[R-sig-ME] Extracting the coefficients from corStruct objects

[Nilsson Fredrik X]

> The reason for the strange parameterization is that the nlme package
> always uses unconstrained parameters for the optimization.  When the
> parameters have natural constraints it is necessary to transform them
> to produce unconstrained versions.
>
It would be interesting to know what the transformation is or how to find out. Is it tanh(x/2)?

Best regards,
Fredrik

>> But here's a solution while waiting from the enlightenment:
>>
>> Look at:
>> intervals(ModelwAR1.lme)$corStruct
>>
>> And then
>>
>> intervals(ModelwAR1.lme)$corStruct[2]
>>
>> is what you look for.
>>
>> Cheers,
>>
>> Fredrik
>>
>> -----Ursprungligt meddelande-----
>> Från: r-sig-mixed-models-bounces at r-project.org [mailto:r-sig-mixed-models-bounces at r-project.org] För Markus Jäntti
>> Skickat: den 9 maj 2007 23:49
>> Till: R-SIG-Mixed-Model
>> Ämne: [R-sig-ME] Extracting the coefficients from corStruct objects
>>
>> Dear All:
>>
>> I have been trying to get the coefficients from corAR1 objects in lme objects in
>> then lme package (3.1-80), but coef(obj$modelStruct$corStruct) returns something
>> rather different that summary(object). I have tried to follow the various
>> summary. methods in the source code, but have kind of gotten lost. Would be
>> greatful if someone could point to a simple solution.
>>
>> It would be helpful to manually extract the AR(1) coefficients, because I am
>> data sets consisting of very large cohorts of individuals. I extract the
>> information I need from each estimated object and store them in small matrices.
>> I have not been able to figure out a way to extract the AR(1) parameters, however.
>>
>> I include example code below, which produces the following results:
>>  > summary(fm2)
>> <snip>
>> Correlation Structure: AR(1)
>>   Formula: ~year | ind
>>   Parameter estimate(s):
>>    Phi
>> 0.424
>> <snip>
>>
>> but
>>
>>  > coef(fm2$modelStruct$corStruct)
>> [1] 0.906
>>
>> Version info for nlme:
>>
>>  > packageDescription("nlme")
>> Package: nlme
>> Version: 3.1-80
>> Date: 2007-03-30
>> Priority: recommended
>> Title: Linear and Nonlinear Mixed Effects Models
>> Author: Jose Pinheiro <Jose.Pinheiro at pharma.novartis.com>, Douglas
>>          Bates <bates at stat.wisc.edu>, Saikat DebRoy
>>          <saikat at stat.wisc.edu>, and Deepayan Sarkar
>>          <deepayan at stat.wisc.edu>
>> Maintainer: R-core <R-core at R-project.org>
>> Description: Fit and compare Gaussian linear and nonlinear
>>          mixed-effects models.
>> Depends: graphics, stats, R (>= 2.3.0)
>> Imports: lattice
>> LazyLoad: yes
>> LazyData: yes
>> License: GPL version 2 or later
>> Packaged: Fri Mar 30 07:37:02 2007; ripley
>> Built: R 2.5.0; i486-pc-linux-gnu; 2007-04-25 03:13:00; unix
>>
>> Example code:
>> <code>
>> library(nlme)
>> ## number of time periods
>> years <- 5
>> ## individuals
>> n <- 100
>>
>> ## the random effects
>> a <- as.vector(sapply(rnorm(n), rep, years))
>>
>> ## indicators for individuals
>> ind <- as.vector(sapply(1:n, rep, years))
>> year <- rep(1:years, n)
>>
>> ## the white noise
>> v <- rnorm(n*years)
>> ## the AR(1) parameter
>> rho <- .5
>> ## generate the AR1 erroer
>> ## randomly draw the starting value of the AR(1) error for each unit
>> u1 <- rnorm(n)
>> ## store values here:
>> u <- numeric(n*years)
>> ## a counter:
>> k <- 1
>> ## there has to be a much more elegant way to do this!
>> for(i in 1:n)
>>    {
>>      u[k] <- u1[i]
>>      for(j in 2:years)
>>        {
>>          k <- k+1
>>          u[k] <- rho*u[k-1] + v[k]
>>        }
>>      k <- k+1
>>    }
>>
>> ## the covariate
>> x <- runif(n*years)
>> ## the intercept and slope
>> b0 <- 10
>> b1 <- 4
>> ## and the dependent variable
>> y <- b0 + b1*x + a + u
>>
>> test.d <- data.frame(y=y, x=x, ind=ind, year=year)
>> rm("x", "y", "ind", "years")
>>
>> fm1 <- lme(y ~ x, random = ~ 1 | ind,
>>             data=test.d)
>> summary(fm1)
>> fm2 <- lme(y ~ x, random = ~ 1 | ind, correlation=corAR1(form=~year),
>>             data=test.d)
>> summary(fm2)
>> coef(fm2$modelStruct$corStruct)
>> </code>
>>
>> --
>> Markus Jantti
>> Abo Akademi University
>> markus.jantti at iki.fi
>> http://www.iki.fi/~mjantti
>>
>> _______________________________________________
>> 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
>>