mmrm 0.3.13's first fit is different than subsequent ones

[luwidmer]

Summary
mmrm 0.3.13 gives a different result on the first call than subsequent ones. Note that with TMB 1.9.15 and mmrm 0.3.12 this does not happen (the result is the same every time).

I've been able to narrow it down to the TMB::config call in the warning about optimize.instantly (https://github.com/openpharma/mmrm/blob/main/R/utils.R#L536 - also, if I manually set optimize.instantly to 1, the warning does not trigger?), so this seems to be a side effect of calling TMB::config to get the configuration even without changing the config? Reproducible example that allows for switching mmrm versions easily below:

mmrm_version <- "0.3.13" # "0.3.12" for comparison
mmrm_libpath <- paste0(".libs-", mmrm_version)

if(!dir.exists(mmrm_libpath)){
  dir.create(mmrm_libpath)
}
.libPaths(c(normalizePath(mmrm_libpath), .libPaths()))

if (!dir.exists(paste0(mmrm_libpath, "/mmrm"))) {
  remotes::install_version("mmrm", version = mmrm_version, type="source", repos="https://stat.ethz.ch/CRAN", dependencies="never")
}

library(mmrm)
library(digest)
options(digits = 22)

res <- character()
cat(as.character(packageVersion("mmrm")), "\n")
for (i in 1:5) {
  set.seed(123)
  mmrmFit <- mmrm(
    formula = FEV1 ~ RACE + SEX + ARMCD * AVISIT + us(AVISIT | USUBJID),
    data = fev_data
  )
  res[i] <- digest(coef(mmrmFit))
  cat("Res ", i, ": ", res[i], "\n")
}

Output for mmrm 0.3.13:

> source("~/.active-rstudio-document")
0.3.13 
Res  1 :  a005171052e35c971e920321dea782f6 
Res  2 :  903104ec4747c253d773e48e6ab96b32 
Res  3 :  903104ec4747c253d773e48e6ab96b32 
Res  4 :  903104ec4747c253d773e48e6ab96b32 
Res  5 :  903104ec4747c253d773e48e6ab96b32

Output for mmrm 0.3.12:

> source("~/.active-rstudio-document")
0.3.12 
Res  1 :  903104ec4747c253d773e48e6ab96b32 
Res  2 :  903104ec4747c253d773e48e6ab96b32 
Res  3 :  903104ec4747c253d773e48e6ab96b32 
Res  4 :  903104ec4747c253d773e48e6ab96b32 
Res  5 :  903104ec4747c253d773e48e6ab96b32 

OS / Environment
Multiple (R 4.4 on ARM64 / MacOS, R 4.3.1 on RHEL 7)

@danielinteractive @kaskr any ideas?

[danielinteractive]

Oh no 😮...

[kaskr]

@luwidmer Surprisingly, this doesn't seem to be a problem with TMB::config? If I out-comment this line:

rlang::warn(msg, .frequency = "once", .frequency_id = "tmb_warn_optimization")

the results become reproducible. Seems like rlang::warn has side effects?

So, why is this warning triggered at all? That's because your are testing for the wrong flag. Should have been:

tmb_config <- TMB::config(DLL = "mmrm")
if ( ! tmb_config$tmbad_deterministic_hash) { ... ## NOTE: NULL on old versions of TMB }

[danielinteractive]

Thanks @kaskr , wow that is really weird that rlang::warn seems to be the problem. I mean then I can also just go with a usual warning() - it will spam the console of the user but I guess that is less problematic than side effects.
Thanks for the pointer to the right flag, I guess I was a bit too naive there when assuming the option has the same name as the list element.

[luwidmer]

@kaskr super interesting, and sorry for blaming TMB, my bad 😢! I did re-run 0.3.13 through RDvalgrind with level 2 valgrind instrumentation and no uninitialized memory is used (also not by rlang) there as far as this instrumentation can detect. Commenting the warning and re-installing also did the trick for me, but I'm a bit at a loss of why rlang::warn would no longer make things reproducible... that's worrying at least 🤔

[luwidmer]

@danielinteractive is it possible a different optimizer is used when this warning is fired? It doesn't show up for me in the above code, even though it should?! I very much doubt the issue is in rlang::warn

[luwidmer]

@danielinteractive yup I think that might be it, changing to warning() gives the following result:

> source("~/.active-rstudio-document")
0.3.13.9000 
Error in refit_multiple_optimizers(fit = fit, control = control) : 
  No optimizer led to a successful model fit. Please try to use a different covariance structure or other covariates.

[danielinteractive]

@luwidmer ah yeah true, the warning will be captured internally and then another optimizer will be chosen, causing the other result. Alright, then I just did a mistake there putting the warning in the TMB fitting function, I guess we will need to move it higher up into mmrm().

[danielinteractive]

tmbad_deterministic_hash

@kaskr For my understanding, we can still now use the TMB tape optimizer but just with deterministic hash and that will produce reproducible results, is that correct? So forcing optimize.instantly = 0 is probably too much because we just disable the whole TMB tape optimizer, right?

[kaskr]

@danielinteractive Yes that is correct!