I had some normally distributed data:
mu <- 3
sigma <- 5
x <- rnorm(1e5, mu, sigma)
I took a kernel density estimate with a fairly high bandwidth:
kernel_density_of_x <- density(x, bw = 5)
Then I differentiated it:
differentiate <- function(x, y)
{
diffOfX <- diff(x)
data.frame(
x = x[-length(x)] + (diffOfX / 2),
dyByDx = diff(y) / diffOfX
)
}
first_derivative <- with(kernel_density_of_x, differentiate(x, y))
This looked just as expected:
library(ggplot2)
(p1 <- ggplot(first_derivative, aes(x, dyByDx)) + geom_line())
When I differentiated again, I expected another smooth curve, but I saw an odd cyclical effect.
second_derivative <- with(first_derivative, differentiate(x, dyByDx))
(p2 <- p1 %+% second_derivative + ylab("d2yByDx2"))
I tried a few different options for the kernel
argument, but the noisiness persisted.
Dropping the bandwidth down to, for example, 0.5
gave a lower frequency noise that dominated the plot (making it nonsense).
Dropping the number of sampling points down from n = 512
to n = 32
stopped the issue, but that causes other problems.
Why does this effect occur? Is it an artifact of the density
function, or have I done something silly?
We can redraw the plot using the probability density function of the normal distribution that x
was generated from to see the shape I expected:
xx <- seq.int(-20, 20, 0.1)
pdf_of_xx <- dnorm(xx, mu, sigma)
first_derivative_of_xx <- differentiate(xx, pdf_of_xx)
second_derivative_of_xx<- with(first_derivative_of_xx, differentiate(x, dyByDx))
ggplot(second_derivative_of_xx, aes(x, dyByDx)) + geom_line()
dnorm
. $\endgroup$R
code; the discrete second derivative operator greatly magnifies the small errors. This has been discussed in another thread which I cannot at the moment find ... When you think about it, smoothness of the derivative (of a kernel density estimate as computed inR
) isn't even defined in this context, because everything in sight is a discrete set of points, anyway. $\endgroup$