jepusto / lmeInfo

Information matrices for fitted nlme::lme and nlme::gls models
https://jepusto.github.io/lmeInfo/
4 stars 2 forks source link

Option to return multiple sigma-sq parameter estimates #40

Closed jepusto closed 2 years ago

jepusto commented 2 years ago

When fitting a gls() or lme() model with weights = varIdent(form = ~ 1 | Grp), the variance parameter estimates are returned as ratios of $\sigma^2$ values relative to a reference level. For example, if Grp has levels A, B, and C, the variance components are parameterized as $\sigma^2_B / \sigma^2_A$ and $\sigma^2_C / \sigma^2_A$. It would be useful to provide an option to return the variance parameter estimates and information matrix in terms of $\sigma^2_A$, $\sigma^2_B$, and $\sigma^2_C$.

The simplest way to implement this would be to add an option separate_variances to the extract_varcomp(), Fisher_info(), and varcomp_vcov() functions. The option would only be relevant for models with variance functions of class varIdent. In extract_varcomp(separate_variances = TRUE), it would just return the separate estimates of $\sigma^2$'s. In Fisher_info(separate_variances = TRUE), it would implement a delta-method transformation to get the correct information matrix. And then in varcomp_vcov(), the option would just need to be passed through to Fisher_info(separate_variances = separate_variances).

jepusto commented 2 years ago

@manchen07 what do you think?

manchen07 commented 2 years ago

@manchen07 what do you think?

Thanks for the explanation! I like this approach. In the above example, is the variance parameter $\delta_1 = 1$(reference), $\delta_2 = \sigma_B/\sigma_A$, and $\delta_3 = \sigma_C/\sigma_A$? To get the variance components with right parameterization, we would calculate $\sigma_A^2 = \sigma^2$, $\sigma_B^2 = \sigma_A^2 \delta_2^2$, and $\sigma_C^2 = \sigma_A^2 \delta_3^2$.

jepusto commented 2 years ago

Yes, I think that’s correct.