教你使⽤MATLAB绘制散点密度图(⼆维核密度)
效果:
原理也很简单,通过matlab⾃带的ksdensity获得⽹格每⼀点密度,通过密度拟合曲⾯,再计算每个数据点对应的概率,并将概率映射到颜⾊即可
为了怕⼤家不到函数这次⼯具函数放到最前⾯
1⼯具函数完整代码
function [CData,h,XMesh,YMesh,ZMesh,colorList]=density2C(X,Y,XList,YList,colorList)
[XMesh,YMesh]=meshgrid(XList,YList);
XYi=[XMesh(:)YMesh(:)];
F=ksdensity([X,Y],XYi);
ZMesh=zeros(size(XMesh));
ZMesh(1:length(F))=F;
h=interp2(XMesh,YMesh,ZMesh,X,Y);
if nargin<5
colorList=[0.270000.3300
0.27000.23000.5100
0.19000.41000.5600
0.12000.56000.5500
0.21000.72000.4700
0.56000.84000.2700
0.99000.91000.1300];
end
colorFunc=colorFuncFactory(colorList);
CData=colorFunc((h-min(h))./(max(h)-min(h)));
colorList=colorFunc(linspace(0,1,100)');
function colorFunc=colorFuncFactory(colorList)
x=(0:size(colorList,1)-1)./(size(colorList,1)-1);
y1=colorList(:,1);y2=colorList(:,2);y3=colorList(:,3);
colorFunc=@(X)[interp1(x,y1,X,'pchip'),interp1(x,y2,X,'pchip'),interp1(x,y3,X,'pchip')];
end
end
2参数说明
输⼊:
X,Y 散点坐标
XList,YList ⽤来构造密度曲⾯⽹格的序列,其实就是把XLim,YLim分成⼩份,例如XList=0:0.1:10
colorList 颜⾊表mx3数组,⽤来构造将⾼度映射到颜⾊函数的数据表
输出:
CData各个点对应颜⾊
h 各个点对应核密度
matlab直方图XMesh,YMesh,ZMesh 核密度曲⾯数据
colorList 插值后更细密的颜⾊表
3使⽤⽅式
假如编写了如下程序:
PntSet1=mvnrnd([23],[10;02],800);
PntSet2=mvnrnd([67],[10;02],800);
PntSet3=mvnrnd([89],[10;01],800);
PntSet=[PntSet1;PntSet2;PntSet3];
scatter(PntSet(:,1),PntSet(:,2),'filled');
结果:
3.1散点赋⾊
将上⾯那段代码改写
PntSet1=mvnrnd([23],[10;02],800);
PntSet2=mvnrnd([67],[10;02],800);
PntSet3=mvnrnd([89],[10;01],800);
PntSet=[PntSet1;PntSet2;PntSet3];
CData=density2C(PntSet(:,1),PntSet(:,2),-2:0.1:15,-2:0.1:15);
scatter(PntSet(:,1),PntSet(:,2),'filled','CData',CData);
3.2等⾼线图
PntSet1=mvnrnd([23],[10;02],800);
PntSet2=mvnrnd([67],[10;02],800);
PntSet3=mvnrnd([89],[10;01],800);
PntSet=[PntSet1;PntSet2;PntSet3];
[~,~,XMesh,YMesh,ZMesh,colorList]=density2C(PntSet(:,1),PntSet(:,2),-2:0.1:12,-2:0.1:12); colormap(colorList)
contourf(XMesh,YMesh,ZMesh,10)
3.3带直⽅图的散点图
colorList=[0.94000.97000.9600
0.89000.93000.9200
0.82000.91000.8800
0.69000.85000.7700
0.59000.78000.6900
0.55000.75000.6500
0.45000.65000.5600
0.40000.58000.4900
0.35000.51000.4200
0.25000.36000.3100
0.13000.17000.1400];
CData=density2C(PntSet(:,1),PntSet(:,2),-2:0.1:15,-2:0.1:15,colorList); set(gcf,'Color',[111]);
%主分布图
ax1=axes('Parent',gcf);hold(ax1,'on')
scatter(ax1,PntSet(:,1),PntSet(:,2),'filled','CData',CData);
ax1.Position=[0.1,0.1,0.6,0.6];
% X轴直⽅图
ax2=axes('Parent',gcf);hold(ax2,'on')
histogram(ax2,PntSet(:,1),'FaceColor',[0.780.880.82],...
'EdgeColor','none','FaceAlpha',0.7)
ax2.Position=[0.1,0.75,0.6,0.15];
ax2.YColor='none';
ax2.XTickLabel='';
ax2.TickDir='out';
ax2.XLim=ax1.XLim;
% Y轴直⽅图
ax3=axes('Parent',gcf);hold(ax3,'on')
histogram(ax3,PntSet(:,2),'FaceColor',[0.780.880.82],...
'EdgeColor','none','FaceAlpha',0.7,'Orientation','horizontal')
ax3.Position=[0.75,0.1,0.15,0.6];
ax3.XColor='none';
ax3.YTickLabel='';
ax3.TickDir='out';
ax3.YLim=ax1.YLim;
3.4带直⽅图的等⾼线图
colorList=[0.93000.95000.9700
0.79000.84000.9100
0.65000.73000.8500
0.51000.62000.7900
0.37000.51000.7300
0.27000.41000.6300
0.21000.32000.4900
0.15000.22000.3500
0.09000.13000.2100
0.03000.04000.0700];
[~,~,XMesh,YMesh,ZMesh,colorList]=density2C(PntSet(:,1),PntSet(:,2),-2:0.1:13,-2:0.1:13,colorList);
set(gcf,'Color',[111]);
%主分布图
ax1=axes('Parent',gcf);hold(ax1,'on')
colormap(colorList)
contourf(XMesh,YMesh,ZMesh,10,'EdgeColor','none')
ax1.Position=[0.1,0.1,0.6,0.6];
ax1.TickDir='out';
% X轴直⽅图
ax2=axes('Parent',gcf);hold(ax2,'on')
[f,xi]=ksdensity(PntSet(:,1));
fill([xi,xi(1)],[f,0],[0.340.470.71],'FaceAlpha',...
0.3,'EdgeColor',[0.340.470.71],'LineWidth',1.2)
ax2.Position=[0.1,0.75,0.6,0.15];
ax2.YColor='none';
ax2.XTickLabel='';
ax2.TickDir='out';
ax2.XLim=ax1.XLim;
% Y轴直⽅图
ax3=axes('Parent',gcf);hold(ax3,'on')
[f,yi]=ksdensity(PntSet(:,2));
fill([f,0],[yi,yi(1)],[0.340.470.71],'FaceAlpha',...
0.3,'EdgeColor',[0.340.470.71],'LineWidth',1.2)
ax3.Position=[0.75,0.1,0.15,0.6];
ax3.XColor='none';
ax3.YTickLabel='';
ax3.TickDir='out';
ax3.YLim=ax1.YLim;
4使⽤⽅式扩展–与ggplot修饰器联动
ggplot风格修饰器:(点击图⽚跳转链接)
⽰例1
ax=gca;
ax.XLim=[-113];
ax.YLim=[-113];
ax=ggplotAxes2D(ax);
CData=density2C(PntSet(:,1),PntSet(:,2),0:0.1:15,0:0.1:15);
scatter(PntSet(:,1),PntSet(:,2),'filled','CData',CData);
是不是瞬间有那味了:
⽰例2
PntSet1=mvnrnd([23],[10;02],800);
PntSet2=mvnrnd([67],[10;02],800);
PntSet3=mvnrnd([89],[10;01],800);
PntSet=[PntSet1;PntSet2;PntSet3];
ax=gca;
ax.XLim=[-313];
ax.YLim=[-313];
ax=ggplotAxes2D(ax);
[~,~,XMesh,YMesh,ZMesh,colorList]=density2C(PntSet(:,1),PntSet(:,2),-2:0.1:12,-2:0.1:12); colormap(colorList)
contourf(XMesh,YMesh,ZMesh,10)
版权声明:本站内容均来自互联网,仅供演示用,请勿用于商业和其他非法用途。如果侵犯了您的权益请与我们联系QQ:729038198,我们将在24小时内删除。
发表评论