dat2bernquaf(6, c(2,10)) # median 1/2 of 2 and 10 is 6 (trivial and fast)
if (FALSE) {
set.seed(5135)
lmr <- vec2lmom(c(1000, 400, 0.2, 0.3, 0.045))
par <- lmom2par(lmr, type="wak")
Q <- rlmomco(83, par) # n = 83 and extremely non-Normal data
lgQ <- max(Q) # 5551.052 by theory
dat2bernquaf(median(Q), Q)$f # returns 0.5100523 (nearly 1/2)
dat2bernquaf(lgQ, Q)$f # unable to root
dat2bernquaf(lgQ, Q, bound.type="sd")$f # unable to root
itf <- c(0.5, 0.99999)
f <- dat2bernquaf(lgQ, Q, interval=itf, bound.type="sd")$f
print(f) # F=0.9961118
qlmomco(f, par) # 5045.784 for the estimate F=0.9961118
# If we were not using the maximum and something more near the center of the
# distribution then that estimate would be closer to qlmomco(f, par).
# You might consider lqQ <- qlmomco(0.99, Q) # theoretical 99th percentile and
# let the random seed wander and see the various results.
}
Run the code above in your browser using DataLab