I think the problem lies in the mixture model. First, why do you have two
free parameters p and q. Shouldnt that be p and (1-p)? Secondly, I am not
sure that your data, d2, is compatible with a binary mixture model. It seems
like a sensible binary mixture model cannot be fitted for your data. It may
be over parametrized. A model with just a single beta distribution seems to
do ok.
loglik2 <- function(par, data) {
sh11 <- par[1]
sh12 <- par[2]
sum(dbeta(data,shape1=sh11,shape2=sh12, log=TRUE))
}
optim(par=c(0.1,0.1), fn=loglik2, data=d2, lower=c(0,0),
control=list(fnscale=-1))
Ravi.
----------------------------------------------------------------------------
-------
Ravi Varadhan, Ph.D.
Assistant Professor, The Center on Aging and Health
Division of Geriatric Medicine and Gerontology
Johns Hopkins University
Ph: (410) 502-2619
Fax: (410) 614-9625
Email: rvarad...@jhmi.edu
Webpage:
http://www.jhsph.edu/agingandhealth/People/Faculty_personal_pages/Varadhan.h
tml
----------------------------------------------------------------------------
--------
-----Original Message-----
From: r-help-boun...@r-project.org [mailto:r-help-boun...@r-project.org] On
Behalf Of Nicolas Turenne
Sent: Wednesday, November 04, 2009 12:17 PM
To: r-help@r-project.org
Subject: [R] pb with optimimization and fitdistr
Hello
i try to fit a data series (N below) with a model consisting of a
mixture of two beta distributions
for that i am using fitdistr of package MASS
as follows
> library(MASS)
> N=c(796,3586,4089,3364,2745,1992,1120,432,99,10,0,0)
> d2 = (N-min(N)+0.000001)/(max(N)-min(N)+0.002)
> mixtBeta <- function(x,p,q,sh11,sh12,sh21,sh22)
p*dbeta(x,shape1=sh11,shape2=sh12)+(q)*dbeta(x,shape1=sh21,shape2=sh22)
> ff= fitdistr(d2,mixtBeta ,
start=list(p=0.2,q=0.1,sh11=2,sh12=9,sh21=6,sh22=10))
Il y a eu 50 avis ou plus (utilisez warnings() pour voir les 50 premiers)
> ff
p q sh11 sh12
sh21 sh22
2.265173e+03 2.751748e+03 6.823009e-02 1.230362e-01
7.130326e-01 1.069649e+00
but manually i find another set of parameters
p=0.26;q=0.12; sh11=2.1, sh12=9, sh21=5.8, sh22=9
S=c(0,1,2,3,4,5,6,7,8,9,10,11)/12
N=c(796,3586,4089,3364,2745,1992,1120,432,99,10,0,0)
d2 = (N-min(N)+0.001)/(max(N)-min(N)+0.002)
x = c(1:100)/100
p=0.26;q=0.12; y = p*dbeta( x, 2.1, 9 ) + q*dbeta( x, 5.8, 9 );
plot(S,d2,xlim=c(0,1), ylim=c(0,1))
lines(x, y )
Is there a problem with fitdistr if one uses a function with parameters ?
or a problem of initialization ?
thanks for help to understand.
regards
Nicolas
--
Nicolas Turenne - INRA, Unité Mathématique Informatique et Génome UR1077,
F-78350 Jouy-en-Josas
http://genome.jouy.inra.fr/~turenne/nt.html
______________________________________________
R-help@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.
______________________________________________
R-help@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.
