I need to demonstrate that the probability density function is the derivative of the CDF. Any distribution will do, but I have been trying with the normal. I have got as far as:
set.seed(53)
b <- rnorm(500)
db <- density(b)
plot(db)
Then I can calculate the cumulative probabilities using pnorm(b), but then I don't know how to differentiate, because D() requires an expression rather than pnorm(). Could anyone help, please?
Here's the console scrape from where I demonstrated the near equality (to 5 or 7 decimal places) of the integral of dnorm to pnorm from -Inf to selected values of "x": The Fundamental Theorem of Calculus says that if the integral of a function f(x) is g(x) then f(x) is the derivative of g(x). (Or words to that effect.)
> sapply(c(0,Inf), function(x) integrate(dnorm, lower=-Inf, upper=x))
[,1] [,2]
value 0.5 1
abs.error 4.680562e-05 9.361124e-05
subdivisions 3 3
message "OK" "OK"
call expression expression
> sapply(c(0,Inf), function(x) integrate(dnorm, lower=-Inf, upper=x)$value)
[1] 0.5 1.0
> sapply(seq(-3,3, by=0.5), function(x) integrate(dnorm, lower=-Inf, upper=x)$value)
[1] 0.001349899 0.006209665 0.022750132 0.066807201 0.158655254 0.308537539
[7] 0.500000000 0.691462461 0.841344751 0.933192799 0.977249868 0.993790335
[13] 0.998650102
> pnorm(seq(-3,3, by=0.5)
+ )
[1] 0.001349898 0.006209665 0.022750132 0.066807201 0.158655254 0.308537539
[7] 0.500000000 0.691462461 0.841344746 0.933192799 0.977249868 0.993790335
[13] 0.998650102
I wasn't sure that the D() was "smart" enough to to the symbolic differentiation, but I shouldn't have been so skeptical. This bit of console interaction was done following the examples on the ?deriv help page:
> D(quote(pnorm(x)), "x")
dnorm(x)
Also ... here's something you can get with deriv:
> norm.expr <- expression(pnorm(x))
> deriv(norm.expr, "x")
expression({
.value <- pnorm(x)
.grad <- array(0, c(length(.value), 1L), list(NULL, c("x")))
.grad[, "x"] <- dnorm(x)
attr(.value, "gradient") <- .grad
.value
})