[R-sig-ME] autocorrelation

Mon Mar 31 17:09:44 CEST 2008

On Mon, Mar 31, 2008 at 9:06 AM, Mike Dunbar <mdu at ceh.ac.uk> wrote:
> Am I missing something here? If you need to estimate an autocorrelation parameter at say lag 1, why not make yourself a new column with your response variable lagged by 1 time unit. Then include that as a fixed effect. Clearly there are issues with missing data, but I'm not aware that nlme does anything more than you can do manually. In fact I wonder if this approach is slightly more flexible, as including a random slope for that lagged variable allows it to vary between groups, and I'm not aware that this is allowed using the in-built structures in nlme. Anything more complex, as I'm continually told, there's always Winbugs...

The approach of incorporating a function of the response as a
predictor in a model generally leads you into a quagmire of modeling
the noise term.  The approach in nlme is to "pre-whiten" the response
and the model matrices.  Given a value of the autocorrelation
parameter, the response is transformed by the inverse of the
corresponding correlation matrix and the model matrices are similarly
transformed.  Following that the log-likelihood is evaluated and
optimized with respect to the variance components and the
autocorrelation parameter simultaneously.

>  >>> "Douglas Bates" <bates at stat.wisc.edu> 30/03/2008 22:43:32 >>>
>
>
> On Sat, Mar 29, 2008 at 2:52 PM, Alan Cobo-Lewis <alanc at umit.maine.edu> wrote:
>  > Doug Bates writes on r-sig-mixed-models at r-project.org on Saturday, March 29, 2008 at 7:00 AM -0500 wrote about his planned book on multilevel modelling in R:
>  >
>  >  >I emphasize graphical displays of the data and aspects
>  >  >of the fitted models and inferences based on MCMC samples from the
>  >  >posterior distribution of the model parameters.
>  >
>  >  (n)lme handled correlated error terms, but lme4 does not.
>
>  So if you want a model with correlated error terms (in addition to the
>  correlation induced by the random effects) then you should use the
>  nlme package.
>
>  > Leaving aside the superior algorithms in lme4, this appears to be the major impediment to considering lme4 capabilities as a superset of (n)lme capabilities.
>
>  I don't recall any statements to the effect that the lme4 capabilities
>  would be a superset of the nlme capabilities.  It seems that whoever
>  made that decision should have informed me of it.
>
>  The development of the lme4 package has been generously funded by
>  several grants, the most important of which was an STTR contract that
>  we had for 3 years.  The purpose of that contract was to develop a
>  package that could fit generalized linear mixed models using the
>  Laplace approximation and allowing for crossed or partially crossed
>  grouping factors for the random effects.  The development is currently
>  funded by another grant specifically to provide for fitting models
>  with crossed and partially crossed random effects and with carryover
>  of random effects from one time period to another.
>
>  Neither generalized linear mixed models nor models with crossed or
>  partially crossed random effects can be fit (well without resorting to
>  egregious kludges) with the nlme package.  Even nonlinear mixed models
>  as fit by nlme are sub-optimal compared to the methods in lme4.  (lme4
>  uses direct optimization of the Laplace approximation to the
>  log-likelihood whereas nlme uses an alternating algorithm that Mary
>  Lindstrom and I proposed.)
>
>  My priorities are to fulfill the tasks that I proposed for these
>  grants and to build the best software that I can.  The beauty of open
>  source software is that if your priorities are different, you have
>  full access to the sources and you can modify them to fulfill your
>  objectives.  So I suggest that you
>   - Continue to use the nlme package if you wish to incorporate
>  (additional) correlation structures in models
>   - Design, code, test and document extensions to the lme4 package to
>  do so and then submit these changes as patches
>   - Develop your own package so you can have things done the way you
>  want them to be done.  You already have access to the lme4 sources so
>  a lot of the heavy lifting has been done for you.
>
>  This developing the software is not as easy as it may seem.  There are
>  many trade-offs and, at least for me, it takes a lot of effort to
>  determine even if it is possible to incorporate various extensions
>  harmoniously.  It is possible to model the mean and variance of the
>  conditional distribution of the response separately when that
>  distribution is multivariate normal.  It is not as easy to do so when
>  that distribution is binomial or Poisson or some other distribution
>  for a generalized linear mixed model.  Because the first purpose of
>  the lme4 package was to allow for generalized linear mixed models I
>  did not incorporate (additional) correlation structures and variance
>  functions.  I'm not even sure it could be done consistently for GLMMs
>  but you are welcome to show us how.
>
>  >  But what do I do if I've got, for example, autocorrelated error terms? Is there a way to "trick" lme4 into handling that (perhaps something analogous to the "random effect variance per treatment group in lmer" thread that David Afshartous and I
>  >  participated in)? Is there instead a good argument for ignoring it? It seems like something that would arise in practice in a non-negligible amount of problems in real data. Will the upcoming book give some advice on how to address this?
>  >
>  >  I can produce self-contained reproducible code if necessary, but I don't think it is.
>  >
>  >  thanks
>  >  alan
>  >
>  >  _______________________________________________
>  >  R-sig-mixed-models at r-project.org mailing list
>  >  https://stat.ethz.ch/mailman/listinfo/r-sig-mixed-models
>  >
>
>  _______________________________________________
>  R-sig-mixed-models at r-project.org mailing list
>  https://stat.ethz.ch/mailman/listinfo/r-sig-mixed-models
>
>  --
>  This message (and any attachments) is for the recipient ...{{dropped:6}}
>
>
>
>  _______________________________________________
>  R-sig-mixed-models at r-project.org mailing list
>  https://stat.ethz.ch/mailman/listinfo/r-sig-mixed-models
>