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

[Spencer Graves]

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 at gmail.com <toby909 at 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.
>