How does one do a Wald test on estimates from two variables?

Question

Given a dataset with two variables $X$ and $Y$, with each observation independent of the others, test the null hypothesis $$\mu_X = \mu_Y\\ \sigma^2_X = \sigma^2_Y$$ using a Wald test.

This question comes from a textbook, but my answer isn't going to be handed in. Hence, I'm posting this question here because I'd like my work to be checked. Our class was given a dataset to use for this question, and I use Stata. I'm not sure of the best way to solve this question—the main struggle I had was with figuring out how to find the variance-covariance matrix for the Wald statistic, as I didn't think it was possible to create a MLE random variable for $\widehat{\mu_X} - \widehat{\mu_Y}$ and $\widehat{\sigma_X^2} - \widehat{\sigma^2_Y}$, in order to find their variances and covariances.

To implement the Wald test, I construct the vector of restrictions such that it is predicted to equal $\boldsymbol 0$. In this case, that is $$\left[\begin{matrix} {{\mu_X - \mu_Y}} \\ \sigma^2_X - \sigma^2_Y \end{matrix}\right] = \left[\begin{matrix} 0 \\ 0 \end{matrix}\right] $$ For notational convenience, label this $\boldsymbol v$. Then $\boldsymbol{{\hat v'}}\left(\Sigma_{\boldsymbol{{\hat v}}}\right)^{-1}\boldsymbol{{\hat v}}$ is the Wald Statistic. The covariance matrix is \begin{align*} \Sigma_{\boldsymbol{{\hat v}}} &= (\boldsymbol{{\hat v}} - \mathrm E[\boldsymbol{{\hat v}}])(\boldsymbol{{\hat v}} - \mathrm E[\boldsymbol{{\hat v}}])' \\ &=\left[\begin{matrix} \mathrm{Var}(\widehat{\mu_X} - \widehat{\mu_Y}) & \mathrm{Cov}\left(\widehat{\mu_X} - \widehat{\mu_Y}, \widehat{\sigma^2_X} - \widehat{\sigma^2_Y}\right) \\ \mathrm{Cov}\left(\widehat{\mu_X} - \widehat{\mu_Y}, \widehat{\sigma^2_X} - \widehat{\sigma^2_Y}\right) & \mathrm{Var}\left(\widehat{\sigma^2_X} - \widehat{\sigma^2_Y}\right) \end{matrix}\right] \end{align*}

Let's look at each element of the matrix in turn: \begin{align*} &&\mathrm{Var}\left(\widehat{\mu_X} - \widehat{\mu_Y}\right) &= \mathrm{Var}\left(\frac 1n\sum_{i=1}^n x_i - \frac 1n\sum_{i=1}^n y_i\right) \\ && &= \mathrm{Var}\left(\frac 1n\left(\sum_{i=1}^n(x_i - y_i)\right)\right) \\ && &=\frac1{n^2}\left(\sum_{i=1}^n \mathrm{Var}(x_i) + \sum_{i=1}^n\mathrm{Var}(y_i) - \sum_{i=1}^n 2\mathrm{Cov}(x_i, y_i)\right) \\ && &= \frac1n\left(\sigma^2_X + \sigma^2_Y - 2\mathrm{Cov}(X, Y)\right) \end{align*} where the observations (individuals) are independent of each other, so $\mathrm{Cov}(x_i, y_j) = 0 $ for $i\neq j$. This scalar will be computed in Stata. Easy enough, but trying to use this same method of deriving formulae for the other entries in terms of $\mu_X, \mu_Y, \sigma_X^2$ and $\sigma_Y^2$ seems to be difficult. (In order to do that, I would use the expected value formula for a variance $\mathrm{Var}(W) = \mathrm E[W^2] - \mathrm E[W]^2$, and substitute $x_i = \mu_X + \varepsilon_{X,i}$ (also with $y_i$) to use that fact that $\mathrm E[\varepsilon_{X,i}^2] = \sigma_X^2$ (also with $\varepsilon_{Y, i}^2$).)

I chose a simpler method of generating variables for $(X-\mu_X)^2$ and $(Y-\mu_Y)^2$ in Stata, in order to compute the covariances according to the formulae written below.

The second and third entries are $$% \mathrm{Cov}\left(\widehat{\mu_X} - \widehat{\mu_Y}, \widehat{\sigma^2_X} - \widehat{\sigma^2_Y}\right) = \mathrm{Cov}\left(\widehat{\mu_X}, \widehat{\sigma^2_X}\right) - \mathrm{Cov}\left(\widehat{\mu_Y}, \widehat{\sigma^2_X}\right) -\mathrm{Cov}\left(\widehat{\mu_X}, \widehat{\sigma^2_Y}\right) + \mathrm{Cov}\left(\widehat{\mu_Y}, \widehat{\sigma^2_Y}\right)% $$ Where \begin{align*} \mathrm{Cov}\left(\widehat{\mu_X}, \widehat{\sigma^2_X}\right) &= \mathrm{Cov}\left(\frac 1n\sum_{i=1}^n x_i, \frac 1n\sum_{i=1}^n (x_i - \mu)^2\right) \\ &= \frac 1{n^2}\sum_{i=1}^n\sum_{j=1}^n\mathrm{Cov}\left(x_i, (x_j - \mu_X)^2\right)\\ &\text{where $x_i$ and $(x_j - \mu)^2$ are not correlated when $i\neq j$} \\ &= \frac 1n\mathrm{Cov}\left(X, (X - \mu_X)^2\right) \end{align*} Similarly, \begin{multline*} \mathrm{Cov}\left(\widehat{\mu_X} - \widehat{\mu_Y}, \widehat{\sigma^2_X} - \widehat{\sigma^2_Y}\right) \\ =\frac 1n\bigg(\mathrm{Cov}\left(X, (X - \mu_X)^2\right) - \mathrm{Cov}\left(Y, (X-\mu_X)^2\right) - \mathrm{Cov}\left(X, (Y-\mu_Y)^2\right) + \mathrm{Cov}\left(Y, (Y - \mu_Y)^2\right)\bigg) \end{multline*}

Lastly, $\mathrm{Var}\left(\widehat{\sigma^2_X} - \widehat{\sigma^2_Y}\right)$ \begin{align*} &= \mathrm{Var}\left(\frac 1n\left(\sum_{i=1}^n(x_i - \mu_X)^2 - \sum_{i=1}^n(y_i - \mu_Y)^2\right)\right) \\ &= \frac 1{n^2}\left(\sum_{i=1}^n\mathrm{Var}\left((x_i - \mu_X)^2\right) + \sum_{i=1}^n\mathrm{Var}\left((y_i - \mu_Y)^2\right) - \sum_{i=1}^n\mathrm{Cov}\Big((x_i - \mu_X)^2, (y_i - \mu_Y)^2\Big)\right) \\ &= \frac 1n\bigg(\mathrm{Var}\left((X - \mu_X)^2\right) + \mathrm{Var}\left((Y - \mu_Y)^2\right) - \mathrm{Cov}\Big((X - \mu_X)^2, (Y - \mu_Y)^2\Big)\bigg) \end{align*}

For each of these formulae for the entries of the matrix, substituting $\mu_X = \widehat{\mu_X},\ \mu_Y = \widehat{\mu_Y},\ \sigma_X^2 = \widehat{\sigma^2_X},$ and $\sigma_Y^2 = \widehat{\sigma^2_Y}$ is asymptotically valid. Hence, my Stata code is this.

I thought of using the Fisher Information matrix $I(\boldsymbol{{\hat v}}) \approx \Sigma_{\boldsymbol{{\hat v}}}^{-1}$ according to the Cramer-Rao Lower Bound, but I don't know how to differentiate the bivariate loglikelihood w.r.t. an expression (i.e. $\mu_X - \mu_Y$). I have found the bivariate loglikelihood function under the assumption of normality.

Answer 1

The OP seeks:

1. $\mathrm{Cov}\left(\widehat{\mu_X} - \widehat{\mu_Y}, \widehat{\sigma^2_X} - \widehat{\sigma^2_Y}\right)$

2. $\mathrm{Var}\left(\widehat{\sigma^2_X} - \widehat{\sigma^2_Y}\right)$

Notation

These are moments of moments problems, in particular finding population moments of sample moments, in a bivariate world. Such problems are not generally easy to solve by hand, and the modus operandi for solving them is to work with power sum notation, which in our bivariate world is of form $s_{r,t}$, namely:

$$s_{r,t}=\sum _{i=1}^n X_i^r Y_i^t$$

In this power sum notation, the sample mean estimators are:

• $\widehat{\mu_X} = \frac{s_{1,0}}{n} \quad \quad \quad \widehat{\mu_Y} = \frac{s_{0,1}}{n}$

... and the MLE sample variance estimators are:

• $\widehat{\sigma^2_X}: \quad m_{2,0} = \frac{s_{2,0}}{n}-\frac{s_{1,0}^2}{n^2} \quad \quad \text{and} \quad \quad \widehat{\sigma^2_Y}: \quad m_{0,2} = \frac{s_{0,2}}{n}-\frac{s_{0,1}^2}{n^2}$

Answer 1

Find: $\mathrm{Cov}\left(\widehat{\mu_X} - \widehat{\mu_Y}, \widehat{\sigma^2_X} - \widehat{\sigma^2_Y}\right)$

Since the covariance operator denotes the $\mu_{1,1}$ Central Moment, we can find the covariance using the mathStatica (for Mathematica) package function:

where the solution is expressed in terms of central moments of the population, namely:

• $\mu _{r,s}$ denotes the product central moment:

$$\mu _{r,s}=E\left[(X-E[X]]^r (Y-E[Y])^s\right]$$

For example, $\mu_{1,1} = \text{Cov}(X,Y)$, $\mu_{2,0}= \text{Var}(X)$ and $\mu_{0,2}= \text{Var}(Y)$.

Answer 2

Find: $\mathrm{Var}\left(\widehat{\sigma^2_X} - \widehat{\sigma^2_Y}\right)$

Since the variance operator denotes the $2^{\text{nd}}$ Central Moment, we can find the desired variance using the mathStatica package function:

More detail

There is an extensive discussion of 'moment of moment' problems in Chapter 7 of our book:

• Rose and Smith, "Mathematical Statistics with Mathematica", Springer, NY

A free download of the written chapter is available here:

http://www.mathstatica.com/book/bookcontents.html

which also includes references for further reading. You may also find of interest Chapter 13 of:

• Stuart and Ord (1994), Kendall's Advanced Theory of Statistics

... but it is very hard to solve such problems by hand.

As disclosure, I should note that I am an author of the CentralMomentToCentral function used above.