Example on Maximum likelihood estimation using R

Hi,


You should re-write your function `fr' as follows:


fr <- function(par){
a <- par[1]
b <- par[2]
p <- par[3]
lambda <- par[4]

l <- 0.5*(lambda + b*p + (1-p)*(lambda-b))
l^2 > lambda*b*p
delta <- sqrt(abs(l^2 - b*p*lambda))
mt <- a/p*(1-exp(-l*t)*cosh(delta*t)-(l-b*p)*exp(-t)*sinh(delta*t)/delta)
logl <- sum(diff(y)*log(diff(mt))-diff(mt)-lfactorial(diff(y)))
return(-logl)
}

However, I don't understand what the following fragment is doing in `fr':

l^2 > lambda*b*p

Can you clarfy that?

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: rvaradhan at jhmi.edu

Webpage:  http://www.jhsph.edu/agingandhealth/People/Faculty/Varadhan.html

 

----------------------------------------------------------------------------
--------


-----Original Message-----
From: r-help-bounces at r-project.org [mailto:r-help-bounces at r-project.org]
On
Behalf Of toh
Sent: Thursday, September 04, 2008 9:15 PM
To: r-help at r-project.org
Subject: Re: [R] Maximum likelihood estimation


Yes I'm trying to optimize the parameters a, b, p and lambda where a > 0,
b

0 and 0 < p < 1. I attached the error message that I got when I run mle. 

t <- c(1:90)
y <-
c(5,10,15,20,26,34,36,43,47,49,80,84,108,157,171,183,191,200,204,211,2
17,226,230,

+
234,236,240,243,252,254,259,263,264,268,271,277,290,309,324,331,346,367,375,
381,
+
401,411,414,417,425,430,431,433,435,437,444,446,446,448,451,453,460,463,463,
464,
+
464,465,465,465,466,467,467,467,468,469,469,469,469,470,472,472,473,473,473,
473,
+ 473,473,473,475,475,475,475)

library(stats4)
fr <- function(a, b, p, lambda){

+ l <- 0.5*(lambda + b*p + (1-p)*(lambda-b))
+ l^2 > lambda*b*p
+ delta <- sqrt(abs(l^2 - b*p*lambda))
+ mt <- 
+ a/p*(1-exp(-l*t)*cosh(delta*t)-(l-b*p)*exp(-t)*sinh(delta*t)/delta)
+ logl <- sum(diff(y)*log(diff(mt))-diff(mt)-lfactorial(diff(y)))
+ return(-logl)
+ }

mle(start=list(a=12,b=0.01,p=0.5,lambda=0.01),fr, method="L-BFGS-B",

+ lower = c(0.002, 0.002, 0.002, 0.002), upper = c(Inf, Inf, 0.999,
Inf),control=list(fnscale=-1))
Error in optim(start, f, method = method, hessian = TRUE, ...) : 
  non-finite finite-difference value [3]












Prof Brian Ripley wrote:

From ?optim

       fn: A function to be minimized (or maximized), with first
           argument the vector of parameters over which minimization is
           to take place.  It should return a scalar result.

I think you intended to optimize over c(a,b,p, lambda), so you need to 
specify them as a single vector.

You may be making life unnecessarily hard for yourself: see function 
mle() in package stats4.

Showing your code without a verbal description of what you are doing 
nor the error message you got is less helpful than we need.

On Wed, 3 Sep 2008, toh wrote:

Hi R-experts,
I'm new to R in mle. I tried to do the following but just couldn't 
get it right. Hope someone can point out the mistakes. thanks a lot.

t <- c(1:90)
y <-
c(5,10,15,20,26,34,36,43,47,49,80,84,108,157,171,183,191,200,204,211,
217,226,230,

234,236,240,243,252,254,259,263,264,268,271,277,290,309,324,331,346,3
67,375,381,

401,411,414,417,425,430,431,433,435,437,444,446,446,448,451,453,460,4
63,463,464,

464,465,465,465,466,467,467,467,468,469,469,469,469,470,472,472,473,473,473,
473,

	473,473,473,475,475,475,475)
fr <- function(a, b, p, lambda){
l <- 0.5*(lambda + b*p + (1-p)*(lambda-b))
l^2 > lambda*b*p
delta <- sqrt(abs(l^2 - b*p*lambda))
mt <- 
a/p*(1-exp(-l*t)*cosh(delta*t)-(l-b*p)*exp(-t)*sinh(delta*t)/delta)
logl <- sum(diff(y)*log(diff(mt))-diff(mt)-lfactorial(diff(y)))
return(-logl)
}
optim(c(15,0.01,0.5,0.01),fr, method="L-BFGS-B", lower = c(0.002, 
0.002, 0.002, 0.002), upper = c(Inf, Inf, 0.999,
Inf),control=list(fnscale=-1))

--
View this message in context:
http://www.nabble.com/Maximum-likelihood-estimation-tp19304249p193042
49.html Sent from the R help mailing list archive at Nabble.com.

-- 
Brian D. Ripley,                  ripley at stats.ox.ac.uk
Professor of Applied Statistics,  http://www.stats.ox.ac.uk/~ripley/
University of Oxford,             Tel:  +44 1865 272861 (self)
1 South Parks Road,                     +44 1865 272866 (PA)
Oxford OX1 3TG, UK                Fax:  +44 1865 272595

--
View this message in context:
http://www.nabble.com/Maximum-likelihood-estimation-tp19304249p19323140.html
Sent from the R help mailing list archive at Nabble.com.