[proofplan]
We first record that the random variables involved are bounded, so all expectations and second moments used below are finite. The mean follows from linearity of expectation and the fact that each [Rademacher random variable](/page/Rademacher%20Random%20Variable) has mean zero. For the variance, we expand $S^2$ as a finite double sum: the diagonal terms contribute $a_i^2$ because $\mathbb E[\varepsilon_i^2]=1$, while the off-diagonal terms vanish because independence gives $\mathbb E[\varepsilon_i\varepsilon_j]=\mathbb E[\varepsilon_i]\mathbb E[\varepsilon_j]=0$ for $i\ne j$.
[/proofplan]
[step:Verify integrability and compute the mean by linearity of expectation]
Since each $\varepsilon_i$ is Rademacher, it takes values in $\{-1,1\}$ almost surely. Therefore $|\varepsilon_i|=1$ almost surely for every $i\in\{1,\ldots,n\}$, and the finite sum $S$ is integrable and square-integrable. By [citetheorem:10110], for every $i\in\{1,\ldots,n\}$,
\begin{align*}
\mathbb E[\varepsilon_i]=0.
\end{align*}
Using finite linearity of expectation and the fact that each coefficient $a_i$ is deterministic, we obtain
\begin{align*}
\mathbb E[S]=\mathbb E\left[\sum_{i=1}^n a_i\varepsilon_i\right].
\end{align*}
Hence
\begin{align*}
\mathbb E[S]=\sum_{i=1}^n a_i\mathbb E[\varepsilon_i].
\end{align*}
Substituting $\mathbb E[\varepsilon_i]=0$ for every $i$ gives
\begin{align*}
\mathbb E[S]=0.
\end{align*}
[/step]
[step:Show that all off-diagonal products have expectation zero]
Let $i,j\in\{1,\ldots,n\}$ with $i\ne j$. Since $\varepsilon_i$ and $\varepsilon_j$ are independent Rademacher random variables, their generated $\sigma$-algebras are independent. In particular, for each $r,s\in\{-1,1\}$,
\begin{align*}
\mathbb P(\varepsilon_i=r,\varepsilon_j=s)=\mathbb P(\varepsilon_i=r)\mathbb P(\varepsilon_j=s).
\end{align*}
Because $\varepsilon_i\varepsilon_j$ takes values in $\{-1,1\}$, its expectation is the finite sum
\begin{align*}
\mathbb E[\varepsilon_i\varepsilon_j]=\sum_{r\in\{-1,1\}}\sum_{s\in\{-1,1\}}rs\,\mathbb P(\varepsilon_i=r,\varepsilon_j=s).
\end{align*}
Using independence in this finite sum gives
\begin{align*}
\mathbb E[\varepsilon_i\varepsilon_j]=\sum_{r\in\{-1,1\}}\sum_{s\in\{-1,1\}}rs\,\mathbb P(\varepsilon_i=r)\mathbb P(\varepsilon_j=s).
\end{align*}
Factoring the product of finite sums yields
\begin{align*}
\mathbb E[\varepsilon_i\varepsilon_j]=\left(\sum_{r\in\{-1,1\}}r\,\mathbb P(\varepsilon_i=r)\right)\left(\sum_{s\in\{-1,1\}}s\,\mathbb P(\varepsilon_j=s)\right).
\end{align*}
The two factors are $\mathbb E[\varepsilon_i]$ and $\mathbb E[\varepsilon_j]$, respectively, so [citetheorem:10110] gives
\begin{align*}
\mathbb E[\varepsilon_i\varepsilon_j]=0.
\end{align*}
[guided]
Fix distinct indices $i,j\in\{1,\ldots,n\}$. The point of this step is to justify precisely why the cross terms in the square of $S$ disappear. Since $\varepsilon_i$ and $\varepsilon_j$ are independent random variables, the events $\{\varepsilon_i=r\}$ and $\{\varepsilon_j=s\}$ are independent for every $r,s\in\{-1,1\}$. Thus
\begin{align*}
\mathbb P(\varepsilon_i=r,\varepsilon_j=s)=\mathbb P(\varepsilon_i=r)\mathbb P(\varepsilon_j=s).
\end{align*}
The product $\varepsilon_i\varepsilon_j$ is a bounded real-valued [random variable](/page/Random%20Variable) taking values in $\{-1,1\}$, so its expectation is computed by summing over the four possible pairs of values:
\begin{align*}
\mathbb E[\varepsilon_i\varepsilon_j]=\sum_{r\in\{-1,1\}}\sum_{s\in\{-1,1\}}rs\,\mathbb P(\varepsilon_i=r,\varepsilon_j=s).
\end{align*}
Substituting the independence factorisation gives
\begin{align*}
\mathbb E[\varepsilon_i\varepsilon_j]=\sum_{r\in\{-1,1\}}\sum_{s\in\{-1,1\}}rs\,\mathbb P(\varepsilon_i=r)\mathbb P(\varepsilon_j=s).
\end{align*}
This double sum separates because it is a finite sum of products:
\begin{align*}
\mathbb E[\varepsilon_i\varepsilon_j]=\left(\sum_{r\in\{-1,1\}}r\,\mathbb P(\varepsilon_i=r)\right)\left(\sum_{s\in\{-1,1\}}s\,\mathbb P(\varepsilon_j=s)\right).
\end{align*}
The first parenthesis is $\mathbb E[\varepsilon_i]$, and the second parenthesis is $\mathbb E[\varepsilon_j]$. By [citetheorem:10110], each Rademacher random variable has mean zero. Therefore
\begin{align*}
\mathbb E[\varepsilon_i\varepsilon_j]=\mathbb E[\varepsilon_i]\mathbb E[\varepsilon_j]=0.
\end{align*}
This proves that every off-diagonal product has expectation zero.
[/guided]
[/step]
[step:Expand the second moment and identify the diagonal contribution]
Since $S$ is square-integrable, its second moment is finite. Expanding the square of the finite sum gives
\begin{align*}
S^2=\left(\sum_{i=1}^n a_i\varepsilon_i\right)^2.
\end{align*}
Thus
\begin{align*}
S^2=\sum_{i=1}^n\sum_{j=1}^n a_i a_j\varepsilon_i\varepsilon_j.
\end{align*}
Taking expectations and using finite linearity of expectation, we get
\begin{align*}
\mathbb E[S^2]=\sum_{i=1}^n\sum_{j=1}^n a_i a_j\mathbb E[\varepsilon_i\varepsilon_j].
\end{align*}
For $i=j$, [citetheorem:10110] gives
\begin{align*}
\mathbb E[\varepsilon_i^2]=1.
\end{align*}
For $i\ne j$, the preceding step gives
\begin{align*}
\mathbb E[\varepsilon_i\varepsilon_j]=0.
\end{align*}
Therefore only the diagonal terms remain:
\begin{align*}
\mathbb E[S^2]=\sum_{i=1}^n a_i^2.
\end{align*}
[/step]
[step:Convert the second moment computation into the variance formula]
By the definition of variance for a square-integrable real-valued random variable,
\begin{align*}
\operatorname{Var}(S)=\mathbb E[(S-\mathbb E[S])^2].
\end{align*}
Equivalently,
\begin{align*}
\operatorname{Var}(S)=\mathbb E[S^2]-(\mathbb E[S])^2.
\end{align*}
The first step proved $\mathbb E[S]=0$, and the preceding step proved
\begin{align*}
\mathbb E[S^2]=\sum_{i=1}^n a_i^2.
\end{align*}
Substituting these two identities yields
\begin{align*}
\operatorname{Var}(S)=\sum_{i=1}^n a_i^2.
\end{align*}
Together with $\mathbb E[S]=0$, this is the desired conclusion.
[/step]
