############################################### ## Author: Joshua M. Tebbs ## Date: 13 January 2019 ## Update: 18 February 2019 ## STAT 512 course notes: R Code ## Chapter 8 ############################################### # Figure 8.1 # Page 75 # Left y = seq(10,30,0.01) plot(y,dnorm(y,20,2.5),type="l",lty=1,bty="n",xaxt="n",yaxt="n",xlab="",ylab="",ylim=c(-0.01,0.16),cex.lab=1.25) abline(h=0) points(x=20,y=0,pch=19,cex=1.25) text(20,-0.01,expression(theta),cex=1.5) text(30,-0.0075,expression(hat(theta)),cex=1.5) # Right y = seq(10,30,0.01) plot(y,dnorm(y,20,2.5),type="l",lty=1,bty="n",xaxt="n",yaxt="n",xlab="",ylab="",ylim=c(-0.01,0.16),cex.lab=1.25) abline(h=0) points(x=17,y=0,pch=19,cex=1.25) text(17,-0.01,expression(theta),cex=1.5) text(30,-0.0075,expression(hat(theta)),cex=1.5) # Figure 8.2 # Page 76 y = seq(0,1,0.01) pdf = rep(1,length(y)) plot(y,pdf,type="l",xlab="y",ylab="PDF",xaxt='n',yaxt='n',xlim=c(0,1.25),ylim=c(0,1.25),cex.lab=1.25) axis(side=1,at=c(0,1),labels=expression(0,theta)) axis(side=2,at=c(0,1),labels=expression(0,1/theta)) lines(c(1,1.5),c(0,0),lty=1) # Figure 8.3 # Page 78 # Left y = seq(10,30,0.01) plot(y,dnorm(y,20,2.5),type="l",lty=1,bty="n",xaxt="n",yaxt="n",xlab="",ylab="",ylim=c(-0.01,0.16),cex.lab=1.25) abline(h=0) points(x=20,y=0,pch=19,cex=1.25) text(20,-0.01,expression(theta),cex=1.5) text(30,-0.0075,expression(hat(theta)[1]),cex=1.5) # Right y = seq(10,30,0.01) plot(y,dnorm(y,20,0.5),type="l",lty=1,bty="n",xaxt="n",yaxt="n",xlab="",ylab="",ylim=c(-0.04,0.8),cex.lab=1.25) abline(h=0) points(x=20,y=0,pch=19,cex=1.25) text(20,-0.04,expression(theta),cex=1.5) text(30,-0.035,expression(hat(theta)[2]),cex=1.5) # Figure 8.4 # Page 80 n = 25 p = seq(0,1,0.001) mse.1 = p*(1-p)/n mse.2 = (n*p*(1-p))/((n+4)^2) + ((n*p+2)/(n+4)-p)^2 plot(p,mse.1,type="l",lty=1,xlab="p",ylab="MSE",xlim=c(0,1),cex.lab=1.25) abline(h=0) lines(p,mse.2,lty=4) legend(0.37,0.004,lty = c(1,4),c(expression(MSE(hat(p)[1])),expression(MSE(hat(p)[2])))) # Figure 8.5 # Page 85 y = seq(10,30,0.01) plot(y,dnorm(y,20,2.5),type="l",bty="n",lty=1,xaxt="n",yaxt="n",xlab="",ylab="",ylim=c(-0.01,0.16),cex.lab=1.25) abline(h=0) points(x=20,y=0,pch=19,cex=1.25) text(20,-0.01,expression(mu),cex=1.5) text(30,-0.0075,expression(bar(y)),cex=1.5) x = seq(25,30,0.01) y = dnorm(x,20,2.5) polygon(c(25,x,30),c(0,y,0),col="lightblue") x = seq(10,15,0.01) y = dnorm(x,20,2.5) polygon(c(10,x,15),c(0,y,0),col="lightblue") text(20,0.07,0.95,cex=1.25) # Figure 8.6 # Page 87 y = seq(10,30,0.01) plot(y,dnorm(y,20,2.5),type="l",bty="n",lty=1,xaxt="n",yaxt="n",xlab="",ylab="",ylim=c(-0.01,0.16),cex.lab=1.25) abline(h=0) points(x=20,y=0,pch=19,cex=1.25) text(20,-0.01,expression(p),cex=1.5) text(30,-0.0075,expression(hat(p)),cex=1.5) x = seq(25,30,0.01) y = dnorm(x,20,2.5) polygon(c(25,x,30),c(0,y,0),col="lightblue") x = seq(10,15,0.01) y = dnorm(x,20,2.5) polygon(c(10,x,15),c(0,y,0),col="lightblue") text(20,0.07,0.95,cex=1.25) # Figure 8.7 # Page 91 birth.weights = read.table("C:\\Users\\tebbs.MATHSTAT\\Documents\\texfiles\\Classes\\USC\\stat512\\s19\\data\\birthweights.txt",header=TRUE) birth.weights = birth.weights[,1] hist(birth.weights,xlab="Birth weight (in grams)",ylab="Count",main="",col="lightblue",ylim=c(0,120)) # Figure 8.8 # Page 93 x = seq(-5,5,0.001) pdf = dt(x,10) plot(x,pdf,type="l",lty=1,xlab="",xaxt="n",yaxt="n",bty="n",ylab="",ylim=c(0,0.4)) abline(h=0) x = seq(-5,qt(0.025,10),0.001) y = dt(x,10) polygon(c(-5,x,qt(0.025,10)),c(0,y,0),col="lightblue") points(x=qt(0.025,10),y=0,pch=19,cex=1.5) x = seq(qt(0.975,10),5,0.001) y = dt(x,10) polygon(c(qt(0.975,10),x,5),c(0,y,0),col="lightblue") points(x=qt(0.975,10),y=0,pch=19,cex=1.5) text(-0.025,0.1,expression(1-alpha),cex=1.25) text(-3.5,0.05,expression(alpha/2),cex=1.25) text(3.5,0.05,expression(alpha/2),cex=1.25) # Figure 8.9 # Page 94 x = seq(0,30,0.001) pdf = dchisq(x,10) plot(x,pdf,type="l",lty=1,xlab="",ylab="",xaxt="n",yaxt="n",bty="n",ylim=c(0,0.1)) abline(h=0) x = seq(0,qchisq(0.025,10),0.001) y = dchisq(x,10) polygon(c(0,x,qchisq(0.025,10)),c(0,y,0),col="lightblue") points(x=qchisq(0.025,10),y=0,pch=19,cex=1.5) x = seq(qchisq(0.975,10),30,0.001) y = dchisq(x,10) polygon(c(qchisq(0.975,10),x,30),c(0,y,0),col="lightblue") points(x=qchisq(0.975,10),y=0,pch=19,cex=1.5) text(9.8,0.02,expression(1-alpha),cex=1.25) text(0.5,0.01,expression(alpha/2),cex=1.25) text(25,0.0075,expression(alpha/2),cex=1.25) # Figure 8.10 # Page 97 x = seq(0,1,0.001) pdf = dbeta(x,2.5,1) plot(x,pdf,type="l",lty=1,xlab="",ylab="",xaxt="n",yaxt="n",bty="n",xlim=c(0,1.05)) abline(h=0) x = seq(0,qbeta(0.025,2.5,1),0.001) y = dbeta(x,2.5,1) polygon(c(0,x,qbeta(0.025,2.5,1)),c(0,y,0),col="lightblue") points(x=qbeta(0.025,2.5,1),y=0,pch=19,cex=1.5) x = seq(qbeta(0.975,2.5,1),1,0.001) y = dbeta(x,2.5,1) polygon(c(qbeta(0.975,2.5,1),x,1),c(0,y,0),col="lightblue") points(x=qbeta(0.975,2.5,1),y=0,pch=19,cex=1.5) text(0.7,0.6,expression(1-alpha),cex=1.25) text(0.1,0.25,expression(alpha/2),cex=1.25) text(1.05,0.25,expression(alpha/2),cex=1.25) # Figure 8.11 # Page 103 x = seq(-5,5,0.001) pdf = dt(x,10) plot(x,pdf,type="l",lty=1,xlab="",xaxt="n",yaxt="n",bty="n",ylab="",ylim=c(0,0.4)) abline(h=0) x = seq(-5,qt(0.025,10),0.001) y = dt(x,10) polygon(c(-5,x,qt(0.025,10)),c(0,y,0),col="lightblue") points(x=qt(0.025,10),y=0,pch=19,cex=1.5) x = seq(qt(0.975,10),5,0.001) y = dt(x,10) polygon(c(qt(0.975,10),x,5),c(0,y,0),col="lightblue") points(x=qt(0.975,10),y=0,pch=19,cex=1.5) text(-0.025,0.1,expression(1-alpha),cex=1.25) text(-3.5,0.05,expression(alpha/2),cex=1.25) text(3.5,0.05,expression(alpha/2),cex=1.25) # Figure 8.12 # Page 105 x = seq(0,30,0.001) pdf = dchisq(x,10) plot(x,pdf,type="l",lty=1,xlab="",ylab="",xaxt="n",yaxt="n",bty="n",ylim=c(0,0.1)) abline(h=0) x = seq(0,qchisq(0.025,10),0.001) y = dchisq(x,10) polygon(c(0,x,qchisq(0.025,10)),c(0,y,0),col="lightblue") points(x=qchisq(0.025,10),y=0,pch=19,cex=1.5) x = seq(qchisq(0.975,10),30,0.001) y = dchisq(x,10) polygon(c(qchisq(0.975,10),x,30),c(0,y,0),col="lightblue") points(x=qchisq(0.975,10),y=0,pch=19,cex=1.5) text(9.8,0.02,expression(1-alpha),cex=1.25) text(0.5,0.01,expression(alpha/2),cex=1.25) text(25,0.0075,expression(alpha/2),cex=1.25) # Figure 8.13 # Page 107 x = seq(-5,5,0.001) pdf = dt(x,10) plot(x,pdf,type="l",lty=1,xlab="",xaxt="n",yaxt="n",bty="n",ylab="",ylim=c(0,0.4)) abline(h=0) x = seq(-5,qt(0.025,10),0.001) y = dt(x,10) polygon(c(-5,x,qt(0.025,10)),c(0,y,0),col="lightblue") points(x=qt(0.025,10),y=0,pch=19,cex=1.5) x = seq(qt(0.975,10),5,0.001) y = dt(x,10) polygon(c(qt(0.975,10),x,5),c(0,y,0),col="lightblue") points(x=qt(0.975,10),y=0,pch=19,cex=1.5) text(-0.025,0.1,expression(1-alpha),cex=1.25) text(-3.5,0.05,expression(alpha/2),cex=1.25) text(3.5,0.05,expression(alpha/2),cex=1.25) # Example 8.15 # Pages 108-109 loc.1 = c(5030,4980,13700,11910,10730,8130,11400,26850,860,17660,2200,22800,4250,1130,15040,1690) loc.2 = c(2800,2810,4670,1330,6890,3320,7720,1230,7030,2130,7330,2190) t.test(loc.1,loc.2,conf.level=0.95,var.equal=TRUE)\$conf.int t.test(loc.1,loc.2,conf.level=0.95,var.equal=FALSE)\$conf.int # Figure 8.14 # Page 110 x = seq(0,30,0.001) pdf = dchisq(x,10) plot(x,pdf,type="l",lty=1,xlab="",ylab="",xaxt="n",yaxt="n",bty="n",ylim=c(0,0.1)) abline(h=0) x = seq(0,qchisq(0.025,10),0.001) y = dchisq(x,10) polygon(c(0,x,qchisq(0.025,10)),c(0,y,0),col="lightblue") points(x=qchisq(0.025,10),y=0,pch=19,cex=1.5) x = seq(qchisq(0.975,10),30,0.001) y = dchisq(x,10) polygon(c(qchisq(0.975,10),x,30),c(0,y,0),col="lightblue") points(x=qchisq(0.975,10),y=0,pch=19,cex=1.5) text(9.8,0.02,expression(1-alpha),cex=1.25) text(0.5,0.01,expression(alpha/2),cex=1.25) text(25,0.0075,expression(alpha/2),cex=1.25)