Estimating the residual variance by using a varFunc in nlme

Question

I'm using a linear model to predict a dependent variable (the score to a test) starting by the age of a person. As shown in the figure below, the variability decreases when the age increases (data are simulated!).

Since this characteristic trend of the variance, I'm using the Generalized Least Squares method implemented in R in the function gls() of the package nlme, by using the varExp variance function:

library("nlme")
vf <- varExp(form=~age)
fit <- gls(score~age, data=dataset, weights=vf)

The variance function provides a good distribution of residuals.

Now, I need to estimate, for some value of age, the corresponding value of the residual variance. However, the function predict() extracts only the fitted values, not the variances. There is a method to estimated the residual variance of each age, starting by the model?

Answer 1

Looking at the help page for varExp, you can see that the function is defined as s2(v) = exp(2* t * v). Where t is the variance function coefficient and v is the covariate (in your case age).

fit$modelStruct$varStruct should give you the parameter estimate for your variance coefficient. You can plot the variance against age by using plot(y= exp( 2* fit$modelStruct$varStruct * dataset$age),x= dataset$age)

Of course you can enter any value for age into exp( 2* fit$modelStruct$varStruct * age) to get a variance estimate for that particular age. Hope this answers your question.