tmpar <- par(mfrow=c(1,2),mar=c(4.1,4.1,0,0))

attach(damage)
mplot(x,y,group=as.character(score),xlab="(a) logit(damage)",ylab="")
tmp <- order(x)
lines(x[tmp],predict(damage.x)[tmp],lty=2)
se.bar(0.5,93,std.dev(damage.x),cap="SD")

plot(score,y,xlim=c(1,5),xaxt="n",xlab="(b) score",ylab="yield")
axis(1,at=1:5)
tmp <- order(score)
lines(score[tmp],predict(damage.s)[tmp],lty=2)
for (i in 1:4)
	lines(i+c(-.4,.4),rep(damage.mean[i],2))
se.bar(3.5,93,std.dev(damage.score),cap="SD")
detach()

par(tmpar)

#damage.level <- matrix(c(
#1,0,5,
#2,5,25,
#3,25,75,
#4,75,95,
#5,95,100),
#5,3,byrow=T,dimnames=list(NULL,c("score","lo","hi")))
#
#plot(c(0,100),c(1,5),type="n",
#	xlab="(c) damage range (%)",ylab="score")
#for (i in 1:5){
#	lines(damage.level[i,c("lo","hi")],rep(damage.level[i,"score"],2))
#}
#
#mplot(100*damage$p,resid(damage.s),group=as.character(damage$score),
#	xlim=c(0,100),xlab="(d) actual damage (%)",ylab="score residual")
#abline(h=0,lty=2)
#abline(h=c(-1,1)*std.dev(damage.score),lty=3)