Está es una la solución para sombrear el área bajo una curva entre \(\pm\) una desviación estándar (DE).

Primero generamos la gráfica para la densidad de edadt2 y añadimos líneas para ubicar la media, mediana y \(\pm\) una DE.

edadt2 <- c(20, 22, 23, 22, 24, 25, 26)
plot(density(edadt2), xlab = "Edad", ylab = "Densidad", main = "Densidad para edadt2")
abline(v = mean(edadt2), lty = 2, col = "darkgreen", lwd = 2)
abline(v = median(edadt2), col = "darkred", lwd = 2)
text(24, 0.17, "Media", col = "darkgreen")
text(22, 0.17, "Mediana", col = "darkred")
sdmas <- mean(edadt2) + sd(edadt2)  #Valor para +1 DE
sdmenos <- mean(edadt2) - sd(edadt2)  #Valor para -1 DE
abline(v = sdmas)  #Añadimos la línea de +1 DE
abline(v = sdmenos)  #Añadimos la línea de -1 DE

Tenemos los valores de sdmenos y sdmas que corresponden a la media menos una DE y la media más una DE, necesitamos la densidad de edadt2 en un objeto que nombraremos densidad.

densidad <- density(edadt2)

Ahora recortamos la densidad que está entre sdmenos y sdmas.

minimo <- min(which(densidad$x >= sdmenos))
maximo <- max(which(densidad$x <= sdmas))

Por último damos las coordenadas en la función polygon() tanto x como y deben tener al menos cuatro valores para el polígono, para este caso en particular estamos incluyendo un continuo de puntos de la densidad y señalando los valores mínimos y máximos para x y y respectivamente.

plot(density(edadt2), xlab = "Edad", ylab = "Densidad", main = "Densidad para edadt2")
abline(v = mean(edadt2), lty = 2, col = "darkgreen", lwd = 2)
abline(v = median(edadt2), col = "darkred", lwd = 2)
text(24, 0.17, "Media", col = "darkgreen")
text(22, 0.17, "Mediana", col = "darkred")
abline(v = sdmas)
abline(v = sdmenos)
abline(h = 0)
with(densidad, polygon(x = c(x[c(minimo, minimo:maximo, maximo)]), y = c(0, 
    y[minimo:maximo], 0), col = rgb(0, 206/255, 209/255, 0.5)))