干货|变邻域搜索(VNS)算法求解Max-Mean Dispersion Problem(附...
Part
1
Max-Mean Dispersion Problem
对MDP和MMDP还是一头雾水？不要着急，今天就和小编一起解决这三个问题—
什么是MDP和MMDP?
它们有什么用?
怎样解决这两个问题?
1.1
Maximum Diversity Problem
要讨论Max-Mean Dispersion Problem,就要首先了解Maximum Diversity Problem (MDP) 。
MDP是一种用来度量元素中差异的问题，通俗来讲，就是要从一个集合中选择一个子集合，使得子集合中的元素差异最大化。
举个例子，对于生态系统的平衡问题，生态系统的存续与否就在于多样性。假如我们拥有任意两个生物之间的差异值，通过找到具有最大多样性的子集，我们就能方便地建立起可行的生态系统。
又比如说在农业育种中，往往需要在子代中挑选出具有理想多样性的种群，问题就又归结到了在子代中找到最大差异化个体的问题上了。
文章开头的表情包，其实质也是一个MDP。在风险投资中，往往要通过找到差异最大的几个市场来进行投资，避免牵一发而动全身的情况。
例子都是多样性在生活中发挥作用的表现，那么我们应该如何最大化这种多样性呢？这就是MDP要解决的问题了。
更多的应用见下图：
1.2
MDP的数学描述
考虑一个元素的集合
,索引集
。每个元素包含着r个属性，我们可以将一个元素用向量
表示。我们的目标就是从n个元素中选择m个元素，最大化所选择的元素的多样性。为了度量元素之间的多样性，我们定义值
来代表元素i,j之间的距离。这个距离有多种算法，如欧几里得距离，曼哈顿距离等。在这里我们使用最为常用的欧几里得距离。 
由此，我们可以定义一组元素的多样性为：
根据这些约定，我们可以通过引入0-1变量的方法，将问题表达为：
1.3
Max-Mean Dispersion Problem
有了对MDP问题的认识，下面我们来看看MMDP。
和MDP要最大化子集的多样性不同，MMDP问题需要最大化子集多样性的平均值。在这里值得注意的一点是，MDP中子集的大小是固定的，是问题给出的。而MMDP中，子集数量的多少需要自己确定。当子集数量确定后，MMDP就转化为了MDP。
还是有些云里雾里？更通俗的讲，假如说问题是在10个人中挑出差异最大的5个人，这自然是一个MDP，但假如说问题是在10个人中挑出几个人，使这几个人的差异最大呢？这时自然不能考虑差异值的和，而是需要考虑差异值的和的平均值，即MMDP了。
我们用一个简单的例子来具体解释MDP和MMDP：
假设给出4个元素A，B，C，D，给出4个元素的距离矩阵如下图：
假如要求是从4个元素中选择3个元素，使它们之间的差异最大，这就是一个MDP。假设选择元素A，B，C，则目标函数的值为1+2+4 = 7.
假如要求是从4个元素中选择任意个元素，使他们之间的平均差异最大，这就是一个MMDP。同样假设选择元素A，B，C，目标函数的值就变为（1+2+4）/ 3 = 7/3。
Part
2
变邻域搜索(VNS)算法再回顾
在之前的推文干货 | 变邻域搜索算法(Variable Neighborhood Search,VNS)超详细一看就懂中，已经对VNS算法有了详细的介绍。
在干货 | 变邻域搜索算法(VNS)求解TSP（附C++详细代码及注释）
中，给出了VNS算法解决问题的实例。
在这里，我们简要地复习下VNS算法的基本内容，详细内容参见以上推文。
2.1
VNS算法介绍
VNS算法的基本思想是在搜索过程中系统地改变邻域结构集来拓展搜索过程，获得局部最优解，再基于此局部最优解重新系统地改变邻域结构集拓展搜索范围找到另一个局部最优解的过程。基本的流程如下图：
正如Hansen在论文Variable neighborhood search Principles and applications一文中提到的，VNS算法本质上还是一种跳出局部最优解的算法。和禁忌搜索与模拟退火算法不同，其算法并不遵循一定的"轨迹"，而是通过shaking动作来跳出当前局部最优解，在不同的邻域中找到其他局部最优解，当且仅当该解优于当前解时进行移动。假如邻域结构可以覆盖整个可行解集，则算法可以找到全局最优解。
Part
3
具体算法介绍
3.1
初始解生成
对于初始解，我们使用贪心的方法来构造。最开始将所有元素都视为已选择，计算出每一元素被移除后，该解目标函数值的提高，不断地移除能提高最多的元素，不断循环，直到不再有元素被移除时目标函数值提高为止。
3.2
邻域动作
我们定义三种邻域动作：
Exchange:从被选择的元素的集合中随机选择元素i，从不被选择的元素的集合中随机选择元素j，交换i,j。
Insert:从不被选择的元素中随机选择元素i，将其从不被选择的元素的集合中移除，并加入到被选择的元素的集合中。
Remove: 从被选择的元素的集合中随机选择元素i，将其从被选择的元素的集合中移除，并加入到不被选择的元素的集合中。
3.3
具体流程
shake函数：我们定义shake函数接受参数k，随机从选择的元素的集合和不被选择的元素的集合中选择k个元素，并交换他们。
通过我们在3.2中定义的邻域动作进行进行搜索，具体流程如下图：
Part
4
代码分享
这里我们分享两份代码，第一份是某位西班牙大佬分享的代码，另一种是小编在他的基础上改编而来的代码，这里展示的是小编实现的版本。在https://github.com/alexfp95/Max-Mean-Dispersion-Problem/tree/master/src中，可以获得西班牙大佬写的版本。
具体实现如下：
1import java.io.*;
2import java.util.ArrayList;
3import java.util.LinkedList;
4import java.util.Queue;
5import java.util.Random;
6
7class Solution //解的类
8{
9 ArrayList select_set = new ArrayList<>();//存放点的集合
10 ArrayList unselect_set = new ArrayList<>();//没选择的点
11 double value;
12 double getValue()
13 {
14 double ans = 0;
15 ArrayList new_set = new ArrayList<>();
16 new_set.addAll(this.select_set);
17 while(new_set.size()!=0)
18 {
19 int i = new_set.get(0);
20 new_set.remove(0);
21 for(Integer j:new_set)
22 {
23 ans+=main.cost_mariax[i][j];
24 }
25 }
26 ans = ans / this.select_set.size();
27 return ans; // 返回目标函数值
28 }
29 double improve_value(int i)//计算返回将这个点转到另一个集合后，value值的改变
30 {
31 double ans;
32 double last_ans = this.value;
33 Integer I = Integer.valueOf(i);
34 Solution new_solution = new Solution();
35 new_solution.select_set.addAll(this.select_set);
36 if(this.select_set.contains(i))
37 {
38
39 new_solution.select_set.remove(I);
40 new_solution.unselect_set.add(I);
41 ans = new_solution.getValue() - last_ans;
42 }
43 else
44 {
45 new_solution.select_set.remove(I);
46 new_solution.unselect_set.add(I);
47 ans = new_solution.getValue() - last_ans;
48
49 }
50 return ans;
51 }
52
53
54}
55abstract class Strategy//策略类，存放算法
56{
57
58 static Solution ini_solution()//初始化一个解，采用贪婪的思想，最开始所有解都在select_set中，随后逐渐减少，每次计算去除点的价值，由大到小，不再有改进
59 {
60 Solution solution = new Solution();
61 for(int i=1;i<=main.CODE_NUMBER;i++)
62 solution.select_set.add(i);//开始时所有解都在select_set中
63 solution.value = solution.getValue();//获得当前解
64 double best = 0;
65 int id = 0;
66 while(true) {
67 best = 0;
68 for (int i : solution.select_set) {
69 //System.out.println(solution.improve_value(i));
70 if (best < solution.improve_value(i)) {
71 best = solution.improve_value(i);
72 id = i;
73 }
74 //System.out.println(solution.select_set.size());
75 }
76 if(best == 0)
77 break;
78 solution.select_set.remove(Integer.valueOf(id));//不断改进
79 solution.unselect_set.add(Integer.valueOf(id));
80 solution.value = solution.getValue();
81 // System.out.println(solution.value);
82
83 }
84 solution.value = solution.getValue();
85 return solution;
86 }
87
88 static Solution exchange(Solution solution)//第一种邻域算子：交换i,j st i在solution中，j不在
89 {
90 Random r = new Random();
91 int i = r.nextInt(solution.select_set.size()-1);
92 int j = r.nextInt(solution.unselect_set.size()-1);//在i,j中随机找两点;
93 Integer mid_one = solution.select_set.get(i);
94 Integer mid_two = solution.unselect_set.get(j);
95 solution.select_set.remove(i);
96 solution.unselect_set.remove(j);
97 solution.unselect_set.add(Integer.valueOf(mid_one));
98 solution.select_set.add(Integer.valueOf(mid_two));
99 solution.value = solution.getValue();
100 return solution;
101 }
102 static Solution insert(Solution solution)//第二种邻域算子：将j从没选择的集合中加入到solution中
103 {
104 Random r = new Random();
105 int j = r.nextInt(solution.unselect_set.size()-1);
106 int mid = solution.unselect_set.get(j);
107 solution.unselect_set.remove(j);
108 solution.select_set.add(Integer.valueOf(mid));
109 solution.value = solution.getValue();
110 return solution;
111 }
112 static Solution remove(Solution solution)//第三种邻域算子：将j从选择的集合中删除
113 {
114 Random r = new Random();
115 int j = r.nextInt(solution.select_set.size()-1);
116 int mid = solution.select_set.get(j);
117 solution.unselect_set.add(Integer.valueOf(mid));
118 solution.value = solution.getValue();
119 return solution;
120 }
121 public static Solution shake(Solution solution,int k)//shake动作，跳出局部最优
122 {
123 for(int i=1;i<=k;i++)
124 {
125 solution = exchange(solution);
126 }
127 return solution;
128 }
129 public static Solution Local_Search(Solution solution)//邻域搜索，获得局部最优
130 {
131 for(int i=1;i<=100;i++)
132 {
133 Solution new_solution = new Solution();
134 new_solution.select_set.addAll(solution.select_set);
135 new_solution.unselect_set.addAll(solution.unselect_set);
136 new_solution.value = solution.value;
137 if(solution.unselect_set.size()==0)
138 {
139 new_solution = Strategy.remove(solution);
140 }
141 else if(solution.select_set.size()==0)
142 {
143 new_solution = Strategy.remove(solution);
144 }
145 else {
146 Random r = new Random();
147 double t = r.nextDouble();
148 if(t>0.6) {
149 new_solution = Strategy.exchange(new_solution);
150 }
151 else if(t>0.3)
152 {
153 new_solution = Strategy.remove(new_solution);
154
155 }
156 else
157 {
158 new_solution = Strategy.insert(new_solution);
159 }
160
161 }
162 if(solution.valuevalue) {
163 solution = new_solution;
164 System.out.println(new_solution.value);
165 }
166 }
167 return solution;
168 }
169 public static Solution V_N_Search(Solution solution)//变邻域搜索
170 {
171 int k = 1;
172 {
173 while(k<=solution.select_set.size())//按照shaking的定义进行shake，不断搜索直到k==被选择的集合的元素个数
174 {
175 System.out.println(k);
176 Solution new_solution = new Solution();
177 new_solution.select_set.addAll(solution.select_set);
178 new_solution.unselect_set.addAll(solution.unselect_set);
179 new_solution.value = solution.value;
180 new_solution = shake(new_solution,k);
181 new_solution = Local_Search(new_solution);
182 if(solution.valuevalue)
183 {
184 solution = new_solution;
185 k=1;
186 }
187 else{
188 k++;
189 }
190 System.out.println(solution.value);
191
192 }
193 }
194 return solution;
195 }
196}
197public class main {
198 static double[][] cost_mariax;//距离矩阵
199 static int CODE_NUMBER;
200 public static void main(String[] args)
201 {
202 CODE_NUMBER = 150;
203 cost_mariax = new double[CODE_NUMBER+1][CODE_NUMBER+1];//初始化
204 for(int i=1;i<=CODE_NUMBER;i++)
205 try {
206 FileReader fr = new FileReader("MDPI1_150.txt");//读入数据
207 BufferedReader br = new BufferedReader(fr);
208 String line = " ";
209 while(true)
210 {
211 line = br.readLine();
212
213 if(line.equals("end"))break;
214 String data[] = line.split("\t");
215 cost_mariax[Integer.valueOf(data[0])][Integer.valueOf(data[1])] = Double.valueOf(data[2]);
216 cost_mariax[Integer.valueOf(data[1])][Integer.valueOf(data[0])] = Double.valueOf(data[2]);
217
218 }
219 }
220 catch(IOException e)
221 {
222 e.printStackTrace();
223 }
224
225
226 for(int i=1;i<=CODE_NUMBER;i++) //初始化
227 for(int j=1;j<=CODE_NUMBER;j++)
228 {
229 if(i == j)
230 {
231 cost_mariax[i][j] = 0;
232 continue;
233 }
234
235
236 }
237 Solution solution = Strategy.ini_solution(); // 初始化解;
239 solution = Strategy.V_N_Search(solution);//VNS搜索
240 System.out.println(solution.value);
241 System.out.println("当前解为"+solution.value);
242 System.out.println("被选择的点集的大小为" + solution.select_set.size());
243
244
245 }
246}
结果如图：
欲下载本文相关代码，测试样例，相关论文，请移步留言区
参考文献：
[1]Martí, Rafael, and Fernando Sandoya. “GRASP and Path Relinking for the Equitable Dispersion Problem.” Computers & Operations Research 40.12 (2013): 3091–3099. Crossref. Web.
[2] Hansen, Pierre , and N. Mladenovi . "Variable neighborhood search: Principles and applications." European Journal of Operational Research 130.3(2001):449-467.
[3]Hansen, Pierre, N. Mladenović, and D. Perez-Britos. "Variable Neighborhood Decomposition Search." Journal of Heuristics 7.4(2001):335-350.
-The End-
文案/代码/排版 ：李博骁
指导老师 / 秦虎 华中科技大学管理学院
如对代码有疑问，可联系小编，无偿提供服务。
李博骁（华中科技大学管理学院本科二年级、1642940680@qq.com)