Center of gaussian mixed distribution area

6 ビュー (過去 30 日間)
Aishwarya Rao
Aishwarya Rao 2023 年 3 月 27 日
コメント済み: Aishwarya Rao 2023 年 4 月 11 日
Hi! I want to calculate the center of a gaussian mixed distribution. The code goes as follows:
dd = importdata('gmdata.txt'); %%%% import data
xx = dd(:,1);
yy = dd(:,2);
obj = gmdistribution.fit([xx,yy],2); %%%% fit a distribution
figure;
scatter(xx,yy,10,'.') %%%% Scatter plot with points of size 10
hold on
gmPDF = @(x,y) arrayfun(@(x0,y0) pdf(obj,[x0 y0]),x,y);
hh = fcontour(gmPDF,[-8 6]);
I am getting the above distribution as output. I want to recreate the value of GM distribution centers as given in "Tanabe, Hiroko, Keisuke Fujii, and Motoki Kouzaki. "Intermittent muscle activity in the feedback loop of postural control system during natural quiet standing." Scientific reports 7.1 (2017): 10631." The sentences read as:
Kindly help me with the code. I guess the center they are mentioning as the center of the fcontour rings, but not sure.
  1 件のコメント
Torsten
Torsten 2023 年 3 月 27 日
編集済み: Torsten 2023 年 3 月 27 日
Your data ressemble a usual multivariate Gaussian, don't they ?

サインインしてコメントする。

回答 (1 件)

the cyclist
the cyclist 2023 年 3 月 27 日
編集済み: the cyclist 2023 年 3 月 27 日
I downloaded and looked over the paper. I am also not certain what they mean by the "center" of the mixed gaussian model.
But, I think a pretty reasonable guess is they mean the peak of the joint PDF. (I think is also what you mean by the "center of the contour rings").
You can use the fminsearch function to find this peak:
dd = importdata('gmdata.txt'); %%%% import data
xx = dd(:,1);
yy = dd(:,2);
obj = gmdistribution.fit([xx,yy],2); %%%% fit a distribution
figure;
scatter(xx,yy,10,'.') %%%% Scatter plot with points of size 10
hold on
gmPDF = @(x,y) arrayfun(@(x0,y0) pdf(obj,[x0 y0]),x,y);
hh = fcontour(gmPDF,[-8 6]);
% Define the function to be minimized (using the negative because we actually want the maximum)
f = @(x) -obj.pdf(x);
% Initial guess as to where the maximum is.
% (This is a relatively poor guess, but the result looks OK anyway)
xy_init = [0 0];
% Use fminsearch to find the maximum
xy_max = fminsearch(f, xy_init)
xy_max = 1×2
-1.9393 -5.2685
Eye-balling the plot, this looks accurate.
  2 件のコメント
Torsten
Torsten 2023 年 3 月 27 日
Yes, it's one of the obj.mu (where both of them look pretty much the same).
Aishwarya Rao
Aishwarya Rao 2023 年 4 月 11 日
Thank you @the cyclist. The code worked!

サインインしてコメントする。

カテゴリ

Help Center および File ExchangeParticle & Nuclear Physics についてさらに検索

製品


リリース

R2023a

Community Treasure Hunt

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

Start Hunting!

Translated by