Ernest Adrogué
2009-Nov-09 15:56 UTC
[R] unexpected results involving the skellam distribution
Hi all, I am new to R. I made this script that plots the deviation of a skellam-distributed random variable with respect to the skellam distribution. I would expect to get random errors, but the plot sistematically shows a non-random pattern (first a peak and then a low). I don't know how to explain this. Can somebody enlighten me?? require(skellam) n <- 5000 lam1 <- 5.45 lam2 <- 2.78 p1 <- rpois(n, lam1) p2 <- rpois(n, lam2) rho <- cor(p1,p2) z <- hist(p1-p2, plot=FALSE) mu1 <- lam1 - (rho*sqrt(lam1*lam2)) mu2 <- lam2 - (rho*sqrt(lam1*lam2)) x <- z$breaks[1:length(z$breaks)-1] y <- dskellam(x, mu1, mu2) plot(x,z$density-y) -- Cheers. Ernest