Can I add 95% confidence ellipses around groups of data in a pca plot (biplot)?

Question

Heidi Hirsh 2020-10-1

0
链接

此问题的直接链接

https://ww2.mathworks.cn/matlabcentral/answers/603379-can-i-add-95-confidence-ellipses-around-groups-of-data-in-a-pca-plot-biplot

评论： Adam Danz 2020-10-4

I am plotting the results of principle component analysis using biplot. I am wondering if there is a way (in matlab) to add confidence ellipses around the groups of data. Maybe I have to use something instead of biplot?

In the example below I would want to draw an ellipse around the data for each site (green and blue points)

Here's how I am currently plotting the data:

[coeff,score,latent,tsquared,explained] = pca(X);
f1=figure(1) 
h= biplot(coeff(:,1:2),'scores',score(:,1:2),'color','k','marker','.','markersize',17,'varlabels',...
       {'DO','O2sat','pH','Temp','Sal','Depth','PAR','|Velocity|'});  
 
%color by site
hID = get(h,'tag'); %identify handle
hPt = h(strcmp(hID,'obsmarker')); %isolate handles to scatter points
grp = findgroups(site);  
grpID = 1:max(grp); 
 
clrMap = winter(length(unique(grp)));   
for i = 1:max(grp)
    set(hPt(grp==i), 'Color', clrMap(i,:), 'DisplayName', sprintf('MSP%d', grpID(i)))
end
 
title('Full Deployment (12h)');
set(gca,'fontsize',18)
xlabel('Component 1 (55.26%)') % how can I add percent ('%s %', explained(1))')
ylabel('Component 2 (21.87%)')
[~, unqIdx] = unique(grp);
legend(hPt(unqIdx))

1 个评论
显示 -1更早的评论隐藏 -1更早的评论

Heidi Hirsh 2020-10-1

Here's a random example I found online for what I want to do:

请先登录，再进行评论。

请先登录，再回答此问题。

Answer 1

Adam Danz 2020-10-1

1
链接

此回答的直接链接

https://ww2.mathworks.cn/matlabcentral/answers/603379-can-i-add-95-confidence-ellipses-around-groups-of-data-in-a-pca-plot-biplot#answer_503758

编辑：Adam Danz 2020-10-4

在 MATLAB Online 中打开

data.mat

If you know the center of the clusters and the CI's along the x and y axes for each cluster, you can use this function to plot the ellipses

https://www.mathworks.com/matlabcentral/answers/476867-how-can-i-draw-std-ellipses-in-a-scatter-plot#answer_388409

A comment below that answer points to an alternative solution as well.

Also check the file exchange.

Addendum

The block of code below implements the ellipses outlined in the first link above with your data which is attached. I made additional changes from your version to clean some stuff up a bit.

Stars mark the center of each cluster using the mean (consider using the median instead since the data are not normally distributed). Major and minor axes of the ellipses represent the 95% CI using the percentile method which is a good method given that the data are not normally distributed.

load('data.mat')
[coeff,score,latent,tsquared,explained] = pca(X);
f1=figure(); 
h= biplot(coeff(:,1:2),'scores',score(:,1:2),'color','k','marker','.','markersize',17,'varlabels',...
    {'DO','O2sat','pH','Temp','Sal','Depth','PAR','|Velocity|'});
hold on
%color by site
hID = get(h,'tag'); %identify handle
hPt = h(strcmp(hID,'obsmarker')); %isolate handles to scatter points
[grp, grpID] = findgroups(site);
clrMap = winter(numel(grpID));
p = 95; % CI level
for i = 1:max(grpID)
    set(hPt(grp==i), 'Color', clrMap(i,:), 'DisplayName', sprintf('MSP%d', grpID(i)))
    
    % Compute centers (means)
    allX = arrayfun(@(hh)hh.XData(1), hPt(grp==i)); 
    allY = arrayfun(@(hh)hh.YData(1), hPt(grp==i)); 
    centers(1) = mean(allX);  % x mean
    centers(2) = mean(allY);  % y mean
    
    % Plot centers, do they make sense?
    plot(centers(1), centers(2), 'rp', 'MarkerFaceColor', clrMap(i,:), 'MarkerSize', 20, 'LineWidth', 1)
    
    % Compute 95% CI using percentile method
    CIx = prctile(allX, [(100-p)/2, p+(100-p)/2]); % x CI [left, right]
    CIy = prctile(allY, [(100-p)/2, p+(100-p)/2]); % y CI [lower, upper]
    CIrng(1) = CIx(2)-CIx(1); % CI range (x)
    CIrng(2) = CIy(2)-CIy(1); % CI range (y)
    
    % Draw ellipses
    llc = [CIx(1), CIy(1)]; % (x,y) lower left corners
    rectangle('Position',[llc,CIrng],'Curvature',[1,1], 'EdgeColor', clrMap(i,:));
    
end
title('Full Deployment (12h)');
set(gca,'fontsize',18)
xlabel('Component 1 (55.26%)') % how can I add percent ('%s %', explained(1))')
ylabel('Component 2 (21.87%)')
[~, unqIdx] = unique(grp);
legend(hPt(unqIdx))

5 个评论
显示 3更早的评论隐藏 3更早的评论

Heidi Hirsh 2020-10-3

编辑：Adam Danz 2020-10-4

在 MATLAB Online 中打开

ellipse attempt.png

I tried adding the following code. I think I got close to a solution but for some reason everything is huge. Can anyone spot what I am doing wrong?

%divide scores by site (first have of observations are site 1 and second half are site 2)
score1 = score(1:34,1:2);
x1 = score1(:,1); y1 = score1(1,:); 
score2 = score(35:end,1:2);
x2 = score2(:,1); y2 = score2(1,:); 
 
%find centers for data from each site
meanx1 = mean(x1);
meany1 = mean(y1);
meanx2 = mean(x2);
meany2 = mean(y2);
 
%confidence intervals
% x is a vector, matrix, or any numeric array of data. NaNs are ignored.
% p is the confidence level (ie, 95 for 95% CI)
% The output is 1x2 vector showing the [lower,upper] interval values.
p=95; 
%find upper and lower bounds for 95% confidence in x and y directions
% CIFcn_x1 = @(x1,p)prctile(x1,abs([0,100]-(100-p)/2));
CIFcn_x1 = prctile(x1,abs([0,100]-(100-p)/2));
CIFcn_y1 = prctile(y1,abs([0,100]-(100-p)/2));
CIFcn_x2 = prctile(x2,abs([0,100]-(100-p)/2));
CIFcn_y2 = prctile(y2,abs([0,100]-(100-p)/2));
 
CI_x1 = abs((CIFcn_x1(1)-CIFcn_x1(2))/2); 
CI_y1 = abs((CIFcn_y1(1)-CIFcn_y1(2))/2); 
CI_x2 = abs((CIFcn_x2(1)-CIFcn_x2(2))/2); 
CI_y2 = abs((CIFcn_y2(1)-CIFcn_y2(2))/2); 
 
mean1 = [meanx1,meany1]
std1 = [CI_x1,CI_y1]
mean2 = [meanx2,meany2]
std2 = [CI_x2,CI_y2]
 
%add ellipses 
E = plotEllipses(mean1,std1); 
E.FaceColor = [0 0 1 .3]; 
E.EdgeColor = clrMap(1,:); 
E.LineWidth = 2; 
 
E = plotEllipses(mean2,std2); 
E.FaceColor = [0 1 0 .3]; 
E.EdgeColor = clrMap(2,:); 
E.LineWidth = 2; 

Heidi Hirsh 2020-10-4

Thank you so much Adam! This is so helpful!

Adam Danz 2020-10-4

Glad I could help.

请先登录，再进行评论。

Can I add 95% confidence ellipses around groups of data in a pca plot (biplot)?

1 个评论
显示 -1更早的评论隐藏 -1更早的评论

采纳的回答

5 个评论
显示 3更早的评论隐藏 3更早的评论

更多回答（0 个）

另请参阅

类别

标签

Community Treasure Hunt

Can I add 95% confidence ellipses around groups of data in a pca plot (biplot)?

1 个评论 显示 -1更早的评论隐藏 -1更早的评论

采纳的回答

5 个评论 显示 3更早的评论隐藏 3更早的评论

更多回答（0 个）

另请参阅

类别

标签

Community Treasure Hunt

1 个评论
显示 -1更早的评论隐藏 -1更早的评论

5 个评论
显示 3更早的评论隐藏 3更早的评论