Model Sum of Squares from ANOVA table Background: I am studying a course on statistical experiments using the textbook by Douglas Montgomery on the analysis of experiments. This is an introductory course and so I am relatively new to the study of these types of models.
Question: Consider the model $Y_{i,j,k} = \mu + R_i + C_j + T_k + \epsilon_{i,j,k}$ where $R_i$ is the impact of from row $i$, $C_j$ is the impact from column $j$ and $T_k$ is the impact from treatment $k$ where there are exactly $n$ rows, $n$ columns, and $n$ treatments. We know that the row sum of squares $S^2_{\text{row}} = \sum _ {i=1}^n \frac{Y_{i \cdot \cdot}^2}{n} - \frac{Y_{\cdot \cdot \cdot}^2}{n^2}$ and similarly for $S^2_{\text{column}}$ and $S^2_{\text{treatment}}$ we have the same form but instead we sum over $j$ and $k$ respectively from $1$ to $n$ (we take the formulae for the sums of squares directly from the corresponding ANOVA table).
I have read in the textbook (without justification) the formula for the expected value $ \color{red}{\mathbb{E}(S^2_{\text{row}}\space)}$ and (again) the same applies for the columns and treatments by replacing $R_i^2$ with $C_j^2$ (for the columns) and $T_k^2$ (for the treatments).
$\text{How do we find the formula for the expectation?}$
Note: we are making the standard modelling assumptions for these types of models where the error terms are independent and normally distributed with mean $0$ and variance $\sigma^2$.

The model you described is LSD.
Here $$ S^2_R = n\sum_i\left(\bar{y}_{i\cdot\cdot}- \bar{y}_{\cdot\cdot\cdot} \right)^2.$$
Then
\begin{align}\mathbb E\left[S^2_R\right] &= \mathbb E\left[n\sum_i\left(R_i + \bar{\varepsilon}_{i\cdot\cdot} - \bar{\varepsilon}_{\cdot\cdot\cdot}\right)^2\right]\\ &= n\left[\sum_iR_i^2+\mathbb E\left\{\sum_i(\bar{\varepsilon}_{i\cdot\cdot} - \bar{\varepsilon}_{\cdot\cdot\cdot})^2\right\} + 2\underbrace{\sum_i\alpha_i\mathbb E(\bar{\varepsilon}_{i\cdot\cdot} - \bar{\varepsilon}_{\cdot\cdot\cdot})}_{=0}\right]\\ &= n\left[\sum_iR_i^2+\mathbb E\left\{\sum_i\bar{\varepsilon}^2_{i\cdot\cdot} - n\cdot \bar{\varepsilon}^2_{\cdot\cdot\cdot}\right\} \right] \\ &= n\left[\sum_iR_i^2 + \sum_i\frac{\sigma^2} n - \frac{\sigma^2}n\right] \\ &= n\sum_i R_i^2 + (n-1)\sigma^2.\end{align}
