MATLAB | 绘图复刻 | 分层聚类分析图:树状图+热图

Posted slandarer

tags:

篇首语:本文由小常识网(cha138.com)小编为大家整理,主要介绍了MATLAB | 绘图复刻 | 分层聚类分析图:树状图+热图相关的知识,希望对你有一定的参考价值。

好久不见啊,今天时绘图复刻第三期,这期画的比较难应该文章也不会太短。。。
前段时间发现公众号SCIPainter发布了一期名为《如何对基因和蛋白质的表达丰度进行相关性分析》,其中有一幅图很好看:

于是我也复刻了一下。由于原文章没有提供数据,我这里随便捏造了点数据进行的绘图,以下是绘制效果:

大家光看图就能看出能展示哪些信息,明显行列相似度大的被调整的挨在一起,且能看出明显颜色分界处也正好是树状图主要分支的分界。

对了,重点需要安装Statistics and Machine Learning Toolbox即统计与机器学习工具箱!!!

0 注

本来发现Bioinformatics Toolbox工具箱就有分层聚类分析的绘制函数clustergram,但是画出来属实太丑了而且配色啊树状图啊统统很难调整:

嗯。。。自己写呗,直接开肝。

1 数据

这里假设X,Y数据为相同数据,计算出的相关性矩阵就是对称矩阵:

X=randn(20,20)+[(linspace(-1,2.5,20)').*ones(1,8),(linspace(.5,-.7,20)').*ones(1,5),(linspace(.9,-.2,20)').*ones(1,7)];
Y=X;
% 计算数据相关性系数
Data=corr(X,Y);

% 输入行列名称
colName='FREM2','ALDH9A1','RBL1','AP2A2','HNRNPK','ATP1A1','ARPC3','SMG5','RPS27A',...
          'RAB8A','SPARC','DDX3X','EEF1D','EEF1B2','RPS11','RPL13','RPL34','GCN1','FGG','CCT3';
rowName=colName;

如果是对称矩阵最终效果:

当然也可以X,Y数据为不同数据:

rng(1)
Y=randn(20,20)+[(linspace(-1,2.5,20)').*ones(1,8),(linspace(.5,-.7,20)').*ones(1,5),(linspace(.9,-.2,20)').*ones(1,7)];
X=randn(20,15)+[(linspace(-1,2.5,20)').*ones(1,4),(linspace(.5,-.7,20)').*ones(1,5),(linspace(.9,-.2,20)').*ones(1,6)];
% 计算数据相关性系数
Data=corr(X,Y);

% 输入行列名称
colName='FREM2','ALDH9A1','RBL1','AP2A2','HNRNPK','ATP1A1','ARPC3','SMG5','RPS27A',...
          'RAB8A','SPARC','DDX3X','EEF1D','EEF1B2','RPS11','RPL13','RPL34','GCN1','FGG','CCT3';
rowName='A1','A2','A3','A4','A5','A6','A7','A8','A9','A10','A11','A12','A13','A14','A15';

如果非对称矩阵最终效果:

2 坐标区域构建

创建坐标区域并调整各个坐标区域的位置Position,中间的热图,两个树状图以及colorbar都分别创建一个坐标区域:

其中比较少见的属性YAxisLocation是用来设置Y轴是在0点处还是坐标区域左侧或者右侧,本期需要设置为右侧,uistack用来调整图形对象层级关系,需要调整一部分坐标区域谁压着谁。

% =========================================================================
% 获取数据矩阵大小
[rows,cols]=size(Data);

fig=figure('Position',[500,200,800,750],'Name','slandarer');


% 坐标区域创建 =============================================================
% 热图坐标区域
placeMat=zeros(7,8);placeMat(2:7,2:7)=1;
axMain=subplot(7,8,find(placeMat'));
axMain.XLim=[1,cols]+[-.5,.5];
axMain.YLim=[1,rows]+[-.5,.5];
axMain.YAxisLocation='right';
axMain.YDir='reverse';
axMain.XTick=1:cols;
axMain.YTick=1:rows; 
hold on
% 树状图坐标区域
axTree1=subplot(7,8,(1:6).*8+1);
axTree1.Position(3)=axTree1.Position(3)+axTree1.Position(1)*4/5;
axTree1.Position(1)=axTree1.Position(1)/5;
axTree1.Position(3)=(axMain.Position(1)-axTree1.Position(1)-axTree1.Position(3))*4/5+axTree1.Position(3);
hold on
axTree2=subplot(7,8,2:7);
axTree2.Position(2)=axMain.Position(2)+axMain.Position(4)+(axTree2.Position(2)-axMain.Position(2)-axMain.Position(4))/5;
axTree2.Position(4)=axTree2.Position(4)+(1-axTree2.Position(2)-axTree2.Position(4))*4/5;
hold on
% colorbar坐标区域
axBar=subplot(7,8,(1:7).*8);
axBar.Position(1)=1/2;
axBar.Position(3)=.92/2;
axBar.Position(2)=axMain.Position(2)+axMain.Position(4)/2;
axBar.Position(4)=axMain.Position(2)+axMain.Position(4)-axBar.Position(2);
axBar.Color='none';
axBar.XColor='none';
axBar.YColor='none';
uistack(axBar,'bottom')
CM=colorbar;
hold on

3 树状图绘制

这里使用linkage获取集聚分层聚类树,再使用dendrogram函数将其绘制,dendrogram函数的Orientation属性可以调整树的方向:

tree1=linkage(Data,'average');
[treeHdl1,~,order1]=dendrogram(tree1,0,'Orientation','left');

设置为黑色并加粗:

% 设置树状图颜色
set(treeHdl1,'Color',[0,0,0]);
set(treeHdl1,'LineWidth',1);


但是注意到此时两个图不在同一figure:

此时我们就想到之前写过一篇《MATLAB | 如何复制figure图窗任意axes的全部信息?》中给了一个复制axes的函数,将其进行微调:

function axbag=copyAxes(fig,k,newAx)
% @author : slandarer
% 公众号  : slandarer随笔
% 知乎    : slandarer
% 
% 此段代码解析详见公众号 slandarer随笔 文章:
%《MATLAB | 如何复制figure图窗任意axes的全部信息?》
% https://mp.weixin.qq.com/s/3i8C78pv6Ok1cmEZYPMyWg
classList(length(fig.Children))=true;
for n=1:length(fig.Children)
    classList(n)=isa(fig.Children(n),'matlab.graphics.axis.Axes');
end
isaaxes=find(classList);
oriAx=fig.Children(isaaxes(end-k+1));
if isaaxes(end-k+1)-1<1||isa(fig.Children(isaaxes(end-k+1)-1),'matlab.graphics.axis.Axes')
    oriLgd=[];
else
    oriLgd=fig.Children(isaaxes(end-k+1)-1);
end
axbag=copyobj([oriAx,oriLgd],newAx.Parent);
axbag(1).Position=newAx.Position;
delete(newAx)
end 

调用该函数复制axes并修饰:

tempFig=treeHdl1(1).Parent.Parent;
% 坐标区域二次修饰
axTree1=copyAxes(tempFig,1,axTree1);
axTree1.XColor='none';
axTree1.YColor='none';
axTree1.YDir='reverse';
axTree1.XTick=[];
axTree1.YTick=[];
axTree1.YLim=[1,rows]+[-.5,.5];
delete(tempFig)

另一个树状图代码几乎一模一样:

% 绘制上方树状图
tree2=linkage(Data.','average');
[treeHdl2,~,order2]=dendrogram(tree2,0,'Orientation','top');
set(treeHdl2,'Color',[0,0,0]);
set(treeHdl2,'LineWidth',1);
tempFig=treeHdl2(1).Parent.Parent;
axTree2=copyAxes(tempFig,1,axTree2);
axTree2.XColor='none';
axTree2.YColor='none';
axTree2.XTick=[];
axTree2.YTick=[];
axTree2.XLim=[1,cols]+[-.5,.5];

4 热图绘制

需要根据linkagedendrogram的聚类结果调整相关系数矩阵行列顺序以及X,Y轴标签顺序:

% 绘制中央热图
Data=Data(order1,:);
Data=Data(:,order2);
imagesc(axMain,Data)
axMain.XTickLabel=colName(order2);
axMain.YTickLabel=rowName(order1);

5 绘制白线

% 绘制白线
LineX=repmat([[1,cols]+[-.5,.5],nan],[rows+1,1]).';
LineY=repmat((.5:1:(rows+.5)).',[1,3]).';
plot(axMain,LineX(:),LineY(:),'Color','w','LineWidth',1);
LineY=repmat([[1,rows]+[-.5,.5],nan],[cols+1,1]).';
LineX=repmat((.5:1:(cols+.5)).',[1,3]).';
plot(axMain,LineX(:),LineY(:),'Color','w','LineWidth',1);

6 调整colorbar范围和配色

范围使用clim调整,配色的话MATLAB自带的配色均可:

colormap(pink)
clim([min(min(Data)),max(max(Data))])

当然自己找一些数据进行插值也可以,这里提供了六种配色,可自行添加,原理就是自己找一些nx3大小的颜色列表,然后插值到上百行,颜色就渐变了:

% 调整colorbar
baseCM以上是关于MATLAB | 绘图复刻 | 分层聚类分析图:树状图+热图的主要内容,如果未能解决你的问题,请参考以下文章

MATLAB | 绘图复刻 | 和弦图+颜色修改+标签旋转

MATLAB | 绘图复刻 | 和弦图+颜色修改+标签旋转

MATLAB | 绘图复刻 | 分组环形热图

MATLAB | 绘图复刻 | 分组环形热图

MATLAB | 绘图复刻 | 热图+差异气泡图

MATLAB | 绘图复刻 | 折线图+误差棒+柱状图+散点抖动+灰色背景+图片叠加