How do I use fitgmdist with a set of data that is not in a histogram format
6 visualizaciones (últimos 30 días)
Mostrar comentarios más antiguos
Max Mortensen
el 18 de En. de 2021
Comentada: Max Mortensen
el 20 de En. de 2021
This set of data below gives me a trimodal distributions and I would like to run a statistical analysis on it. When I run fitgmdist it takes the data and turns it into a histogram and then tells me mixing prop and mean. Is there a way to run my data as it looks without Matlab turning it into a histogram and causing a large shift between 0-1?
X = [0.91,3.55,1.20,2.35,1.58,0.75,6.65,10.22,15.30,19.38,6.98,7.79,3.43,3.68,2.04,1.93,1.03,1.14,2.16,5.65,0.48,0.75,0.28,0.20,0.15,0.20,0.11,0.07,0.05];
bar(X, 'stacked')
figure(1)
3 comentarios
Adam Danz
el 19 de En. de 2021
With so few data and only 1 "repetition" of data for each x-value, curve fitting will be difficult.
An alternative would be to split the bars into 3 groups, one for each mode. That could be done with differentiation and smoothing or perhaps with findpeaks. Then fit each mode; fitdist could come in handy for that.
Respuesta aceptada
Adam Danz
el 19 de En. de 2021
Editada: Adam Danz
el 20 de En. de 2021
You could fit the distribution to a tri-modal mixture of gaussians but very little data and a very narrow peaks (1-2 bars), you'll need to set bounds to the parameters. Even then you won't get great fits.
The demo below walks through the use of fmincon and fits the data to a tri-modal mixture of gaussians. Its loosely based on a fitting proceedure I recently used to fit a bi-modal distribution (with much more data) and you'll likely need to tweek it. The results are highly sensitive to the bounds. See inline comments for details.
% Original data to fit
X = [0.91,3.55,1.20,2.35,1.58,0.75,6.65,10.22,15.30,19.38,6.98,7.79,3.43,3.68,2.04,1.93,1.03,1.14,2.16,5.65,0.48,0.75,0.28,0.20,0.15,0.20,0.11,0.07,0.05];
xvals = 1:numel(X);
% Define a tri-modal gaussian with 3 independent peaks with 11 parameters:
% b(1) vertical off
% b(2) slope
% b(3) amp1
% b(4) mu1
% b(5) std1
% b(6) amp2
% b(7) mu2
% b(8) std2
% b(9) amp3
% b(10) mu3
% b(11) std3
mixGausMdl = @(b,x) b(1) + b(2) * x(:, 1) + b(3) * exp( -(x(:, 1) - b(4)).^2/(2*b(5)^2)) + b(6) * exp(-(x(:, 1) - b(7)).^2/(2*b(8)^2)) + b(9) * exp(-(x(:, 1) - b(10)).^2/(2*b(11)^2));
% Set initial guesses, lower and upper bounds to each parameter
initGuess = [0, 0.1, 3, 2, 1.06, 19, 10, 5.88 4 20 1];
lowbounds = [0 -inf 0 0 .1 0 0 0.1 0 20 .1];
highbounds = [0 inf 10 10 10 inf 25 10 6 30 10];
% Poisson log-likelihood; fmincon finds the minimum of this function. Thus,
% will find the max of the poisson loglikelihood.
poissonLogLikelihood = @(r,lambda) r(:)'*log(lambda(:)) - sum(lambda);
% loss functions for fitting
nllFcn = @(params) -poissonLogLikelihood(X, mixGausMdl(params, xvals(:)));
% Do fitting
opts = optimset('LargeScale', 'off', 'MaxIter', 50000, 'Display', 'off');
[paramEst, nll, exitflag] = fmincon(nllFcn, initGuess, [], [], [], [], lowbounds, highbounds,[],opts);
% show results
T = table((1:numel(paramEst))',paramEst(:), initGuess(:), lowbounds(:), highbounds(:), ...
'VariableNames', {'idx','Estimates','initGuess','lowBound','HiBound'}, ...
'RowNames',{'VertOff','slope','amp1','mu1','std1','amp2','mu2','std2','amp3','mu3','std3'});
disp(T)
% Plot results
clf()
bar(xvals, X, 'FaceColor',[.8 .8 .8])
hold on
% Plot fit
xx = linspace(xvals(1),xvals(end),100);
yy = mixGausMdl(paramEst,xx(:));
plot(xx,yy,'r-','LineWidth',5)
% Show each guassian separately
gaus = @(x,mu,sig,amp,vo)amp*exp(-(((x-mu).^2)/(2*sig.^2)))+vo;
y1 = gaus(xx, paramEst(4), paramEst(5), paramEst(3), paramEst(1));
y2 = gaus(xx, paramEst(7), paramEst(8), paramEst(6), paramEst(1));
y3 = gaus(xx, paramEst(10), paramEst(11), paramEst(9), paramEst(1));
plot(xx,y1, 'c:', 'LineWidth',3)
plot(xx,y2, 'k:', 'LineWidth',3)
plot(xx,y3, 'b:', 'LineWidth',3)
Más respuestas (0)
Ver también
Community Treasure Hunt
Find the treasures in MATLAB Central and discover how the community can help you!
Start Hunting!