Superponer curva normal al histograma en R
Me las arreglé para encontrar en línea cómo superponer una curva normal a un histograma en R, pero me gustaría conservar el eje y de "frecuencia" normal de un histograma. Vea dos segmentos de código a continuación y observe cómo en el segundo, el eje y se reemplaza por "densidad". ¿Cómo puedo mantener ese eje y como "frecuencia", como está en el primer gráfico?
COMO BONIFICACIÓN: también me gustaría marcar las regiones SD (hasta 3 SD) en la curva de densidad. ¿Cómo puedo hacer esto? Lo intenté abline
, pero la línea se extiende hasta la parte superior del gráfico y se ve fea.
g = d$mydata
hist(g)
g = d$mydata
m<-mean(g)
std<-sqrt(var(g))
hist(g, density=20, breaks=20, prob=TRUE,
xlab="x-variable", ylim=c(0, 2),
main="normal curve over histogram")
curve(dnorm(x, mean=m, sd=std),
col="darkblue", lwd=2, add=TRUE, yaxt="n")
Vea cómo en la imagen de arriba, el eje y es "densidad". Me gustaría que eso fuera "frecuencia".
Aquí hay una manera fácil y agradable que encontré:
h <- hist(g, breaks = 10, density = 10,
col = "lightgray", xlab = "Accuracy", main = "Overall")
xfit <- seq(min(g), max(g), length = 40)
yfit <- dnorm(xfit, mean = mean(g), sd = sd(g))
yfit <- yfit * diff(h$mids[1:2]) * length(g)
lines(xfit, yfit, col = "black", lwd = 2)
Necesita encontrar el multiplicador correcto para convertir la densidad (una curva estimada donde el área debajo de la curva es 1) en recuentos. Esto se puede calcular fácilmente a partir del hist
objeto.
myhist <- hist(mtcars$mpg)
multiplier <- myhist$counts / myhist$density
mydensity <- density(mtcars$mpg)
mydensity$y <- mydensity$y * multiplier[1]
plot(myhist)
lines(mydensity)
Una versión más completa, con una densidad normal y líneas en cada desviación estándar alejada de la media (incluida la media):
myhist <- hist(mtcars$mpg)
multiplier <- myhist$counts / myhist$density
mydensity <- density(mtcars$mpg)
mydensity$y <- mydensity$y * multiplier[1]
plot(myhist)
lines(mydensity)
myx <- seq(min(mtcars$mpg), max(mtcars$mpg), length.out= 100)
mymean <- mean(mtcars$mpg)
mysd <- sd(mtcars$mpg)
normal <- dnorm(x = myx, mean = mymean, sd = mysd)
lines(myx, normal * multiplier[1], col = "blue", lwd = 2)
sd_x <- seq(mymean - 3 * mysd, mymean + 3 * mysd, by = mysd)
sd_y <- dnorm(x = sd_x, mean = mymean, sd = mysd) * multiplier[1]
segments(x0 = sd_x, y0= 0, x1 = sd_x, y1 = sd_y, col = "firebrick4", lwd = 2)