EDIT 2023-04-17
The answer below is wrong because $R \Sigma R$ is not lower-triangular in the general case. I leave it here for posterity. Down-votes may help ensure people know not to trust it.
Original post (which is incorrect)
Thanks for the nudges to explain more about my issue. The overarching problem is a statistical issue regarding multivariate normal distributions. Suppose random variable $y$ is a vector of length $n$ with a multivariate normal distribution with mean vector $\mu$ and positive-definite symmetric covariance matrix $\Sigma$. We denote this by:
$$
y \sim \text{MVN}(\mu, \Sigma)
$$
My goal is to efficiently work with the marginal distribution of a subset of the components of $y$:
$$
S R y \sim \text{MVN} ( SR\mu, SR\Sigma R^TS^T)
$$
In this problem, we assume we already know the Cholesky factorization $\Sigma = L L^T$. In order to quickly evaluate the density of the marginal distribution above, we want a shortcut to get the Cholesky factor of $SR\Sigma R^TS^T$ without having to run the standard factorization algorithm from scratch. This is because multi_normal_cholesky_lpdf() is much faster than multi_normal_lpdf() in Stan. Fortunately, a shortcut exists, as I prove below.
Here, $R$ is an $n \times n$ permutation matrix. In my case, $Ry$ is a column vector that puts the desired elements of the subset of $y$ at the front of the vector, while preserving the relative order of the desired components. To explain more about what I mean, suppose $y = [a \ \ b \ \ c \ \ d]^T$ and I want the marginal distribution of $[b \ \ c]^T$. Then, $Ry = [b \ \ c \ \ a \ \ d]^T$, with $b$ and $c$ at the front and $b$ before $c$.
Also, $S$ is an $m \times n$ matrix, where $m < n$ is the cardinality of the desired subvector $S R y$. Here, $S = [I_m \ \ 0_{m \times (n - m)}]$, where $I_m$ is the $m \times m$ identity matrix and $0_{m \times (n - m)}$ is the appropriately sized matrix of zeroes.
Now, let $L$ be the lower-triangular Cholesky factor of $\Sigma$ (with $\Sigma = L L^T$). We know $R L R^T$ is lower triangular because $R$ permutes the rows of $L$ in the exact same way as $R^T$ permutes the columns of $L$. (The previous sentence is definitely wrong, and it is straightforward to come up with small counterexamples.) In addition, $R^T = R^{-1}$ because $R$ is a permutation matrix and thus orthonormal. Thus:
$$
R L R^T(R L R^T)^T = R L R^T R L^T R^T = R L L^T R^T = R \Sigma R^T
$$
That means $R L R^T$ is the Cholesky factor of $R \Sigma R^T$. Finally, the matrix $S$ is constructed such that $S R \Sigma R^T S^T$ is the upper-left $m \times m$ block of $R \Sigma R^T$ and $S R L R^T S^T$ is the upper-left $m \times m$ block of $R \Sigma R^T$. By the lemma below, this makes $S R L R^T S^T$ the lower-triangular Cholesky factor of $S R \Sigma R^T S^T$. Best of all, $S R L R^T S^T$ is super fast to compute because it is just a subset of the rows and columns of $L$.
Lemma
Define $\Sigma$ to be an $n \times n$ positive definite symmetric matrix with lower-triangular covariance matrix $L$. Write these matrices in block notation:
$$
\Sigma = \begin{bmatrix}
A & B \\
C & D
\end{bmatrix} \qquad
L = \begin{bmatrix}
X & O \\
Y & Z
\end{bmatrix}
$$
where $A$ and $X$ are $m \times m$ blocks $(m < n)$ and $O$ is the appropriately sized matrix of zeroes. Then, $X$ is the $m \times m$ lower-triangular Cholesky factor of $m \times m$ matrix $A$.
Proof of the Lemma
By construction, matrix $X$ is already lower-triangular. From matrix multiplication and the definition of the Cholesky factorization, we find:
$$
\begin{aligned}
\Sigma &= L L^T \\
\begin{bmatrix}
A & B \\
C & D
\end{bmatrix} &= \begin{bmatrix}
X & O \\
Y & Z
\end{bmatrix} \begin{bmatrix}
X^T & Y^T \\
O^T & Z^T
\end{bmatrix} = \begin{bmatrix}
XX^T + OO^T & XY^T + OZ^T \\
YX^T + ZO^T & YY^T + ZZ^T
\end{bmatrix} = \begin{bmatrix}
XX^T & XY^T + OZ^T \\
YX^T + ZO^T & YY^T + ZZ^T
\end{bmatrix}
\end{aligned}
$$
Thus, $A = XX^T$.
