I'm stuck with the following integral: $I=\int ... \int exp\Big(-\frac{1}{2} \sum \limits_{t=1}^{n} x_{t}^T{\Sigma_{x}}^{-1} x_{t}+\sum \limits_{t=1}^{n} x_{t}^T{\Sigma_{x}}^{-1} z_{t} -\frac{1}{2} \sum \limits_{t=1}^{n} z_{t}^T{\Sigma_{x}}^{-1} z_{t}-\frac{1}{2} a_{1}^T \Big({\Sigma_{a}}^{-1}+n {\Sigma_{z}}^{-1}\Big)a_{1}+ (z_{n}-z{1})^T{\Sigma_{z}}^{-1} a_{1}-\frac{1}{2} \sum \limits_{t=1}^{n} \Delta z_{t}^T{\Sigma_{z}}^{-1}\Delta z_{t} \Big) da_{1}dz_{1} dz_{2}...dz_{n}$
Because $\sum \limits_{t=1}^{n} x_{t}^T{\Sigma_{x}}^{-1} z_{t}=Tr(({\Sigma_{x}}^{-1}X)^T Z)$ and $\sum \limits_{t=1}^{n} z_{t}^T{\Sigma_{x}}^{-1} x_{t}=Tr(Z^T{\Sigma_{x}}^{-1}Z)$, we can use the matrix cookbook formula for Multivariate Gaussians (eq 346, version 2012) to get: \ $\int ... \int exp \Big(-\frac{1}{2} Tr(({\Sigma_{x}}^{-1}X)^T Z)+Tr(Z^T{\Sigma_{x}}^{-1}Z) \Big)dZ= {\sqrt{det(2\pi\Sigma_{x})}}^{n} exp \Big(-\frac{1}{2}Tr(X^T{\Sigma_{x}}^{-1}X) \Big)$\
But then, what do we do with the fact that we have $z_{1}$ and $z_{n}$ mixed with the $a_{1}$ term, as well as $\sum \limits_{t=1}^{n} \Delta z_{t}^T{\Sigma_{z}}^{-1}\Delta z_{t} $ ???? \
If we expand it, we get: $\sum \limits_{t=1}^{n} \Delta z_{t}^T{\Sigma_{z}}^{-1}\Delta z_{t}=\sum \limits_{t=1}^{n} z_{t}^T{\Sigma_{z}}^{-1} z_{t}+2\sum \limits_{t=1}^{n} z_{t}^T{\Sigma_{z}}^{-1} z_{t-1}+\sum \limits_{t=1}^{n} z_{t-1}^T{\Sigma_{z}}^{-1} z_{t-1} $\
How can one deal with the cross terms $z_{t}z_{t-1} $? Is there a way to integrate these terms together with the previous ones and make better use of the trace to solve the integral $I$? Any help would be very much appreciated.