(This article was originally published at Analysis with Programming, and syndicated at StatsBlogs.)
Problem
- Suppose that we have two independent random samples $X_1,\cdots, X_n$ are exponential($\theta$), and $Y_1,\cdots, Y_m$ are exponential($\mu$).
- Find the LRT (Likelihood Ratio Test) of $H_0:\theta=\mu$ versus $H_1:\theta\neq\mu$.
- Show that the test in part (a) can be based on the statistic $$ T=\frac{\sum X_i}{\sum X_i+\sum Y_i}. $$
- Find the distribution of $T$ when $H_0$ is true.
Solution
- The Likelihood Ratio Test is given by $$ \lambda(\mathbf{x},\mathbf{y}) = \frac{\displaystyle\sup_{\theta = \mu,\mu>0}\mathrm{P}(\mathbf{x},\mathbf{y}|\theta,\mu)}{\displaystyle\sup_{\theta > 0,\mu>0}\mathrm{P}(\mathbf{x}, \mathbf{y}|\theta,\mu)}, $$ where the denominator is evaluated as follows: $$ \sup_{\theta > 0,\mu>0}\mathrm{P}(\mathbf{x}, \mathbf{y}|\theta,\mu)= \sup_{\theta > 0}\mathrm{P}(\mathbf{x}|\theta)\sup_{\mu > 0}\mathrm{P}(\mathbf{y}|\mu),\quad\text{by independence.} $$ So that, $$ \begin{aligned} \sup_{\theta > 0}\mathrm{P}(\mathbf{x}|\theta)&=\sup_{\theta>0}\prod_{i=1}^{n}\frac{1}{\theta}\exp\left[-\frac{x_i}{\theta}\right]=\sup_{\theta>0}\frac{1}{\theta^n}\exp\left[-\frac{\sum_{i=1}^{n}x_i}{\theta}\right]\\ &=\frac{1}{\bar{x}^n}\exp\left[-\frac{\sum_{i=1}^{n}x_i}{\bar{x}}\right]=\frac{1}{\bar{x}^n}\exp[-n], \end{aligned} $$ since $\bar{x}$, or the sample mean is the MLE of $\theta$. Also, $$ \begin{aligned} \sup_{\mu > 0}\mathrm{P}(\mathbf{y}|\mu)&=\sup_{\mu>0}\prod_{j=1}^{m}\frac{1}{\mu}\exp\left[-\frac{y_j}{\mu}\right]=\sup_{\mu>0}\frac{1}{\mu^m}\exp\left[-\frac{\sum_{j=1}^{m}y_j}{\mu}\right]\\ &=\frac{1}{\bar{y}^m}\exp\left[-\frac{\sum_{j=1}^{m}y_j}{\bar{y}}\right]=\frac{1}{\bar{y}^m}\exp[-m]. \end{aligned} $$ Now the numerator is evaluated as follows, $$ \begin{aligned} \sup_{\theta = \mu,\mu>0}\mathrm{P}(\mathbf{x},\mathbf{y}|\theta,\mu)&=\sup_{\theta=\mu,\mu>0}\mathrm{P}(\mathbf{x}|\theta)\mathrm{P}(\mathbf{y}|\mu),\quad\text{by independence.}\\ &=\sup_{\theta=\mu,\mu>0}\prod_{i=1}^{n}\frac{1}{\theta}\exp\left[-\frac{x_i}{\theta}\right]\prod_{j=1}^{m}\frac{1}{\mu}\exp\left[-\frac{y_j}{\mu}\right]\\ &=\sup_{\theta=\mu,\mu>0}\frac{1}{\theta^n}\exp\left[-\frac{\sum_{i=1}^nx_i}{\theta}\right]\frac{1}{\mu^m}\exp\left[-\frac{\sum_{j=1}^m y_j}{\mu}\right]\\ &=\sup_{\mu>0}\frac{1}{\mu^n}\exp\left[-\frac{\sum_{i=1}^nx_i}{\mu}\right]\frac{1}{\mu^m}\exp\left[-\frac{\sum_{j=1}^m y_j}{\mu}\right]\\ &=\sup_{\mu>0}\frac{1}{\mu^{n+m}}\exp\left\{-\frac{1}{\mu}\left[\sum_{i=1}^nx_i+\sum_{j=1}^m y_j\right]\right\} \end{aligned} $$ Note that $\mu$ is a nuisance parameter, and so we will also maximize this over its domain. And to do that we take the log-likeihood function first, $$ \begin{aligned} \ell(\mu|\mathbf{x},\mathbf{y})&=-\log(\mu^{n+m})-\frac{1}{\mu}\left[\sum_{i=1}^nx_i+\sum_{j=1}^m y_j\right]\\ &=-(n+m)\log(\mu)-\frac{1}{\mu}\left[\sum_{i=1}^nx_i+\sum_{j=1}^m y_j\right]. \end{aligned} $$ Taking the derivative with respect to $\mu$, gives us $$ \frac{\operatorname{d}}{\operatorname{d}\mu}\ell(\mu|\mathbf{x},\mathbf{y})=-(n+m)\frac{1}{\mu}+\frac{1}{\mu^2}\left[\sum_{i=1}^nx_i+\sum_{j=1}^m y_j\right], $$ equate this to zero to obtain the stationary point, $$ \begin{aligned} -(n+m)\frac{1}{\mu}+\frac{1}{\mu^2}\left[\sum_{i=1}^nx_i+\sum_{j=1}^m y_j\right]&=0\\ -(n+m)\mu+\left[\sum_{i=1}^nx_i+\sum_{j=1}^m y_j\right]&=0\\ \mu&=\frac{1}{n+m}\left[\sum_{i=1}^nx_i+\sum_{j=1}^m y_j\right]. \end{aligned} $$ To verify if this is the MLE, we take the second derivative test for the log-likelihood function, $$ \frac{\operatorname{d}^2}{\operatorname{d}\mu^2}\ell(\mu|\mathbf{x},\mathbf{y})=(n+m)\frac{1}{\mu^2}-\frac{2}{\mu^3}\left[\sum_{i=1}^nx_i+\sum_{j=1}^m y_j\right]<0, $$ since $\frac{1}{\mu^2}<\frac{2}{\mu^3}$, implying $\hat{\mu}=\displaystyle\frac{1}{n+m}\left[\sum_{i=1}^nx_i+\sum_{j=1}^m y_j\right]$ is the MLE of $\mu$. Thus the LRT, $\lambda(\mathbf{x},\mathbf{y})$ would be, $$ \begin{aligned} \lambda(\mathbf{x},\mathbf{y})&=\frac{\sup_{\mu>0}\displaystyle\frac{1}{\mu^{n+m}}\exp\left\{-\frac{1}{\mu}\left[\sum_{i=1}^nx_i+\sum_{j=1}^m y_j\right]\right\}}{\displaystyle\frac{1}{\bar{x}^n}\frac{1}{\bar{y}^m}\exp[-(n+m)]}\\ &=\left(\frac{1}{\frac{1}{{(n+m)}^{n+m}}\left[\sum_{i=1}^nx_i+\sum_{j=1}^m y_j\right]^{n+m}}\times\right.\\ &\qquad\left.\exp\left\{-\frac{1}{\frac{1}{n+m}\left[\sum_{i=1}^nx_i+\sum_{j=1}^m y_j\right]}\left[\sum_{i=1}^nx_i+\sum_{j=1}^m y_j\right]\right\}\right)\bigg/\\ &\qquad\qquad\qquad\displaystyle\frac{1}{\bar{x}^n}\frac{1}{\bar{y}^m}\exp[-(n+m)] \end{aligned} $$ $$ \begin{aligned} &=\frac{\displaystyle\frac{1}{\displaystyle\frac{1}{{(n+m)}^{n+m}}\left[\displaystyle\sum_{i=1}^nx_i+\sum_{j=1}^m y_j\right]^{n+m}}\times\exp[-(n+m)]}{\displaystyle\frac{1}{\bar{x}^n}\frac{1}{\bar{y}^m}\exp[-(n+m)]}\\[.3cm] &=\frac{\displaystyle \bar{x}^n \bar{y}^m}{\displaystyle\frac{1}{{(n+m)}^{n+m}}\left[\sum_{i=1}^nx_i+\sum_{j=1}^m y_j\right]^{n+m}}. \end{aligned} $$ And we say that $H_0$ is rejected if $\lambda(\mathbf{x},\mathbf{y})\leq c$.
- If we do some algebra on the LRT in part (a), we obtain the following: $$ \begin{aligned} \lambda(\mathbf{x},\mathbf{y})&=\frac{\displaystyle \bar{x}^n \bar{y}^m}{\displaystyle\frac{1}{{(n+m)}^{n+m}}\left[\sum_{i=1}^nx_i+\sum_{j=1}^m y_j\right]^{n+m}}\\ &=\frac{\displaystyle\frac{1}{n^n}\left(\sum_{i=1}^{n}x_i\right)^{n}\frac{1}{m^{m}}\left(\sum_{j=1}^{m}y_j\right)^{m}}{\displaystyle\frac{1}{(n+m)^{n+m}}\left[\sum_{i=1}^{n}x_i+\sum_{j=1}^{m}y_j\right]^{n+m}}\\ &=\frac{\displaystyle (n+m)^{n+m}\left(\sum_{i=1}^{n}x_i\right)^{n}\left(\sum_{j=1}^{m}y_j\right)^{m}}{\displaystyle n^{n}m^{m}\left[\sum_{i=1}^{n}x_i+\sum_{j=1}^{m}y_j\right]^{n+m}}\\ &=\frac{(n+m)^{n+m}}{n^nm^{m}}\left[\frac{\displaystyle \sum_{j=1}^{m}y_j}{\displaystyle\sum_{i=1}^{n}x_i+\sum_{j=1}^{m}y_j}\right]^{m}\left[\frac{\displaystyle \sum_{i=1}^{n}x_i}{\displaystyle\sum_{i=1}^{n}x_i+\sum_{j=1}^{m}y_j}\right]^{n}\\ &=\frac{(n+m)^{n+m}}{n^nm^{m}}\left[1-\frac{\displaystyle \sum_{i=1}^{n}x_i}{\displaystyle\sum_{i=1}^{n}x_i+\sum_{j=1}^{m}y_j}\right]^{m}\left[\frac{\displaystyle \sum_{i=1}^{n}x_i}{\displaystyle\sum_{i=1}^{n}x_i+\sum_{j=1}^{m}y_j}\right]^{n}\\ &=\frac{(n+m)^{n+m}}{n^nm^{m}}\left[1-T\right]^{m}\left[T\right]^{n}. \end{aligned} $$ Hence, the LRT can be based on the statistic $T$.
- The distribution of $\sum X_i$ is obtain using the MGF (Moment Generating Function) technique, that is $$ \begin{aligned} \mathrm{M}_{\Sigma X_i}(t)&=\mathrm{E}\exp[t\Sigma X_i]=\mathrm{E}\exp[tX_1 +\cdots + tX_n]\\ &=\mathrm{E}\exp[tX_1]\times\cdots\times\mathrm{E}\exp[tX_n],\quad\text{by independence.}\\ &=\frac{1}{1-\theta t}\times\cdots\times\frac{1}{1-\theta t}\\ &=\left(\frac{1}{1-\theta t}\right)^{n}=\text{MGF of gamma}(n,\theta). \end{aligned} $$ Now, when $H_0$ is true then $\sum X_i$ is gamma($m,\theta$). For brevity, let $X=\sum_{i=1}^{n} X_i$ and $Y=\sum_{j=1}^{m}Y_j$. The joint distribution of $X$ and $Y$ is given below, $$ f_{XY}(x, y)=\frac{1}{\Gamma (n)\theta^{n}}x^{n-1}\exp[-x/\theta]\times\frac{1}{\Gamma (m)\theta^{m}}y^{m-1}\exp[-y/\theta]. $$ Let $U=\frac{X}{X+Y}$ and $V=X+Y$, then the support of $(X,Y)$ is $\mathcal{A}=\left\{(x,y)\in \mathbb{R}^{+}\times \mathbb{R}^{+}\right\}$. Since the transformations $U$ and $V$ is one-to-one and onto, then $\mathcal{B}=\left\{(u,v)\in [0,1]\times \mathbb{R}^{+}\right\}$. Consider the following transformations $$ u=g_{1}(x,y)=\frac{x}{x+y}\quad\text{and}\quad v=g_{2}(x,y)=x+y. $$ Then, \begin{equation} \label{eq:bvt1} u=\frac{x}{x+y}\Rightarrow x=\frac{uy}{1-u} \end{equation} and \begin{equation} \label{eq:bvt2} v=x+y\Rightarrow y = v-x. \end{equation} Substitute Equation (\ref{eq:bvt2}) to Equation (\ref{eq:bvt1}), then $$ \begin{aligned} x=\frac{u(v-x)}{1-u}&\Rightarrow x(1-u)=u(v-x)\\ x-ux=uv-ux&\Rightarrow x=uv=h_{1}(u,v). \end{aligned} $$ Substitute $x$ above to Equation (\ref{eq:bvt2}) to obtain, $$y=v(1-u)=h_2(u,v).$$ And the Jacobian matrix is, $$ \mathbf{J}=\bigg| \begin{array}{cc} v&u\\[.2cm] -v&1-u \end{array} \bigg|=v(1-u)+uv=v. $$ So that, $$ \begin{aligned} f_{UV}(u,v)&=f_{XY}(h_1(u,v),h_2(u,v))\lvert \mathbf{J}\rvert=f_{XY}(uv,v(1-u))\lvert v\rvert\\ &=\frac{1}{\Gamma (n)\theta^{n}}(uv)^{n-1}\exp[-uv/\theta]\times\\ &\quad\;\frac{1}{\Gamma (m)\theta^{m}}(v(1-u))^{m-1}\exp[-v(1-u)/\theta]v\\ &=\frac{1}{\Gamma (n)\theta^{n}}(uv)^{n-1}\exp[-uv/\theta]\times\\ &\quad\;\frac{1}{\Gamma (m)\theta^{m}}(v(1-u))^{m-1}\exp[-v/\theta]\exp[uv/\theta]v\\ &=\frac{1}{\Gamma (n)\theta^{n}}u^{n-1}v^{n-1}\times\frac{1}{\Gamma (m)\theta^{m}}v^{m-1}(1-u)^{m-1}\exp[-v/\theta]v\\ &=\frac{1}{\Gamma (n)}\underbrace{u^{n-1}(1-u)^{m-1}}_{\text{Beta}(n,m)\text{ kernel}}\frac{1}{\Gamma (m)\theta^{m+n}}v^{m-1}v^{n-1}\exp[-v/\theta]v\\ &=\frac{\Gamma(m)\Gamma(m+n)}{\Gamma(m)\Gamma(m+n)}\frac{u^{n-1}(1-u)^{m-1}}{\Gamma (n)}\times\\ &\quad\;\frac{1}{\Gamma (m)\theta^{m+n}}v^{m-1}v^{n}\exp[-v/\theta]\\ &=\underbrace{\frac{\Gamma(m+n)}{\Gamma (n)\Gamma(m)}u^{n-1}(1-u)^{m-1}}_{\text{Beta}(n,m)}\times\\ &\quad\;\underbrace{\frac{1}{\Gamma(m+n)\theta^{m+n}}v^{m+n-1}\exp[-v/\theta]}_{\text{Gamma}(m+n,\theta)}. \end{aligned} $$ So that the marginal density of $U=\displaystyle\frac{\sum X_i}{\sum X_i +\sum Y_i}$ is Beta($n,m$).
Reference
Please comment on the article here: Analysis with Programming
