Problem 1

Say \(Y\) is bivariate normal with mean vector \(\pmb{\mu}=\begin{pmatrix} 1 \\ 2\end{pmatrix}\) and covariance matrix \(\pmb{\Sigma}=\begin{pmatrix} 2 & 2 \\ 2 & 4\end{pmatrix}\). Find \(P(X<Y)\) and \(P\left(X<1\vert Y=0.5\right)\).

First note that \(P(X<Y)=P(X-Y<0)\). Now E

\[ \begin{aligned} &E[X-Y] = E[X]-E[Y]=1-2=-1\\ &var(X-Y) =var(X)+var(Y)-2cov(X,Y)=2+4-2\times 2=2 \end{aligned} \] and so \(X-Y\sim N(-1,1)\), and \(P(X<Y)=0.76\)

pnorm(0, -1, sqrt(2))
## [1] 0.7602499

By (5.2.13) \(X\vert Y=0.5\) has a normal distribution with mean

\[\theta=\mu_x+\Sigma_{xy}\Sigma_{yy}^{-1}(y-\mu_y)=1+2(4)^{-1}(0.5-2)=0.25\] and variance

\[\tau^2=\Sigma_{xx}-\Sigma_{xy}\Sigma_{yy}^{-1}\Sigma_{yx}=2-2(4)^{-1}2=1\]

and so \(P\left(X<1\vert Y=0.5\right)=0.77\)

pnorm(1, 0.25, 1)
## [1] 0.7733726

Let’s check:

set.seed(1111)
library(mvtnorm)
x=rmvnorm(1e6, mean=c(1, 2), sigma=matrix(c(2,2,2,4), 2, 2))
round(sum(x[,1]<x[,2])/1e6, 3)
## [1] 0.761
z=x[abs(x[,2]-0.5)<0.01, 1]
c(mean(z), var(z), length(z))
## [1]    0.2746113    1.0209370 3114.0000000
round(sum(z<1)/length(z), 3)
## [1] 0.773

Problem 2

  1. Say \(\pmb{X}\sim N(\begin{pmatrix}1\\1\end{pmatrix}, \begin{pmatrix}2&-1\\-1& 2\end{pmatrix}\). Find \(var(X_1^2+X_2^2)\).

Note

\[x^2+y^2 = \begin{pmatrix}x&y\end{pmatrix}\pmb{I}\begin{pmatrix}x\\y\end{pmatrix}\]

so by theorem (5.3.8)

\[var(X_1^2+X_2^2)=var(\pmb{X'IX})=2tr(\pmb{\Sigma}^2)+4\pmb{\mu'\Sigma\mu}\]

\[ \begin{aligned} &\pmb{\Sigma}^2 = \begin{pmatrix}2&-1\\-1& 2\end{pmatrix} \begin{pmatrix}2&-1\\-1& 2\end{pmatrix}= \begin{pmatrix}5&-4\\-4& 5\end{pmatrix}\\ &2tr(\pmb{\Sigma}^2) = 20\\ &\pmb{\mu'\Sigma\mu} = \begin{pmatrix}1&1\end{pmatrix}\begin{pmatrix}2&-1\\-1& 2\end{pmatrix}\begin{pmatrix}1\\1\end{pmatrix}=\\ &\begin{pmatrix}1&1\end{pmatrix}\begin{pmatrix}1\\1\end{pmatrix}=2\\ &var(X_1^2+X_2^2)=20+4\times 2=28 \end{aligned} \]

Let’s check:

library(mvtnorm)
x=rmvnorm(1e5, c(1,1), matrix(c(2,-1,-1,2), 2, 2))
round(var(x[,1]^2+x[,2]^2), 1)
## [1] 28.4
  1. \(X=\begin{pmatrix}X_1&...&X_n\end{pmatrix}\sim N(\pmb{0}, \pmb{\Sigma})\) with \(var(X_i)=2\) and \(cov(X_i,X_j)=1\) if \(\vert i-j\vert=1\), 0 otherwise. Find \(var(\sum_{i=1}^n X_i^2)\).

\[ \pmb{\Sigma} = \begin{pmatrix} 2 & 1& 0 & 0 &..&0 \\ 1 & 2& 1 & 0&..&0 \\ 0 & 1& 2 & 1 &..&0 \\ \vdots & &&& &\vdots \\ 0 & 0 & 0 &..& 1&2 \\ \end{pmatrix} \] now

\[ \pmb{\Sigma}^2 = \begin{pmatrix} 5 & 4& 1 & 0 &0&..&0 \\ 4 & 6& 4 & 1 &0 &..&0 \\ 1 & 4& 6 & 4 & 1&..&0 \\ \vdots & &&& &\vdots \\ 0 & 0 & ... &1&4& 6&4 \\ 0 & 0 & ... & 0&1& 4&5 \\ \end{pmatrix}\\ tr(\pmb{\Sigma}^2) = 2\times5+(n-2)6 = 6n-2 \]

and so

\[var(\sum_{i=1}^n X_i^2)=2tr(\pmb{\Sigma}^2)+4\pmb{\mu'\Sigma\mu}=12n-4\] Let’s check:

A=function(n) {
  A=2*diag(n)
  for(i in 2:n) A[i-1,i]=1
  for(i in 2:n) A[i,i-1]=1
  A
}
for(n in 3:7) {
  x=rmvnorm(1e5, sigma=A(n))
  print(c(n,round(var(apply(x^2,1,sum)), 1), 12*n-4))
}
## [1]  3.0 32.1 32.0
## [1]  4.0 43.8 44.0
## [1]  5 56 56
## [1]  6.0 67.9 68.0
## [1]  7.0 79.8 80.0

Problem 3

Say \(\pmb{X}\sim N(\begin{pmatrix}1\\1\end{pmatrix}, \begin{pmatrix}2&-1\\-1& 2\end{pmatrix}\). Find the \(cov(X_1+X_2, X_1^2+X_2^2)\).

By theorem (5.3.12) we have

\[cov(\pmb{BX,X'AX})=2\pmb{B\Sigma A\mu}\]

here we have \(\pmb{B}=\begin{pmatrix}1&1\end{pmatrix}\) and \(\pmb{A}=\pmb{I}\), so

\[cov(\pmb{BX,X'AX})=2\pmb{B\Sigma A\mu}=\\ 2\begin{pmatrix}1&1\end{pmatrix}\begin{pmatrix}2&-1\\-1& 2\end{pmatrix}\begin{pmatrix}1\\1\end{pmatrix}=\\ 2\begin{pmatrix}1&1\end{pmatrix}\begin{pmatrix}1\\1\end{pmatrix}=4\]

library(mvtnorm)
x=rmvnorm(1e5, c(1,1), matrix(c(2,-1,-1,2), 2, 2))
round(cov(x[,1]+x[,2],x[,1]^2+x[,2]^2), 1)
## [1] 4

Problem 4

Say we have a \(\pmb{y}=(y_1\text{ .. }y_n)'\) sample from \(N(\mu,\sigma^2)\). We want to test

\[H_0:\mu=\mu_0\text{ vs }H_a:\mu\ne\mu_0\] The usual test uses the test statistic \(T=\sqrt{n}\frac{\bar{y}-\mu_0}{s}\) and rejects \(H_0\) if \(|T|>t_{\alpha/2,n-1}\). Find the power of this test if in fact \(\mu=\mu_1\). As a numeric example use \(\mu_0=1,\mu_1=1.5,\sigma=2,n=100,\alpha=0.05\).

  1. directly

  2. using the noncentral t distribution

\[ \begin{aligned} &P(|T|>t_{\alpha/2,n-1}\vert\mu=\mu_1) = \\ &P(|\sqrt{n}\frac{\bar{y}-\mu_0}{s}|>t_{\alpha/2,n-1}\vert\mu=\mu_1)=\\ &1-P(-t_{\alpha/2,n-1}<\sqrt{n}\frac{\bar{y}-\mu_0}{s}<t_{\alpha/2,n-1}\vert\mu=\mu_1)=\\ &1-P(-t_{\alpha/2,n-1}<\sqrt{n}\frac{\bar{y}-\mu_1+(\mu_1-\mu_0)}{s}<t_{\alpha/2,n-1}\vert\mu=\mu_1)=\\ &1-P(-t_{\alpha/2,n-1}<\sqrt{n}\frac{\bar{y}-\mu_1}{s}+\sqrt{n}\frac{\mu_1-\mu_0}{s}<t_{\alpha/2,n-1}\vert\mu=\mu_1)=\\ &1-P(-t_{\alpha/2,n-1}-\sqrt{n}\frac{\mu_1-\mu_0}{s}<\sqrt{n}\frac{\bar{y}-\mu_1}{s}<t_{\alpha/2,n-1}-\sqrt{n}\frac{\mu_1-\mu_0}{s}\vert\mu=\mu_1)=\\ &1-\left(pt(t_{\alpha/2,n-1}-\sqrt{n}\frac{\mu_1-\mu_0}{s})-pt(-t_{\alpha/2,n-1}-\sqrt{n}\frac{\mu_1-\mu_0}{s})\right) \end{aligned} \]

mu0=1;mu1=1.5;sigma=2;n=100;alpha=0.05
1-(pt(qt(1-alpha/2,n-1)-sqrt(n)*(mu1-mu0)/sigma,n-1)-pt(-qt(1-alpha/2,n-1)-sqrt(n)*(mu1-mu0)/sigma,n-1))
## [1] 0.6964319

Under the alternative distribution \(y_i/\sigma\sim N(\mu_1/\sigma,1)\), and so

\(\sqrt{n}\bar{y}\sim N(\sqrt{n}(\mu_1-\mu_0)/\sigma,1)\) and so \(\sqrt{n}\bar{y}/s\) has a non-central t distribution with n-1 degrees of freedom and noncentrality parameter \(\delta=\sqrt{n}(\mu_1-\mu_0)/\sigma\). Therefore

\[ \begin{aligned} &P(|T|>t_{\alpha/2,n-1}\vert\mu=\mu_1) = \\ &1-\left(pt(t_{\alpha/2,n-1},n-1,\delta)-pt(-t_{\alpha/2,n-1},n-1,\delta)\right) \end{aligned} \]

mu0=1;mu1=1.5;sigma=2;n=100;alpha=0.05
delta=sqrt(n)*(mu1-mu0)/sigma
delta
## [1] 2.5
1-(pt(qt(1-alpha/2,n-1),n-1,delta)-pt(-qt(1-alpha/2,n-1),n-1, delta))
## [1] 0.6969803

Let’s check with a simulation:

n=100;mu0=1;mu1=1.5;sigma=2;alpha=0.05;B=10000
crit=qt(1-alpha/2, n-1)
counter=0
for(i in 1:B) {
  x=rnorm(n,mu1,sigma)
  TS=sqrt(n)*(mean(x)-mu0)/sd(x)
  if(abs(TS)>crit) counter=counter+1
}
counter/B
## [1] 0.6979