I have some question.
Let say we have sample with size $3$ from $N(\mu,\sigma)$ with known $\sigma$ How can I calculate the Variance of the estimator $W(X_1,X_2,X_3)$?
What I've tried to do:
$$ \theta = \text{$p^{\text{th}}$ quantile of $X$ distribution}, \qquad p=\frac{k}{n+1} $$
For $f(\theta) \ne 0$:
$$ E[X_k] \approx \theta, \qquad \operatorname{Var}[X_k]\approx\frac{p(1-p)}{(n+2)(f(\theta))^2}$$
For the normal distribution: $p=0.5$, $\theta=\phi^{-1}(p)=\mu$, and
$$ \operatorname{Var[X_2]} \approx\frac{0.25}{5\Bigl(\frac{1}{\sqrt{2\pi}\sigma}\exp\Bigl(-\frac{(\mu-\mu)^2}{2\sigma^2}\Bigr)\Bigr)^2} =\frac{2\pi\sigma^2}{20} \approx31.4 \quad !?$$
This is my R code, here I get something around 45.30
N=10000
n=3
mu=75
sigma=10
W2 = c(0,N)
c_var = c()
for (i in 1:N) {
X=rnorm(n, mean = mu, sd = sigma)
W2[i] = median(X)
c_var = c(c_var, var(W2[1:i]))
}
plot(c_var, type = "l", lwd = 2, col = "blue", ylab = "Oberved Var - Normal",
xlab = "Number of Obs", main='Variance of median sample of Normal distribution')
print(var(W2))
Writing $Z_i = (X_i - \mu)/\sigma$ and noting that the order statistics satisfy $X_{(i)} = \mu + \sigma Z_{(i)}$, it suffices to study the variance of $Z_{(2)}$. Now, it is well-known that the PDF of $Z_{(2)}$ is given by
$$ f_{Z_{(2)}}(x) = 6\Phi(x)(1-\Phi(x))\phi(x), $$
where $\Phi$ and $\phi$ are the CDF and PDF of the standard normal distribution, respectively. Now by noting that
$$ f_{-Z_{(2)}}(x) = f_{Z_{(2)}}(-x) = 6\Phi(-x)(1-\Phi(-x))\phi(-x) = 6(1-\Phi(x))\Phi(x)\phi(x) = f_{Z_{(2)}}(x), $$
it follows that the distribution of $Z_{(2)}$ is symmetric about $0$ and hence $ \mathbf{E}[Z_{(2)}] = 0 $. So, the variance of $Z_{(2)}$ is the same as $\mathbf{E}[Z_{(2)}^2]$. However, by utilizing the formula
$$ \int x \phi(x) \, \mathrm{d}x =- \phi(x) + \mathsf{C}, $$
which can be easily proved by u-substitution, we get
\begin{align*} \mathbf{Var}(Z_{(2)}) &= \mathbf{E}[Z_{(2)}^2] \\ &= \int_{-\infty}^{\infty} 6x^2\Phi(x)(1-\Phi(x))\phi(x) \, \mathrm{d}x \\ &= \int_{-\infty}^{\infty} 6 \bigl( x\Phi(x)(1-\Phi(x)) \bigr)'\phi(x) \, \mathrm{d}x \tag{$\because$ IbP} \\ &= \int_{-\infty}^{\infty} \bigl( \underbrace{6 \Phi(x)(1-\Phi(x)) \phi(x)}_{=f_{Z_{(2)}}(x)} + \underbrace{6x\phi(x)^2}_{\text{odd}} - 12x\phi(x)^2\Phi(x) \bigr) \, \mathrm{d}x \\ &= 1 + 0 - \int_{-\infty}^{\infty} 12x\phi(x)^2\Phi(x) \, \mathrm{d}x. \end{align*}
So it suffices to compute the last integral:
\begin{align*} &\int_{-\infty}^{\infty} 12x\phi(x)^2\Phi(x) \, \mathrm{d}x \\ &= \frac{6}{\sqrt{2}\pi^{3/2}} \int_{-\infty}^{\infty} \int_{-\infty}^{\infty} xe^{-x^2-y^2/2} \mathbf{1}_{\{ y \leq x \}} \, \mathrm{d}x \mathrm{d}y \\ &= \frac{6}{\pi^{3/2}} \int_{-\infty}^{\infty} \int_{-\infty}^{\infty} xe^{-x^2-y^2} \mathbf{1}_{\{ \sqrt{2}y \leq x \}} \, \mathrm{d}x \mathrm{d}y \tag{$y \mapsto \sqrt{2}y$} \\ &= \frac{6}{\pi^{3/2}} \biggl( \int_{\arctan(1/\sqrt{2})-\pi}^{\arctan(1/\sqrt{2})} \cos \theta \, \mathrm{d}\theta \biggr) \biggl( \int_{0}^{\infty} r^2 e^{-r^2} \, \mathrm{d}r \biggr)\\ &= \frac{6}{\pi^{3/2}} \biggl( \frac{2}{\sqrt{3}} \biggr) \biggl( \frac{\sqrt{\pi}}{4} \biggr) \\ &= \frac{\sqrt{3}}{\pi}. \end{align*}
Therefore
$$ \mathbf{Var}(Z_{(2)}) = 1 - \frac{\sqrt{3}}{\pi} \approx 0.448671. $$
The following is a histogram for the $10^8$ simulated medians of $(Z_1, Z_2, Z_3)$ together with the summary statistics:
Here is the Mathematica code used to generate the above histogram: