In the nlme package, I understand that to extract variance components in the random effects there's a VarCorr function. However, I recently also came across how it is done in another snippet of package
VC <- lapply(pdMatrix(lmeObject$modelStruct$reStruct), "*", lmeObject$sigma^2)[]
where it multiplies the matrix of the random effects with the within-subject standard variance.
May I know why is this done?