Re: Drift in plotting of confidence region ellipsoid

1 回表示 (過去 30 日間)
Ravi mutturi
Ravi mutturi 2016 年 5 月 12 日
I am trying to plot the confidence region ellipsoid using following code for a (3x3) covariance matrix into 3 seperate 2-dimensional plots with midpoints [16.4083;96.9713;9.4083]. However my confidence ellipsoids are having drift from the above midpoints. I am assuming chisquare distribution with 95%confidence region. Kindly help.
covar_matrix=[ 88.5333 -33.6000 -5.3333
-33.6000 15.4424 2.6667
-5.3333 2.6667 0.4848];
cov_m1m2 = [88.5333 -33.6000
-33.6000 15.4424];
cov_m2m3 = [15.4424 2.6667
2.6667 0.4848];
cov_m1m3 = [88.5333 -5.3333
-5.3333 0.4848];
mL2m1m2 = [16.4083
96.9713];
mL2m2m3 = [96.9713
9.4083];
mL2m1m3 = [16.4083
9.4083];
m1=-50:2.44:50;
m2=80:1:120;
m3=7:0.122:12;
chi_1 = zeros(41,41);
chi_2 = zeros(41,41);
chi_3 = zeros(41,41);
for j=1:1:41
for i=1:1:41
p_5 = m1(i);
p_6 = m2(j);
mtrue =[p_5 p_6]';
p_4 = mtrue - mL2m1m2;
epp = (p_4'/cov_m1m2)*p_4;
chi_1(i,j) = chi2pdf(epp,2);
if chi_1(i,j) > 5.9915
chi_1(i,j) = 5.9915;
end
%chi_1(i,j) = chi2inv(0.95,chi_a);
clear p4 p5 p6 epp
end
end
clear i j
for j=1:1:41
for i=1:1:41
p_5 = m2(i);
p_6 = m3(j);
mtrue =[p_5 p_6]';
p_4 = mtrue - mL2m2m3;
epp = (p_4'/cov_m2m3)*p_4;
chi_2(i,j) = chi2pdf(epp,2);
if chi_2(i,j) > 5.9915
chi_2(i,j) = 5.9915;
end
%chi_2(i,j) = chi2inv(0.95,chi_a);
clear p4 p5 p6 epp
end
end
clear i j
for j=1:1:41
for i=1:1:41
p_5 = m1(i);
p_6 = m3(j);
mtrue =[p_5 p_6]';
p_4 = mtrue - mL2m1m3;
epp = (p_4'/cov_m1m3)*p_4;
chi_3(i,j) = chi2pdf(epp,2);
if chi_3(i,j) > 5.9915
chi_3(i,j) = 5.9915;
end
clear p4 p5 p6 epp
end
end
clear i j
chi_1;
x_cord_1 = m1;
x_cord_2 = m2;
x_cord_3 = m3;
[X_1,X_2] = meshgrid(x_cord_1,x_cord_2);
figure
hold on
Z_1 = chi_1;
subplot(3,1,1);
contour(X_1,X_2,Z_1)
[X_2,X_3] = meshgrid(x_cord_2,x_cord_3);
Z_2 = chi_2;
subplot(3,1,2);
contour(X_2,X_3,Z_2)
[X_1,X_3] = meshgrid(x_cord_1,x_cord_3);
Z_3 = chi_3;
subplot(3,1,3);
contour(X_1,X_3,Z_3)

回答 (0 件)

カテゴリ

Help Center および File ExchangeSurface and Mesh Plots についてさらに検索

Community Treasure Hunt

Find the treasures in MATLAB Central and discover how the community can help you!

Start Hunting!

Translated by