This is more or less a conceptual question — I will write it out using math instead of Stan code.

I have a multilevel model where group-level y_i \in \mathbb{R}^n are modeled as a multivariate normal distribution. The twist is that for each y_i, some elements are parameters and some are data, depending on i. Estimating a centered parametrization is straightforward, but I run into convergence issues.

I am trying to convert to a non-centered parametrization. Let L be a Cholesky correlation factor, and \sigma a vector of standard deviations. Then

would allow me to go back and forth between standard normal IID \varepsilon's and y. But given that some elements of each y_i are data, I am constrained in the elements of \varepsilon's.

This is how I conceptualize this. Let

for easier notation. Assume that y_i and \varepsilon_i can be partitioned as *p*arameters and *d*ata, so

Then I parametrize in terms of \varepsilon_p for each i, and have

I calculate the log likelihood based on [\varepsilon_{i,p}, \varepsilon_{i,d}] as being standard normals, and my log Jacobian determinant correction is

which thankfully is the same for all i.

All I have to deal with is a heinous amount of bookeeping. I am asking for feedback on whether the above would be the recommended approach, or I missed an obvious clean solution.