Multivariate lognormal distribution in Stan

General
1

Hi,

Currently, I am trying to fit a lognormal model for one of the parameters in the Stan program. For one-dimension case, the easy way is just using the built-in Stan function as:

\beta \sim lognormal(\mu, \sigma)

However, currently my parameter \beta is a vector of length 3, where \beta = (\beta_1, \beta_2, \beta_3) and I would like to jointly model the three components as multivariate lognormal distribution with some mean vector and covariance matrix (maybe treat them as hyperparameter in the model). May I ask that does Stan provide such built-in functions for multivariate lognormal distribution?

Thanks so much for your help!

2

The easiest way to do this would be to take the logarithm of your data in the transformed data block. Tutorial here.

3

Usually when we include multivariate random effects in a model, we are modeling covariances between multiple vectors, not multiple elements of a single vector. If \beta_1, \beta_2, and \beta_3 are all elements of a single random-effect vector, then there is just the variance and no additional covariance. If, on the other hand, \beta_1, \beta_2, and \beta_3 all represent single-element random effect vectors, then how do you expect to be able to estimate the underlying means, variances, and covariances? If the means, variances, and covariances are tightly constrained a priori then this could be a useful approach, but that’s not what you seem to be driving at when you write

4

We don’t have a multivariate lognormal built in, but it’s easy to code by inverse transforming a multivariate normal.

data {
  cov_matrix[3] Sigma;
}
transformed data {
  cholesky_factor_cov[3] L_Sigma = cholesky_decompose(Sigma);
}
parameters {
  vector[3] log_beta;  
}
transformed parameters {
  vector[3] beta = exp(log_beta);
}
model {
  log_beta ~ multi_normal_cholesky(mu, L_Sigma);
}

This implies that beta has a multi_lognormal(mu, Sigma) distribution, where Sigma = L_Sigma * L_Sigma'. With this code, you can use beta just as you would had it been given a multivariate lognormal directly.

Multivariate logistic normal works the same way, only with inv_logit() in place of exp().

[edit: added he covariance matrix as data and its cholesky factor as a transform—you could also take this to be a parameter, for which you might consider the LKJ prior on a cholesky factor of correlation matrices scaled by a vector of scales with their own independent prior]

4 Likes
5

Hi Bob,

Thanks so much for your help!