I’m not sure how to construct this log abs determinant of the jacobian. Let’s assume, I have a symmetric and positive definite matrix X. Where each element of X is constructed from functions of parameters.
\begin{equation*}
X_{n,n} =
\begin{pmatrix}
f_{11}(a_{1,1}) & f_{12}(a_{1,2}) & \cdots & f_{1n}(a_{1,n}) \\
f_{21}(a_{2,1} )& f_{22}(a_{2,2}) & \cdots & f_{2n}(a_{2,n}) \\
\vdots & \vdots & \ddots & \vdots \\
f_{n1}(a_{n,1} ) & f_{n2}(a_{n,2} )& \cdots & f_{nn}(a_{n,n} )
\end{pmatrix}
\end{equation*}
I want to work on the Cholesky form because it’s more stable and much faster to do the next things. My question is if I get the jacobian from the change of variables into the Cholesky factor do I then have to apply the chain rule to that since I have functions on each element in X? Or can I just use the result below?
that is from
http://web.mit.edu/18.325/www/handouts/handout2.pdf
I found this https://rpubs.com/bbolker/cholblocks, which suggests I can just construct this matrix with the derivatives of the elements needed. Anyway, just seeing if there’s a simple way to construct this.
3 Likes
According to the “exercise” I posted do I only need to know the derivatives of the elements on the diagonal of the cholesky factor to compute the determinant of J(L)? I emphasized only because the diagonal of the Cholesky matrix contains the sum of the j^{th} row elements (ie j, 1:j-1).
I’m doing all this stuff and it’s getting hard to keep track of what the end result Jacobian should be. So I have an n \times n correlation matrix where each entry is from a vector size n(n-1)/2 parameter in the parameters block drawn from a normal distribution with mean = 0 and sigma chosen to get the desired shape of the marginal distribution of the correlations.
In the transformed parameters block, I transform the parameters vector so that each element is between -1, 1 using -1 + 2 / (\exp(x + 1). Then I’m constructing the correlation matrix that obey upper and lower bounds of the entries by transforming into spherical coordinates (see this question on stackexchange or the paper by Pinheiro and Bates). I want to add the jacobian transform, where I’m using L in a multi_normal_cholesky.
I found that the determinant jacobian for spherical coordinates is given as
\frac{\partial (x_1,\ldots,x_n)}{\partial(r,\theta_1,\ldots,\theta_n)} = r^{n-1}(\sin\theta_1)^{n-2}\cdots(\sin\theta_{n-2})
where r will be 1 in my case, from http://web.mit.edu/18.325/www/handouts/handout2.pdf.
Maybe an easier question is, how would I know if I did the determinant of the Jacobian right?
If those f functions are univariate, then the change of variables from \mathbf{a} to \mathbf{X} has a diagonal Jacobian matrix that is n \times n. Then there is another change of variables to the Cholesky factor, which has a triangular Jacobian that is also n \times n but the diagonal simplifies to the expression in Edelman’s lecture notes.
To check if you did things right, numerically differentiate the log-kernel with respect to the elements of \mathbf{a}.
2 Likes
Is it possible to take the C++ code generated from the stan functions blocks and use the jacobian() function in the stan math library to numerically differentiate?
Not very easily but there is most of an example in the vignette
1 Like
