Error reported when summary(fit)
Closed this issue · 4 comments
Dear mmrm developers,
I have been heavily use this package recently and really appreciate your contributions towards this package!
Here is an error reported when I called for coefficients table using summary() after fitting the model. Is there any solutions towards this issue? Many thanks!
Error message:
> summary(fit)
Error in h_df_1d_list(est = est, var = var, v_num = v_num, v_denom = v_denom) :
Assertion on 'v_denom' failed: Element 1 is not >= 2.22507e-308.
Outputting directly the fitted model:
mmrm fit
Formula: change ~ base + treatArm * visitNum + us(visitNum | patientID)
Data: cut_info (used 377 observations from 99 subjects with maximum 7 timepoints)
Covariance: unstructured (28 variance parameters)
Method: REML
Deviance: -1827.339
Coefficients:
(Intercept) base treatArm2 visitNumWEEK 16 visitNumWEEK 20
0.086617083 -0.018586266 0.026598894 0.006935757 0.010958797
visitNumWEEK 24 visitNumWEEK 28 visitNumWEEK 4 visitNumWEEK 8 treatArm2:visitNumWEEK 16
0.005375097 -0.063040454 0.006567467 0.009339331 -0.008895166
treatArm2:visitNumWEEK 20 treatArm2:visitNumWEEK 24 treatArm2:visitNumWEEK 28 treatArm2:visitNumWEEK 4 treatArm2:visitNumWEEK 8
-0.019292566 -0.013914643 0.159960538 -0.012124116 -0.012901674
Model Inference Optimization:
Converged with code 0 and message: convergence: rel_reduction_of_f <= factr*epsmch
R session info
R version 4.0.3 (2020-10-10)
Platform: x86_64-pc-linux-gnu (64-bit)
Running under: Red Hat Enterprise Linux
Matrix products: default
BLAS/LAPACK: /usr/lib64/libopenblas-r0.2.20.so
locale:
[1] LC_CTYPE=en_US.UTF-8 LC_NUMERIC=C LC_TIME=en_US.UTF-8 LC_COLLATE=en_US.UTF-8
[5] LC_MONETARY=en_US.UTF-8 LC_MESSAGES=en_US.UTF-8 LC_PAPER=en_US.UTF-8 LC_NAME=C
[9] LC_ADDRESS=C LC_TELEPHONE=C LC_MEASUREMENT=en_US.UTF-8 LC_IDENTIFICATION=C
attached base packages:
[1] stats graphics grDevices utils datasets methods base
other attached packages:
[1] ggpubr_0.4.0 forcats_0.5.0 stringr_1.4.0 dplyr_1.0.2 purrr_0.3.4 readr_1.4.0 tidyr_1.1.2 tibble_3.0.4
[9] ggplot2_3.3.2 tidyverse_1.3.0
loaded via a namespace (and not attached):
[1] nlme_3.1-149 fs_1.5.0 lubridate_1.7.9 httr_1.4.2 numDeriv_2016.8-1.1 TMB_1.9.4
[7] tools_4.0.3 backports_1.1.10 R6_2.5.0 DBI_1.1.0 colorspace_1.4-1 withr_2.3.0
[13] tidyselect_1.1.0 emmeans_1.5.2-1 curl_4.3 compiler_4.0.3 cli_2.1.0 rvest_0.3.6
[19] xml2_1.3.2 sandwich_3.0-0 labeling_0.4.2 scales_1.1.1 checkmate_2.0.0 mvtnorm_1.1-1
[25] digest_0.6.27 foreign_0.8-80 rmarkdown_2.5 mmrm_0.2.2 rio_0.5.16 pkgconfig_2.0.3
[31] htmltools_0.5.0 dbplyr_1.4.4 rlang_0.4.8 readxl_1.3.1 rstudioapi_0.11 farver_2.0.3
[37] generics_0.1.0 zoo_1.8-8 jsonlite_1.7.1 zip_2.1.1 car_3.0-10 magrittr_1.5
[43] Matrix_1.2-18 Rcpp_1.0.5 munsell_0.5.0 fansi_0.4.1 abind_1.4-5 lifecycle_0.2.0
[49] stringi_1.5.3 multcomp_1.4-14 yaml_2.2.1 carData_3.0-4 gbRd_0.4-11 MASS_7.3-53
[55] grid_4.0.3 blob_1.2.1 parallel_4.0.3 crayon_1.3.4 lattice_0.20-41 cowplot_1.1.0
[61] haven_2.3.1 splines_4.0.3 hms_0.5.3 knitr_1.30 pillar_1.4.6 estimability_1.3
[67] ggsignif_0.6.0 codetools_0.2-16 reprex_0.3.0 glue_1.4.2 evaluate_0.14 data.table_1.13.2
[73] modelr_0.1.8 vctrs_0.3.4 Rdpack_2.0 cellranger_1.1.0 gtable_0.3.0 assertthat_0.2.1
[79] xfun_0.19 openxlsx_4.2.3 rbibutils_1.3 xtable_1.8-4 broom_0.7.2 coda_0.19-4
[85] rstatix_0.6.0 survival_3.2-7 TH.data_1.0-10 ellipsis_0.3.1
Thanks @Dominic7227 , could you add a reprex for this?
As a workaround you could try to use the development version of the package and use another d.f. adjustment method (because the error is specific for Satterthwaite)
Sure @danielinteractive, here is the demo data and the lines I used.
fit <- mmrm::mmrm(
formula = change ~ base + treatArm*visitNum + us(visitNum | patientID),
data = cut_info,
control = mmrm::mmrm_control(
optimizer = c( "BFGS", "CG", "SANN", "Brent", "Nelder-Mead", "L-BFGS-B","nlminb"),
n_cores=8
)
)
summary(fit)
Thanks @Dominic7227 !
I looked at the example. It is interesting that on my machine and with the latest development version of the package the model does not fit successfully in the first place, which already points to an overparametrized covariance structure. For example, if I switch it to adh, toeph or ar1h they all work.
The specific error you got is also a symptom of too many parameters in the covariance structure. Apologies that the error message is very cryptic, but I think this is the first edge case we have seen this assertion being triggered.
You can try to add method = "Residual" after installing the development version of the package (e.g. from R-Universe most conveniently, see https://openpharma.github.io/mmrm/main/#development, and then see what you get as summary output - probably you'll see which covariance entries cannot be estimated e.g.
closing as stale