Elliot Tanis
April 20, 2006
Exercise 5.3-12  (Exercise_5_3-12.mws)

Exercise 5.3-12

(a) Simulations of T  when sampling from N (1, 1).

>    randomize():
for k from 1 to 500 do
  X := NormalS(1,1,7):
  TT[k] := evalf((Mean(X) - 1)/(StDev(X)/sqrt(7))):
od:
T := [seq(TT[k], k = 1 .. 500)]:
Min(T),Max(T);
tbar := Mean(T);
mu := 0;
svar := Variance(T);
var := 6/4;

-4.018484602, 4.169998454

tbar := .6379135082e-1

mu := 0

svar := 1.214785714

var := 3/2

>    xtics := [seq(-5 + 0.5*k, k = 0 .. 20)]:
ytics := [seq(0.05*k, k = 1 .. 9)]:
P1 := plot(TPDF(6,x), x = -5 .. 5, color=black, thickness=2, xtickmarks=xtics, ytickmarks=ytics, labels=[``,``]):
P2 := HistogramFill(T, -5.5 .. 5.5, 22):
display({P1, P2});

[Maple Plot]

>    quantiles := [seq(TP(6, k/501), k = 1 .. 500)]:

>    T := sort(T):

>    xtics := [seq(-5.5 + k*.5, k = 0 .. 22)]:
ytics := xtics:
P1 := plot([[-5,-5],[5,5]], x = -5.5 .. 5.5, y = -5.5 .. 5.5, color=black, thickness=2, xtickmarks=xtics, ytickmarks=ytics, labels=[``,``]):
P2 := ScatPlotCirc(T, quantiles):
display({P1, P2});

[Maple Plot]

(b) Now sample from an exponential distribution with mean 1 and observe values of " T ".

>    for k from 1 to 500 do
  X := ExponentialS(1,7):
  TT[k] := evalf((Mean(X) - 1)/(StDev(X)/sqrt(7))):
od:
T := [seq(TT[k], k = 1 .. 500)]:
Min(T),Max(T);
tbar := Mean(T);
mu := 0;
svar := Variance(T);
var := 6/4;

-19.81023059, 3.295568446

tbar := -.5882264744

mu := 0

svar := 3.763433258

var := 3/2

>    xtics := [seq(-7 + 0.5*k, k = 0 .. 28)]:
ytics := [seq(0.05*k, k = 1 .. 9)]:
P1 := plot(TPDF(6,x), x = -7 .. 7, color=black, thickness=2, xtickmarks=xtics, ytickmarks=ytics, labels=[``,``]):
P2 := HistogramFill(T, -7 .. 7, 28):
display({P1, P2});

`WARNING: There are`, 7, ` data points out of plot range.`

[Maple Plot]

If you have already done the previous simulation, you do not have to again calculate the quantiles for the t-distribution with with 6 degrees of freedom. If this is the first simulation, remove # and run the following line.

>    #quantiles := [seq(TP(6, k/501), k = 1 .. 500)]:

>    T := sort(T):

>    xtics := [seq(-10 + k, k = 0 .. 18)]:
ytics := xtics:
P1 := plot([[-7,-7],[7,7]], x = -7 .. 7, y = -7 .. 7,color=black, thickness=2, xtickmarks=xtics, ytickmarks=ytics, labels=[``,``]):
P2 := ScatPlotCirc(T, quantiles):
display({P1, P2});

[Maple Plot]

(c) Now try sampling from the Cauchy distribution.  The following procedure can be used.

>    CauchyS(7);

[-.4430582587, .8814225711, 2.606831677, -1.240217020, .4413527455e-2, .1575246206, .4338170866e-1]

>    for k from 1 to 500 do
  X := CauchyS(7):
  TT[k] := evalf(Mean(X)/(StDev(X)/sqrt(7))):
od:
T := [seq(TT[k], k = 1 .. 500)]:
Min(T),Max(T);
tbar := Mean(T);
mu := 0;
svar := Variance(T);
var := 6/4;

-2.886241409, 3.676701498

tbar := .4985784196e-1

mu := 0

svar := 1.181171956

var := 3/2

>    xtics := [seq(-7 + 0.5*k, k = 0 .. 28)]:
ytics := [seq(0.05*k, k = 1 .. 9)]:
P1 := plot(TPDF(6,x), x = -4 .. 4, color=black, thickness=2, xtickmarks=xtics, ytickmarks=ytics, labels=[``,``]):
P2 := HistogramFill(T, -4 .. 4, 16):
display({P1, P2});

[Maple Plot]

>    T := sort(T):

>    xtics := [seq(-5 + k*.5, k = 0 .. 20)]:
ytics := xtics:
P1 := plot([[-4,-4],[4,4]], x = -4 .. 4, y = -4 .. 4, color=black, thickness=2, xtickmarks=xtics, ytickmarks=ytics, labels=[``,``]):
P2 := ScatPlotCirc(T, quantiles):
display({P1, P2});

[Maple Plot]

>   

Return to Menu