Plotting a two peak model as a fitted line on a probplot.

5 次查看（过去 30 天）
Andrew Feenan
Andrew Feenan 2022-10-5
评论： dpb 2022-10-7
Hi,
I am trying to plot a two peak distribution model on a probplot.
For example, Fit Custom Distributions - MATLAB & Simulink Example - MathWorks United Kingdom
and then using this example from MathWorks: Probability plots - MATLAB probplot - MathWorks United Kingdom
However, as my distribution uses a custom mle function is doesn't seem to work as desired. This is the example code I am working with.
%example from fist link above
clear; clc;
rng(10) % For reproducibility
x3 = [trnd(20,1,50) trnd(4,1,100)+3];
histogram(x3)
pdf_normmixture = @(x3,p,mu1,mu2,sigma1,sigma2) ...
p*normpdf(x3,mu1,sigma1) + (1-p)*normpdf(x3,mu2,sigma2);
pStart = .5;
muStart = quantile(x3,[.25 .75]);
sigmaStart = sqrt(var(x3) - .25*diff(muStart).^2);
start = [pStart muStart sigmaStart sigmaStart];
lb = [0 -Inf -Inf 0 0];
ub = [1 Inf Inf Inf Inf];
options = statset('MaxIter',300,'MaxFunEvals',600);
paramEsts = mle(x3,'pdf',pdf_normmixture,'Start',start, ...
'LowerBound',lb,'UpperBound',ub,'Options',options);
xgrid = linspace(1.1*min(x3),1.1*max(x3),200);
pdfgrid = pdf_normmixture(xgrid, ...
paramEsts(1),paramEsts(2),paramEsts(3),paramEsts(4),paramEsts(5));
figure(1)
histogram(x3,'Normalization','pdf')
hold on
plot(xgrid,pdfgrid,'-')
hold off
xlabel('x3')
ylabel('Probability Density')
legend('Sample Data','Fitted pdf','Location','best')
%end of first example
%example from second link above
figure(2)
probplot(x3,'noref')
h = probplot(gca,pdf_normmixture,paramEsts);
h.Color = 'r';
h.LineStyle = '-';
%end of second example
Is there a way to display something like this on a probplot?
  3 个评论
显示 1更早的评论
Andrew Feenan
Andrew Feenan 2022-10-7
Thanks for your response.
The distribution fitter app isn't particularly suitable for my use case but thanks for looking into it.
It' frustrating this isn’t possible on a probplot like it is on a histogram etc.
dpb
dpb 2022-10-7
编辑：dpb 2022-10-7
Agreed. This penchant of having multiple disparate objects that are similar but not the same is a pain. It seems the Statistics Toolbox has been particularly bad about introducing stuff before it's fully ripe (the ill-fated datatset comes to mind, too.)
Looks like you'll have to build the probablllity plot for the custom distribution by evaluating it directly without the benefit of the prepared probplot.
I've not looked inside it so see just how much one might be able to cop as a starting point or if it is exceedingly complex making use of the zillion undocumented behind-the-sceenes pieces.

请先登录，再进行评论。

请先登录，再回答此问题。

采纳的回答

dpb
dpb 2022-10-7
编辑：dpb 2022-10-7
Ran in:
Actually, when I had a few minutes to think about it, it's not too hard; you just have to define the associated CDF function to go w with your PDF and use it instead...
%example from fist link above
rng(10) % For reproducibility
x3 = [trnd(20,1,50) trnd(4,1,100)+3];
%histogram(x3)
pdf_normmixture = @(x3,p,mu1,mu2,sigma1,sigma2) ...
p*normpdf(x3,mu1,sigma1) + (1-p)*normpdf(x3,mu2,sigma2);
pStart = .5;
muStart = quantile(x3,[.25 .75]);
sigmaStart = sqrt(var(x3) - .25*diff(muStart).^2);
start = [pStart muStart sigmaStart sigmaStart];
lb = [0 -Inf -Inf 0 0];
ub = [1 Inf Inf Inf Inf];
options = statset('MaxIter',300,'MaxFunEvals',600);
paramEsts = mle(x3,'pdf',pdf_normmixture,'Start',start, ...
'LowerBound',lb,'UpperBound',ub,'Options',options);
xgrid = linspace(1.1*min(x3),1.1*max(x3),200);
pdfgrid = pdf_normmixture(xgrid, ...
paramEsts(1),paramEsts(2),paramEsts(3),paramEsts(4),paramEsts(5));
probplot(x3,'noref')
% define a cdf function to match...
cdf_normmixture = @(x3,p,mu1,mu2,sigma1,sigma2) ...
p*normcdf(x3,mu1,sigma1) + (1-p)*normcdf(x3,mu2,sigma2);
h = probplot(gca,cdf_normmixture,paramEsts);
h.Color = 'r';
h.LineStyle = '-';
title('Probability plot for Mixed Normal Distribution')
  2 个评论
Andrew Feenan
Andrew Feenan 2022-10-7
Thank you very much. You have saved me so much time. I was in the middle of building the probability plot for my distribution function.
dpb
dpb 2022-10-7
Glad to help -- I knew that was what needed to be able to add; just wasn't thinking clearly how simple the CDF function was to write until came to me in a realization was over-complicating things...

请先登录，再进行评论。

更多回答（0 个）

类别

在 Help Center 和 File Exchange 中查找有关 Uniform Distribution (Continuous) 的更多信息

标签

产品


版本

R2022b

Community Treasure Hunt

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

Start Hunting!

Translated by