[R] mle
Peter Dalgaard
P.Dalgaard at biostat.ku.dk
Tue Jan 8 17:51:50 CET 2008
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
More information about the R-help
mailing list