Zero Gaussian partial correlation is conditional independence

Partition $X=(Y,Z)$ with $Y=(X_i,X_j)$ and $Z=X_S$; if $S=\mathrm{\varnothing }$, omit the $Z$ block. For nonempty $S$, ${\mathrm{\Sigma }}_{S,S}$ is positive definite because every nonzero vector supported on $S$ has positive $\mathrm{\Sigma }$-quadratic form, so it is invertible. Put $A={\mathrm{\Sigma }}_{Y,S}{\mathrm{\Sigma }}_{S,S}^{-1}$ and $U=Y-{\mu }_Y-A(Z-{\mu }_Z)$. Then $\mathrm{Cov}(U,Z)=0$ and $\mathrm{Cov}(U)=K_{ij\mid S}$. For every nonzero $u\in {\mathbb{R}}^2$, the quadratic form of $\mathrm{\Sigma }$ at $(u,-{\mathrm{\Sigma }}_{S,S}^{-1}{\mathrm{\Sigma }}_{S,Y}u)$ equals $u^T{K}_{ij\mid S}u>0$. Thus $K_{ij\mid S}$ is positive definite; for empty $S$ this follows directly from the principal $2\times 2$ covariance block.

The vector $(U,Z)$ is a linear image of $X$, hence is jointly Gaussian by [F1]. Its covariance is block diagonal with blocks $K_{ij\mid S}$ and ${\mathrm{\Sigma }}_{S,S}$. By [F2], the joint density factors as the product of the two Gaussian densities, since the inverse and determinant of a block-diagonal covariance split by blocks. Integrating over measurable rectangles shows $U$ and $Z$ are independent by [F5]. It follows that, conditionally on $Z$, $Y$ is the translate ${\mu }_Y+A(Z-{\mu }_Z)+U$ and has conditional covariance $K_{ij\mid S}$.

If $(K_{ij\mid S}{)}_{12}=0$, the positive-definite bivariate Gaussian density of $U$ factors into its one-coordinate densities. Thus the coordinates of $U$ are independent. Since $U$ is independent of $Z$, the conditional law of $X_i,X_j$ given $Z$ is the product of its conditional marginals. Integrating bounded measurable functions against this product conditional law gives the identity in [F3], hence conditional independence.

Conversely assume the conditional independence in the Statement. Apply [F3] to the sine and cosine functions of each coordinate and combine the four real identities to obtain factorization of the conditional characteristic function at $(1,1)$. By step 2.1, the joint conditional characteristic function is the product of the marginal factors times $\exp (-(K_{ij\mid S}{)}_{12})$. All other factors are nonzero; therefore $\exp (-(K_{ij\mid S}{)}_{12})=1$. The covariance entry is real, so $(K_{ij\mid S}{)}_{12}=0$.

The diagonal entries of $K_{ij\mid S}$ are positive by step 1.1. Its off-diagonal entry is zero exactly when the normalized entry ${\rho }_{ij\mid S}$ is zero, proving both implications.