Confidence Interval for Pairs of Normal Data Suppose we have $\alpha_1,..\alpha_n,\beta_1,..\beta_n$~$N(0,1)$.
Define $x_i = u_i + \sigma\alpha_i$, $y_i = u_i + \sigma\beta_i$, for $(\mu_1,..,\mu_n) \in \mathbb{R}^n$, $\sigma > 0$.
Consider the model $(x_1,y_1),..,(x_n,y_n)$, for it we get the likelihood function which will yield:
$\mu_i* = \frac{x_i+y_i}{2}$ as the MLE estimator for $\mu_i$, and $\sigma* = \frac{1}{4n}\sum^{23}_{i=1}(x_i-y_i)^2$.
It follows the $\mu_i*$ is unbiased for $\mu_i$, and $2\sigma*$ for $\sigma$.
I want to build a $0.99$ confidence interval for $\mu_i$, supposing $n=23$;
I'm unsure which is the correct way to calculate this - 
either $[\mu_i* \pm \frac{\sqrt{\sigma*}}{\sqrt{23}}t_{22,.975}]$
or $[\mu_i* \pm \frac{\sqrt{\sigma*}}{\sqrt{2}}t_{22,.975}]$. I think the latter is correct since we would base this CI on the Cental Limit Theorem using the mean of $x_i$ and $y_i$.
 A: 
Suppose we have $\alpha_1,..\alpha_n,\beta_1,..\beta_n$~$N(0,1)$.... Define $x_i = u_i + \sigma\alpha_i$, $y_i = u_i + \sigma\beta_i$, for $(\mu_1,..,\mu_n) \in \mathbb{R}^n$, $\sigma > 0$... $\sigma* = \frac{1}{4n}\sum^{23}_{i=1}(x_i-y_i)^2$... I want to build a $0.99$ confidence interval for $\mu_i$, supposing $n=23$;

When considering the confidence interval for $\mu_i$, there are two things that need to be taken into account:


*

*The estimate of $\sigma$ might be off.

*For the true $\sigma$, the outcome of $x_i - y_i$ might be off.
Let's consider the first item. $ x_i - y_i  \sim \mathcal{N}(0, 2\sigma^2)$. As such, 
$$
\frac{1}{2n} \sum_{i = 1}^n \left( x_i - y_i \right)^2
$$
is a known-mean unbiased estimator for $\sigma^2$ (why $4n$ in your question?). However, $
\frac{1}{2n} \sum_{i = 1}^n \left( x_i - y_i \right)^2
$ is itself a random variable, and so has variance. Using the properties of chi-square distributions (see this question), an $\epsilon$ confidence interval for $2 \sigma^2$ is
$$
\left(\frac{S^2 n}{\chi^2_{1-\epsilon2}},
 \frac{S^2 n}{\chi^2_{\epsilon/2}}\right). 
$$
Obviously, the worst case is the right boundary of this interval. If you plug in this value for $2 \sigma^2$, you can calculate a $\delta$ interval for the normal distribution using the usual method.
However, note that you need to take both uncertainties into account. One way of upper-bounding the probability of error would be to upper-bound it by taking the probability of error as $(1 - \epsilon) ( 1 - \delta)$. If you're looking for a 99% CI, you need to use $\epsilon, \delta$, s.t. 
$(1 - \epsilon) ( 1 - \delta) \geq 0.99$. One way of doing so (not necessarily the optimal) would be taking $\epsilon = \delta = 1 - \sqrt{0.99}$.
