37363 - Lab/Tutorial 4 Solutions

Scott Alexander

Q1

(a)

The MGF of (Subscript[X, k]^2) is
(Subscript[M, (Subscript[X, k]^2)])(u) := 𝔼[(e^(u (Subscript[X, k]^2)))] = (∫_(-∞)^∞)(e^(u (x^2)))(1/(Sqrt[2π]))(e^(-(x^2)/2))dx = (∫_(-∞)^∞)(1/(Sqrt[2π]))(e^(-(1-2u)(x^2)/2))dx
  = (Sqrt[(1/(1-2u))])(∫_(-∞)^∞)(1/(Sqrt[2π (1/(1-2u))]))(e^(-((x^2)/(2(1/(1-2u))))))dx
  = (Sqrt[(1/(1-2u))])
as     (1/(Sqrt[2π (1/(1-2u))]))(e^(-((x^2)/(2(1/(1-2u))))))     is density of N(0,(1/(1-2u))) so long as u<(1/2).

That is,   (Subscript[M, (Subscript[X, k]^2)])(u)=(((1-2u))^(-1/2))    with    u<(1/2).

(b)

The MGF of Y is
(Subscript[M, (Subscript[Y, n])])(u) := 𝔼[(e^(u (Subscript[Y, n])))] = 𝔼[(e^(u (∑_(k=1)^n)(Subscript[X, k]^2)))] = 𝔼[(∏_(k=1)^n)(e^(u (Subscript[X, k]^2)))]
  =(∏_(k=1)^n)𝔼[(e^(u (Subscript[X, k]^2)))]    ((Subscript[X, k]) are independent from each other)
  = (∏_(k=1)^n)(Subscript[M, (Subscript[X, k]^2)])(u) = (∏_(k=1)^n)(((1-2u))^(-1/2))    (using result from (a))
  =(((1-2u))^(-n/2))    with    u<(1/2).

Calculations

Clear[n];
MomentGeneratingFunction[ChiSquareDistribution[n],u]

(c)

To apply CLT we need to write 
(Subscript[Z, n]) = c (((1/n)(∑_(k=1)^n)(Subscript[X, k]^2)-𝔼[(Subscript[X, k]^2)])/((Sqrt[var((Subscript[X, k]^2))])/(Sqrt[n])))
where c is a constant to be determined.

First note that
𝔼[(Subscript[Y, 1])] = 𝔼[(Subscript[X, 1]^2)] = var((Subscript[X, 1]^2))+((𝔼[(Subscript[X, 1])])^2) = 1+0 = 1
and
var((Subscript[Y, 1])) = var((Subscript[X, 1]^2)) = 𝔼[(Subscript[X, 1]^4)]-((𝔼[(Subscript[X, 1]^2)])^2) = 3-1 = 2.

Then
(Subscript[Z, n]) = (((Subscript[Y, n])-n)/(8(Sqrt[3n]))) = (((∑_(k=1)^n)(Subscript[X, k]^2)-n)/(8(Sqrt[3n]))) = (((1/n)((∑_(k=1)^n)(Subscript[X, k]^2)-n))/((1/n)8(Sqrt[3n])))
  = (((1/n)(∑_(k=1)^n)(Subscript[X, k]^2)-1)/(8(Sqrt[3])/(Sqrt[n]))) = ((Sqrt[2])/(8(Sqrt[3])))(((1/n)(∑_(k=1)^n)(Subscript[X, k]^2)-1)/((Sqrt[2])/(Sqrt[n])))
  = ((Sqrt[2])/(8(Sqrt[3])))(((1/n)(∑_(k=1)^n)(Subscript[X, k]^2)-𝔼[(Subscript[X, k]^2)])/((Sqrt[var((Subscript[X, k]^2))])/(Sqrt[n]))).

But by CLT  Chapter 3 Notes
(lim_(n->∞))(((1/n)(∑_(k=1)^n)(Subscript[X, k]^2)-𝔼[(Subscript[X, k]^2)])/((Sqrt[var((Subscript[X, k]^2))])/(Sqrt[n]))) ~ N(0,1)
so
(lim_(n->∞)) (Subscript[Z, n]) = (lim_(n->∞)) ((Sqrt[2])/(8(Sqrt[3])))(((1/n)(∑_(k=1)^n)(Subscript[X, k]^2)-𝔼[(Subscript[X, k]^2)])/((Sqrt[var((Subscript[X, k]^2))])/(Sqrt[n])))
  = ((Sqrt[2])/(8(Sqrt[3])))(lim_(n->∞)) (((1/n)(∑_(k=1)^n)(Subscript[X, k]^2)-𝔼[(Subscript[X, k]^2)])/((Sqrt[var((Subscript[X, k]^2))])/(Sqrt[n]))) ~ N(0,(((((Sqrt[2])/(8(Sqrt[3])))))^2)=(1/96)).

Calculations

Mean[ChiSquareDistribution[1]]

Variance[ChiSquareDistribution[1]]

Expectation[(x^2),x\[Distributed]NormalDistribution[]]

Expectation[(x^4),x\[Distributed]NormalDistribution[]]-((Expectation[(x^2),x\[Distributed]NormalDistribution[]])^2)

(((((Sqrt[2])/(8(Sqrt[3])))))^2)

(d)

From (c) we know that for large n we have the approximation
(Subscript[Z, n]) = (((Subscript[Y, n])-n)/(8(Sqrt[3n])))~N(0,(1/96))
or after re-arranging the approximation
(Subscript[Y, n]) = n+8(Sqrt[3n])(Subscript[Z, n])~N(n,2n).

To see this note that (Subscript[Y, n]) is approximately Guassian (multiple of an approximate Gaussian RV) with
𝔼[(Subscript[Y, n])] = 𝔼[n+8(Sqrt[3n])(Subscript[Z, n])] = n+8(Sqrt[3n])𝔼[(Subscript[Z, n])]  = n
and
var((Subscript[Y, n])) = var(n+8(Sqrt[3n])(Subscript[Z, n])) = (((8(Sqrt[3n])))^2)var((Subscript[Z, n])) = (((8(Sqrt[3n])))^2)(1/96) = 2n.

Using Mathematica

pApprox[n_]:=CDF[NormalDistribution[n,Sqrt[2n]],(((n^2)+1)/n)];
pExact[n_]:=CDF[ChiSquareDistribution[n],(((n^2)+1)/n)];
Table[pApprox[n]/pExact[n],{n,{(10^4),(10^5)}}]//N

Using R

See R code file.

Q2

(a)

Let X~U(1,8) which has density
(Subscript[f, X])(x) = (1/7)I(1≤x≤8) = ({{{, ({{1/7,, 1≤x≤8}, {0,, otherwise}})}}).

Setting
g(x) = 21x (e^(-x/2))
we have
G = (∫_1^8)3x (e^(-x/2))dx = (∫_(-∞)^∞)21x (e^(-x/2))(1/7)I(1≤x≤8)dx = (∫_(-∞)^∞)g(x)(Subscript[f, X])(x)dx
  = 𝔼[g(X)].

Calculations

PDF[UniformDistribution[{1,8}],x]

(b)

First note that
∫x (e^(-x/2))dx = -2∫x (d/(dx))((e^(-x/2)))dx = -2(x (e^(-x/2))-∫(e^(-x/2))(d/(dx))xdx)    (integration by parts)
  = -2(x (e^(-x/2))-∫(e^(-x/2))dx) = -2(x (e^(-x/2))+2(e^(-x/2)))
  = -2(x+2)(e^(-x/2)).

Then the expectation
G = 𝔼[g(X)] = (∫_(-∞)^∞)g(x)(Subscript[f, X])(x)dx = (∫_(-∞)^∞)21x (e^(-x/2))(1/7)I(1≤x≤8)dx
  = (∫_1^8)3x (e^(-x/2))dx = 3(Subscript[(-2)[(x+2)(e^(-x/2))], 1]^8)    (using result above)
  = -6(10(e^(-4))-3(e^(-1/2)))
  = -60(e^(-4))+18(e^(-1/2))
  ≈9.81861354150335`.

Calculations

∫3x (E^(-x/2))\[DifferentialD]x//FullSimplify

(∫_1^8)3x (E^(-x/2))\[DifferentialD]x//Expand

Expectation[21x (E^(-x/2)),x\[Distributed]UniformDistribution[{1,8}]]//Expand

(c)

First note that
∫(x^2) (e^(-x))dx = -∫(x^2) (d/(dx))((e^(-x)))dx = -((x^2) (e^(-x))-∫(e^(-x))(d/(dx))(x^2)dx)    (integration by parts)
  = -((x^2) (e^(-x))-2∫x (e^(-x))dx) = -((x^2) (e^(-x))+2∫x (d/(dx))(e^(-x))dx)
  = -((x^2) (e^(-x))+2(x (e^(-x))-∫(e^(-x)) (d/(dx))xdx))    (integration by parts)
  = -((x^2) (e^(-x))+2(x (e^(-x))-∫(e^(-x)) dx)) = -((x^2) (e^(-x))+2(x (e^(-x))+(e^(-x))))
  = -((x^2) +2x +2)(e^(-x)).

Then the second moment
𝔼[(g^2)(X)] = (∫_(-∞)^∞)(g^2)(x)(Subscript[f, X])(x)dx = (∫_(-∞)^∞)(((21x (e^(-x/2))))^2)(1/7)I(1≤x≤8)dx
  = (∫_1^8)63(x^2) (e^(-x))dx = -(Subscript[63[((x^2) +2x +2)(e^(-x))], 1]^8)    (using result above)
  = -63(64 +16 +2)(e^(-8))+63(1 +2 +2)(e^(-1))
  = -5166(e^(-8))+315(e^(-1)).

From this weobtain the variance as
(σ^2)(g) = var(g(X)) = 𝔼[(g^2)(X)]-((𝔼[g(X)])^2)
  = -5166(e^(-8))+315(e^(-1))-(((-60(e^(-4))+18(e^(-1/2))))^2)
  = -8766(e^(-8))+2160(e^(-9/2))-9(e^(-1))
  ≈ 23.62518294825382`.

Calculations

∫(x^2) (E^(-x))\[DifferentialD]x

(∫_1^8)63(x^2) (E^(-x))\[DifferentialD]x//Expand

-63(64 +16 +2)(E^(-8))+63(1 +2 +2)(E^(-1))//Expand

-5166(E^(-8))+315(E^(-1))-(((-60(E^(-4))+18(E^(-1/2))))^2)//Expand

Expectation[(((21x (E^(-x/2))))^2),x\[Distributed]UniformDistribution[{1,8}]]-((Expectation[21x (E^(-x/2)),x\[Distributed]UniformDistribution[{1,8}]])^2)//Expand

(d)

n=10^5;
g[x_]:=21x (E^(-x/2));
expG=Expectation[g[x],x\[Distributed]UniformDistribution[{1,8}]]; (* 𝔼[g(X)] *)
expG2=Expectation[((g[x])^2),x\[Distributed]UniformDistribution[{1,8}]]; (* 𝔼[(g^2)(X)] *)
varG=expG2-(expG^2); (* var(g(X)) *)
rv=RandomReal[UniformDistribution[{1,8}],n];
rv=g[rv];
Gn=Mean[rv];
sdGn=(Sqrt[varG/n]);    (* standard error of estimate using population variance *)
sdGnSample=(Sqrt[Variance[rv]/n]);    (* standard error of estimate using sample variance *)
{Gn,sdGn,sdGnSample}//N

(e)

See R code file.

(f)

The approximation
(((Subscript[G, n])-G)/(Sqrt[var((Subscript[G, n]))]))~N(0,1)    (distribution holds approximately)
so if (Subscript[z, 0.975]) is the 0.975 quantile from N(0,1) then
0.975 = P((((Subscript[G, n])-G)/(Sqrt[var((Subscript[G, n]))]))≤(Subscript[z, 0.975]))
which by symmatery of the N(0,1) distribution means
0.95 = P(|(((Subscript[G, n])-G)/(Sqrt[var((Subscript[G, n]))]))|≤(Subscript[z, 0.975]))
  = P(-(Subscript[z, 0.975])≤(((Subscript[G, n])-G)/(Sqrt[var((Subscript[G, n]))]))≤(Subscript[z, 0.975]))
  = P(-(Subscript[z, 0.975])(Sqrt[var((Subscript[G, n]))])≤(Subscript[G, n])-G≤(Subscript[z, 0.975])(Sqrt[var((Subscript[G, n]))]))
  = P(-(Subscript[G, n])-(Subscript[z, 0.975])(Sqrt[var((Subscript[G, n]))])≤-G≤-(Subscript[G, n])+(Subscript[z, 0.975])(Sqrt[var((Subscript[G, n]))]))
  = P((Subscript[G, n])+(Subscript[z, 0.975])(Sqrt[var((Subscript[G, n]))])≥G≥(Subscript[G, n])-(Subscript[z, 0.975])(Sqrt[var((Subscript[G, n]))]))
  = P((Subscript[G, n])-(Subscript[z, 0.975])(Sqrt[var((Subscript[G, n]))])≤G≤(Subscript[G, n])+(Subscript[z, 0.975])(Sqrt[var((Subscript[G, n]))])).

So the approximate 95% CI is
(Subscript[G, n])-(Subscript[z, 0.975])(Sqrt[var((Subscript[G, n]))]) ≤ G ≤ (Subscript[G, n])+(Subscript[z, 0.975])(Sqrt[var((Subscript[G, n]))]).

Calculations

z975=InverseCDF[NormalDistribution[],0.975];
lb=Gn-z975*sdGn;(* LB using population std err *)
ub=Gn+z975*sdGn;(* UB using population std err *)
lbSample=Gn-z975*sdGnSample;(* LB using sample std err *)
ubSample=Gn+z975*sdGnSample;(* UB using sample std err *)
{lb,ub}    
{lbSample,ubSample}

(g)

See R code file.

(Extra Work)

Let X~Exp(1/2) which has density
(Subscript[f, X])(x) = (1/2)(e^(-x/2))I(x≥0) = ({{{, ({{(1/2)(e^(-x/2)),, x≥0}, {0,, otherwise}})}}).

Setting
g(x) = 6x 𝕀(1≤x≤8)
we have
G = (∫_1^8)3x (e^(-x/2))dx = (∫_(-∞)^∞)6x 𝕀(1≤x≤8)(1/2)(e^(-x/2))I(x≥0)dx = (∫_(-∞)^∞)g(x)(Subscript[f, X])(x)dx
  = 𝔼[g(X)].

Calculations

PDF[ExponentialDistribution[1/2],x]

(* Monte Carlo *)

n=(10^5);
g[x_]:=If[1≤x≤8,6x,0];
expG=Expectation[g[x],x\[Distributed]ExponentialDistribution[1/2]]; (* 𝔼[g(X)] *)
expG2=Expectation[((g[x])^2),x\[Distributed]ExponentialDistribution[1/2]]; (* 𝔼[(g^2)(X)] *)
varG=expG2-(expG^2); (* var(g(X)) *)
rv=RandomReal[ExponentialDistribution[1/2],n];
rv=Map[g,rv];
Gn=Mean[rv];
sdGn=(Sqrt[varG/n]);    (* standard error of estimate using population variance *)
sdGnSample=(Sqrt[Variance[rv]/n]);    (* standard error of estimate using sample variance *)
{Gn,sdGn,sdGnSample}//N

(* Confidence Interval *)
z975=InverseCDF[NormalDistribution[],0.975];
lb=Gn-z975*sdGn;(* LB using population std err *)
ub=Gn+z975*sdGn;(* UB using population std err *)
lbSample=Gn-z975*sdGnSample;(* LB using sample std err *)
ubSample=Gn+z975*sdGnSample;(* UB using sample std err *)
{lb,ub}    
{lbSample,ubSample}