nx <- 10
ny <- 8
shp <- 3
smpx <- rgamma(nx, shape = shp) - shp
smpy <- rgamma(ny, shape = shp)
mx <- mean(smpx)
my <- mean(smpy)
m <- numeric(6)
for (j in 2:6) {
m[j] <- mean(c((smpx - mx)^j, (smpy - my)^j))
}
uM6pool(m[2], m[3], m[4], m[6], nx, ny)
Run the code above in your browser using DataLab