##Paul Johnson
##2010-01-29
##

N <- 400
x<- 50+ 50*rnorm(N)
y <- 450 + 1.4 *x + 100 * rnorm(N)

## plot to the on-screen device (whatever is default on your system")

plot(x,y,xlab="School Spending", ylab="Standardized Test Scores")
abline(lm(y~x))


## Dump that into an EPS file

dev.copy(postscript, file="school-1.eps", horizontal=F,onefile=F)
dev.off()

## Same as shortcut "convenience" program

dev.copy2eps(file="school-2.eps")

##Caution, no need for "dev.off" because it kills the graph on screen.
##dev.off()

## Specify Sizes, experiment

dev.copy(postscript, file="school-3.eps", onefile=F, horizontal=F, height=2, width=5, pointsize=12, family="Times")
dev.off()

## Very interestin behavior. It scrunched the graph, but not the label words!

##Most of the time, the output from dev.copy is adequate.


postscript(file="school-4.eps", onefile=F, horizontal=F, height=5, width=5, pointsize=12, family="Times")

plot(x,y,xlab="School Spending", ylab="Standardized Test Scores")
abline(lm(y~x))
dev.off()

### Let's see how well it copies to pdf
plot(x,y,xlab="School Spending", ylab="Standardized Test Scores")
abline(lm(y~x))
dev.copy(pdf, file="school-1.pdf", onefile=F, height=5, width=5, pointsize=12, family="Times", paper="special")
dev.off()


pdf(file="school-2.pdf", onefile=F, height=5, width=5, pointsize=12, family="Times",paper="special")

plot(x,y,xlab="School Spending", ylab="Standardized Test Scores")
abline(lm(y~x))
dev.off()


dev.copy(xfig,file="school.fig")
dev.off()

dev.copy(svg, file="school.svg")
dev.off()





# Suppose you needed to scan through 10 datasets and
# create scatterplots for each one in order to compare
# the relationship between two variables.  In the
# "point and click" era, you would do a lot of clicking.
# Instead, in R, you can master the commands to save
# graphs into files and get the job done.

for (i in 1:10){
  beta <- rnorm(1)
  x<- 50+ 50*rnorm(N)
  y <- 450 + beta *x + 100 * rnorm(N)
  
  pdf(file=paste("schsp",i,".pdf",sep=""), onefile=F, height=5.5, width=5.5,paper="special", pointsize=12, family="Times")
  plot(x,y,xlab="School Spending", ylab="Standardized Test Scores")
  abline(lm(y~x))
  dev.off()
}


### Boring to re-type options every time. Run this once per program
pdf.options(onefile=F, height=5.5, width=5.5,paper="special", pointsize=12, family="Times")

for (i in 1:10){
  beta <- rnorm(1)
  x<- 50+ 50*rnorm(N)
  y <- 450 + beta *x + 100 * rnorm(N)
  myfilename = paste("schsp",i,".pdf",sep="")
  pdf(file = myfilename)
  plot(x,y,xlab="School Spending", ylab="Standardized Test Scores")
  abline(lm(y~x))
  dev.off()
}


### Similarly, I want encapsulated postscript output.
### You get EPS output from a postscript device by adding the options
### onefile=F and paper="special" (that creates a "bounding box")
ps.options(onefile=F, height=5.5, width=5.5,paper="special", pointsize=12, family="Times", horizontal=F)

level2beta <- 0.4
for (i in 1:10){
  beta <- level2beta + rnorm(1)
  x<- 50+ 50*rnorm(N)
  y <- 450 + beta *x + 100 * rnorm(N)
  
  postscript(file=paste("schsp",i,".eps",sep=""))
  plot(x,y,xlab="School Spending", ylab="Standardized Test Scores")
  mod <- lm(y~x)
  abline(mod)
  estbeta <- round(coef(mod),2)
  eqtext <- bquote(widehat(Scores[i]) == .(estbeta[1]) + .(estbeta[2]) * Spending[i])
  myyposition <- ifelse(estbeta[2]>0, min(y), max(y))
  text(max(x), myyposition, eqtext, pos=2)
  abline(450, level2beta, col="blue",lty=4)
  mylegendposition <- ifelse(estbeta[2]>0, "topleft","bottomleft")

  legend(mylegendposition, legend=c("subsample estimate","level 2 beta"), col=c("black","blue"), lty=c(1,4))
  
  dev.off()
}