\begin{Schunk}
\begin{Sinput}
 par(mfrow=c(2,2))
 N <- 2000
 x <- seq(0.0000001, 2, length.out = 4000)
 y1 <- rgamma(N, shape = 1/3, scale = 3)
 h1 <- hist(log(y1), main = "shape = 1/3", breaks = 50,  prob = TRUE,
        #     xlim = rockchalk::magRange(log(y1), c(1.25,1.25)),
             xlab = expression(log(delta[i])))
 y2 <- rgamma(N, shape = 2/3, scale = 3/2)
 hist(log(y2), main = "shape = 2/3", breaks = 50, prob = TRUE,
       # xlim = rockchalk::magRange(log(y1), c(1.25,1.25)),
             xlab = expression(log(delta[i])))
 y3 <- rgamma(N, shape = 1/1, scale = 1)
 hist(log(y3), main="shape = 1", breaks = 50, prob = TRUE,
        xlab = expression(log(delta[i])))
 y4 <- rgamma(x, shape = 3/2, scale = 2/3)
 hist(log(y4), main="shape = 1.5", breaks = 50, prob = TRUE,
  xlab = expression(log(delta[i])))
 par(mfrow=c(1,1))
\end{Sinput}
\end{Schunk}
