PS 8 Solutions 

11.8 

>	read "/home/fac/little/public_html/ProbStat0910/MSP.map";

 

 

Warning, the use of _seed is deprecated.  Please consider using one of the alternatives listed on the _seed help page.
	(1)

 

>	XList:=[39.0,37.5,22.2,17.5,.64,.45,2.62,2.36,32.0,.77];

 

(2)

 

>	YList:=[23.0,22.3,9.4,9.7,.15,.28,.75,.51,28.,.39];

 

(3)

 

>	Xbar:=Mean(XList);

 

(4)

 

>	Ybar:=Mean(YList);

 

(5)

 

>	Sxx:=add((XList[i]-Xbar)^2,i=1..nops(XList));

 

(6)

 

>	Sxy:=add((XList[i]-Xbar)*(YList[i]-Ybar),i=1..nops(XList));

 

(7)

 

>	Syy:=add((YList[i]-Ybar)^2,i=1..nops(YList));

 

(8)

 

>	beta[1]:=Sxy/Sxx;

 

(9)

 

>	beta[0]:=Ybar-beta[1]*Xbar;

 

(10)

 

>	with(plots):

 

>	PP:=plot([seq([XList[i],YList[i]],i=1..nops(XList))],style=point,symbol=diamond,color=blue):

 

>	LP:=plot(beta[0]+beta[1]*x,x=0..40,color=red):

 

>	display(PP,LP);

 

Estimate of value for x = 12: 

>	beta[0]+beta[1]*12;

 

(11)

 

11.12 

a) 

>	XList:=[0.0,4,14,32,52];

 

(12)

 

>	YList:=[19.8,16.5,12.8,8.1,7.5];

 

(13)

 

>	Xbar:=Mean(XList);

 

(14)

 

>	Ybar:=Mean(YList);

 

(15)

 

>	Sxx:=add((XList[i]-Xbar)^2,i=1..nops(XList));

 

(16)

 

>	Sxy:=add((XList[i]-Xbar)*(YList[i]-Ybar),i=1..nops(XList));

 

(17)

 

>	Syy:=add((YList[i]-Ybar)^2,i=1..nops(YList));

 

(18)

 

>	beta[1]:=Sxy/Sxx;

 

(19)

 

>	beta[0]:=Ybar-beta[1]*Xbar;

 

(20)

 

b) 

>	with(plots):

 

>	PP:=plot([seq([XList[i],YList[i]],i=1..nops(XList))],style=point,symbol=diamond,color=blue):

 

>	LP:=plot(beta[0]+beta[1]*x,x=0..55,color=red):

 

>	display(PP,LP);

 

c)  Estimate for x = 20: 

>	beta[0]+20*beta[1];

 

(21)

 

11.13 

>	XList:=[7.23,8.53,9.82,10.26,8.96,12.27,10.28,4.45,1.78,4.0,3.3,4.3,0.8,0.5];

 

(22)

 

>	YList:=[190,160,134,129,172,197,167,239,542,372,245,376,454,410];

 

(23)

 

>	Xbar:=add(XList[i],i=1..nops(XList))/nops(XList);

 

(24)

 

>	Ybar:=add(evalf(YList[i]),i=1..nops(YList))/nops(YList);

 

(25)

 

>	Sxx:=add((XList[i]-Xbar)^2,i=1..nops(XList));

 

(26)

 

>	Sxy:=add((XList[i]-Xbar)*(YList[i]-Ybar),i=1..nops(XList));

 

(27)

 

>	Syy:=add((YList[i]-Ybar)^2,i=1..nops(YList));

 

(28)

 

>	beta[1]:=Sxy/Sxx;

 

(29)

 

>	beta[0]:=Ybar-beta[1]*Xbar;

 

(30)

 

So the line is:   y = x 

b)  The plot of the data points with the regression line: 

>	with(plots):

 

>	PP:=plot([seq([XList[i],YList[i]],i=1..nops(XList))],style=point,symbol=diamond,color=blue):

 

>	LP:=plot(beta[0]+beta[1]*x,x=0..14,color=red):

 

>	display(PP,LP);

 

11.24 

We must test the null hypothesis   :    against the alternative   

>	SSE:=Syy-beta[1]*Sxy;

 

(31)

 

>	S:=sqrt(SSE/(nops(XList)-2));

 

(32)

 

>	t:=beta[1]/(S*sqrt(1/Sxx));

 

(33)

 

>	nops(XList)-2;

 

(34)

 

The  p-value of the test is < 2(.005)  since     so we reject the  

null hypothesis if .   A more accurate estimate of the p-value is: 

>	p:=2*(1-TCDF(12,5.775));

 

(35)

 

>

 

11.26 

a) 

>	XList:=[29.4,39.2,49.0,58.8,68.6,78.4];

 

(36)

 

>	YList:=[4.25,5.25,6.5,7.85,8.75,10.0];

 

(37)

 

>	Xbar:=add(XList[i],i=1..nops(XList))/nops(XList);

 

(38)

 

>	Ybar:=add(evalf(YList[i]),i=1..nops(YList))/nops(YList);

 

(39)

 

>	Sxx:=add((XList[i]-Xbar)^2,i=1..nops(XList));

 

(40)

 

>	Sxy:=add((XList[i]-Xbar)*(YList[i]-Ybar),i=1..nops(XList));

 

(41)

 

>	Syy:=add((YList[i]-Ybar)^2,i=1..nops(YList));

 

(42)

 

>	beta[1]:=Sxy/Sxx;

 

(43)

 

>	beta[0]:=Ybar-beta[1]*Xbar;

 

(44)

 

b) 

>	S:=sqrt((Syy-beta[1]*Sxy)/4);

 

(45)

 

>	c11:=1/Sxx;

 

(46)

 

From table in text,    

>	beta[1]-2.776*S*sqrt(c11),beta[1]+2.776*S*sqrt(c11);

 

(47)

 

c) 

>	c00:=add(XList[i]^2,i=1..6)/(6*Sxx);

 

(48)

 

>	t:=beta[0]/(S*sqrt(c00));

 

(49)

 

>	2*(1 - TCDF(4,t));

 

(50)

 

p = attained sig. level is = is between 2(.01) = .02 and 2(.005) = .01 

since  t_(.01)(4) = 3.747 and t_(.005)(4) = 4.604.  At , 

we would reject the null hypothesis and conclude  . 

(11.27 on other sheets) 

11.69 

a) 

>	XList:=[-7.,-5.,-3.,-1.,1.,3.,5.,7.];

 

(51)

 

>	YList:=[18.5,22.6,27.2,31.2,33.0,44.9,49.4,35.0];

 

(52)

 

>	Xbar:=Mean(XList);

 

(53)

 

>	Ybar:=Mean(YList);

 

(54)

 

>	Sxx:=add((XList[i]-Xbar)^2,i=1..nops(XList));

 

(55)

 

>	Sxy:=add((XList[i]-Xbar)*(YList[i]-Ybar),i=1..nops(XList));

 

(56)

 

>	Syy:=add((YList[i]-Ybar)^2,i=1..nops(YList));

 

(57)

 

>	beta[1]:=Sxy/Sxx;

 

(58)

 

>	beta[0]:=Ybar-beta[1]*Xbar;

 

(59)

 

Using matrix formuation.  The model  Y =  

>	with(linalg):

 

>	X:=transpose(stackmatrix([seq(1,i=1..nops(XList))],[XList]));

 

(60)

 

>	M:=multiply(transpose(X),X);

 

(61)

 

>	linsolve(M,multiply(transpose(X),YList));

 

(62)

 

The model  Y =  

>	X2:=transpose(stackmatrix([seq(1,i=1..nops(XList))],[XList],[seq(XList[i]^2,i=1..nops(XList))]));

 

(63)

 

>	M2:=multiply(transpose(X2),X2);

 

(64)

 

>	linsolve(M2,multiply(transpose(X2),YList));

 

(65)

 

>	with(plots):

 

>	PP:=plot([seq([XList[i],YList[i]],i=1..nops(XList))],style=point,symbol=diamond,color=blue):

 

>	LP:=plot(32.725+1.812*x,x=-7..7,color=red):

 

>

 

>	display(PP,LP,QP);

 

>

 

In this case , it is more or less clear from the graphs that the point [7,35]  from  

the year 2003 is the ``culprit'' here.  If we excluded that point, the model Y = would fit very well and the quadratic model Y = would not be 

much of an improvement.  With that point, neither model is all that good for  (!) 

>