《Geo Word Clouds》paper implementation

Oct 18, 2021 30 min read

Geo Word Clouds 论文复现（2020.8）

准备阶段

地图准备：用 echarts.js 绘制山东省地图，并修改图片大小为 800×500 并保存，作为词云背景
初步确定 Geo Map 范围实验初始阶段的想法是用山东的17个地级市（这样做结果更加直观）作为输入的标签数据，来测试词云对地理信息的反映程度，所以先用 Matlab App 中的 Image Viewer 大致测取各地级市的范围（单位为像素点），结果如下

City	min_x	max_x	min_y	max_y
威海		617		742		111		209
烟台		471		675		81		226
青岛		468		594		162		342
潍坊		354		509		141		324
日照		391		480		290		401
东营		360		448		42		183
淄博		311		380		147		301
临沂		285		439		275		479
滨州		276		364		26		210
莱芜		281		330		235		290
济南		186		310		115		284
泰安		170		335		240		321
济宁		153		302		299		474
枣庄		235		315		375		473
德州		144		301		79		246
聊城		108		209		177		319
菏泽		66		187		313		460

最后根据各地市的面积比例以及限定的坐标范围来生成随机数据
确定实验流程

实验阶段

这一部分只简要介绍实验流程，实现的具体细节见代码实现部分

初始随机数据集在地图上的分布情况如下首先使用 Python 中的 wordcloud 库对生成的随机数据进行可视化，结果如下：可以看出单词的分布比较随机杂乱，没有包含任何地理信息，且单词的大小几乎相同，难以读出频率信息在不聚类的情况下，只加入地理信息因素，绘制的词云图如下：可以看出，虽然现在的词云图已反映出了地理信息，但在地图中为每个点放置一个单词，很难找到出现频率高的单词，同时相邻近的单词颜色比较相近，难以区分，视觉效果不佳

加入原论文中的颜色分配算法，结果如下：

之后对初始数据集进行 K-means 聚类，继续绘制词云图，结果如下：其中各单词的位置与实际地理位置相吻合，同时，单词的大小也在一定程度上反映了该地市的面积（由于这里只着重测试地理位置，聚类较少，所以地图形状空缺较大，用户可以自定义填充单词进行形状填充）叠加图层观察可视化结果，可以看出，各聚类中心大部分被单词所覆盖，虽有个别点出现偏差，但也在可接受的范围之内。接下来再对另一个的数据集进行测试：数据大小为 9935 ，执行 K-means 聚类后生成 126 个聚类运行结果如下：再对几个出现频率较高的单词进行叠加图层比较，结果如下：

Edith

单词覆盖基本准确 Alice

左边的聚类覆盖基本吻合，右边的聚类由于单词 Nick 在聚类中心，同时如果在聚类中心放置的话，由于 Alice 的字体较大，无法放置，所以偏移到了左边 Bettina

左边的聚类基本吻合，右边的聚类由于与 Gasol 冲突，进行了左移调整可视化衡量指标的结果将在下一部分进行描述

代码实现部分

生成初始数据集

为了方便实验，这里的坐标单位均为像素，实际应用的话，只需将经纬度坐标按一定比例转换即可。输入山东省各地市的地理位置信息以及自己所指定的标签名称，来生成包含地理位置坐标的标签数据集

使用 C++ 中的 random_device 随机数引擎生成服从正态分布的随机数据： RandomData.cpp

<div class="highlight highlight-source-c++ position-relative overflow-auto" data-snippet-clipboard-copy-content="#include
using namespace std;

//Geo Info
struct Region{
int min_x,max_x,min_y,max_y,area;
void Output(){
cout<<min_x<<" "<<max_x<<" "<<min_y<<" "<<max_y<<" "<<area<<endl; } }geo[20]; tags info struct tags{ string name; int num,sz; idx[10]; void output(){ cout<<num<<" "<<name<<" "; for(int i="1;i<=sz;i++){" cout<<idx[i]<<" cout<
>city; in>>Geo[i].min_x>>Geo[i].max_x>>Geo[i].min_y>>Geo[i].max_y>>Geo[i].area; } in.close(); /* for(int i=1;i<=17;i++){ Geo[i].Output(); } */ freopen(“TagsInfoVer2.txt”,”r”,stdin); freopen(“CityTags3.txt”,”w”,stdout); // Total Points int Count=0; //Input Tags Region vector
a; string str; stringstream ss; bool flag=false; for(int i=1;i<=60;i++){ str=””; getline(cin,str); ss.clear(); ss<

>tags[i].num; ss>>tags[i].name; while(ss){ int t; ss>>t; a.push_back(t); } a.pop_back(); tags[i].sz=a.size(); for(int j=0;j
<a.size();j++){ tags[i].idx[j+1]="a[j];" } a.clear(); count+="tags[i].num*tags[i].sz;" tags[i].output(); generate data for(int i="1;i<=60;i++){" j="1;j<=tags[i].sz;j++){" frequency int num="((double)tags[i].num/(double)Count)*10000;" cout<<num<<endl; random_device rd; mt19937_64 eng(rd()); coordinate uniform_int_distribution
distrx(Geo[tags[i].idx[j]].min_x, Geo[tags[i].idx[j]].max_x); uniform_int_distribution

distry(Geo[tags[i].idx[j]].min_y, Geo[tags[i].idx[j]].max_y); // Tags x y frequency for(int k=0;k
<num;k++){ cout<<tags[i].name<<" "<<distrx(eng)<<" "<<distry(eng)<<" "<<tags[i].idx[j]<

#include <bits/stdc++.h>
using namespace std;

//Geo Info
struct Region{
	int min_x,max_x,min_y,max_y,area;
	void Output(){
		cout<
        <<
        " "<
        
         <<
         " "<
         
          <<
          " "<
          
           <<
           " "<
           <
           <endl; } }geo[
            20];


            //Tags Info

            struct 
            Tags{
	string name;
	
            int num,sz;
	
            int idx[
            10];
	
            void 
            Output(){
		cout<
            
             <<
             " "<
             
              <<
              " "; 
              for(
              int i=
              1;i<=sz;i++){ cout<
              
               <<
               " "; } cout<
               <endl; } }tags[
                100]; 
                int 
                main() { string city,TagName; 
                int min_x,max_x,min_y,max_y,area; ifstream in; in.
                open(
                "SDCities.txt"); 
                //Input Geo Info 
                for(
                int i=
                1;i<=
                17;i++){ in>>city; in>>Geo[i].
                min_x>>Geo[i].
                max_x>>Geo[i].
                min_y>>Geo[i].
                max_y>>Geo[i].
                area; } in.
                close(); 
                /* 
                 for(int i=1;i<=17;i++){ 
                 Geo[i].Output(); 
                 } 
                 */ 
                freopen(
                "TagsInfoVer2.txt",
                "r",stdin); 
                freopen(
                "CityTags3.txt",
                "w",stdout); 
                // Total Points 
                int Count=
                0; 
                //Input Tags Region  vector<
                int> a; string str; stringstream ss; 
                bool flag=
                false; 
                for(
                int i=
                1;i<=
                60;i++){ str=
                ""; 
                getline(cin,str); ss.
                clear(); ss<
                
                 >tags[i].
                 num; ss>>tags[i].
                 name; 
                 while(ss){ 
                 int t; ss>>t; a.
                 push_back(t); } a.
                 pop_back(); tags[i].
                 sz=a.
                 size(); 
                 for(
                 int j=
                 0;j
                 
                  size();j++){ tags[i].
                  idx[j+
                  1]=a[j]; } a.
                  clear(); Count+=tags[i].
                  num*tags[i].
                  sz; 
                  //tags[i].Output(); } 
                  //Generate Data 
                  for(
                  int i=
                  1;i<=
                  60;i++){ 
                  for(
                  int j=
                  1;j<=tags[i].
                  sz;j++){ 
                  // Frequency 
                  int num=((
                  double)tags[i].
                  num/(
                  double)Count)*
                  10000; 
                  //cout<
                   
                    <
                    
                    random_device rd; mt19937_64 
                  eng(
                  rd()); 
                  // Coordinate uniform_int_distribution<
                  unsigned 
                  long 
                  long> 
                  distrx(Geo[tags[i].
                  idx[j]].
                  min_x, Geo[tags[i].
                  idx[j]].
                  max_x); uniform_int_distribution<
                  unsigned 
                  long 
                  long> 
                  distry(Geo[tags[i].
                  idx[j]].
                  min_y, Geo[tags[i].
                  idx[j]].
                  max_y); 
                  // Tags x y frequency 
                  for(
                  int k=
                  0;k
                  <num;k++){ cout<<tags[i].
                   name<<
                   " "<<
                   distrx(eng)<<
                   " "<<
                   distry(eng)<<
                   " "<
                   
                    idx[j]<
                    <endl; } 
                     return 
                     0; }

数据预处理

K-means 聚类

原文中只简单的提及了使用 K-means 算法对标签进行聚类，但没有具体说明实现的细节，而 K 值的选取对最后的聚类结果有着较大的影响。由于我们要在程序中对每一个标签都要执行 K-means 算法，所以，利用肘部法则，根据变化趋势选取 K 值是不可行的，或者通过不同 K 值轮廓系数的对比来选取合适的 K 值，但这样的话要对同一标签进行多次聚类后，再根据结果进行比较，时间开销偏大，这里选择的是根据地理位置进行聚类。即根据标签点所落在的区域来确定 K 值，比如某一标签落在了（威海，济南）那么我们选取 K 值为 2 即可，但节省时间的同时，这样也存在一些问题，比如某一标签落在了（济南，泰安），由于这两个区域距离很近，其实聚为一类结果可能是最佳的，现在解决这个问题的想法是：如果是做成交互式的自动可视化的应用的话，可以在绘制词云前，首先展示我们的聚类结果，然后用户可以根据聚类结果对某一标签进行调整。

K_means_Region_Cluster.m

45
Dir = [Dir;0];
else
Dir = [Dir;1];
end
sz = sum(idx==cluster);
tmpX = sort(X(idx==cluster,:));
min_x = tmpX(ceil(0.05*sz),1);
max_x = tmpX(ceil(0.95*sz),1);
min_y = tmpX(ceil(0.05*sz),2);
max_y = tmpX(ceil(0.95*sz),2);
Range = [Range;min_x,max_x,min_y,max_y];
C = [C;D(cluster,1) D(cluster,2) sz cnt];
cnt = cnt + 1;
x = [];
y = [];
Tags = [Tags;string(pre)];
end
end
x = [x A(i,1)];
y = [y A(i,2)];
Region = [Region A(i,3)];
pre=cur;
end
X = [x’ y’];
Region = unique(Region);
Rsz = size(Region);
K = Rsz(2);
[idx,D] = kmeans(X,K);
for cluster = 1:K
CurCluster = [X(idx==cluster,1) X(idx==cluster,2)];
CurCluster_Size = size(CurCluster);
tmp = zeros(CurCluster_Size(1),1);
tmp = tmp + cnt;
cluster_points = [cluster_points;X(idx==cluster,1) X(idx==cluster,2) tmp];
direction = abs(PCA_Rotation(CurCluster’));
if direction > 45
Dir = [Dir;0];
else
Dir = [Dir;1];
end
sz = sum(idx==cluster);
Range = [Range;min(X(idx==cluster,1)),max(X(idx==cluster,1)),min(X(idx==cluster,2)),max(X(idx==cluster,2))];
C = [C;D(cluster,1) D(cluster,2) sz cnt];
cnt = cnt + 1;
x = [];
y = [];
Tags = [Tags;string(cur)];
end
Map = [];
GeoSize = size(Geo);
InfoSize = size(C);
for i = 1:InfoSize(1)
for j = 1:GeoSize(1)
cx = C(i,1);
cy = C(i,2);
if cx >= Geo(j,1) && cx = Geo(j,3) && cy

clc,clear
[A,B] = xlsread('CityTags.xlsx');
[Geo,GeoName] = xlsread('Region.xlsx');
Geo = Geo(:,1:4);
num = length(A);
pre = B(1);
x = [];
y = [];
idx = [];
C = [];
D = [];
Region = [];
Tags = [];
Range = [];
% 单词方向 ：水平 = 1，竖直 = 0
Dir = [];
% 第 cnt 个聚类
cnt = 1;
% 为输入信息 打上聚类标签
cluster_points = [];
for i = 1:num
    cur = B(i);
    if ~isequal(pre,cur)
        X = [x' y'];
        Region = unique(Region);
        Rsz = size(Region);
        Region = [];
        K = Rsz(2);
        [idx,D] = kmeans(X,K);
        for cluster = 1:K
            CurCluster = [X(idx==cluster,1) X(idx==cluster,2)];
            CurCluster_Size = size(CurCluster);
            tmp = zeros(CurCluster_Size(1),1);
            tmp = tmp + cnt;
            cluster_points = [cluster_points;X(idx==cluster,1) X(idx==cluster,2) tmp];
            direction = abs(PCA_Rotation(CurCluster'));
            if direction > 45
                Dir = [Dir;0];
            else
                Dir = [Dir;1];
            end
            sz = sum(idx==cluster);
            tmpX = sort(X(idx==cluster,:));
            min_x = tmpX(ceil(0.05*sz),1);
            max_x = tmpX(ceil(0.95*sz),1);
            min_y = tmpX(ceil(0.05*sz),2);
            max_y = tmpX(ceil(0.95*sz),2);
            Range = [Range;min_x,max_x,min_y,max_y];
            C = [C;D(cluster,1) D(cluster,2) sz cnt];
            cnt = cnt + 1;
            x = [];
            y = [];
            Tags = [Tags;string(pre)];
        end
    end
    x = [x A(i,1)];
    y = [y A(i,2)];
    Region = [Region A(i,3)];
    pre=cur;
end
X = [x' y'];
Region = unique(Region);
Rsz = size(Region);
K = Rsz(2);
[idx,D] = kmeans(X,K);
for cluster = 1:K
    CurCluster = [X(idx==cluster,1) X(idx==cluster,2)];
    CurCluster_Size = size(CurCluster);
    tmp = zeros(CurCluster_Size(1),1);
    tmp = tmp + cnt;
    cluster_points = [cluster_points;X(idx==cluster,1) X(idx==cluster,2) tmp];
    direction = abs(PCA_Rotation(CurCluster'));
    if direction > 45
        Dir = [Dir;0];
    else
        Dir = [Dir;1];
    end
    sz = sum(idx==cluster);
    Range = [Range;min(X(idx==cluster,1)),max(X(idx==cluster,1)),min(X(idx==cluster,2)),max(X(idx==cluster,2))];
    C = [C;D(cluster,1) D(cluster,2) sz cnt];
    cnt = cnt + 1;
    x = [];
    y = [];
    Tags = [Tags;string(cur)];
end
Map = [];
GeoSize = size(Geo);
InfoSize = size(C);
for i = 1:InfoSize(1)
    for j = 1:GeoSize(1)
        cx = C(i,1);
        cy = C(i,2);
        if cx >= Geo(j,1) && cx <= Geo(j,2) && cy >= Geo(j,3) && cy <= Geo(j,4)
            Map = [Map;string(GeoName(j))];
            break;
        end
    end
end
xlswrite('CityTagsCluster1.xlsx', [Tags,C,Dir Map,Range], 'A1')
xlswrite('CityTagsWithIndex1.xlsx',cluster_points,'A1')

PCA 确定聚类主方向 论文中 Placement Algorithm 的第三步需要确定单词的初始方向，方法如下所示这里采用PCA算法来确定聚类的主方向：参考对协方差矩阵进行特征分解，按照特征值从大到小的顺序，有特征矩阵 V，其每一行对应一个特征向量，解方程 V * T = (1 0)’ 向量 T 经过 V 的投影之后，在主方向上为 1，在垂直主方向上为 0。故 T 即指示了主方向。

PCA_Rotation.m

<div class="highlight highlight-source-python position-relative overflow-auto" data-snippet-clipboard-copy-content="function [Dir,res] = PCA(Data)
% 以[-pi/2,pi/2]之间的角度pi表示方向
% 1.去除均值
[buf K] = size(Data);
miu = mean(Data')';
for k=1:K
Data(:,k) = Data(:,k)-miu;
end;
sigma = zeros(2,2);
% 2.计算协方差
for k=1:K
x = Data(:,k);
sigma = sigma+x*x';
end;
sigma = sigma/K;
% 3.特征分解
[V,D] = eig(sigma);
if (D(1,1)

function [Dir,res] = PCA(Data)
% 以[-pi/2,pi/2]之间的角度pi表示方向
% 1.去除均值
[buf K] = size(Data);
miu = mean(Data')';
for k=1:K
    Data(:,k) = Data(:,k)-miu; 
end;
sigma = zeros(2,2);
% 2.计算协方差
for k=1:K
    x = Data(:,k);
    sigma = sigma+x*x';
end;
sigma = sigma/K;  
% 3.特征分解
[V,D] = eig(sigma);
if (D(1,1)<D(2,2)) % 把较大的特征值对应的向量挪到第一行
     buf = V(1,:);
     V(1,:) = V(2,:);
     V(2,:) = buf;
end;
% 4.求解主方向向量
res = inv(V)*[1;0];
Dir = atan(res(2)/res(1))/pi*180;

核心放置算法

Step 1 : Load Data

# Load Geo Map
dataf = pd.read_excel('GeoRegion.xlsx')
wordmap = dataf.values
sz = wordmap.shape
geomap = dict()
for i in range(0,sz[0]):
    geomap[wordmap[i][0]] =[wordmap[i][1],wordmap[i][2],wordmap[i][3],wordmap[i][4]]

# Load Data Set
# Load Cluster Info
df = pd.read_excel('CityTagsClusterVer4.xlsx')
words = df.values
# Load Tags labeled by cluster index
df2 = pd.read_excel('CityTagsWithIndex2.xlsx')
points = df2.values
cluster_size = words.shape[0] + 10
cluster_points = [[] for i in range(cluster_size)]

for i in range(0, points.shape[0]):
    cluster_points[points[i][2]].append((points[i][0],points[i][1]))

Step 2 : Cluster 由 K_means_Region_Cluster.m 实现

WordCluster 类声明

other.num

def __str__(self):
return ‘(‘ + self.word + ‘, ‘ + str(self.x) + ‘, ‘ + str(self.y) + ‘, ‘ + str(self.num) + ‘, ‘+str(self.direction)+’)’

“>

class WordCluster(object):
    def __init__(self,word,x,y,num,direction,region,min_x,max_x,min_y,max_y,idx):
        self.word = word
        self.x = x
        self.y = y
        self.num = num
        self.direction = direction
        self.region = region
        self.min_x = min_x
        self.max_x = max_x
        self.min_y = min_y
        self.max_y = max_y
        self.idx = idx

    def __lt__(self, other):
        return self

    def __lt__(self, other):
        return self.num > other.num

    def __str__(self):
        return '(' + self.word + ', ' + str(self.x) + ', ' + str(self.y) + ', ' + str(self.num) + ', '+str(self.direction)+')'

Step 3 : Assign Attributes For Each Cluster
- (1) Choose the Font Size 求解 $A_m$

# compute the area of M in pixel : Am
img = cv2.imread("SD.png")
gray = cv2.cvtColor(img, cv2.COLOR_BGRA2GRAY)
h, w = gray.shape[:2]
m = np.reshape(gray, [1, w * h])
mean = m.sum() / (w * h)
ret, binary = cv2.threshold(gray, mean, 255, cv2.THRESH_BINARY)
Am = len(binary[binary == 0])

求解 $Font Size$：

font_size = int((cmath.sqrt((freq * Am)/n).real))

(2) Choose the Rotation 由 PCA_Rotation.m 实现

Step 4 : Greedy Placement 论文中提到，如果某个单词不能被放在它准确的位置，那么就缩小它，重新加入队列，关于放置的方法没有提及（只是说了尽可能覆盖聚类中更多的点) 下面就依照我自己对问题的理解来实现放置算法：关于放置的具体位置我尝试了质心周围矩形范围内放置、按区域放置、按聚类结果放置，
第一种方法，由于重点一直关注在质心，所以为了实现放置可能会出现很多冲突的单词由于被缩放多次直至小于 minimal size 被抛弃，导致很多单词没有被绘制，实验结果也证实了这种猜测，由于存在多次迭代，运行时间较长，且被绘制在词云上的单词数减少了
第二种方法，规定范围放置，比如要放置的单词位于青岛，那么就在青岛这个范围内，随便找取一个空闲位置放置，通过实验，这种方法虽然运行速度很快，但由于限定的范围太过宽泛，地理偏差较大
第三种方法，在观察上述两种方法的词云与聚类中心点叠加图层后的结果后受到启发，取了二者的折中，将限定范围缩小，缩小为聚类所占据的范围，为了避免一些聚类边缘点的干扰，首先对聚类的坐标点进行排序，然后取 5%到95% 的数据作为聚类的范围，然后在寻找位置时，在所有可行的坐标中取距离聚类中心点距离最短的一组作为放置位置，放置算法实现代码如下：

<div class="highlight highlight-source-python position-relative overflow-auto" data-snippet-clipboard-copy-content="def sample_position_v2(self, size_x, size_y, b_x, b_y, e_x, e_y, center_x, center_y):
dis_to_center = 99999999
ans_x = 0
ans_y = 0
for y in range(max(1, b_y), min(e_y – size_y, self.height – size_y)):
for x in range(max(1,b_x), min(e_x – size_x,self.width – size_x)):
area = self.integral[y – 1, x – 1] + self.integral[y + size_y – 1, x + size_x – 1]
area -= self.integral[y – 1, x + size_x – 1] + self.integral[y + size_y – 1, x – 1]
if not area:
tmpdis = Euclidean_distance(x + size_x/2 , center_x, y + size_y/2, center_y).real
if tmpdis

def sample_position_v2(self, size_x, size_y, b_x, b_y, e_x, e_y, center_x, center_y):
    dis_to_center = 99999999
    ans_x = 0
    ans_y = 0
    for y in range(max(1, b_y), min(e_y - size_y, self.height - size_y)):
        for x in range(max(1,b_x), min(e_x - size_x,self.width - size_x)):
            area = self.integral[y - 1, x - 1] + self.integral[y + size_y - 1, x + size_x - 1]
            area -= self.integral[y - 1, x + size_x - 1] + self.integral[y + size_y - 1, x - 1]
            if not area:
                tmpdis = Euclidean_distance(x + size_x/2 , center_x, y + size_y/2, center_y).real
                if tmpdis < dis_to_center:
                    dis_to_center = tmpdis
                    ans_x = x
                    ans_y = y
    if dis_to_center == 99999999:
        return None
    else:
        return ans_y, ans_x