#------------------------
# This generates scatter plots  used in class7-prep

makeCorPdf = FALSE  #SET TO TRUE TO OUTPUT pdf files
dotRadius = .5
runOneCor = function(n, overlap, npts, col, makePdf)
{
    #X = sum of n uniform(0,1)
    #Y = sum of n uniform(0,1)
    #overlap = number in common between X and Y
    ncol = 2*n - overlap
    tmp = runif(npts*ncol,0,1)
    data = matrix(tmp, nrow=npts,ncol=ncol)  #raw data
  
    xcols = 1:n                         #x data = first n columns
    ycols = (n-overlap+1):(2*n-overlap) #y data overlaps x data in overlap columns
    if (n == 1)  #special case, only one column so rowSums barfs
    {
        x = data[,xcols]
        y = data[,ycols]
    }
    else
    {
      x = rowSums(data[,xcols])
      y = rowSums(data[,ycols])
    }
    theoryCor = overlap/n       #Theoretical correlation
    sampleCor = cor(x,y)        #Sample correlation
    s = sprintf("(%d, %d) cor=%.2f, sample_cor=%.2f", n, overlap, theoryCor, sampleCor)

    if (makePdf)
    {
      pdfwidth = 4
      pdfheight = 3.8
      pdfname = sprintf("../pdftex/img/class7-cor_%d_%d.pdf", n,overlap)
      pdf(file=pdfname, width=pdfwidth, height=pdfheight)
    }
    plot(x,y, col=col, cex=dotRadius, main=s, xlab="")
    title(xlab='x',mgp=c(2,1,0))
    abline(v= mean(x))
    abline(h=mean(y)) 
    if (makePdf)
      dev.off()
}
npts = 1000
runOneCor(1,0,npts,rgb(.2,.4,.9),makeCorPdf)
runOneCor(2,1,npts,'red',makeCorPdf)
runOneCor(5,1,npts,'orange',makeCorPdf)
runOneCor(5,3,npts,'cyan',makeCorPdf)
runOneCor(10,5,npts,'purple',makeCorPdf)
runOneCor(10,8,npts,'magenta',makeCorPdf)

#------------------------
# Bivariate normal scatter plots

makeBivarPdf = FALSE  #SET TO TRUE TO OUTPUT pdf files
dotRadius = .5
runOneBivariateNormal = function(rho, npts, col, makePdf)
{
  #draw (u,v) independent standard normal
  u = rnorm(npts,0,1) 
  v = rnorm(npts,0,1)
  # manipulate to bivariate normal with means m1,m2, stdev s1, s2, cor rho
  m1=0; m2=0; s1=1; s2=1;
  x = m1 + s1*u
  y = m2 + s2*(rho*u + sqrt(1 - rho^2)*v)
  
  sampleCor = cor(x,y)            #Sample correlation
  s = sprintf("rho=%.2f",rho)
  
  if (makePdf)
  {
    pdfwidth = 4
    pdfheight = 3.8
    f = round(10*rho)
    pdfname = sprintf("../pdftex/img/class7-bivarnorm_%d.pdf", f)
    pdf(file=pdfname, width=pdfwidth, height=pdfheight)
  }
  plot(x,y, col=col, cex=dotRadius, main=s, xlab="")
  title(xlab='x',mgp=c(2,1,0))
  abline(v= mean(x))
  abline(h=mean(y)) 
  if (makePdf)
    dev.off()
}
npts = 1000
runOneBivariateNormal(0.0,npts,rgb(.2,.4,.9),makeBivarPdf)
runOneBivariateNormal(0.3,npts,'red',makeBivarPdf)
runOneBivariateNormal(0.7,npts,'orange',makeBivarPdf)
runOneBivariateNormal(1.0,npts,'cyan',makeBivarPdf)
runOneBivariateNormal(-0.5,npts,'purple',makeBivarPdf)
runOneBivariateNormal(-0.9,npts,'magenta',makeBivarPdf)

#-------------------------------