my_p = 0.25; y = rexp(nn <- 1000)
(myexp = qeexp(my_p))
sum(myexp - y[y <= myexp]) / sum(abs(myexp - y)) # Should be my_p
par(mfrow=c(2,1))
yy = seq(-0, 4, len = nn)
plot(yy, deexp(yy), col = "blue", ylim = 0:1, xlab = "y", ylab = "g(y)",
type = "l", main = "g(y) for Exp(1); dotted green is f(y) = dexp(y)")
lines(yy, dexp(yy), col="darkgreen", lty="dotted", lwd=2) # 'original'
plot(yy, peexp(yy), type = "l", col = "blue", ylim = 0:1,
xlab = "y", ylab = "G(y)", main = "G(y) for Exp(1)")
abline(v = 1, h = 0.5, col = "red", lty = "dashed")
lines(yy, pexp(yy), col = "darkgreen", lty = "dotted", lwd = 2)
Run the code above in your browser using DataLab