###################################################
### chunk number 1: startup
###################################################
    rm(list=ls())
    set.seed(439)
    options(width=55,digits=4,continue=" ")
    par(ask=FALSE, mfrow=c(1,1))
    setwd("c:/myd2/JSM Workshop")


###################################################
### chunk number 2: Uniform
###################################################
  n <- 1000
  u <- runif(n)         #vector of UNIF(0,1)
  v <- runif(n, -1, 1)  #vector of UNIF(-1,1)


###################################################
### chunk number 3: Uniform
###################################################
  summary(cbind(u, v))  #descriptive stats for both


###################################################
### chunk number 4: cumsum
###################################################
  z <- runif(5, -1, 1)
  x <- cumsum(c(0, z)) #cumulative sum
  z
  x


###################################################
### chunk number 5: tinyRW
###################################################
  plot(x, type="l")


###################################################
### chunk number 6: RW
###################################################
  n <- 1000
  z <- runif(n, -1, 1)
  x <- cumsum(c(0, z))
  plot(x, type="l")


###################################################
### chunk number 7: Normal
###################################################
  n <- 1000
  y <- rnorm(n, 2)  #vector of N(2,1)


###################################################
### chunk number 8: norm6
###################################################
  head(y)  #the first few entries of y


###################################################
### chunk number 9: Normalplot
###################################################
    hist(y, prob = TRUE, breaks = "scott",
        ylim = c(0, 1/sqrt(2*pi)))
    curve(dnorm(x,2), add=TRUE)


###################################################
### chunk number 10: perm
###################################################
  sample(letters[1:5])


###################################################
### chunk number 11: sample
###################################################
  sample(1:20, size=10)


###################################################
### chunk number 12: srs
###################################################
  x <- sample(state.region, size=100, replace=TRUE)
  cbind(table(state.region)/50 ,table(x)/100)


###################################################
### chunk number 13: biased.coin
###################################################
  x <- sample(c("H","T"), size=1000, replace=TRUE,
      prob=c(.4, .6))
  table(x)


###################################################
### chunk number 14: biased.coin2
###################################################
  table(sample(c(1,0), size=1000, replace=TRUE,
      prob=c(.4, .6)))


###################################################
### chunk number 15: binom.coin
###################################################
    y <- rbinom(1000, size=1, prob=.4)
    table(y)


###################################################
### chunk number 16: binom.coin
###################################################
    mean(y)


###################################################
### chunk number 17: propTRUE
###################################################
    x <- runif(10)
    round(x, 3)
    x < .25
    mean(x < .25)


###################################################
### chunk number 18: forloop
###################################################
  m <- 10000; n <- 200
  xbar1 <- numeric(m)    #empty vector length m
  for (i in 1:m) {
      x <- rnorm(n)
      xbar1[i] <- mean(x)
  }


###################################################
### chunk number 19: replicate
###################################################
    m <- 10000; n <- 200
    xbar2 <- replicate(m, mean(rnorm(n)))


###################################################
### chunk number 20: CLT
###################################################
  par(mfrow=c(2,2))
  hist(xbar1, prob = TRUE, main="for loop")
    curve(dnorm(x, 0, 1/sqrt(n)), add=TRUE)
  hist(xbar2, prob = TRUE, main="replicate")
    curve(dnorm(x, 0, 1/sqrt(n)), add=TRUE)
  par(mfrow=c(1,1))


###################################################
### chunk number 21: replicate2
###################################################
  m <- 1000; n <- 200
  stats <- replicate(m, expr={
      x <- rnorm(n)
      c(mean(x), var(x))
      })
  dim(stats)


###################################################
### chunk number 22: rep2plot
###################################################
  xbar <- stats[1,]   #means in row 1
  v <- stats[2,]      #vars in row 2
  cor(xbar, v)
  plot(xbar, v, main="N(0,1) data", cex=.5)


###################################################
### chunk number 23: rep2xplot
###################################################
  stats <- replicate(m, expr={
      x <- rlnorm(n)
      c(mean(x), var(x))
      })
  xbar <- stats[1,]
  v <- stats[2,]
  cor(xbar, v)


###################################################
### chunk number 24: rep2x
###################################################
  plot(xbar, v, main="LogN(0,1) data", cex=.5)


###################################################
### chunk number 25: scott
###################################################
  nclass.scott(xbar1)
  nclass.Sturges(xbar1)  #hist default


###################################################
### chunk number 26: abbr
###################################################
x <- sapply(state.name, FUN=abbreviate)
head(x)  #display the first few states


###################################################
### chunk number 27: payoff
###################################################
print(payoff <- matrix(rgeom(16, prob=.1), 4, 4))


###################################################
### chunk number 28: minmax
###################################################
  cMax <- apply(payoff, 2, max)
  rMin <- apply(payoff, 1, min)
  cMax
  rMin
  list(minmax=min(cMax), maxmin=max(rMin))


###################################################
### chunk number 29: normalmixture
###################################################
  p <- .75
  m <- 10000
  u <- runif(m)
  mu <- ifelse(u <= p, 0, 3)
  x1 <- rnorm(m, mean=mu, sd=1)


###################################################
### chunk number 30: mu
###################################################
  head(mu)


###################################################
### chunk number 31: normalmixture2
###################################################
  mu <- sample(c(0,3), size=m, replace=TRUE,
          prob=c(p, 1-p))
  x2 <- rnorm(m, mean=mu, sd=1)


###################################################
### chunk number 32: mu
###################################################
  head(mu)


###################################################
### chunk number 33: normalmixtureplot
###################################################
  par(mfrow=c(1,2))
  hist(x1, breaks="scott", prob=TRUE)
  curve(p*dnorm(x)+(1-p)*dnorm(x,3), add=TRUE)
  hist(x2, breaks="scott", prob=TRUE)
  curve(p*dnorm(x)+(1-p)*dnorm(x,3), add=TRUE)
  par(mfrow=c(1,1))