I am trying to find the gradient and hessian of a matrix-t-distribution representing a generalized linear dynamic linear model. At each time step k, we get a marginal t-distribution whose log is given by - \begin{equation} \propto \frac{-(v+p)}{2}\log \left (1+\frac{1}{v}(\eta_{k}-f_{k})^{T}X^{-1}(\eta_{k}-f_{k})\right) \end{equation} where $X^{-1} = c\Xi_{k-1}$. $\eta$ and $f$ are $p$-dimensional column vectors, $v$ is a constant, and $\Xi$ is a $p \times p$ dimensional symmetric matrix. And $k$ here represents one-time point in the time series.
I consider the log part in the equation above and it's a little simplified version as the term $S$ given by the following equation - \begin{equation} S = \frac{1}{v} \left [(\eta_{k} - f_{k})^{T} \ X^{-1} \ (\eta_{k} - f_{k}) \right ] \end{equation}
Furthermore, these matrices are related to each other by the following recursive equations: \begin{equation} A_{k} = G_{k}M_{k-1} \end{equation} \begin{equation} f_{k}^T = F_{k}^TA_{k} \end{equation} \begin{equation} e_{k}^T = \eta_{k}^T - f_{k}^T \end{equation} \begin{equation} M_{k} = A_{k} + F_{k}e_{k}^T \end{equation} \begin{equation} \Xi_{k} = \Xi_{k-1} + e_{k}e_{k}^T \end{equation}
where $G$ is a $q \times q$ matrix, $M$ and $A$ are $q \times p$ matrix, $F$ is a $q \times 1$ vector and $\eta$ is a $1 \times p$ vector. I was able to calculate the gradient and got this: \begin{equation} \frac{dS}{d\eta_{k}} = \frac{2}{v_{k-1}}(\eta_{k}^{T}-f_{k}^{T})X^{-1} \end{equation}
But I am having trouble finding the hessian and would really appreciate your help with it.
SOMETHING I TRIED TO DO:
I thought the hessian's off-diagonal elements can be calculated using: \begin{equation} \frac{dS}{d\eta_{k}d\eta_{k-1}} \end{equation} After performing some calculations, I got this - \begin{equation} dS = [d(\eta_{k}^{T}-f_{k}^{T})X^{-1} + (\eta_{k}^{T}-f_{k}^{T})d(X^{-1})]d\eta_{k} \end{equation}
\begin{equation} dS = [(0 - F_{k}^TG_{k}d(M_{k-1}))X^{-1} + (\eta_{k}^{T}-f_{k}^{T})(-X^{-1}dXX^{-1})]d\eta_{k} \end{equation}
\begin{equation} dS = [-F_{k}^TG_{k}F_{k-1}d\eta_{k-1}^TX^{-1} + (\eta_{k}^{T}-f_{k}^{T})(-X^{-1}dXX^{-1})]d\eta_{k} \end{equation}
Now we need to find $dX$, \begin{equation} dX = c\left [ d\eta_{k-1}(\eta_{k-1}^T - f_{k-1}^T) + (\eta_{k-1} - f_{k-1})d\eta_{k-1}^T\right] \end{equation}
But can't further simplify the expression.