最近帮导师搬砖看了篇论文,记录一下。
a scalable approach for general correlation clustering
该论文有两个创新点:
1.聚类不用预定义最终聚类的个数k
2.据说效率很高迭代一次大约是O(n2)的,所以适合大规模的数据
该文章将待聚类的对象看成一个个点,然后对象之间的关系分为两类,一种是“attract”,用一条标为“+”的边来表示,一种是“repel”,用一条标位“-”的边来表示
同时标为“+”和“-”的边有一个权值,表示 相近|相异 的程度。用原文的话来说就是一个n×n的矩阵,即:
$$W\begin{cases}>0 &\text{: u and v attract each other}\\<0 &\text{: u and v repel each other} \\=0 &\text{: the relation between u and v is uncertain}\end{cases}$$
该文章于是给出了一个算法,只要输入一个矩阵W,就能得到一个最终的聚类方案。
首先用一个n维的向量C表示聚类方案,C(i)表示点i分在哪个类。
那么误差可以表述为:
1.W[u,v] > 0 即u,v相近但是C[u]!=C[v],即没有在一个类中
2.W[u,v] < 0 即u,v相斥但是C[u]==C[v], 即分在一个类中了
形式化的表示为:
$$ E = \sum\limits_{e(u,v)=+} W[u,v][C[u] \neq C[v]] + \sum\limits_{e(u,v)=-} W[u,v][C[u] = C[v]]$$
然后我们就要求出一个C 使得 E达到最小值,
然后可以引入一个k×n的L矩阵来表示分类的具体情况,L[i,u]表示第i类中是否含有点u:
$$L[i,u]\begin{cases}1 & u \in i\\-1 & u \not\in i\end{cases}$$
可以将目标变形为:
$$\min\limits_{L,C,k} \sum\limits_{W[u,v]L[C[u],v] < 0} W[u,v]L[C[u],v]$$
我们可以发现,向量C和矩阵L很像EM算法中的隐藏变量和其现有估计值。(C和L其实是一个东西,知道了L就知道C,知道了C就知道了L)
于是该论文指出,可以用一种EM算法的思想来迭代计算出最后的C。
每次迭代的过程即:
$$C^{t+1}[u] = \arg\limits_{C[u]\in \{1,2,...,k\}} \min \sum\limits_{W[u,v]L^t[C[u],v] < 0} -W[u,v]L^t[C[u],v]$$
于是得出算法1:(伪EM算法)
在这个算法中,每次迭代都有可能有一些类没有包含任何的点所以要将这些空类给删除然后重新编号。
同时你可能注意到这个k参数一开始是没有给出来的。k是在迭代过程中逐渐稳定下来的,而且k只会减少不会增加。
一开始你可以将初始的L设置为一个n×n的矩阵,表示一共有n个类且每个类中只包含一个点,但是这样效率不高,
于是论文中给出了一个特别的初始化方法。
首先将所有的点按照含有标记为+的边的个数降序排序,得到一个list,然后从第一个开始向最后一个扫描。
如果扫描到当前点u,那么所有与u有标为“+”边的点将和u分为一类,且这些点将从该list中删除。这样知道list为空。
这个初始化算法,想法就是把这个k在一开始就压缩的小一点。不过在数据量很小的时候可能会出现一点问题,由于k是不会减小的
所以要是一开始初始化的时候就把所有的点分为一个类中,那么后续的迭代也就无法进行了,所以在点少的时候还是初始化为
n×n的矩阵,每个类一个点比较好。
由次得出算法2:(初始化L,大数据量时的优化)
用n,a,r表示点数,标为+的边数,标为-的边数。
算法二的时间复杂度O(a + nlogn)
算法一每次迭代的复杂度是O(a + r + n)
我简单的实现了一下,贴上一个C++代码,实现的不是很好,意思算是差不多。
1 //################################ 2 //# Author : septicmk 3 //# HomePage : www.cnblogs.com/septic 4 //# Email : septic.mk@gmial.com 5 //# Date : 2014/5/27 6 //# FileName : Cluster 7 //################################ 8 #include <iostream> 9 #include <cstdio> 10 #include <cstring> 11 #include <algorithm> 12 #include <vector> 13 #include <map> 14 15 #define vec vector<int> 16 using namespace std; 17 18 class ClusterTool{ 19 20 //初始化根据+的度数排序的数据结构 21 struct point{ 22 int pdegree; 23 int id; 24 point(int u,int val){ 25 id = u; 26 pdegree = val; 27 } 28 bool operator < (const point &hs)const{ 29 return pdegree > hs.pdegree; 30 } 31 }; 32 33 //判断两次迭代结果是否有差异 34 bool check(vec C_old,vec C_new){ 35 for(int i = 0 ; i < C_old.size() ; i ++){ 36 if(C_old[i] != C_new[i]) return true; 37 } 38 return false; 39 } 40 41 //枚举u所在的类,S.T. 误差值最小 42 int optimal(int u,vector<vec> W,vector<vec> L){ 43 int k = L.size(); 44 int n = W.size(); 45 int tmp,best; 46 for(int i = 0 ; i < k ; i ++){ 47 int sum = 0; 48 for(int j = 0 ; j < n ; j ++){ 49 sum += -W[u][j]*L[i][j]; 50 } 51 if(i == 0){ 52 tmp = sum; 53 best = 0; 54 }else{ 55 if(tmp > sum){ 56 tmp = sum; 57 best = i; 58 } 59 } 60 } 61 return best; 62 } 63 64 65 //由C得到L 66 vector<vec> update(vec C){ 67 map<int,int> m; 68 m.clear(); 69 int ind = 0; 70 for(int i = 0 ; i < C.size() ; i ++){ 71 if(m.find(C[i]) == m.end()){ 72 m.insert(pair<int,int>(C[i],ind)); 73 C[i] = ind; 74 ind++; 75 }else{ 76 C[i] = m[C[i]]; 77 } 78 } 79 80 vector<vec> L; 81 L.clear(); 82 vec l; 83 for(int i = 0 ; i < ind ; i ++){ 84 l.clear(); 85 for(int i = 0 ; i < C.size() ; i ++) 86 l.push_back(-1); 87 L.push_back(l); 88 } 89 90 for(int i = 0 ; i < C.size() ; i ++){ 91 L[C[i]][i] = 1; 92 } 93 return L; 94 } 95 96 //由L得到C 97 vec update(vector<vec> L){ 98 vec C; 99 C.clear(); 100 for(int j = 0; j < L[0].size() ; j++){ 101 for(int i = 0 ; i < L.size() ; i ++){ 102 if(L[i][j] == 1){ 103 C.push_back(i); 104 break; 105 } 106 } 107 } 108 return C; 109 } 110 111 //迭代过程 112 vec Pseudo_EM(vector<vec> W,vector<vec> L,bool mode){ 113 int t = 0; 114 vec C_cur; 115 vec C_nxt; 116 vec C_test = update(L); 117 for(int i = 0 ; i < C_test.size() ; i ++) 118 cout << C_test[i] << " "; 119 puts(""); 120 do{ 121 C_cur = update(L); 122 C_nxt.clear(); 123 for(int u = 0 ; u < W.size() ; u ++){ 124 C_nxt.push_back(optimal(u,W,L)); 125 if(mode){ 126 L[C_cur[u]][u] = -1; 127 L[C_nxt[u]][u] = 1; 128 } 129 } 130 L = update(C_nxt); 131 t++; 132 }while(check(C_cur,C_nxt)); 133 return C_cur; 134 } 135 136 //初始化L0 137 vector<vec> initialize(vector<vec> W){ 138 vector<vec> L; 139 L.clear(); 140 vector<vec> G; 141 vector<point> P; 142 vec use; 143 int n = W.size(); 144 for(int u = 0 ; u < W.size() ; u ++){ 145 vec g;g.clear(); 146 int cnt = 0; 147 for(int v = 0 ; v < W[u].size() ; v ++){ 148 int t = W[u][v] > 0 ? 1: -1; 149 if(t == 1) cnt++; 150 if(u == v) t = 1; 151 g.push_back(t); 152 } 153 P.push_back(point(u,cnt)); 154 use.push_back(0); 155 G.push_back(g); 156 } 157 158 sort(P.begin(),P.end()); 159 160 int rest = n; 161 int i = 0,u; 162 vec l; 163 while(rest){ 164 for(;i<P.size();i++){ 165 if(use[P[i].id] == 0){ 166 u = P[i].id; 167 //cout << i << endl; 168 break; 169 } 170 } 171 l.clear(); 172 for(int j = 0 ; j < n ; j ++){ 173 l.push_back(-1); 174 } 175 176 for(int j = 0 ; j < P.size() ; j ++){ 177 int v = P[j].id; 178 if(!use[v] && G[u][v] == 1){ 179 l[v] = 1; 180 use[v] = 1; 181 rest--; 182 } 183 } 184 L.push_back(l); 185 } 186 return L; 187 } 188 189 //对外的方法,输入W得到C 190 public: vec cluster(vector<vec> W){ 191 //vector<vec> L0 = initialize(W); 192 int n = W.size(); 193 vector<vec> L0; 194 L0.clear(); 195 vec l; 196 for(int i = 0 ; i < n ; i ++){ 197 l.clear(); 198 for(int j = 0 ; j < n; j ++){ 199 l.push_back(-1); 200 if(i == j) l.push_back(1); 201 } 202 L0.push_back(l); 203 } 204 puts("initizlized"); 205 vec C = Pseudo_EM(W,L0,true); 206 return C; 207 } 208 }; 209 210 int main(){ 211 int n; 212 cin >> n; 213 vector<vec> W; 214 215 for(int i = 0 ; i < n ; i ++){ 216 vec w; 217 w.clear(); 218 int tmp; 219 for(int j = 0 ; j < n ; j ++){ 220 cin >> tmp; 221 w.push_back(tmp); 222 } 223 W.push_back(w); 224 } 225 ClusterTool *CT = new ClusterTool(); 226 vec C = CT->cluster(W); 227 for(int i = 0 ; i < C.size() ; i ++){ 228 if(i == 0) cout << C[i]; 229 else cout << " " << C[i]; 230 }puts(""); 231 return 0; 232 }