Difficulty showing Cochran's theorem for the chi-squared distribution due to dependence of variables
We first note that, if $\mathbf{1} = (1, \ldots, 1)$ and $\mathrm{X} = (X_1,\ldots,X_n)$ are regarded as column vectors, then $S$ is a quadratic form in $X$:
$$ S = \sum_{k=1}^{n} (X_k - \overline{X})^2 = \sum_{k=1}^{n} X_k^2 - \frac{1}{n} (X^{\mathsf{T}} \mathbf{1})(\mathbf{1}^{\mathsf{T}} X) = X^{\mathsf{T}} P X, $$
where $P = I_n - \frac{1}{n}\mathbf{1}\mathbf{1}^{\mathsf{T}}$ is the orthogonal projection onto the orthogonal complement $\operatorname{span}(\mathbf{1})^{\perp}$. So, the moment generating function of $S$ is given by
\begin{align*} \mathbf{E}[e^{\lambda S}] = \int_{\mathbb{R}^n} e^{\lambda x^{\mathsf{T}}Px} \frac{e^{-|x|^2/2}}{(2\pi)^{n/2}} \, \mathrm{d}x = \frac{1}{(2\pi)^{n/2}} \int_{\mathbb{R}^n} e^{-x^{\mathsf{T}}\left( \frac{1}{2}I_n - \lambda P \right)x} \, \mathrm{d}x. \end{align*}
Now if $\Sigma$ is an $n\times n$ positive-definite matrix, then it is well-known that
$$ \int_{\mathbb{R}^n} e^{-x^{\mathsf{T}}\Sigma x} \, \mathrm{d}x = \frac{\pi^{n/2}}{\sqrt{\det \Sigma}}. $$
Since the eigenvalues of $\frac{1}{2}I_n - \lambda P$ are $\frac{1}{2}, \frac{1}{2}-\lambda, \ldots, \frac{1}{2}-\lambda$, it follows that $\frac{1}{2}I_n - \lambda P$ is positive-definite for $\lambda < \frac{1}{2}$ and
\begin{align*} \mathbf{E}[e^{\lambda S}] &= \frac{1}{2^{n/2}\sqrt{\det\left( \frac{1}{2}I_n - \lambda P \right)}} = \frac{1}{\sqrt{\det\left( I_n - 2 \lambda P \right)}} = \frac{1}{(1-2\lambda)^{(n-1)/2}}. \end{align*}
This is precisely the m.g.f. of $\chi^2_{n-1}$, and therefore $S \sim \chi^2_{n-1}$.