Priors to get Positive-definite covariance matrix

Modeling
1

Hello, I need to draw the structural model from a trivariate normal distribution. I have three variables: y, x, z. I know that \sigma_z = 1. This is expressed by indicating scale_new[3] = 1. I also know that the correlation between the error of x and the error of z is 0, i.e., \rho_{23} = 0 . This is defined by cov_mat[2,3] = 0 and cov_mat[3, 2] = 0. The thing is now: I am getting divergent transitions, and most probably because my variance-covariance matrix is not positive definite. What kind of priors do i need to define separately for the correlations (i.e. \rho_{13} and \rho_{12}) in order to get positive-definite matrix? What goes wrong with beta priors for the correlations? I attach the code:

data {
  int<lower=0> N;
  int<lower=1> NX; // dimension of endogenous var//i.e. categories 
  vector[N] y; //outcome vector
  vector[N] x; //endogenous variable
  vector[N] z; // instrument
  int sim_data;
}

transformed data {
  vector[3] Y[N];
  for(i in 1:N) {
    Y[i][1] = y[i];
    Y[i][2] = x[i];
    Y[i][3] = z[i];
  }
}

parameters {
  ordered[NX] beta_x;
  vector[2] beta_null;
  real alpha_z;
  vector<lower=0>[3] scale;
  corr_matrix[3] Omega;
}

transformed parameters {
  vector[3] mu[N];
  matrix[3, 3] cov_mat;
  vector[3] scale_new;
  for(i in 1:N) {
    mu[i][1] = beta_null[1] +  x[i] * beta_x[NX];
    mu[i][2] = beta_null[2] + z[i] * alpha_z;
    mu[i][3] = 0;
  }
  cov_mat = Omega;
  cov_mat[2, 3] = 0;
  cov_mat[3, 2] = 0;
  scale_new = scale;
  scale_new[3] = 1;
}

model {
  
  //priors
  beta_x ~ normal(0,1);
  beta_null[1] ~ normal(13,3); 
  beta_null[2] ~ normal(5,3);
  alpha_z ~ normal(0,1);
  scale_new[1] ~ inv_gamma(2,2);
  scale_new[2] ~ inv_gamma(2,2);
  cov_mat[1, 2] ~ beta(5, 5);
  cov_mat[1, 3] ~ beta(2, 8);
  cov_mat[2, 1] ~ beta(5, 5);
  cov_mat[3, 1] ~ beta(2, 8);
  
  
  if(sim_data==0) {
    Y ~ multi_normal(mu, quad_form_diag(cov_mat, scale_new));
  }
}
2

@bgoodri

3

Hey, I’m not an expert but trying to get more into answering questions on here to challenge my own understanding. Do you not want cov_mat to be symmetric?

cov_mat[1, 2] ~ beta(5, 5);
cov_mat[1, 3] ~ beta(2, 8);
cov_mat[2, 1] ~ beta(5, 5);
cov_mat[3, 1] ~ beta(2, 8);

You’ve placed the same priors on the off diagonal elements but they will be independent draws, maybe you want something more like this:

cov_mat[1, 2] ~ beta(5, 5);
cov_mat[1, 3] ~ beta(2, 8);
cov_mat[2, 1]  = cov_mat[1, 2];
cov_mat[3, 1]  = cov_mat[1, 3];

If the matrix is still not positive definite then it might be worth adding a nugget term to your exact values to give the linear algebra some extra wiggle room.

What goes wrong with beta priors for the correlations?

I’m also not certain about the beta priors (I’m not super familiar with structural models so I might be way off here). Constraining covariance to be positive and between zero and one are both quite strong assumptions (the latter particularly when you don’t know what the scale is).

At the moment this isn’t a correlation matrix as your first and second variances are free so you are setting beta priors on covariances not on correlations. If you do set the diagonals to 1 then a beta might be more appropriate but under the current formulation you’re assuming positive correlation which may or may not be reasonable depending on your system.

I hope this isn’t unhelpful and maybe gives you a couple of things to try until someone more qualified can respond!

4

good point about symmetry! I changed that. However I still get the same divergent transitions. Can you clarify what you mean here a bit more:

At the moment this isn’t a correlation matrix as your first and second variances 
are free so you are setting beta priors on covariances not on correlations.
 If you do set the diagonals to 1 then a beta might be more appropriate
 but under the current formulation you’re assuming positive correlation 
which may or may not be reasonable depending on your system.

Why is this not a prior for correlations? When I generate corr_matrix[3] Omega, that automatically sets diagonal elements to 1 right? Also I set priors separately for \sigma-s and correlations and then use quad_form_diag, so that should overall end up with a covariance matrix.

I tried to read up on this paper: Generating random correlation matrices based on partial correlations - ScienceDirect which is about generating partial correlations from the beta distribution. Here I change priors to be:

cov_mat[1, 3] = (1-2*rho_13)*5
cov_mat[1, 2] = (1-2*rho_12)*5
where
rho_13 ~ beta(1, 1)
rho_12 ~ beta(1, 1)

This should be on the support (-1, 1). However, still cannot avoid divergent transitions.

5

@bgoodri will definitely have better answers.

Why is this not a prior for correlations? When I generate corr_matrix[3] Omega, that automatically sets diagonal elements to 1 right? Also I set priors separately for σ -s and correlations and then use quad_form_diag, so that should overall end up with a covariance matrix.

Ahhh I see, I think I got confused because cov_mat isn’t constrained in the transformed parameters block but what you say sounds right. My bad.

This should be on the support (-1, 1).

This looks better to me. Sorry it’s not a solution to your divergent transitions!! Have you tried visualising where the transitions are occurring using ShinyStan or similar?

6

ah I think I see what the potential problem could be. I use simulated data where i generate errors of x and y from a multivariate distribution but forget to generate errors of x, y and z together, not only x and y. That might be causing the issues as the model will not be well suited to the data.
I will change and update, thanks for help so far:)

1 Like