Center of gaussian mixed distribution area

3 次查看（过去 30 天）
Aishwarya Rao
Aishwarya Rao 2023-3-27
评论： Aishwarya Rao 2023-4-11
Ran in:
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 个）

the cyclist
the cyclist 2023-3-27
编辑：the cyclist 2023-3-27
Ran in:
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.

类别

在 Help Center 和 File Exchange 中查找有关 Exploration and Visualization 的更多信息

标签

产品


版本

R2023a

Community Treasure Hunt

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

Start Hunting!

Translated by