class = sdme

[ASKurz]

I'm refitting a model with measurement error on the predictor and the response, using the mi() and me() syntax (fit b14.2_mi from here). After several warnings about large Rhat values, discovered the cause was with the sdme_memar_obs parameter. Based on my read of the b14.2_mi$model output, this is a latent SD.

// generated with brms 2.17.0
functions {
}
data {
  int<lower=1> N;  // total number of observations
  vector[N] Y;  // response variable
  // data for measurement-error in the response
  vector<lower=0>[N] noise;
  // information about non-missings
  int<lower=0> Nme;
  int<lower=1> Jme[Nme];
  int<lower=1> K;  // number of population-level effects
  matrix[N, K] X;  // population-level design matrix
  int<lower=1> Ksp;  // number of special effects terms
  // data for noise-free variables
  int<lower=1> Mme_1;  // number of groups
  vector[N] Xn_1;  // noisy values
  vector<lower=0>[N] noise_1;  // measurement noise
  int prior_only;  // should the likelihood be ignored?
}
transformed data {
}
parameters {
  vector[N] Yl;  // latent variable
  vector[K] b;  // population-level effects
  vector[Ksp] bsp;  // special effects coefficients
  real<lower=0> sigma;  // dispersion parameter
  // parameters for noise free variables
  vector[Mme_1] meanme_1;  // latent means
  vector<lower=0>[Mme_1] sdme_1;  // latent SDs
  vector[N] zme_1;  // standardized latent values
}
transformed parameters {
  vector[N] Xme_1;  // actual latent values
  real lprior = 0;  // prior contributions to the log posterior
  // compute actual latent values
  Xme_1 = meanme_1[1] + sdme_1[1] * zme_1;
  lprior += normal_lpdf(b | 0, 10);
  lprior += normal_lpdf(bsp | 0, 10);
  lprior += cauchy_lpdf(sigma | 0, 2.5)
    - 1 * cauchy_lccdf(0 | 0, 2.5);
  lprior += cauchy_lpdf(sdme_1 | 0, 2.5)
    - 1 * cauchy_lccdf(0 | 0, 2.5);
}
model {
  // likelihood including constants
  if (!prior_only) {
    // initialize linear predictor term
    vector[N] mu = X * b;
    for (n in 1:N) {
      // add more terms to the linear predictor
      mu[n] += (bsp[1]) * Xme_1[n];
    }
    target += normal_lpdf(Yl | mu, sigma);
  }
  // priors including constants
  target += lprior;
  target += normal_lpdf(Y[Jme] | Yl[Jme], noise[Jme]);
  target += normal_lpdf(Xn_1 | Xme_1, noise_1);
  target += std_normal_lpdf(zme_1);
}
generated quantities {
}

But when I look at the default priors for the latent SD, they're flat and unbounded.

b14.2_mi$prior

          prior  class            coef group resp dpar nlpar lb ub       source
  normal(0, 10)      b                                                     user
  normal(0, 10)      b               A                             (vectorized)
  normal(0, 10)      b       Intercept                             (vectorized)
  normal(0, 10)      b memar_obsmar_sd                             (vectorized)
         (flat) meanme                                                  default
         (flat) meanme       memar_obs                             (vectorized)
         (flat)   sdme                                                  default
         (flat)   sdme       memar_obs                             (vectorized)
 cauchy(0, 2.5)  sigma                                        0            user

Is this a mistake or am I misunderstanding something? Wouldn't it make more sense to at least default to lb = 0 for a latent SD parameter?

Regardless, it would be nice if there was documentation on priors of class = sdme in the brms reference manual. I didn't see any in the current version (2.17.0).

[paul-buerkner]

Yeah, it should have updated the boundary on the sdme parameter as well. Did you updat the old model (via update) or did you run the model from scratch with the latest brms version? In any case I would need a reproducible example, because when I just check the generated Stan code of my me example, all looks good to me.

I will add the infor about meanme and sdme to the doc.

[paul-buerkner]

Nevermind, I think I would the issue. Will fix it asap.

[ASKurz]

You don't need the example? [I'm half-way through preparing it.]

[paul-buerkner]

No I don't. Sorry for the confusion.

[ASKurz]

No problem. In addition to lb = 0, will your update include a default half-cauchy distribution for priors of class = sdme, or will you keep the prior flat?

[paul-buerkner]

I will keep them flat for now because that is another issue. I may replace me by mi terms completely in brms 3.0 and then will have to tackle the whole prior specification issue of these terms in detail anyway.

[ASKurz]

Good deal. I get the impression the 3.0 update will be a big one.

[paul-buerkner]

The bug should now be fixed. Can you double check?

[ASKurz]

Yep, the version 2.17.3 update now works on my end. Thanks, Paul!