Hi everyone,

I am seeking help on a series of two connected problems. First, how can I extend @bgoodri’s multivariate probit model here to cases where each the ordinal response variable is polychotomous and each category may be of different size (e.g., one variable has 3 categories, another may have 14). Second and connected to this problem, the response vector y is mixed with a combination of continuous and ordinal values. How (or is it possible) can I model this complex data vector?

I have fit a multivariate normal (and Student T) model to up to 6 continuous response variables (see model below). The next step is how do I combine these continuous variables with the ordinal variables that all make up a complex response vector. Let x be a scalar independent variable and **Y** be a vector of response variables. J indexes the ordinal variables, and K the continuous variables, so Y is a vector of size J + K. I assume the overall response vector is distributedas y_n ~ N (f(g,h), \sigma(x)) where h is a function for the continuous mean and g for that of an ordinal variable and \sigma(x) a function for the standard deviation. The continuous mean function is h_k(x) = a_k*x^r_k + b_k. The ordinal is g_k(x)=x^r_k. Here a = 1 and b=0 for identifiability. For the sd I assume linearity where \sigma(x) = s_k[1+\kappa_k*x]. For each ordinal variable j, only m categories are observable and where I need help is mapping this ordinal response to a latent continuous response z via any number of threshold parameters \tau where \tau_0 = -\infty and \tau_{m+1}= \infty.

Apologies if this is a bit confusing. Essentially, I would like to combine continuous and ordinal response variables in a principled manner. My domain knowledge suggests the mean is related to an offset power law and the noise is heteroskedastic and scales linearly with x. I provide the continuous model below, which is tested on simulated and real data and returns good results. I am seeking help to extend this to the ordinal + continuous sense. I am thinking it maybe something like this. But, I could be totally wrong.

Thank you for all the help!

```
data{
int N; // # of individuals
int K; // # of responses
vector[K] y[N]; // vector of growth responses per individual
real X[N]; //age
}
parameters{
vector<lower=0>[K] a; // multiplicative constant
vector<lower=0>[K] r; // scaling parameter
vector[K] b; // offset
vector<lower=0>[K] s_scale; //constant noise
vector<lower=0>[K] kappa; // slope/gradient of linear noise function
cholesky_factor_corr[K] L_Omega; // cholesky factor the corr matrix for efficiency
}
transformed parameters{
vector[K] mu[N]; //mean array of size N containing vectors of length K
vector<lower=0>[K] s[N]; //sd array of size N containing vectors of length K
matrix[K,K] Sigma[N]; // cov array of size N with K x K cov matrices
for(i in 1:N){
for(k in 1:K){
mu[i,k] = a[k]*X[i]^(r[k])+b[k]; // mean function
s[i,k] = s_scale[k]*(1 + kappa[k]*X[i]); // sd function
}
}
for(i in 1:N){
Sigma[i] = diag_pre_multiply(s[i],L_Omega); // combining the cholesky corr matrix with the noise
}
}
model{
//priors
a ~ normal(0,10); // half-normal
r ~ normal(0,1); //half-normal
b ~ normal(0,10);
kappa ~ normal(0,1); //half-normal
s_scale ~ cauchy(0,5); //half-cauchy
L_Omega ~ lkj_corr_cholesky(0.5);
for(i in 1:N){
y[i] ~ multi_normal_cholesky(mu[i], Sigma[i]); //likelihood
}
}
generated quantities{
vector[K] Y_mean[N]; //posterior mean
vector[K] Y_pred[N]; //posterior predictive check
corr_matrix[K] Omega; //re-factor corr matrix
cov_matrix[K] Sigma_new[N];//put cov matrix back together
vector[N] log_lik; // log-likelihood for model comparison
vector[K] s_new[N]; // posterior variance
Omega = multiply_lower_tri_self_transpose(L_Omega); //corr matrix
for(i in 1:N){
for(k in 1:K){
Y_mean[i,k] = a[k]*X[i]^r[k] +b[k]; //posterior mean
s_new[i,k] = s_scale[k]*(1+kappa[k]*X[i]);
}
Sigma_new[i] = quad_form_diag(Omega, s_new[i]); //posterior cov matrix
Y_pred[i] = multi_normal_rng(Y_mean[i],Sigma_new[i]); //posterior predictive
}
for(n in 1:N){
log_lik[n] = multi_normal_cholesky_lpdf(y[n]|mu[n],Sigma[n]); // log lik of each obs
}
}
```