[R] mle

[Peter Dalgaard]

Antonio Gasparrini wrote:
>  
> Hello,
>  
> I'm trying to obtain a maximum likelyhood estimate of a gaussian model
> by the MLE command, as I did with a Poisson model:
>  
> x <- rep(1:2,each=500)
> y <- rnorm(length(x), mean=10+3*x, sd=1)
>  
> glm1 <- glm(y ~ x , family=gaussian())
>  
> library(stats4)
> func1 <- function(alfa=10, beta=3, sigma=1)
>     -sum(dnorm(y, mean=alfa+beta*x, sd=sigma), log=FALSE)
> mle(func1, method = "BFGS")
>
> func2 <- function(alfa=10, beta=3, sigma=1)
>    
> -sum((1/sqrt(2*pi*sigma^2))*exp(-0.5*(((y-alfa-beta*x)^2)/sigma^2)))
> mle(func2, method = "BFGS")
>  
> I don't understand why it doesn't work.
> Have you some suggestions?
>  
> Thank you so much for your help
>   
You need the minus LOG likelihood, so in func1, you need log=TRUE, and
it belongs inside dnorm(...).

In func2 you need the similar change, or replace sum() with log(prod()).
(The latter might be somewhat less stable numerically, though)

In either case remember that MLE for sigma is not the usual least
squares estimate.
> Antonio Gasparrini
> Public and Environmental Health Research Unit (PEHRU)
> London School of Hygiene & Tropical Medicine
> Keppel Street, London WC1E 7HT, UK
> Office: 0044 (0)20 79272406
> Mobile: 0044 (0)79 64925523
> www.lshtm.ac.uk/pehru/ 
> antonio.gasparrini at lshtm.ac.uk
>
> 	[[alternative HTML version deleted]]
>
> ______________________________________________
> R-help at 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.
>   


-- 
   O__  ---- Peter Dalgaard             Øster Farimagsgade 5, Entr.B
  c/ /'_ --- Dept. of Biostatistics     PO Box 2099, 1014 Cph. K
 (*) \(*) -- University of Copenhagen   Denmark          Ph:  (+45) 35327918
~~~~~~~~~~ - (p.dalgaard at biostat.ku.dk)                  FAX: (+45) 35327907