2
$\begingroup$

I'm having some trouble with understanding how to calculate the variance of a non-parametric estimator.

The example comes from Wasserman's "All of statistics book"

Let $X_1, \ldots,X_n \sim \text{Uniform}(a,b)$ where $a$ and $b$ are unknown parameters and $a \lt b$.

Let $\tau = \int x dF(x)$. Let $\hat{\tau}$ be the maximum likelihood estimator (MLE) of $\tau$. Let $\tilde{\tau}$ be the nonparametric plug-in estimator of $\tau = \int{xdF(x)}$.

Suppose that $a=1$, $b=2$, and $n=10$. Find the mean-square-error (MSE) of $\hat{\tau}$ by simulation. Find the MSE of $\tilde{\tau}$ analytically. Compare

So I believe the MLE of $\tau$ is $\hat{\tau} = \frac{\hat{a}+\hat{b}}{2}$ where $\hat{a} = min(X_i)$ and $\hat{b} = max(X_i)$ for $i \in\{1,\ldots,n\}$ are the ML-estimators of the parameters since: $$\hat{\tau} = \int_\hat{a}^\hat{b}{\frac{x}{\hat{b}-\hat{a}}} = \frac{\hat{a}+\hat{b}}{2}$$

(I believe I used the correct notation here, but PLEASE let me know if I am not)

Using this equation, you can simulate $N$ datasets of $X_1,\ldots,X_{10} \sim \text{Uniform}(a,b)$ (since $n=10$) and calculate the MSE using the formula:

$$\text{MSE}_{\hat{\tau}} = \frac{1}{N-1}\sum_{j=0}^N{(\hat{\tau}_j - \bar{\hat{\tau}})^2}$$

where $\hat{\tau}_j = \frac{a_j+b_j}{2}$, $a_j = min(X_{ij})$, $b_j = max(X_{ij})$ for independent simulation $j$ of 10 datapoints drawn from $\text{Uniform}(1,3)$. $\bar{\hat{\tau}}$ is the mean of all the estimates of tau from the simulations.

Performing this simulation, you get the following distribution of MSE of $\hat{\tau}$:

import numpy as np
import matplotlib.pyplot as plt
MSE = lambda x: np.sum([np.power(_-np.mean(x),2) for _ in x])/float(len(x)-1)
plt.hist([MSE([tau(np.random.uniform(1,3,10)) for _ in range(10)]) for __ in range(1000)],bins=50)
plt.title('MSE of tau by simulation')

MSE <span class=$\hat{\tau}$">

The mean of the MSEs by MLE is ~0.015

Now, the non-parametric portion of the question is getting me a bit confused, and this is where I'm more unsure if I'm using the correct notation.

I believe that the non-parametric plug-in estimator of $\tau$ is: $$\tilde{\tau} = \int{xd\tilde{F}(x)} = \frac{1}{N}\sum_{i=1}^N{X_i}$$

Now, to get the MSE of $\tilde{\tau}$, I imagine you'd want to find the variance: $$V[\tilde{\tau}] = V\left[\int{xd\tilde{F}(x)}\right] = V\left[\frac{1}{N}\sum_{i=1}^N{X_i}\right] = \frac{V[X_i]}{N} = \frac{\sum_{i=0}^N{(X_i - \frac{1}{N}\sum_{i=1}^N{X_i})^2}}{N(N-1)} = \frac{\sum_{i=0}^N{(X_i - \bar{X_N})^2}}{N(N-1)}$$

Or using the information that $a=1$, $b=3$, and $N=10$ the MSE is $$\frac{V[X_i]}{N} = \frac{(3-1)^2}{12 N} = \frac{4}{120} = 0.0333\ldots$$

Is this method correct? Am I missing some key understanding here?

$\endgroup$
3
$\begingroup$

First of all, the estimator $\text{MSE}_{\hat{\tau}}$ you have given is a consistent estimator of the MSE only if $\text{Bias}(\widehat{\tau},\tau)=0$. To see this, recall that \begin{equation} \text{MSE}_{\tau}(\hat{\tau})=\mathbb{E}_{\tau}\left[(\hat{\tau}-\tau)^2\right]=\text{Var}_{\tau}(\hat{\tau})+\text{Bias}(\widehat{\tau},\tau)^2 \end{equation} Thus using $\tau=\mathbb{E}(X_1)=\frac{a+b}{2}=2$ you might obtain a consistent estimator of $\text{MSE}_{\tau}(\hat{\tau})$ as follows: \begin{align} \widehat{\text{MSE}}_{\tau}(\hat{\tau})&=\frac{1}{N}\sum_{j=1}^{N}(\widehat{\tau}_{j}-\tau)^2\\ &=\frac{1}{N}\sum_{j=1}^{N}(\widehat{\tau}_{j}-2)^2, \end{align} where $N$ corresponds to the number of simulation runs.

In the specific case at hand, both estimators are unbiased. Thus, your estimator $\text{MSE}_{\hat{\tau}}$ would consistently estimate the MSE's of $\tilde{\tau}$ and $\hat{\tau}$. First consider $\tilde{\tau}=\bar{X}=\frac{1}{n}\sum_{i=1}^{n}X_i$. We have \begin{align} \mathbb{E}(\bar{X})&=\frac{1}{n}\sum_{i=1}^{n}\mathbb{E}(X_i)=\mathbb{E}(X_1)=2\qquad\text{and}\\ \text{Var}(\bar{X})&=\frac{1}{n^2}\sum_{i=1}^{n}\text{Var}(X_i)=\frac{1}{n}=\frac{(b-a)^2}{12}\frac{1}{n}=\frac{1}{30} \end{align} the same you computed above. Note that I've assumed $X_1,\ldots,X_n$ to be independent since the exercise says nothing about their dependence structure.

For $\hat{\tau}$, we need to know $\mathbb{E}(\hat{a})$ and $\mathbb{E}(\hat{b})$. See this question on how to compute $\mathbb{E}(\hat{a})$ and this for both formulas. We obtain \begin{equation} \mathbb{E}(\hat{\tau})=\frac{1}{2}\left(\mathbb{E}(\hat{a})+\mathbb{E}(\hat{b})\right)=\frac{1}{2}\left(\frac{b+na}{n+1}+\frac{a+bn}{n+1}\right)=\frac{1}{2}(a+b)=\tau, \end{equation} thus establishing unbiasedness. Since the computation of $\text{Var}(\hat{\tau})$ is more involved, we resort to a simulation study.

#set seed for reproducibility
set.seed(124)
sim <- function(a,b,n,N){
  a <- a
  b <- b
  tau <- (1/2)*(a+b)
  tauhat <- rep(0,N)
  tautilde <- rep(0,N)

  for(j in seq_along(tauhat)){
    samp <- runif(n,min=a,max=b)
    tauhat[j] <- (min(samp)+max(samp))/2
    tautilde[j] <- mean(samp)
  }

  MSEtauhat <- (1/N)*sum((tauhat-tau)^2)
  MSEtautilde <- (1/N)*sum((tautilde-tau)^2)
  vartautilde <-  (1/(N-1))*sum((tautilde-mean(tautilde))^2)
  vartauhat <- (1/(N-1))*sum((tauhat-mean(tauhat))^2)
  results <- list(MSEhat=MSEtauhat,Esthat=tauhat,Varhat=vartauhat,
                  MSEtilde=MSEtautilde,Esttilde=tautilde,Vartilde=vartautilde)
  return(results)
}
res <- sim(1,3,10,10e5)
round(c(res$MSEhat,res$MSEtilde),8)
[1] 0.01516704 0.03331627
round(c(res$Varhat,res$Vartilde),8)
[1] 0.01516705 0.03331624

Obtaining almost the same numbers for $\widehat{\text{Var}}$ and $\widehat{\text{MSE}}$ confirms unbiasedness of both estimators $\hat{\tau}$ and $\tilde{\tau}$. Furthermore, the MSE of $\tilde{\tau}$ is the smaller one.

$\endgroup$

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.