Error only occurring on Amazon linux instance

Posted on Tue, 11/16/2021 - 14:27

ssciarra Joined: 01/08/2021

Forums

I am trying to run simulations on an Amazon linux instance and get the following error when I use `mxAutoStart()` to generate starting values :

library(devtools)

#install_github("sciarraseb/nonlinSims", dependencies = T, force=T)

library(easypackages)

library(nonlinSims)

library(parallel)

library(tidyverse)

library(OpenMx)

library(data.table)

test_data <- read_csv(file = 'test_data.csv') model <- create_logistic_growth_model(data_wide = test_data, model_name = 'test') model <- mxAutoStart(model) Error in solve.default(I - A) : system is computationally singular: reciprocal condition number = 0

Here is the output provided when calling `traceback()`:
12: solve.default(I - A) 11: solve(I - A) 10: genericGetExpected(model[[subname]]$expectation, model, component, defvar.row, subname) 9: genericGetExpected(model[[subname]]$expectation, model, component, defvar.row, subname) 8: mxGetExpected(model, c("covariance", "means", "thresholds"), subname = subname) 7: autoStartDataHelper(model, type = type) 6: mxModel(model, autoStartDataHelper(model, type = type)) 5: omxBuildAutoStartModel(model, type) 4: is(model, "MxModel") 3: warnModelCreatedByOldVersion(model) 2: mxRun(omxBuildAutoStartModel(model, type), silent = TRUE) 1: mxAutoStart(model)

I have provided other pertinent information below (i.e., R version, optimizer used in OpenMx, etc.):
OpenMx version: 2.19.8 [GIT v2.19.8] R version: R version 4.0.2 (2020-06-22) Platform: x86_64-koji-linux-gnu Default optimizer: SLSQP NPSOL-enabled?: No OpenMP-enabled?: Yes

Interestingly, the error appearing on the Amazon instance does not appear when I run the code offline in RStudio (on either Mac or Windows). Here is the output provided by `mxVersion()`in the offline version of R that I am using:
OpenMx version: 2.19.8 [GIT v2.19.8] R version: R version 4.0.5 (2021-03-31) Platform: x86_64-apple-darwin17.0 MacOS: 12.0.1 Default optimizer: SLSQP NPSOL-enabled?: No OpenMP-enabled?: No

I have also attached the data set (`test_data.csv`). I have written `create_logistic_growth_model()` and posted it in a GitHub repository (hopefully you can download the package to use the function). I can provide more information about this function if necessary.

Replied on Tue, 11/16/2021 - 15:47

AdminRobK Joined: 01/24/2014

sessionInfo()

Under Amazon Linux, what do you get from running `sessionInfo()` at the R prompt? I'm primarily curious about the 'Matrix products', 'BLAS', and 'LAPACK'.

Log in or register to post comments

Replied on Tue, 11/16/2021 - 16:53

ssciarra Joined: 01/08/2021

Here is the output from

Here is the output from `SessionInfo()`:

R version 4.0.2 (2020-06-22)

Platform: x86_64-koji-linux-gnu (64-bit)

Running under: Amazon Linux 2
Matrix products: default

BLAS/LAPACK: /usr/lib64/R/lib/libRblas.so
locale:

 [1] LC_CTYPE=en_US.UTF-8       LC_NUMERIC=C               LC_TIME=en_US.UTF-8

 [4] 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

[10] LC_TELEPHONE=C             LC_MEASUREMENT=en_US.UTF-8 LC_IDENTIFICATION=C       
attached base packages:

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

 [1] forcats_0.5.1         stringr_1.4.0         dplyr_1.0.7           purrr_0.3.4

 [5] readr_2.1.0           tidyr_1.1.4           tibble_3.1.6          ggplot2_3.3.5

 [9] tidyverse_1.3.1       nonlinSims_0.0.0.9000 devtools_2.4.2        usethis_2.1.3

[13] data.table_1.14.2     OpenMx_2.19.8

loaded via a namespace (and not attached): [1] Rcpp_1.0.7 lubridate_1.8.0 lattice_0.20-41 prettyunits_1.1.1 [5] ps_1.6.0 assertthat_0.2.1 rprojroot_2.0.2 digest_0.6.28 [9] utf8_1.2.2 cellranger_1.1.0 R6_2.5.1 backports_1.3.0 [13] reprex_2.0.1 httr_1.4.2 pillar_1.6.4 rlang_0.4.12 [17] readxl_1.3.1 rstudioapi_0.13 callr_3.7.0 Matrix_1.2-18 [21] desc_1.4.0 bit_4.0.4 munsell_0.5.0 broom_0.7.10 [25] compiler_4.0.2 modelr_0.1.8 pkgconfig_2.0.3 pkgbuild_1.2.0 [29] tidyselect_1.1.1 fansi_0.5.0 crayon_1.4.2 tzdb_0.2.0 [33] dbplyr_2.1.1 withr_2.4.2 MASS_7.3-51.6 grid_4.0.2 [37] jsonlite_1.7.2 gtable_0.3.0 lifecycle_1.0.1 DBI_1.1.1 [41] magrittr_2.0.1 scales_1.1.1 RcppParallel_5.1.4 vroom_1.5.6 [45] stringi_1.7.5 cli_3.1.0 cachem_1.0.6 fs_1.5.0 [49] remotes_2.4.1 testthat_3.1.0 xml2_1.3.2 ellipsis_0.3.2 [53] generics_0.1.1 vctrs_0.3.8 tools_4.0.2 bit64_4.0.5 [57] glue_1.5.0 hms_1.1.1 processx_3.5.2 pkgload_1.2.3 [61] fastmap_1.1.0 colorspace_2.0-2 sessioninfo_1.2.1 rvest_1.0.2 [65] memoise_2.0.0 haven_2.4.3

Log in or register to post comments

Replied on Tue, 11/16/2021 - 16:27

AdminRobK Joined: 01/24/2014

works for me; NAs

At any rate, I do not reproduce the error on my system. My `mxVersion()` output:



OpenMx version: 2.19.8.59 [GIT v2.19.8-59-g5f6ffe95-dirty]

R version: R version 4.0.2 (2020-06-22)

Platform: x86_64-pc-linux-gnu

Default optimizer: SLSQP

NPSOL-enabled?: Yes

OpenMP-enabled?: Yes

I am using the default BLAS and LAPACK implementations.

As a workaround, try adjusting the initial values your function uses for some of your MxPaths and MxMatrix elements.

Log in or register to post comments

Replied on Thu, 11/18/2021 - 19:41

ssciarra Joined: 01/08/2021

Starting value experimentation + matrix products

I have experimented starting values but had no success (unfortunately). I also posted the output from `SessionInfo()` above for you to look at the matrix products.

Log in or register to post comments

Replied on Fri, 11/19/2021 - 13:19

AdminRobK Joined: 01/24/2014

linear algebra; debugger

I also posted the output from `SessionInfo()` above for you to look at the matrix products.

Yes, I saw. My suspicion was that your instance of Amazon linux was building R with non-default implementations of BLAS/LAPACK, which might have a different tolerance for when a matrix is considered "computationally singular". Evidently, that's not what's going on.

Can you make any headway under Amazon Linux if you don't use `mxAutoStart()`?

You might need to run your script under R's debugger on an Amazon Linux instance.

Log in or register to post comments

Replied on Sat, 11/20/2021 - 00:34

ssciarra Joined: 01/08/2021

Debugging attempt

I did my best to understand where the error may be occurring (using the debugging link you posted). It appears that the error is a result of `NaN`s appearing in the A matrix (which contains the factor loadings from the latent variables to the manifest variables; see output below). The outputs from environments 9–10 are shown below (I can provide more output if necessary).

Error in solve.default(I - A) :

  system is computationally singular: reciprocal condition number = 0

Available environments had calls:

1: mxAutoStart(model)

2: mxRun(omxBuildAutoStartModel(model, type), silent = TRUE)

3: warnModelCreatedByOldVersion(model)

4: is(model, "MxModel")

5: omxBuildAutoStartModel(model, type)

6: mxModel(model, autoStartDataHelper(model, type = type))

7: autoStartDataHelper(model, type = type)

8: mxGetExpected(model, c("covariance", "means", "thresholds"), subname = subname)

9: genericGetExpected(model[[subname]]$expectation, model, component, defvar.row, subname)

10: genericGetExpected(model[[subname]]$expectation, model, component, defvar.row, subname)

11: solve(I - A)

12: solve.default(I - A)
Enter an environment number, or 0 to exit

Selection: 10

Browsing in the environment with call:

   genericGetExpected(model[[subname]]$expectation, model, component, defvar.row, subname)

Called from: debugger.look(ind)

Browse[1]> A

       t1_1 t2_45 t3_90 t4_135 t5_180 t6_225 t7_270 t8_315 t9_360 diff beta gamma

t1_1      0     0     0      0      0      0      0      0      0  NaN  NaN   NaN

t2_45     0     0     0      0      0      0      0      0      0    1  NaN   NaN

t3_90     0     0     0      0      0      0      0      0      0    1  NaN   NaN

t4_135    0     0     0      0      0      0      0      0      0    1  NaN   NaN

t5_180    0     0     0      0      0      0      0      0      0    1  NaN   NaN

t6_225    0     0     0      0      0      0      0      0      0    1  NaN   NaN

t7_270    0     0     0      0      0      0      0      0      0    1  NaN   NaN

t8_315    0     0     0      0      0      0      0      0      0    1  NaN   NaN

t9_360    0     0     0      0      0      0      0      0      0    1  NaN   NaN

diff      0     0     0      0      0      0      0      0      0    0    0     0

beta      0     0     0      0      0      0      0      0      0    0    0     0

gamma     0     0     0      0      0      0      0      0      0    0    0     0
Enter an environment number, or 0 to exit

Selection: 9

Browsing in the environment with call:

   genericGetExpected(model[[subname]]$expectation, model, component, defvar.row, subname)

Called from: debugger.look(ind)

Browse[1]> model$A

FullMatrix 'A' 
$labels

       t1_1 t2_45 t3_90 t4_135 t5_180 t6_225 t7_270 t8_315 t9_360 diff      beta      gamma

t1_1   NA   NA    NA    NA     NA     NA     NA     NA     NA     "Dl[1,1]" "Bl[1,1]" "Gl[1,1]"

t2_45  NA   NA    NA    NA     NA     NA     NA     NA     NA     "Dl[2,1]" "Bl[2,1]" "Gl[2,1]"

t3_90  NA   NA    NA    NA     NA     NA     NA     NA     NA     "Dl[3,1]" "Bl[3,1]" "Gl[3,1]"

t4_135 NA   NA    NA    NA     NA     NA     NA     NA     NA     "Dl[4,1]" "Bl[4,1]" "Gl[4,1]"

t5_180 NA   NA    NA    NA     NA     NA     NA     NA     NA     "Dl[5,1]" "Bl[5,1]" "Gl[5,1]"

t6_225 NA   NA    NA    NA     NA     NA     NA     NA     NA     "Dl[6,1]" "Bl[6,1]" "Gl[6,1]"

t7_270 NA   NA    NA    NA     NA     NA     NA     NA     NA     "Dl[7,1]" "Bl[7,1]" "Gl[7,1]"

t8_315 NA   NA    NA    NA     NA     NA     NA     NA     NA     "Dl[8,1]" "Bl[8,1]" "Gl[8,1]"

t9_360 NA   NA    NA    NA     NA     NA     NA     NA     NA     "Dl[9,1]" "Bl[9,1]" "Gl[9,1]"

diff   NA   NA    NA    NA     NA     NA     NA     NA     NA     NA        NA        NA

beta   NA   NA    NA    NA     NA     NA     NA     NA     NA     NA        NA        NA

gamma  NA   NA    NA    NA     NA     NA     NA     NA     NA     NA        NA        NA       
$values

       t1_1 t2_45 t3_90 t4_135 t5_180 t6_225 t7_270 t8_315 t9_360 diff beta gamma

t1_1      0     0     0      0      0      0      0      0      0    0    0     0

t2_45     0     0     0      0      0      0      0      0      0    0    0     0

t3_90     0     0     0      0      0      0      0      0      0    0    0     0

t4_135    0     0     0      0      0      0      0      0      0    0    0     0

t5_180    0     0     0      0      0      0      0      0      0    0    0     0

t6_225    0     0     0      0      0      0      0      0      0    0    0     0

t7_270    0     0     0      0      0      0      0      0      0    0    0     0

t8_315    0     0     0      0      0      0      0      0      0    0    0     0

t9_360    0     0     0      0      0      0      0      0      0    0    0     0

diff      0     0     0      0      0      0      0      0      0    0    0     0

beta      0     0     0      0      0      0      0      0      0    0    0     0

gamma     0     0     0      0      0      0      0      0      0    0    0     0
$free: No free parameters.
$lbound: No lower bounds assigned.

$ubound: No upper bounds assigned.

At this point, it may be important to mention the model I am running. In short, I am running a three-parameter logistic curve model (with the growth parameters treated as latent variables). I have provided a relevant code chunk (from `create_logistic_growth_model()` and attached a path diagram showing the specific model that I am using. Seven parameters are estimated: fixed- and random-effects parameters for $diff$, $\beta$, and $\gamma$ and one residual component for each time point ($\epsilon$). Note that a structured latent curve model is implemented.
mxModel(model = model_name, type = 'RAM', independent = T, mxData(observed = data_wide, type = 'raw'),


                   manifestVars = manifest_vars,

                   latentVars = latent_vars,
                   #Residual variances; by using one label, they are assumed to all be equal

                   #(homogeneity of variance)

                   mxPath(from = manifest_vars,

                          arrows=2, free=TRUE, values = 1, labels='epsilon',  lbound = 1e-3),
                   #Set manifest means to observed means

                   mxPath(from = 'one', to = manifest_vars, free = FALSE, arrows = 1, values = manifest_means),
                   #Latent variable covariances and variances

                   mxPath(from = latent_vars,

                          connect='unique.pairs', arrows=2,

                          #aa(diff_rand),

                          #ab(cov_diff_beta), ac(cov_diff_gamma),

                          #bb(beta_rand), bc(var_beta_gamma), cc(gamma_rand)

                          free = c(TRUE,

                                   FALSE, FALSE,

                                   TRUE, FALSE, TRUE),

                          values=c(1,

                                   NA,NA,

                                   1, NA, 1),

                          labels=c('diff_rand',

                                   'NA(cov_diff_beta)','NA(cov_diff_gamma)',

                                   'beta_rand', 'NA(var_beta_gamma)', 'gamma_rand'),

                          lbound = c(1e-3,

                                     NA, NA,

                                     2, NA, 1)), 
                   #Latent variable means (linear parameters). Note that the nonlinear parameters of

                   #beta and gamma do not have estimated means

                   mxPath(from = 'one', to = 'diff', free = TRUE, arrows = 1, values = 1,

                          labels = 'diff_fixed', lbound = 0.1, ubound = 2),
                   #Functional constraints

                   mxMatrix(type = 'Full', nrow = length(manifest_vars), ncol = 1, free = TRUE,

                            values = 1, labels = 'diff_fixed', name = 'd', lbound = 0.1, ubound = 2),
                   mxMatrix(type = 'Full', nrow = length(manifest_vars), ncol = 1, free = TRUE,

                            values = 1, labels = 'beta_fixed', name = 'b', lbound = 0.1, ubound = 360),
                   mxMatrix(type = 'Full', nrow = length(manifest_vars), ncol = 1, free = TRUE,

                            labels = 'gamma_fixed', name = 'g', lbound = 0.1, ubound = 360),
                   mxMatrix(type = 'Full', nrow = length(manifest_vars), ncol = 1, free = FALSE,

                            values = measurement_days, name = 'time'),
                   #Algebra specifying first partial derivatives

                   mxAlgebra(expression = 1/(1 + exp((b - time)/g)), name="Dl"),

                   mxAlgebra(expression =  -(d * (exp((b - time)/g) * (1/g))/(1 + exp((b - time)/g))^2),

                   name = 'Bl'),

                   mxAlgebra(expression =  d * (exp((b - time)/g) * ((b - time)/g^2))/(1 + exp((b - time)/g))^2,

                   name = 'Gl'),
                   #Factor loadings; all fixed and, importantly, constrained to change according to their

                   #partial derivatives (i.e., nonlinear functions). Note that ampersands precede % symbols

                   #in post to prevent exiting in HTML code tags (source code does not actually include ampersands).

                   mxPath(from = 'diff', to = manifest_vars, arrows=1, free=FALSE,

                          labels = sprintf(fmt = 'Dl[&%d,1]', 1:(ncol(data_wide)-1))),

                   mxPath(from='beta', to = manifest_vars, arrows=1,  free=FALSE,

                          labels =  sprintf(fmt = 'Bl[&%37;d,1]', 1:(ncol(data_wide)-1))),

                   mxPath(from='gamma', to = manifest_vars, arrows=1,  free=FALSE,

                          labels =  sprintf(fmt = 'Gl[&%d,1]', 1:(ncol(data_wide)-1))),
                   mxFitFunctionML(vector = FALSE)
  )

return(model) }

File attachments

path_diagram.png

Log in or register to post comments

Log in or register to post comments

News

Recent Posts

Error only occurring on Amazon linux instance

sessionInfo()

Here is the output from

works for me; NAs

Starting value experimentation + matrix products

linear algebra; debugger

Debugging attempt

News

Recent Posts