How to calculate the variance of an estimator with simulation in R

Question

So I have a question in the exam of Statistics that I dont know how to do it:

There is $X\sim \text{Poisson}(\theta)$.

I want to estimate $\exp[-\theta] = P(X=0)$.

I know by the maximum likelihood method that a consistent estimator is $\exp[-\bar{X}]$.

Now is asked me to use R commands to calculate an approximation of $\text{Var}(\exp[-\bar{X}])$, with $n=20$.

I am really not getting how I should do this if I don't have any knowledge about the parameter $\theta$. How I should simulate the data?

Answer 1

Here is a short R code comparing the estimators $$\frac{1}{n}\sum_{t=1}^n \mathbb{I}_0(X_t)$$ (in yellow) and $$\exp\left\{-\sum_{t=1}^nX_t\big/n\right\}$$ (in brown) for $n=20$ for a range of values of $\theta$

compest<-function(T=1e3,tmin=.1,tmax=10){

 varz=matrix(0,50,2)
 theta=seq(tmin,tmax,le=50)
 for (t in 1:50){
  sampl=matrix(rpois(n=20*T,lambda=theta[t]),ncol=20)
  varz[t,]=c(var(apply(sampl==0,1,mean)),
             var(exp(-apply(sampl,1,mean))))}
 return(varz)}

This shows that the estimator based on $\bar{X}_n$ is leading to a smaller variance than the one based on the frequency of zero draws.

The continuous curves are the theoretical values of the variances, namely $e^{-\theta}(1-e^{-\theta})/n$ for the Binomial proportion of zero draws and $e^{-2\theta}\theta/n$ for the exponential of the average. (This variance is a delta-method approximation of the exact variance, but the fit is very good!)

Actually, as pointed out by George Henry on my blog, the derivation of the mean and variance of $$\exp\left\{-\sum_{t=1}^nX_t\big/n\right\}$$ is quite manageable: since $n\bar{X}_n$ is a Poisson $\mathscr{P}(n\theta)$ variable \begin{align*} \mathbb{E}[\exp\{-\bar{X}_n\}]&=\sum_{i=0}^\infty \exp\{-i/n\}\frac{(n\theta)^i}{i!}\exp\{-n\theta\}\\ &=\sum_{i=0}^\infty \left(\exp\{-1/n\} n\theta \right)^i \frac{\exp\{-\theta\}}{i!}\\ &=\exp\left\{-n\theta+n\theta\exp\{-1/n\} \right\}\\ &=\exp\left\{-n\theta[1-\exp\{-1/n\}]\right\} \end{align*} and \begin{align*} \mathbb{E}[\exp\{-\bar{X}_n\}^2]&=\sum_{i=0}^\infty \exp\{-2i/n\}\frac{(n\theta)^i}{i!}\exp\{-n\theta\}\\ &=\sum_{i=0}^\infty \left(\exp\{-2/n\} n\theta \right)^i \frac{\exp\{-\theta\}}{i!}\\ &=\exp\left\{-n\theta+n\theta\exp\{-2/n\} \right\}\\ &=\exp\left\{-n\theta[1-\exp\{-2/n\}]\right\} \end{align*} Hence $$\text{var}(\exp\{-\bar{X}_n\})=\exp\left\{-n\theta[1-\exp\{-2/n\}]\right\}-\exp\left\{-2n\theta[1-\exp\{-1/n\}]\right\}$$ As shown on the plot below, the difference with the approximation is hard to spot!