[R] weighted likelihood for lme

[Marco Geraci]

--- Spencer Graves <spencer.graves at pdf.com> wrote:

> 	  Thank you for providing such a marvelous example.
>  I wish I could 
> reward your dilligence with a simple, complete
> answer.  Unfortunately, 
> the best I can offer at the moment is a guess and a
> reference.  First, I 
> believe you are correct in that the "weights"
> argument describes "the 
> within-group heteroscedasticity structure".  To
> specify between-group 
> heterscadisticity, have you considered the
> following:
> 
>   foo <- Orthodont
>   foo$w <- c(rep(1, 5*4), rep(0.5, 22*4))
>  > lme(distance ~ 1, random = ~ w-1|Subject,
> +     method="ML", data = foo)
> Linear mixed-effects model fit by maximum likelihood
>    Data: foo
>    Log-likelihood: -258.8586
>    Fixed: distance ~ 1
> (Intercept)
>      23.8834
> 
> Random effects:
>   Formula: ~w - 1 | Subject
>                 w Residual
> StdDev: 3.370796 2.233126
> 
> Number of Observations: 108
> Number of Groups: 27
> 

Nice! I haven't considered this alternative. As you
said, it's not a complete answer, but, I say, it's a
smart start.

> 	  Second, have you consulted Pinheiro and Bates
> (2000) Mixed-Effects 
> Models in S and S-Plus (Springer)?  If you have not
> already, I encourage 
> you to spend some quality time with that book.  For
> me, this book helped 
> transform "lme" from an inaequately documented and
> unusable black box 
> into a simple, understandable, elegant tool.  I may
> have recommeded it 
> to more people than any other single work over the
> past five years.

I know the book by Pinheiro and Bates and I do really
have to spend some time with it.

Thanks very much,
Marco

> 
> 	  hope this helps.
> 	  spencer graves
> 
> Marco Geraci wrote:
> 
> >     Dear R users,
> >   I'm trying to fit a simple random intercept
> model with a fixed intercept. 
> >     Suppose I want to assign a weight w_i to the
> i-th contribute to the log-likelihood, i.e.
> >    
> >   w_i * logLik_i
> >    
> >   where logLik_i is the log-likelihood for the
> i-th subject.
> >   I want to maximize the likelihood for N subjects
> >    
> >   Sum_i  {w_i * logLik_i}
> >    
> >   Here is a simple example to reproduce
> >    
> >   # require(nlme)
> >     > foo <- Orthodont
> > 
> >>lme(distance ~ 1, random = ~ 1|Subject,
> method="ML", data = foo)
> > 
> > Linear mixed-effects model fit by maximum
> likelihood
> >   Data: foo 
> >   Log-likelihood: -257.7456
> >   Fixed: distance ~ 1 
> > (Intercept) 
> >    24.02315 
> >   Random effects:
> >  Formula: ~1 | Subject
> >         (Intercept) Residual
> > StdDev:    1.888748 2.220312
> >   Number of Observations: 108
> > Number of Groups: 27 
> > 
> > Then I assign arbitrary weights, constant within
> the group. I want to give weight 1 to the first 5
> subjects, and weight 0.5 to the others 22 (4 is the
> number of repeated measurements for each subject)
> > 
> >    
> >   > foo$w <- c(rep(1, 5*4), rep(0.5, 22*4))
> >    
> >   Maybe I am missing something, but I believe that
> >    
> >   > lme(distance ~ 1, random = ~ 1|Subject,
> method="ML", data = foo, weight=~w)
> >    
> >   does not maximize the likelihood Sum_i  {w_i *
> logLik_i}, 
> since 'weight' describes the with-in
> heteroscedasticity structure.
> >   I think I need something like the option
> 'iweight' 
> (importance weight) for the command 'xtreg' of
> Stata.
> >    
> >   Any suggestion for R?
> >    
> >   Thanks in advance,
> >    
> >   Marco Geraci
> >    
> >   >  sessionInfo()
> > R version 2.2.1, 2005-12-20, i386-pc-mingw32 
> >   attached base packages:
> > [1] "methods"   "stats"     "graphics" 
> "grDevices" "utils"    
> > [6] "datasets"  "base"     
> >   other attached packages:
> >     nlme 
> > "3.1-66" 
> > 
> > 
> > 
> > 
> > 
> > 		
> > ---------------------------------
> >  
> >  What are the most popular cars?  Find out at
> Yahoo! Autos
> > 	[[alternative HTML version deleted]]
> > 
> > ______________________________________________
> > R-help at stat.math.ethz.ch mailing list
> > https://stat.ethz.ch/mailman/listinfo/r-help
> > PLEASE do read the posting guide!
> http://www.R-project.org/posting-guide.html
>