LS estimates using spatial covariance

[hqian-lilly]

Summary

When fitting mixed model using spatial covariance, the model fitting appears fine. However, extracting lsmeans estimates with "lsmeans" does not yield all the expected results. Lots of lsmeans and their differences are labelled as not estimable. All these lsmeans and their differences should be estimable.

mydata <- bcva_data %>% mutate(time=VISITN, visit=factor(VISITN))
myfit <- mmrm(formula = BCVA_CHG~BCVA_BL + visit + ARMCD + ARMCD:visit + sp_exp(time | USUBJID), data=mydata, method="Kenward-Roger", vcov="Kenward-Roger-Linear")
mylsmeans <- lsmeans(myfit, ~ARMCD*visit)
summary(mylsmeans)
mydiffs <- summary(pairs(mylsmeans, adjust="none")) 
mydiffs[1:20, ]

R session info

# R -e "utils::sessionInfo()" output goes here
R version 4.3.2 (2023-10-31)
Platform: x86_64-pc-linux-gnu (64-bit)
Running under: Red Hat Enterprise Linux

Matrix products: default
BLAS/LAPACK: /lrlhps/apps/intel/intel-2020/compilers_and_libraries_2020.0.166/linux/mkl/lib/intel64_lin/libmkl_gf_lp64.so;  LAPACK version 3.7.0

locale:
 [1] LC_CTYPE=en_US.UTF-8       LC_NUMERIC=C               LC_TIME=en_US.UTF-8        LC_COLLATE=en_US.UTF-8     LC_MONETARY=en_US.UTF-8    LC_MESSAGES=en_US.UTF-8   
 [7] LC_PAPER=en_US.UTF-8       LC_NAME=C                  LC_ADDRESS=C               LC_TELEPHONE=C             LC_MEASUREMENT=en_US.UTF-8 LC_IDENTIFICATION=C       

time zone: America/Indiana/Indianapolis
tzcode source: system (glibc)

attached base packages:
[1] stats     graphics  grDevices utils     datasets  methods   base     

other attached packages:
 [1] emmeans_1.8.9   mmrm_0.3.11     nortest_1.0-4   ggsci_3.0.0     lubridate_1.9.3 forcats_1.0.0   stringr_1.5.1   dplyr_1.1.4     purrr_1.0.2     readr_2.1.4    
[11] tidyr_1.3.0     tibble_3.2.1    ggplot2_3.4.4   tidyverse_2.0.0

loaded via a namespace (and not attached):
 [1] gtable_0.3.4       TMB_1.9.9          xfun_0.41          lattice_0.22-5     tzdb_0.4.0         vctrs_0.6.4        tools_4.3.2        Rdpack_2.6        
 [9] generics_0.1.3     parallel_4.3.2     sandwich_3.0-2     fansi_1.0.5        highr_0.10         pkgconfig_2.0.3    Matrix_1.6-4       checkmate_2.3.0   
[17] lifecycle_1.0.4    compiler_4.3.2     munsell_0.5.0      codetools_0.2-19   htmltools_0.5.7    yaml_2.3.7         pillar_1.9.0       MASS_7.3-60       
[25] multcomp_1.4-25    nlme_3.1-164       tidyselect_1.2.0   digest_0.6.33      mvtnorm_1.2-4      stringi_1.8.2      splines_4.3.2      fastmap_1.1.1     
[33] grid_4.3.2         colorspace_2.1-0   cli_3.6.1          magrittr_2.0.3     survival_3.5-7     utf8_1.2.4         TH.data_1.1-2      withr_2.5.2       
[41] scales_1.3.0       backports_1.4.1    estimability_1.4.1 timechange_0.2.0   rmarkdown_2.25     zoo_1.8-12         hms_1.1.3          coda_0.19-4       
[49] evaluate_0.23      knitr_1.45         rbibutils_2.2.16   rlang_1.1.2        Rcpp_1.0.11        xtable_1.8-4       glue_1.6.2         rstudioapi_0.15.0 
[57] R6_2.5.1       

OS / Environment

• OS: x86_64-pc-linux-gnu

[danielinteractive]

Thanks Hui-Rong for reporting this problem! We are going to look at this with high priority.

[clarkliming]

root cause of the issue:

1. if visit is not a fixed effect, then it should not be included as a valid choice in emmeans spec. e.g. the formula is y ~ arm + us(visit|id), then visit should not be used in emmeans(fit, ~ visit)
2. however in previous versions we allow for that (to allow users to see the ls means for each visit, although they are the same)
3. for spatial covairance structure, the visit is a numeric value. the coefficient for that value is NA. in addition, X is not full rank (because the numerical visit is not part of the fixed effect) and there is no linear combinations of x to create the L such that $L\beta$ is estimable

current implementation will discard the visit if the covariance structure is spatial, but I think it will be a better choice to gradually deprecate the use of emmeans(fit, ~ visit) if visit is not part of a fixed effect.

[hqian-lilly]

Hi Daniel, Not sure how this works. Does this fix mean that in the updated mmrm package, we will be able to derive the LSmeans values using emmeans? Another question, in SAS proc mixed, there is an option “nobound” for the variance to allow negative variance, is this implemented in mmrm?

Also I am not completely follow Liming’s initial comments, since in my model, I did have visit as a fixed effect in the model, and generated a second numeric variable (time) to model the variace. I noticed that the variance matrix results is displayed a little bit differently when using sp_exp as an option, which might contribute to the NA calculating the LSmeans values?

I like the mmrm package and really appreciate your effort. Please let me know if I can do anything to help.

Best, Hui-Rong

From: Daniel Sabanes Bove @.> Date: Monday, July 15, 2024 at 7:09 AM To: openpharma/mmrm @.> Cc: Hui-Rong Qian @.>, Author @.> Subject: [EXTERNAL] Re: [openpharma/mmrm] LS estimates using spatial covariance (Issue #446) EXTERNAL EMAIL: Use caution before replying, clicking links, and opening attachments.

Closed #446https://urldefense.com/v3/__https:/github.com/openpharma/mmrm/issues/446__;!!A_kSX19cfVk!BAh3MCJqkc43A94EUa38BG0Uw31DL1oA0MVeipjpMKWmaUarG2nYT4n3poTm3gVU_0G5DpW_HHWBfuQGmZ_rEw$ as completed via #450https://urldefense.com/v3/__https:/github.com/openpharma/mmrm/pull/450__;!!A_kSX19cfVk!BAh3MCJqkc43A94EUa38BG0Uw31DL1oA0MVeipjpMKWmaUarG2nYT4n3poTm3gVU_0G5DpW_HHWBfuT3PFDn1Q$.

— Reply to this email directly, view it on GitHubhttps://urldefense.com/v3/__https:/github.com/openpharma/mmrm/issues/446*event-13503872454__;Iw!!A_kSX19cfVk!BAh3MCJqkc43A94EUa38BG0Uw31DL1oA0MVeipjpMKWmaUarG2nYT4n3poTm3gVU_0G5DpW_HHWBfuT8qlarpA$, or unsubscribehttps://urldefense.com/v3/__https:/github.com/notifications/unsubscribe-auth/AFTLDN2UZB22RZ6BDQPH2O3ZMOUXNAVCNFSM6AAAAABJUVRMACVHI2DSMVQWIX3LMV45UABCJFZXG5LFIV3GK3TUJZXXI2LGNFRWC5DJN5XDWMJTGUYDGOBXGI2DKNA__;!!A_kSX19cfVk!BAh3MCJqkc43A94EUa38BG0Uw31DL1oA0MVeipjpMKWmaUarG2nYT4n3poTm3gVU_0G5DpW_HHWBfuQPdT1qaA$. You are receiving this because you authored the thread.Message ID: @.***>

[danielinteractive]

Thanks @hqian-lilly , @clarkliming would you like to have a look?

[clarkliming]

Not sure how this works. Does this fix mean that in the updated mmrm package, we will be able to derive the LSmeans values using emmeans?

Yes the error is fixed. you are able to derive the ls means with emmeans and see update https://github.com/openpharma/mmrm/pull/450/files (there is also a test)

Another question, in SAS proc mixed, there is an option “nobound” for the variance to allow negative variance, is this implemented in mmrm?

I believe this is totally another story, to allow negative variance. The current parameterization does not allow for that.

Also I am not completely follow Liming’s initial comments, since in my model, I did have visit as a fixed effect in the model, and generated a second numeric variable (time) to model the variace. I noticed that the variance matrix results is displayed a little bit differently when using sp_exp as an option, which might contribute to the NA calculating the LSmeans values?

here for spatial covariance structure, I was using "visit" to refer to all the timing variables, that includes the "visit" you used as fixed effect, and also the "numeric variable (time)". The results should be displayed differently because spatial exponential only calculates the distance between each visit (under my terminology) so we only provide the correlation matrix where the distance is 1. This is not related to the bug you mentioned, and this bug is resolved in the latest mmrm package. Please try it out on github version and let me know if everything is working for you!

Thank you again for spotting this issue and feedback!