[R] Difficulty in calculating MLE through NLM

[Ravi Varadhan]

Hi Madan,

You are trying to find the MLE of a binary mixture distribution.  So, there are constraints on the parameters, as you have indicated.  The nlm() function cannot handle constraints. I would recommend one of the following functions (not necessarily in any order), all of which can handle box-constraints:

1. optim(), with method="L-BFGS-B")
2. nlminb()
3. spg() in package "BB"
4. Write an EM algorithm; this automatically imposes all the required constraints, and is guaranteed to converge (to a local maximum).

John Nash and I are writing a package called "optimx", which contains a function called optimx() that integrates all the different optimization functions for smooth nonlinear problems.  Thsi would make it extremely easy for you to do steps (1) - (3) above, using a single call to optimx().  Let me know if you are interested.

Hope this helps,
Ravi.

____________________________________________________________________

Ravi Varadhan, Ph.D.
Assistant Professor,
Division of Geriatric Medicine and Gerontology
School of Medicine
Johns Hopkins University

Ph. (410) 502-2619
email: rvaradhan at jhmi.edu


----- Original Message -----
From: Madan Gopal Kundu <Madan.Kundu at ranbaxy.com>
Date: Wednesday, July 1, 2009 9:25 am
Subject: [R] Difficulty in calculating MLE through NLM
To: r-help <r-help at stat.math.ethz.ch>


> Hi R-friends,
>  
>  Attached is the SAS XPORT file that I have imported into R using 
> following code
>  library(foreign)
>  mydata<-read.xport("C:\\ctf.xpt")
>  print(mydata)
>  
>  I am trying to maximize logL in order to find Maximum Likelihood 
> Estimate (MLE) of 5 parameters (alpha1, beta1, alpha2, beta2, p) using 
> NLM function in R as follows.
>  
>  # Defining Log likelihood - In the function it is noted as logL
>  > library(stats)
>  > loglike1<- function(x)
>  + {
>  + alpha1<-x[1]
>  + beta1<-x[2]
>  + alpha2<-x[3]
>  + beta2<-x[4]
>  + p<-x[5]
>  + n<- mydata[3]
>  + e<-mydata[4]
>  + f1<- ((1+beta1/e)^(-n))*((1+e/beta1)^(-alpha1))*gamma(alpha1+n)/(gamma(n+1)*gamma(alpha1))
>  + f2<- ((1+beta2/e)^(-n))*((1+e/beta2)^(-alpha2))*gamma(alpha2+n)/(gamma(n+1)*gamma(alpha2))
>  + logL=sum(log(p*f1+(1-p)*f2))
>  + logL<- -logL
>  + }
>  
>  # Supplying starting parameter values
>  > theta<-c(.2041,.0582, 1.4150, 1.8380,0.0969)
>  
>  # Calculating MLE using NLM function
>  > result<- nlm(loglike1, theta, hessian=TRUE, print.level=1)
>  
>  Now the problem is, this is not working as there is no improvement in 
> final parameter estimate over starting values and NLM just stops just 
> after 1 iteration with gradient value of all the 5 parameters as zero. 
> I have tried other set of starting values, but then also I am getting 
> final parameter estimates similar to starting values and iteration 
> stops just after one step. When I check for warnings, R displays 
> following kind of warnings:
>  
>   *   NA/Inf replaced by maximum positive value
>   *   value out of range in 'gammafn'
>  
>  Please suggest what I should do. I am expecting the final MLE of 
> alpha1, alpha2, beta1 and beta2 greater than 0 and P should lie 
> between 0 to 1.
>  
>  Thanks & Regards,
>  Madan Gopal Kundu
>  Biostatistician, CDM, MACR, Ranbaxy Labs. Ltd.
>  Tel(O): +91 (0) 1245194045 - Mobile: +91 (0) 9868788406
>  
>  
>  
>  (i) The information contained in this e-mail message is intended only 
> 
>  for the confidential use of the recipient(s) named above. This 
> message 
>  is privileged and confidential. If the reader of this message is not 
> 
>  the intended recipient or an agent responsible for delivering it to 
> the 
>  intended recipient, you are hereby notified that you have received 
> this 
>  document in error and that any review, dissemination, distribution, 
> or 
>  copying of this message is strictly prohibited. If you have received 
> this 
>  communication in error, please notify us immediately by e-mail, and 
> delete 
>  the original message. 
>  
>  (ii) Madan.Kundu at ranbaxy.com confirms that Ranbaxy shall not be 
> responsible if this 
>  email message is used for any indecent, unsolicited or illegal 
> purposes, 
>  which are in violation of any existing laws and the same shall solely 
> be 
>  the responsibility of Madan.Kundu at ranbaxy.com and that Ranbaxy shall 
> at all times be 
>  indemnified of any civil and/ or criminal liabilities or consequences 
> there. 
> ______________________________________________
>  R-help at r-project.org mailing list
>  
>  PLEASE do read the posting guide 
>  and provide commented, minimal, self-contained, reproducible code.