pdf("F:\\Rmisc\\graphs\\ballistic1.pdf") ballistic <- read.table("http://www.stat.ufl.edu/~winner/data/ballistic1.dat",header=F, col.names=c("bulltype","layers","v50","sharp","fsp")) attach(ballistic) v50_2 <- (v50/100)^2 int.mod <- lm(v50_2 ~ layers + sharp + fsp + I(layers*sharp) + I(layers*fsp)) summary(int.mod) anova(int.mod) add.mod <- lm(v50_2 ~ layers + sharp + fsp) summary(add.mod) anova(add.mod) anova(add.mod,int.mod) xp <- 0:50 pred.int1 <- int.mod$coef[1]+int.mod$coef[2]*xp pred.int2 <- int.mod$coef[1]+int.mod$coef[2]*xp+int.mod$coef[3]+int.mod$coef[5]*xp pred.int3 <- int.mod$coef[1]+int.mod$coef[2]*xp+int.mod$coef[4]+int.mod$coef[6]*xp pred.add1 <- add.mod$coef[1]+add.mod$coef[2]*xp pred.add2 <- add.mod$coef[1]+add.mod$coef[2]*xp+add.mod$coef[3] pred.add3 <- add.mod$coef[1]+add.mod$coef[2]*xp+add.mod$coef[4] plot(layers,v50_2,pch=as.numeric(bulltype),xlim=c(0,45),ylim=c(0,1.2*max(v50_2)),main="Interaction Model") lines(xp,pred.int1,lty=1) lines(xp,pred.int2,lty=2) lines(xp,pred.int3,lty=3) legend("topleft",inset=.05,title="Bullet Type", c("Round","Sharp","FSP"), lty=c(1,2,3), pch=c(1,2,3)) plot(layers,v50_2,pch=as.numeric(bulltype),xlim=c(0,45),ylim=c(0,1.2*max(v50_2)),main="Additive Model") lines(xp,pred.add1,lty=1) lines(xp,pred.add2,lty=2) lines(xp,pred.add3,lty=3) legend("topleft",inset=.05,title="Bullet Type", c("Round","Sharp","FSP"), lty=c(1,2,3), pch=c(1,2,3)) dev.off()