Hello. I am having a hard time making my PyMC model estimate work the same way my Stan code does/did. My model is sort of like an ordered logit.
One difference, maybe, is that I am using different priors. I have no idea how to implement the
induced_dirichlet_lpdf prior
defined in this example: Ordinal Regression
I mindlessly copied that code for my model in Stan. It is:
functions {
real induced_dirichlet_lpdf(vector c, vector alpha, real phi) {
int K = num_elements(c) + 1;
vector[K - 1] sigma = inv_logit(phi - c);
vector[K] p;
matrix[K, K] J = rep_matrix(0, K, K);
// Induced ordinal probabilities
p[1] = 1 - sigma[1];
for (k in 2:(K - 1))
p[k] = sigma[k - 1] - sigma[k];
p[K] = sigma[K - 1];
// Baseline column of Jacobian
for (k in 1:K) J[k, 1] = 1;
// Diagonal entries of Jacobian
for (k in 2:K) {
real rho = sigma[k - 1] * (1 - sigma[k - 1]);
J[k, k] = - rho;
J[k - 1, k] = rho;
}
return dirichlet_lpdf(p | alpha)
+ log_determinant(J);
}
}
and invoked with
c ~ induced_dirichlet(rep_vector(1, K), 0);
I would like to have the same prior in PyMC, and LLMs can’t tell me how to write this.
All I’ve been using so far is
cuts = pm.Normal('cuts',
mu= np.linspace(0, K,
nCutpoints),
sigma=sigma_cuts,
transform=pm.distributions.transforms.ordered,
)
which is giving me undesirable fits. (Let’s just say qualitatively different from what Stan gives).
I’m not sure why the priors could be the cause of my problem, but I’d like to make sure.
(1) Are these induced dirichlet priors generally a good choice for an ordered choice model?
(2) Can anyone write my code for me ie translate the relevant part of the Betancourt example linked above? I can’t understand many of the words in his section 2.2 there, let alone port the code.