hist-loess.R

set.seed(110951)
rm(list = ls())
library(MASS)
N <- 1000
ymax <- 300
data1 <- rgamma(N, shape = 3)
data2 <- 8 + rgamma(N, shape =  6, scale = 1/6)
data <- ifelse(runif(N)< 0.8, data1, data2)
max(data)

## [1] 10.24648

truehist(data, breaks = 0:12, prob = FALSE, xlim = c(0,12), ymax = ymax)
x <- (0:11) + 0.5
y <- hist(data, breaks = 0:12, plot = FALSE)$counts
max(y)

## [1] 226

points(x, y)
xnew <- (0:110)/10
rooty <- sqrt(y)
loess.out <- loess(rooty ~ x, span = 3/7)
rootynew <- predict(loess.out, newdata = data.frame(x = xnew))
lines(xnew, rootynew^2, lwd = 2)
yfit <- predict(loess.out, newdata = data.frame(x = x))
points(x, yfit^2, pch = "*", cex = 2)