[This question was formerly called "On Non-normal distributions with zero skewness and zero excess kurtosis" and relabeled to better reflect its focus.]
I am trying to write a little simulation using @Glen_b's very nice answer in this post: Non-normal distributions with zero skewness and zero excess kurtosis?
In particular, I want to show that (what is known in econometrics as) the Jarque-Bera test (which was actually considered earlier in for example Bowman and Shenton (1975). "Omnibus contours for departures from normality based on √b1 and b2". Biometrika 62 (2): 243–250, see @Glen_b's comment) lacks power against non-normal but symmetric distributions without excess kurtosis. The figure shows the distribution of the p-values of the simulation.
For his uniform (left panel) and Poisson (right panel) example I do get a distribution of the p-values that leads to no power beyond size/a conservative test (depending on whether you call JB a test of normality or a test of the two moments), but in the gamma example (middle panel) there does even seem to be some power.
In neither case do I get a uniform distribution of the p-values though, although I (believe to) simulate data under (what I think is) the null - symmetry and no excess kurtosis.
Thoughts on how/why that happens?
CODE:
library(tseries)
library(MASS)
n <- 1e5
a <- sqrt(5+sqrt(24))
b <- (sqrt(13)+1)/2
lambda <- .5
reps <- 1000
JBpval <- matrix(rep(NA,3*reps),ncol=3)
for (i in 1:reps) {
#(a)
u1a <- runif(n/2,-1,1)
u2a <- runif(n/2,-a,a)
#(b)
u1b <- rgamma(n/2,shape = b, scale = 1)
u2b <- -rgamma(n/2,shape = b, scale = 1)
#(c)
u1c <- sqrt(rpois(n/2,lambda = lambda))
u2c <- -sqrt(rpois(n/2,lambda = lambda))
ua <- c(u1a,u2a)
ub <- c(u1b,u2b)
uc <- c(u1c,u2c)
JBpval[i,1] <- jarque.bera.test(ua)$p.value
JBpval[i,2] <- jarque.bera.test(ub)$p.value
JBpval[i,3] <- jarque.bera.test(uc)$p.value
}
par(mfrow=c(1,3))
truehist(JBpval[,1])
truehist(JBpval[,2])
truehist(JBpval[,3])