Azzera filtri
Azzera filtri

How to create a Gaussian mixed model (GMM) with many components?

3 visualizzazioni (ultimi 30 giorni)
Hi, I have a series of 1D data which is in the histogram. I want to use a GMM equivalent model for this data. I have written the following code, can you guide me in the continuation of this code?
clear all ; close all ; clc
warning off ;
load Data
% type 1: section [1,37,78,187,249,288,364,396,451,532,563,607]
X = 1:607 ;
D = D(X) ;
ft = fittype( 'smoothingspline' );
opts = fitoptions( 'Method', 'SmoothingSpline' );
opts.SmoothingParam = 0.015;
[SF, ~] = fit( X', D, ft, opts );
D = SF(X) ;
D = normalize(D-min(D),'norm') ;
%% Get Data from bar
figure ;
B = bar(X,D,'EdgeColor','none') ; hold on
Y = B.YData ;
%% peak clusters
[P data] = findpeaks(-Y+max(Y),'MinPeakWidth',4.5,'MinPeakDistance',10) ; % 'WidthReference','halfheight'
plot(data,P,'vg') ; data = data' ;
indexSep = [1 data' length(X)] ;
K = length(indexSep)-1 ;
index = zeros(size(X)) ;
for i=1:K
index(indexSep(i):indexSep(i+1)) = i ;
end
index = index' ;
%% Show Clusters
G = {'r*','b*','g*','k*','c*','y*','m*','ro','bo','go','ko','co','yo'} ;
for i = 1:K
plot(X(index==i),Y(index==i),[G{i}])
end
%% Does it need weight?
%% fit components
Model = 'a*exp(-((x-b)/c)^2)' ;
for i = 1:K
opts = fitoptions( 'Method', 'NonlinearLeastSquares','MaxIter', 500,'Robust','Bisquare');
needindex = find(maxk(Y(index==i),1)==Y(index==i)) ;
opts.StartPoint = [max(Y(index==i)) indexSep(i)+X(needindex(1)) (max(X(index==i))-min(X(index==i)))/2];
opts.Lower = [max(Y(index==i))-5 opts.StartPoint(2)-10 0] ;
opts.Upper = [max(Y(index==i))+5 opts.StartPoint(2)+10 Inf] ;
f = fit(X(index==i)',Y(index==i)',Model,opts) ;
AICoeff = coeffnames(f) ;
AIValue = cellstr(num2str(coeffvalues(f)')) ;
fs = replace(Model,AICoeff,AIValue) ;
fs = replace(fs,{'--','-+','+-','++','*','^'},{'+','-','-','+','.*','.^'}) ;
eval(['F' num2str(i) ' = fs ;'])
eval(['f' num2str(i) ' = inline(F' num2str(i) ') ;'])
% need fit with full sum(FK)
end
%% Show
newx = -400:1200 ;
for i = 1:K
eval(['newy(' num2str(i) ',:) = f' num2str(i) '(newx) ;'])
plot(newx,newy(i,:),'LineWidth',1,'Color','k')
end
axis([X(1)-10 X(end)+10 0 max(Y)])
%% In this section I want to see an equivalent GMM that matches this data
options = statset('MaxIter', 500);
gm = gmdistribution.fit(Y',K-1, 'Options', options);
gm = fitgmdist(Y',K) ;
y = pdf(gm,newx');
plot(newx,y,'LineWidth',2,'Color','r')
Indeed, I need a way to calculate GMM coefficients(correct gm.mu,gm.Sigma,\pi_{k}) and match each component.
  1 Commento
Mohamed Hassan
Mohamed Hassan il 22 Giu 2022
add this part to your code:
ComponentMeans = gm.mu;
ComponentCovariance = gm.Sigma;
MixtureProportions = gm.ComponentProportion
I have a question please. Could I apply your code to an image? I wanna show the gaussian distributions on the histogram.

Accedi per commentare.

Risposte (0)

Community Treasure Hunt

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

Start Hunting!

Translated by