Sampling not done because 'y is not positive definite'

[Jeroen CORSTJENS]

Hi,

I'm fitting a non-linear mixed effects model with 25 random effects and 5 group-level predictors. 12 of these random effects are interactions, so x and z are random effects as well as the interaction xz. I do not specify any priors myself since I have no knowledge about the probability distributions of the various variables.

I have a data set of 4000 observations (200 groups with 20 observations per group). All continuous predictor variables were mean centred. The reciprocal transformation of the response variable (multiplied by large constant) was taken to ameliorate heteroscedasticity in the model.

When fitting the model using the R package 'brms' (which relies on rstan) I get the following error:

SAMPLING FOR MODEL 'gaussian(identity) brms-model' NOW (CHAIN 1).

[1] "Error in sampler$call_sampler(args_list[[i]]) : "                                        

[2] "  model14a86bff6439_file14a877266e59_namespace::write_array: y is not positive definite."

[1] "error occurred during calling the sampler; sampling not done"

When using 'rstanarm' (with QR=FALSE and QR=TRUE) I get the following error

fixed-effect model matrix is rank deficient so dropping 128 columns / coefficients

Both errors suggests possible collinearity issues, while the latter error suggests the presence of collinearity in my fixed effects design matrix, meaning some columns can be expressed as linear combinations of other columns. However, trying to identify these colums with the function findLinearCombos(Fixed Effects Matrix) (package caret) did not show any linear dependencies. Executing the same function on the random effects design matrix Zt did show many linear dependencies and suggests me to remove, for example, "b[Greedy problem:1]"  and "b[Regret2 problem:1]" from the model.

When I remove 3 random effects, sampling begins, but the model does not converge, with Rhats > 1.1. Doubling the number of observations (200 groups with 40 observations per group) does not help either. I have tried fitting a model that considers the interactions effects as fixed effects instead of random effects and this works fine, so the issue lies with considering these interactions as random effects.

My question: Is it simply impossible given the data set I use to find a stable fit for a model that considers (this many) 2-way interactions as random effects? Is there anything else I could try to resolve the 'not positive definite' issue that prevents sampling to start.

p.s. I also tried to fit the fit model using rstan directly. However, the stancode provided by the brm-model gave the following error in rstan

Error in basename(fname) : path too long

In addition: Warning message:

In normalizePath(path.expand(path), winslash, mustWork) :

  path[1]="// generated with brms 1.5.0

functions { 

} 

...

 vector[N] Z_1_25; 

  int<lower=1> NC_1; 

  int prior_only;  // should the likeliho [... truncated]

Is there a maximum code length for rstan?

[Andrew Gelman]

Hi, just a quick comment which is that you could clean up your code a huge amount by using arrays:  for example:

vector[N] Z_1[25];

and so forth.

Also use unit scales.  It can't be right to have priors with scale 2776.

Finally, it should help to put priors on the parameters b.

A

--
You received this message because you are subscribed to the Google Groups "Stan users mailing list" group.
To unsubscribe from this group and stop receiving emails from it, send an email to stan-users+...@googlegroups.com.
To post to this group, send email to stan-...@googlegroups.com.
For more options, visit https://groups.google.com/d/optout.
<Stan_code_full_model.stan>