[R] how to avoid NaN in optim()

Ravi Varadhan rvaradhan at jhmi.edu
Thu Sep 30 21:04:28 CEST 2010


You also need the constrain that par[1] + par[2] < 1 in order to avoid NaNs.

You can do this using the `projectLinear' argument in  `spg'.

library(BB)
?spg


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: Ravi Varadhan <rvaradhan at jhmi.edu>
Date: Thursday, September 30, 2010 2:54 pm
Subject: Re: [R] how to avoid NaN in optim()
To: arindam fadikar <arindam.fadikar at gmail.com>
Cc: r-help at r-project.org


> Change the `else NA' to `else Inf' in your loglik function  and then 
> run the following:
>  
>  
>  library(BB)
>  
>  p0 <- runif(2)
>  
>  spg(p0, lik (176,182 , 60 ,17) , lower=0,  upper=1)
>  
>  
>  This will give you correct results.
>  
>  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: arindam fadikar <arindam.fadikar at gmail.com>
>  Date: Thursday, September 30, 2010 2:17 pm
>  Subject: [R] how to avoid NaN in optim()
>  To: r-help at r-project.org
>  
>  
>  > hi ,
>  >  
>  >  lik <- function(nO, nA, nB, nAB){
>  >  
>  >  loglik <- function(par)
>  >    {
>  >        p=par[1]
>  >        q=par[2]
>  >        r <- 1 - p - q
>  >  
>  >        if (c(p,q,r) > rep(0,3) && c(p,q,r) < rep(1,3) )
>  >  
>  >          {
>  >            -(2 * nO * log (r) + nA * log (p^2 + 2 * p * r)
>  >              + nB * log (q^2 + 2 * q * r)
>  >              + nAB * (log(2) +log(p) +log(q)))
>  >         }
>  >        else
>  >          NA
>  >    }
>  >  
>  >  loglik
>  >  
>  >  }
>  >  
>  >  
>  >  i want to maximize this likelihood function over the range (0,0,0) 
> to
>  >  (1,1,1). so i tried
>  >  
>  >  optim ( c(0.3,0.3), lik ( 176,182 , 60 ,17) , method = "CG")
>  >  
>  >  and obtained the following :
>  >  
>  >  > optim(c(0.3,0.3), fn, method="CG")
>  >  $par
>  >  [1] 0.26444187 0.09316946
>  >  
>  >  $value
>  >  [1] 492.5353
>  >  
>  >  $counts
>  >  function gradient
>  >       130       19
>  >  
>  >  $convergence
>  >  [1] 0
>  >  
>  >  $message
>  >  NULL
>  >  
>  >  Warning messages:
>  >  1: In log(q^2 + 2 * q * r) : NaNs produced
>  >  2: In log(q) : NaNs produced
>  >  3: In log(q^2 + 2 * q * r) : NaNs produced
>  >  4: In log(q) : NaNs produced
>  >  
>  >  
>  >  please help ...
>  >  
>  >  
>  >  -- 
>  >  Arindam Fadikar
>  >  M.Stat
>  >  Indian Statistical Institute.
>  >  New Delhi, India
>  >  
>  >  	[[alternative HTML version deleted]]
>  >  
>  >  ______________________________________________
>  >  R-help at r-project.org mailing list
>  >  
>  >  PLEASE do read the posting guide 
>  >  and provide commented, minimal, self-contained, reproducible code.
>  
>  ______________________________________________
>  R-help at r-project.org mailing list
>  
>  PLEASE do read the posting guide 
>  and provide commented, minimal, self-contained, reproducible code.



More information about the R-help mailing list