center, major and minor axis of ellipsoid
6 views (last 30 days)
Show older comments
Here is my code to plot the ellipsoids od 50,80,95,99% confidence intervals from my data (x1_bladder, y1_rectum), and it works, but i don not know how to plot the major and minor axis of this ellipsoid, can anyone help me? Thanks in advance
% my Data
Ith_MinT = [x1_bladder];
Can_MinT = [y1_rectum];
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%%% Calculate moments of the data sets %%%
% Observed
Ith_mean = mean(Ith_MinT); % x mean
Can_mean= mean(Can_MinT); % y mean
CV = cov(Ith_MinT,Can_MinT); % the angle of ellipse is determined by the cov of data
[Evec, Eval]=eig(CV); % Eigen values and vectors of covariance matrix
%%% Plot observed multivariate contours %%
% Observed data
xCenter = Ith_mean; % ellipses centered at sample averages, center at x,y data mean
yCenter = Can_mean;
theta = 0 : 0.01 : 2*pi; % angles used for plotting ellipses 0:360 degree
% compute angle for rotation of ellipse
% rotation angle will be angle between x axis and first eigenvector
x_vec= [1 ; 0]; % vector along x-axis
cosrotation =dot(x_vec,Evec(:,1))/(norm(x_vec)*norm(Evec(:,1)));
rotation =pi/2-acos(cosrotation); % rotation angle
R = [sin(rotation) cos(rotation); ...
-cos(rotation) sin(rotation)]; % create a rotation matrix
% create chi squared vector
chisq = [1.368 3.2188 4.605 5.991 9.210]; % percentiles of chi^2 dist df=2% the scale of the ellipsoid depends on confidence interval chosen
% chisq = [1.368 3.2188 4.605 5.991];
% size ellipses for each quantile
for i = 1:length(chisq)
% calculate the radius of the ellipse
xRadius(i)=(chisq(i)*Eval(1,1))^.5; % primary
yRadius(i)=(chisq(i)*Eval(2,2))^.5; % secondary
% lines for plotting ellipse
x{i} = xRadius(i)* cos(theta);
y{i} = yRadius(i) * sin(theta);
% rotate ellipse
rotated_Coords{i} = R*[x{i} ; y{i}];
% center ellipse
x_plot{i}=rotated_Coords{i}(1,:)'+xCenter;
y_plot{i}=rotated_Coords{i}(2,:)'+yCenter;
end
% Set up plot
figure
set(gca,'Title',text('String','Probability Ellipses of M1 full weighted of patient plan',...
'FontAngle', 'italic', 'FontWeight', 'bold'),...
'xlabel',text('String', '$\mathbf{M1 full weighted bladder}$', 'Interpreter', 'latex'),...
'ylabel',text('String', '$\mathbf{M1 full weighted rectum}$', 'Interpreter', 'latex'))
hold on
% Plot data points
plot(Ith_MinT,Can_MinT,'o');% scattered data plot
% Plot contours
for j = 1:length(chisq)
plot(x_plot{j},y_plot{j})
end
legend('Data points','50%', '80%', '90%', '95%', '99%')
hold on
0 Comments
Answers (2)
KSSV
on 12 Apr 2019
% my Data
Ith_MinT = rand(1,100) ;
Can_MinT = rand(1,100) ;
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%%% Calculate moments of the data sets %%%
% Observed
Ith_mean = mean(Ith_MinT); % x mean
Can_mean= mean(Can_MinT); % y mean
CV = cov(Ith_MinT,Can_MinT); % the angle of ellipse is determined by the cov of data
[Evec, Eval]=eig(CV); % Eigen values and vectors of covariance matrix
%%% Plot observed multivariate contours %%
% Observed data
xCenter = Ith_mean; % ellipses centered at sample averages, center at x,y data mean
yCenter = Can_mean;
theta = 0 : 0.01 : 2*pi; % angles used for plotting ellipses 0:360 degree
% compute angle for rotation of ellipse
% rotation angle will be angle between x axis and first eigenvector
x_vec= [1 ; 0]; % vector along x-axis
cosrotation =dot(x_vec,Evec(:,1))/(norm(x_vec)*norm(Evec(:,1)));
rotation =pi/2-acos(cosrotation); % rotation angle
R = [sin(rotation) cos(rotation); ...
-cos(rotation) sin(rotation)]; % create a rotation matrix
% create chi squared vector
chisq = [1.368 3.2188 4.605 5.991 9.210]; % percentiles of chi^2 dist df=2% the scale of the ellipsoid depends on confidence interval chosen
% chisq = [1.368 3.2188 4.605 5.991];
% size ellipses for each quantile
for i = 1:length(chisq)
% calculate the radius of the ellipse
xRadius(i)=(chisq(i)*Eval(1,1))^.5; % primary
yRadius(i)=(chisq(i)*Eval(2,2))^.5; % secondary
% lines for plotting ellipse
x{i} = xRadius(i)* cos(theta);
y{i} = yRadius(i) * sin(theta);
% rotate ellipse
rotated_Coords{i} = R*[x{i} ; y{i}];
% center ellipse
x_plot{i}=rotated_Coords{i}(1,:)'+xCenter;
y_plot{i}=rotated_Coords{i}(2,:)'+yCenter;
end
% Set up plot
figure
set(gca,'Title',text('String','Probability Ellipses of M1 full weighted of patient plan',...
'FontAngle', 'italic', 'FontWeight', 'bold'),...
'xlabel',text('String', '$\mathbf{M1 full weighted bladder}$', 'Interpreter', 'latex'),...
'ylabel',text('String', '$\mathbf{M1 full weighted rectum}$', 'Interpreter', 'latex'))
hold on
% Plot data points
plot(Ith_MinT,Can_MinT,'o');% scattered data plot
% Plot contours
for j = 1:length(chisq)
P = [x_plot{j} y_plot{j}] ;
plot(x_plot{j},y_plot{j})
% Get center
C = mean(P) ;
plot(C(1),C(2),'*')
% Get distances
d = sqrt(sum((C-[P(:,1) P(:,2)]).^2,2)) ;
M1 = [C(1)-min(d),C(2) ; C(1)+min(d),C(2)] ;
plot(M1(:,1),M1(:,2),'r')
M2 = [C(1),C(2)-max(d) ; C(1) ,C(2)+max(d)] ;
plot(M2(:,1),M2(:,2),'r')
end
legend('Data points','50%', '80%', '90%', '95%', '99%')
hold on
Waqar Khan
on 18 Mar 2021
Hello, I need to make ellipse using two points such as X=0.098, and Y=0.765. how can i find center from these X and Y and after that draw an ellipse. Please need help.
3 Comments
Waqar Khan
on 18 Mar 2021
Thank you for correcting me here, How can i do it from one point, is it possible?
See Also
Categories
Find more on Contour Plots in Help Center and File Exchange
Community Treasure Hunt
Find the treasures in MATLAB Central and discover how the community can help you!
Start Hunting!