Parametric mixture model of survival time data in matlab
显示 更早的评论
Good day every one, I am new to matlab and working on parametric mixture model of survival time data. I would like to fit a mixture of weibull and gamma distributions to estimate the mixing probabilities and the parameters of the distributions via Expectation Maximization Algorithm EM. If any one how has experience in modeling survival mixtures to help and put me through.
回答(2 个)
Tom Lane
2012-4-22
Do you need to use the EM algorithm? The mle function may be able to handle a problem like this:
x = [gamrnd(2,2,50,1); wblrnd(20,2,100,1)];
f = @(x,p,a,b,c,d) p*gampdf(x,a,b) + (1-p)*wblpdf(x,c,d);
mle(x,'pdf',f,'start',[.5 2 2 10 2])
2 个评论
nasrin
2015-8-21
Hi Tom,
I am using mle function to fit mixture of gamma distributions. so I do:
f = @(x,p,a,b,c,d) p*gampdf(x,a,b) + (1-p)*gampdf(x,c,d); mle(x,'pdf',f,'start',[.5 2 2 10 1])
but I get this error:
The PDF function returned negative or zero values
I think the problem is with choice of starting values. Is there any way to infer them from my data? Thank you :)
okoth ochola
2022-7-28
Hi Tom, kindly explain the meaning of codes in words so that one can easily apply the same in other data sets, for example what does line one give, why have you used the commands and are the meaning of the values inside the brackets? What do they do. Kindly explain if you dont mind. Take for example gamrnd(2,2,50,1); what does it mean?
x = [gamrnd(2,2,50,1); wblrnd(20,2,100,1)];
f = @(x,p,a,b,c,d) p*gampdf(x,a,b) + (1-p)*wblpdf(x,c,d);
mle(x,'pdf',f,'start',[.5 2 2 10 2])
Russ Adheaux
2012-9-17
0 个投票
I have a related question. I have a a 1D array, x: 0<x<5; 300 data points; sorted ascending. I want to fit Rayleigh distribution in the range 0<x<=xo and 2-parameter weibull in the range xo<x<5. Then I want to plot the CDF of these two distributions on a single (Rayleigh probability plot) plot. Tom's answer gave some clue, but I don't seem to get it working. Thanks.
3 个评论
Tom Lane
2012-9-18
Maybe you could clarify some more. It sounds like you want to estimate a Rayleigh parameter and two Weibull parameters at least. Do you know xo or do you need to estimate that? I guess you want to truncate the Rayleigh to [0,xo] and the Weibull to [xo,5]. Do you also expect to estimate the amount of probability to allocate to the Rayleigh vs. the Weibull? Then is it true that you want a probability plot against an untruncated Rayleigh?
Russ Adheaux
2012-9-19
Sorry for the missing info...
Yes, I know xo. (Once I get this working I will try to estimate it)
Yes, I want to truncate Rayleigh to [0,xo] and the 2-parameter Weibull to [xo,5].
Because I am specifying xo, I thought I am specifying probability allocation, in a way.
I would like to see 3 things on my Rayleigh Probability Plot 1)data 2) Straight line indicating Rayleigh distrib 3) Another bent line from xo onwards showing 2-parameter fitted Weibull distrib.
I expect my data to lie on Rayleigh line up to xo and then on the Weibull line from xo onwards.
This is trickier. Here's an example that works sometimes. First generate samples from two truncated distributions:
a = raylrnd(1,100,1);
b = wblrnd(4,10,100,1);
x0 = 3;
x = [a(a<x0); b(b>x0 & b<5)];
[N,X] = hist(x);
bar(X,N/((X(2)-X(1))*length(x)),'hist')
shg
Define truncated Rayliegh and Weibull distributions and a mixture.
rpdf = @(x,p1) (x<x0) .* raylpdf(x,p1)/max(.001,raylcdf(x0,p1));
wpdf = @(x,p2,p3) (x>x0 & x<5) .* wblpdf(x,p2,p3) / max(.001,(wblcdf(5,p2,p3)-wblcdf(x0,p2,p3)));
mypdf = @(x,p1,p2,p3,p4) (exp(p4)/(1+exp(p4)))*rpdf(x,abs(p1)) + (1/(1+exp(p4)))*wpdf(x,abs(p2),abs(p3));
Use mle to fit the mixture.
opt = statset('mlecustom');
opt = statset(opt,'FunValCheck','off','MaxIter',1e6);
p = mle(x,'pdf',mypdf,'start',[1 4 10 0])
Superimpose on the data, along with the separate distributions for reference.
xx = linspace(0,5);
line(xx,rpdf(xx,p(1)),'Color','g')
line(xx,wpdf(xx,p(2),p(3)),'Color','c')
line(xx,mypdf(xx,p(1),p(2),p(3),p(4)),'Color','r')
legend('hist(x)','Rayleigh','Weibull','Combined')
Note that I parametrized the mixture using p4 rather than working with the mixture proportion directly, because p4 doesn't need to be constrained to the interval [0,1].
[Note that I edited an earlier version of this post where I neglected to restrict the Rayleigh and Weibull functions to their range.]
类别
在 帮助中心 和 File Exchange 中查找有关 Weibull Distribution 的更多信息
Community Treasure Hunt
Find the treasures in MATLAB Central and discover how the community can help you!
Start Hunting!