[R] 'varFunc' class with additive variance? (was: can I get same results using lme and gls?)

From: Spencer Graves <spencer.graves_at_pdf.com>
Date: Mon, 21 May 2007 10:35:40 -0700

Hi, Doug and others:

      What might be the best tools for modeling an additive variance structure for residuals, something like the following:

           var(resid) = s0.2*(1+var.pred) + daysFromTraining*var(process migration per day),

      where var.pred = relative variance of prediction error = something roughly like crossprod(x, solve(crossprod(X), x)),

      and var(process migration per day) = the variance of a random walk of some process characteristic.

      My data are residuals from predictions from models produced by a data mining algorithm with various choices for training and test sets. I've been using 'gls' to fit models using 'varFunc' classes. However, 'varComb' models the relative standard deviation as a product of components from different sources. I'm thinking of creating 'varSumSq' functions by modifying your 'varComb' code to model an additive (not multiplicative) variance structure.

      Might there be other tools for modeling an additive variance structure? Alternatively, does it sound sensible to create 'varSumSq' functions similar to 'varComb'?

      Thanks,
      Spencer Graves    

Douglas Bates wrote:
> On 5/20/07, toby909@gmail.com <toby909@gmail.com> wrote:
>
>
>> I was wondering how to get the same results with gls and lme. In my lme, the
>> design matrix for the random effects is (should be) a identity matrix and
>> therefore G should add up with R to produce the R matrix that gls would report
>> (V=ZGZ'+R). Added complexity is that I have 3 levels, so I have R, G and say H
>> (V=WHW'+ZGZ'+R). The lme is giving me the correct results, I am having trouble
>> finding the right corresponding specification for the gls.
>>
>
> Thanks for including a reproducible example. However, I'm a bit at a
> loss as to why you would want to try to create a gls model that fits a
> mixed-effects model that has random effects for intercept and slope at
> two nested levels. I don't think that corCompSymm will do what you
> want but, to tell the truth, I have difficulty in thinking of the
> model in that form. I much prefer the mixed-effects form.

>
>
>
>> Thanks for your help.
>>
>> Toby
>>
>>
>> dtaa =
>> read.table("http://www.ats.ucla.edu/stat/mplus/examples/ma_snijders/mlbook1.dat",
>> sep=",")
>> dta1 = reshape(dtaa, list(c("V10","V12")), "score", direction="long",
>> drop=c("V2","V3","V4","V5","V6","V7","V8","V9","V11","V13","V14","V15","V16","V17","V18","V19","V20","V21","V22","V23","V24","V25"))
>> colnames(dta1)[1] = "schoolNR"
>> dta2 = dta1[order(dta1$id),]
>> head(dta2)
>> timef = factor(dta2$time)
>>
>> summary(mdl1l <- lme(score~timef-1, dta2, ~timef-1|schoolNR/id,,,,"ML"))
>> summary(mdl1g <- gls(score~timef-1, dta2, corCompSymm(, ~timef|schoolNR/id),
>> varIdent(, ~1|id*timef),,"ML"))
>>
>> ______________________________________________
>> 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
>> and provide commented, minimal, self-contained, reproducible code.
>>
>>
>
> ______________________________________________
> 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
> and provide commented, minimal, self-contained, reproducible code.
>



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 and provide commented, minimal, self-contained, reproducible code. Received on Mon 21 May 2007 - 17:42:50 GMT

Archive maintained by Robert King, hosted by the discipline of statistics at the University of Newcastle, Australia.
Archive generated by hypermail 2.2.0, at Mon 21 May 2007 - 18:31:37 GMT.

Mailing list information is available at https://stat.ethz.ch/mailman/listinfo/r-help. Please read the posting guide before posting to the list.