# variance of sample median of normal distribution

Does anyone tell me how to find the variance of sample median when random samples are normally distributed?? I know that I need to use the order statistics, but actually I don't understand how to apply it.

I will make some illustrations with R. First, using results discussed at Central limit theorem for sample medians the asymptotic distribution of sample median of $N$ iid observations from a normal distribution is itself a normal distribution centered at the theoretical median $\mu$ and with variance $\frac{2\pi \sigma^2}{4 N}$ where $\sigma^2$ is the variance of the parent distribution. The sample median is the central order statistic, and in the case $N=2s+1$ is odd, it is the order statistic of order $s+1$, which have a distribution we can calculate exactly rather easily (theory given in the cited post and elsewhere at this site). I will use some R function I detail at the end of this post to compare the exact and approximate distributions of the sample median. We use the standard normal for parent distribution.

First with sample size $N=11$: The approximation is already quite good! Let us look at the case $N=101$: Here the approximation is virtually exact. Some R functions for distributions of order statistics:

porder  <-  function(x, N=1, r=1, basename="norm", ..., log.p=FALSE) {
pfun  <-  get(paste("p", basename, sep=""))
stopifnot(r  <=  N)
prob  <-  pfun(x, ...)
retval  <-  pbeta(prob, r, N-r+1, log.p=log.p)
retval
}

dorder <-  function(x, N=1, r=1, basename="norm", ..., log=FALSE)  {
pfun  <-  get(paste("p", basename, sep=""))
dfun  <-  get(paste("d", basename, sep=""))
stopifnot(r <= N)
logdens  <-  -lbeta(r, N-r+1)+(r-1)*pfun(x, ..., log.p=TRUE) +(N-r)*pfun(x, ..., lower.tail=FALSE, log.p=TRUE) + dfun(x, ..., log=TRUE)
retval  <-  if (log) logdens else exp(logdens)
retval
}

qorder  <-  function(p, N=1, r=1, basename="norm",  ..., log.p=FALSE) {
pfun  <-  get(paste("p", basename, sep=""))
qfun  <-  get(paste("q", basename, sep=""))
stopifnot( r <= N)
retval  <-  qfun(qbeta( p, r, N-r+1, log.p=log.p), ..., log.p=FALSE)
retval
}

rorder  <-  function(n, N=1, r=1, basename="norm", ...)  {
qfun  <-  get(paste("q", basename, sep=""))
qfun( rbeta(n, r, N-r+1), ... )
}