Matlab数据可视化（5）：二维数据可视化 I

以下介绍数据的二维可视化。

一. 二维散点图

（源代码：scatter.m）

我们用著名的Iris数据集（Fisher, 1936）作为绘图实例。Iris数据集包含3种鸢尾花的150个样本数据，每个数据都有4个属性（花萼和花瓣的长度及宽度）。

1) 基本散点图

我们用其中两个属性值作为X和Y轴，另一个属性值表示点的大小。（图1）

[plain] view plaincopy

1. %%  加载数据集  
2. [attrib className] = xlsread('iris.xlsx');  
3.   
4. %% 绘制基本的散点图  
5. figure('units','normalized','Position',[0.2359    0.3009    0.4094    0.6037]);  
6. scatter(attrib(:,1),attrib(:,2),10*attrib(:,3),[0 1 0],'filled','Marker','^');  
7. set(gca,'Fontsize',12);  
8. title({'Iris数据集包含150个数据，每个数据含4个属性',...  
9.     '第32个属性值扩大后用来决定标志的大小',...  
10.     '标志使用自定义的大小和样式'});  
11. xlabel('属性1'); ylabel('属性2');  
12. box on;  
13. set(gcf,'color',[1 1 1],'paperpositionmode','auto');  

图 1

2) 散点矩阵图

plotmatrix可以将每两个属性组合的散点图以矩阵的形式绘制出来。通过观察如图2，我们知道，3个各类的花朵不能通过某两个属性就可以加以区分。

[plain] view plaincopy

1. %% 绘制散点矩阵图  
2. % The output parameters have  
3. % A matrix of handles to the objects created in H,  
4. % A matrix of handles to the individual subaxes in AX,  
5. % A handle to a big (invisible) axes that frames the subaxes in BigAx,  
6. % A matrix of handles for the histogram plots in P.  
7. % BigAx is left as the current axes so that a subsequent title, xlabel, or ylabel command is centered with respect to the matrix of axes.  
8. figure('units','normalized','Position',[ 0.2359    0.3009    0.4094    0.6037]);  
9. [H,AX,BigAx,P] = plotmatrix(attrib,'r.');  
10. attribName = {[char(10) '花萼长度'],[char(10) '花萼宽度'],[char(10) '花瓣长度'],[char(10) '花瓣宽度']};  
11. % 添加标注  
12. for i = 1:4  
13.     set(get(AX(i,1),'ylabel'),'string',['属性' num2str(i) attribName{i}]);  
14.     set(get(AX(4,i),'xlabel'),'string',['属性' num2str(i) attribName{i}]);  
15. end  
16. set(get(BigAx,'title'),'String',{'散点矩阵图 (数据含3类Iris花)', ...  
17.     '如图所示这3类不能通过某2个属性进行区分'},'Fontsize',14);  
18. set(gcf,'color',[1 1 1],'paperpositionmode','auto');  

图 2

3) 带直方图的散点图

我们可以在绘制两个属性的散点图的同时，绘制出各属性分布的直方图。（图3）

[plain] view plaincopy

1. %% 在轴线上画出分布图  
2. figure('units','normalized','Position',[0.2589    0.3296    0.3859    0.5750]);  
3. H = scatterhist(attrib(:,2),attrib(:,3),'nbins',50,'direction','out');  
4. mainDataAxes = H(1);  
5. xhistAxes = H(2);  
6. yhistAxes = H(3);  
7. set(get(mainDataAxes,'title'),'String','有一个类别基于如下2个属性是可分的','Fontsize',12);  
8. set(get(mainDataAxes,'xlabel'),'String','x (属性2, 花萼宽度)','Fontsize',14);  
9. set(get(mainDataAxes,'ylabel'),'String','y (属性3, 花瓣长度)','Fontsize',14);  
10. set(gcf,'color',[1 1 1],'paperpositionmode','auto');  

图 3

二.  用热图表示数据点分布的密度

（源代码：scattersmooth.m）

当数据量很大时，散点图之间会有显著的重叠，这时密度图往往能够更加清晰的反映出数据的分布情况。我们以用如下方法生成的1000个数据点为例进行说明。

[plain] view plaincopy

1. %% 生成数据  
2. z = [repmat([1 2],1000,1) + randn(1000,2)*[1 .5; 0 1.32];...  
3.      repmat([9 1],1000,1) + randn(1000,2)*[1.4 .2; 0 0.98];...  
4.      repmat([4 8],1000,1) + randn(1000,2)*[1 .7; 0  0.71];];  

如果直接绘制散点图（如图4），我们会发现象非常严重的重叠问题。 

[plain] view plaincopy

1. %% 绘制原始数据  
2. figure('units','normalized','position',[ 0.4458    0.6296    0.1995    0.2759]);  
3. plot(z(:,1),z(:,2),'.');  
4. title({'原始数据','散点视图有显著的重叠'});  
5. set(gcf,'color',[1 1 1],'paperpositionmode','auto');  

图 4

我们通过将数据进行抽象，绘制密度图来解决上述问题。之前，我们介绍过热图（heat map）的绘制，当时用到的是imagesc和pcolor，而这里我们利用绘制二维曲面的surf命令，接合适当的视角来完成。（图5）

[plain] view plaincopy

1. xx = z(:,1);  
2. yy = z(:,2);  
3. gridSize = 100;  
4. % 选择colormap  
5. colormap(summer);  
6.   
7. % 建立网格  
8. x=linspace(min(xx),max(xx),gridSize);  
9. y=linspace(min(yy),max(yy),gridSize);  
10. gridEval = zeros(length(x)-1,length(y)-1);  
11.   
12. % 计算每个网格中点的频数  
13. for cnt_x=1:length(x)-1  
14.     for cnt_y=1:length(y)-1  
15.         x_ind=intersect(find(xx>x(cnt_x)),find(xx<=x(cnt_x+1)));                                                      
16.         xy_ind=intersect(find(yy(x_ind)>y(cnt_y)), find(yy(x_ind)<=y(cnt_y+1)));       
17.         gridEval(cnt_y, cnt_x)=length(xy_ind);  
18.     end  
19. end  
20.   
21. % surface函数绘制热图  
22. surf((x(1:end-1)+ x(2:end))/2,(y(1:end-1)+y(2:end))/2,gridEval); view(2);   
23. shading interp;  hold on;  
24. axis([min(xx),max(xx) min(yy),max(yy)]);  
25.   
26. % 添加标注  
27. title(['密度图, 网格大小: ' num2str(gridSize) ' x ' num2str(gridSize) ' 个网格'],'Fontsize',14);  
28. xlabel('x','Fontsize',14);  
29. ylabel('y','Fontsize',14);  
30. axis tight;  
31. h1 = gca; % 保存句柄，以便后面添加边框  
32.   
33. % 添加颜色条  
34. h=colorbar;  
35. axes(h);ylabel('密度, 每个网格的点数','Fontsize',14);  
36. set(gcf,'color',[1 1 1],'paperpositionmode','auto');  
37.   
38. % 添加黑色的边框  
39. axes(h1);  
40. line(get(gca,'xlim'),repmat(min(get(gca,'ylim')),1,2),'color',[0 0 0],'linewidth',1);  
41. line(get(gca,'xlim'),repmat(max(get(gca,'ylim')),1,2),'color',[0 0 0],'linewidth',2);  
42. line(repmat(min(get(gca,'xlim')),1,2),get(gca,'ylim'),'color',[0 0 0],'linewidth',2);  
43. line(repmat(max(get(gca,'xlim')),1,2),get(gca,'ylim'),'color',[0 0 0],'linewidth',1);  

图 5

直接计算密度得到的密度可能会显得不太自然，我们可以第个网格点处对所有数据进行加权，这样便可以成到平滑的作用。（图6）

[plain] view plaincopy

1. xx = z(:,1);  
2. yy = z(:,2);  
3. sigma = 0.1;  
4. gridSize = 100;  
5.   
6. % 选择colormap  
7. colormap(summer);  
8.   
9. x=linspace(min(xx),max(xx),gridSize);  
10. y=linspace(min(yy),max(yy),gridSize);  
11. gridEval = zeros(length(x)-1,length(y)-1);  
12.   
13. % 计算每个点处的高斯函数  
14. for i = 1:length(x)-1  
15.     for j = 1:length(y)-1  
16.     %calculate a Gaussian function on the grid with each point in the center and add them up  
17.         gridEval(j,i) = gridEval(j,i) + sum(exp(-(((x(i)-xx).^2)./(2*sigma.^2) + ((y(j)-yy).^2)./(2*sigma.^2))));         
18.     end  
19. end  
20.   
21. % 绘制热图  
22. surf((x(1:end-1)+ x(2:end))/2,(y(1:end-1)+y(2:end))/2,gridEval); view(2); shading interp;  
23. axis([min(xx),max(xx) min(yy),max(yy)]);  
24. h1 = gca; % 保存句柄，以便后面添加边框  
25.   
26. % 添加标注  
27. title(['平滑散点图, \sigma = ' num2str(sigma) ', 网格数： ' num2str(gridSize) ' x' num2str(gridSize)],'Fontsize',14);  
28. xlabel('x','Fontsize',14);  
29. ylabel('y','Fontsize',14);  
30.   
31. % 添加颜色条  
32. h=colorbar;  
33. axes(h);ylabel('Intensity','Fontsize',14);  
34.   
35. % 添加黑色的边框  
36. axes(h1);  
37. line(get(gca,'xlim'),repmat(min(get(gca,'ylim')),1,2),'color',[0 0 0],'linewidth',1);  
38. line(get(gca,'xlim'),repmat(max(get(gca,'ylim')),1,2),'color',[0 0 0],'linewidth',2);  
39. line(repmat(min(get(gca,'xlim')),1,2),get(gca,'ylim'),'color',[0 0 0],'linewidth',2);  
40. line(repmat(max(get(gca,'xlim')),1,2),get(gca,'ylim'),'color',[0 0 0],'linewidth',1);  
41.   
42. set(gcf,'color',[1 1 1],'paperpositionmode','auto');  

图 6

三. 双向误差图（bidirectional error bar）

（源代码：bidirerrorbars.m）

对某些数据而言，X轴和Y轴上的数据都有一定的波动（误差）范围，有时我们需要同时在图上反映出这些波动信息。这称为双向误差图。我们曾通过errorbar命令绘制过Y轴方向上的误差图，在此基础上，我们绘制双向误差图。

（图7）

[plain] view plaincopy

1. %% 加载数据  
2. load flatPlateBoundaryLayerData  
3.   
4. xx = laminarFlow(:,2);  
5. yy = laminarFlow(:,3);  
6.   
7. % 划定网格  
8. xg=linspace(min(xx),max(xx),6);  
9. yg=linspace(min(yy),max(yy),6);  
10.   
11. % 计算各个网格中点的频数，以及它们的均值和方差  
12. for cnt_x=1:length(xg)-1  
13.     x_ind=intersect(find(xx>xg(cnt_x)),find(xx<=xg(cnt_x+1)));    
14.     x(cnt_x)=mean(xx(x_ind));  
15.     e_x(cnt_x)=std(xx(x_ind));  
16. end  
17. for cnt_y=1:length(yg)-1  
18.     y_ind=intersect(find(yy>yg(cnt_y)),find(yy<=yg(cnt_y+1)));    
19.     y(cnt_y)=mean(yy(y_ind));  
20.     e_y(cnt_y)=std(yy(y_ind));  
21. end  
22.   
23.    
24. figure('units','normalized','Position',[0.0750    0.5157    0.5703    0.3889]);  
25. axes('Position',[0.0676    0.1100    0.8803    0.8150]);  
26.   
27. %% 绘制双向误差图  
28. h = biDirErrBar(x,y,e_x,e_y);  
29.   
30. %% 标注  
31. set(get(h,'title'),'string','Laminar Flow的双向误差图','Fontsize',15);  
32. set(get(h,'xlabel'),'string','位置测度','Fontsize',15);  
33. set(get(h,'ylabel'),'string','速度测度','Fontsize',15);  
34. set(gcf,'color',[1 1 1],'paperpositionmode','auto');  

图 7

其中，自定义函数biDirErrBar在Matlab命令errorbar基础上，添加X轴方向上的误差。定义如下：

[plain] view plaincopy

1. function h = biDirErrBar(x,y,e_x,e_y)  
2.   
3. % bidirErrBar(x,y,e_x,e_y) plots bi directional error bars and returns the  
4. % handle to the axes for annotation  
5. % x, y are the data vectors with n elements each. e_x and e_y are the error bars to be  
6. % positioned at each point in xi,yi,   
7.   
8. errorbar(x,y,e_y); hold on;  
9. plot(x,y,'Color',[0 0 0]);  
10. x_lower=x-e_x; x_upper=x+e_x;   
11. y_lower=y-.0001*y; y_upper=y+.0001*y;  
12. line([x_lower; x_upper],[y; y],'Color',[0 0 1]);  
13.   
14. hold on; line([x_lower; x_lower], [y_lower; y_upper],'Color',[0 0 1]);  
15. hold on; line([x_upper; x_upper], [y_lower; y_upper],'Color',[0 0 1]);  
16.   
17. h = gca;  

四. 二维关系图

我们知道，在关系矩阵中，1表示连接，0表示非去接。进一步，非零整数可以表示两个结点间的连接数。我们可以利用热图表示出这种扩展的关系。以维克多·雨果《悲剧世界》中人物同时登场的次数数量为实例，我们绘制二维的关系图。（图8）

[plain] view plaincopy

1. %% 加载数据  
2. load characterCoOccurences  
3.   
4. %% 设置图像和坐标轴   
5. figure('units','normalized','Position',[ 0.1990    0.1324    0.4854    0.7731]);  
6. mainAx = axes('position',[ 0.1361    0.0143    0.8042    0.8038]);  
7.   
8. % 设置colormap。由于关系图是稀疏的，因此这里我们反转colormap的顺序  
9. m = colormap(copper);  
10. m = m(end:-1:1,:);  
11. colormap(m);  
12.   
13. % 用surf创建热图  
14. r=surf(lesMiserables);view(2);  
15. % 将edgecolor设置为半透明  
16. set(r,'edgealpha',0.2);  
17. set(gca,'clim',[min(lesMiserables(:)) max(lesMiserables(:))]);  
18.   
19. % 在上方添加颜色条  
20. h=colorbar('northoutside');   
21.   
22. % 添加刻度标签  
23. set(mainAx,'xAxisLocation','top','xtick',0:78,'xticklabel',{' ' LABELS{:} ' '},...  
24.          'ytick',0:78,'yticklabel',{' ' LABELS{:} ' '},'ticklength',[0 0],'fontsize',8);  
25. axis tight;  
26. rotateXLabels(gca,90);  
27.   
28. % 调整颜色条和坐标轴的位置，因为他们会受琶刻度标签旋转的影响  
29. set(h,'position',[0.1006    0.9409    0.8047    0.0128]);  
30. set(get(h,'title'),'String','维克多·雨果《悲剧世界》中人物同时出场次数','Fontsize',14);  
31. set(mainAx,'position',[ 0.1361    0.0143    0.8042    0.8038]);  
32.   
33. set(gcf,'color',[1 1 1],'paperpositionmode','auto');  
34. box on;  

图 8

五. 绘制系统树图（dendrogram）

(源代码：dendrocluster.m)

图 9

我们首先利用pdist命令得到基因之间的两两距离。然后，linkage命令通过聚类的方法，得到层级树形层级结构。层级结构的可视化则由dendrogram命令完成。

[plain] view plaincopy

1. %% 加载数据  
2. load 14cancer.mat  
3. data = [Xtrain(find(ytrainLabels==9),genesSet); Xtest(find(ytestLabels==9),genesSet)];  
4.   
5. %% 癌症样本及100个感兴趣的基因绘制系统树图  
6. Z_genes =   linkage(pdist(data'));  
7. Z_samples = linkage(pdist(data));  
8.   
9. %% 设定图像位置  
10. figure('units','normalized','Position',[0.5641    0.2407    0.3807    0.6426]);  
11. mainPanel = axes('Position',[.25 .08 .69 .69]);  
12. leftPanel = axes('Position',[.08 .08 .17 .69]);  
13. topPanel =  axes('Position',[.25 .77 .69 .21]);  
14.   
15. %% 较低的值颜色较浅  
16. m = colormap(pink); m = m(end:-1:1,:);  
17. colormap(m);  
18.   
19. %% 绘制系统树图 - 展示所有基因的节点（如果没有参数“0”，则默认值为最多30个节点）  
20. axes(leftPanel);h = dendrogram(Z_samples,'orient','left'); set(h,'color',[0.1179         0         0],'linewidth',2);  
21. axes(topPanel); h = dendrogram(Z_genes,0);set(h,'color',[0.1179         0         0],'linewidth',2);  
22.   
23. %% 获取系统树图和热图数据生成的次序  
24. Z_samples_order = str2num(get(leftPanel,'yticklabel'));  
25. Z_genes_order = str2num(get(topPanel,'xticklabel'));  
26. axes(mainPanel);  
27. surf(data(Z_samples_order,Z_genes_order),'edgecolor',[.8 .8 .8]);view(2);  
28. set(mainPanel,'Xticklabel',[],'yticklabel',[]);  
29.   
30. %% 对齐X轴和Y轴  
31. set(leftPanel,'ylim',[1 size(data,1)],'Visible','Off');  
32. set(topPanel,'xlim',[1 size(data,2)],'Visible','Off');  
33. axes(mainPanel);axis([1 size(data,2) 1 size(data,1)]);  
34.   
35. %% 添加标注  
36. axes(mainPanel); xlabel('30个不同的基因','Fontsize',14);   
37. colorbar('Location','northoutside','Position',[ 0.0584    0.8761    0.3082    0.0238]);   
38. annotation('textbox',[.5 .87 .4 .1],'String',{'基因表达水平', '白血病'},'Linestyle','none','fontsize',14);  
39.   
40. %% 在坐标轴不可见的情况下，显示标签  
41. set(leftPanel,'yaxislocation','left');  
42. set(get(leftPanel,'ylabel'),'string','样本','Fontsize',14);  
43. set(findall(leftPanel, 'type', 'text'), 'visible', 'on');  
44.   
45. set(gcf,'color',[1 1 1],'paperpositionmode','auto');