#### R tips ######################################################### ## R tips: how to use R helps? # want to know how to use function "runif" ?runif # want to know how to generate random numbers from normal distribution help.search("normal distribution") ?Normal ######################################################### ######################################################### ## R tips: how to print or save R graphs? # if use menu operations: [1] click the R graph window; # [2] (menu)File --> print or File --> Save as # if use R command lines: x <- rnorm(1000, mean=1, sd=2) postscript("normal.ps",paper="special",width=8,height=7,horizontal=F) par(mfrow=c(1,1)) hist(x, freq=F) curve(dnorm(x, mean=1, sd=2), lwd=2, lty=2, add=T) dev.off() ######################################################### ### Example: Table 8.1-1 on page 440 ## read data from a file, if data file is in the working folder data = read.table("table811.dat", head=T) x = data\$Weight y = data\$GPM ## plot Figure 8.1-2 par(mfrow=c(1,1)) plot(x, y, xlab="Weight (in 1000 lb)", ylab="GPM (gallongs/100 miles)", main="Figure 8.1-2: Scatter plot of gpm against weight") # print graph into ps file postscript("figure812.ps",paper="special",width=8,height=7,horizontal=F) par(mfrow=c(1,1)) plot(x, y, xlab="Weight (in 1000 lb)", ylab="GPM (gallongs/100 miles)", main="Figure 8.1-2: Scatter plot of gpm against weight") dev.off() # print graph into eps file postscript("figure812.eps",paper="special",width=8,height=7,horizontal=T) par(mfrow=c(1,1)) plot(x, y, xlab="Weight (in 1000 lb)", ylab="GPM (gallongs/100 miles)", xlim=c(0,5), ylim=c(0,8), main="Figure 8.1-2: Scatter plot of gpm against weight") dev.off() ## estimate coefficients (beta0, beta1) xbar = mean(x) ybar = mean(y) xybar = mean(x*y) x2bar = mean(x*x) beta1 = (xybar - xbar*ybar)/(x2bar - xbar^2) beta0 = ybar - beta1*xbar beta1 #[1] 1.638996 beta0 #[1] -0.3630888 ## use R function "lm" directly lm(y~x) # see more output of "lm" fit1 = lm(y~x) summary(fit1) # plot output of "lm" par(mfrow=c(2,2)) plot(fit1) ## plot residual plot directly par(mfrow=c(1,1)) plot(fit1\$fitted, fit1\$residual, xlab="Fitted values", ylab="Residuals", main="Residuals vs. Fitted") abline(0, 0, lty=2)