DPQ (version 0.4-3)

dtWV: Noncentral t Distribution Density by W.V.

Description

Compute the density function \(f(x)\) of the t distribution with df degrees of freedom and non-centrality parameter ncp, according to Wolfgang Viechtbauer's proposal in 2002.

Usage

dtWV(x, df, ncp = 0, log = FALSE)

Arguments

x

numeric vector.

df

degrees of freedom (\(> 0\), maybe non-integer). df = Inf is allowed.

ncp

non-centrality parameter \(\delta\); If omitted, use the central t distribution.

log

logical; if TRUE, \(log(f(x))\) is returned instead of \(f(x)\).

Value

numeric vector of density values, properly recycled in (x, df, ncp).

Details

The formula used is “asymptotic”: Resnikoff and Lieberman (1957), p.1 and p.25ff, proposed to use recursive polynomials for (integer !) degrees of freedom \(f = 1,2,\dots, 20\), and then, for df\( = f > 20\), use the asymptotic approximation which Wolfgang Viechtbauer proposed as a first version of a non-central t density for R (when dt() did not yet have an ncp argument).

References

Resnikoff, George J. and Lieberman, Gerald J. (1957) Tables of the non-central t-distribution; Technical report no. 32 (LIE ONR 32), April 1, 1957; Applied Math. and Stat. Lab., Stanford University. https://statistics.stanford.edu/research/tables-non-central-t-distribution-density-function-cumulative-distribution-function-and

See Also

dt, R's (C level) implementation of the (non-central) t density; dntJKBf, for Johnson et al.'s summation formula approximation.

Examples

# NOT RUN {
tt <- seq(0, 10, len = 21)
ncp <- seq(0, 6, len = 31)
dt3R  <- outer(tt, ncp, dt  , df = 3)
dt3WV <- outer(tt, ncp, dtWV, df = 3)
all.equal(dt3R, dt3WV) # rel.err 0.00063
dt25R  <- outer(tt, ncp, dt  , df = 25)
dt25WV <- outer(tt, ncp, dtWV, df = 25)
all.equal(dt25R, dt25WV) # rel.err 1.1e-5

x <- -10:700
fx  <- dt  (x, df = 22, ncp =100)
lfx <- dt  (x, df = 22, ncp =100, log=TRUE)
lfV <- dtWV(x, df = 22, ncp =100, log=TRUE)

head(lfx, 20) # shows that R's dt(*, log=TRUE) implementation is "quite suboptimal"

## graphics
opa <- par(no.readonly=TRUE)
par(mar=.1+c(5,4,4,3), mgp = c(2, .8,0))
plot(fx ~ x, type="l")
par(new=TRUE) ; cc <- c("red", adjustcolor("orange", 0.4))
plot(lfx ~ x, type = "o", pch=".", col=cc[1], cex=2, ann=FALSE, yaxt="n")
sfsmisc::eaxis(4, col=cc[1], col.axis=cc[1], small.args = list(col=cc[1]))
lines(x, lfV, col=cc[2], lwd=3)
dtt1 <- "      dt"; dtt2 <- "(x, df=22, ncp=100"; dttL <- paste0(dtt2,", log=TRUE)")
legend("right", c(paste0(dtt1,dtt2,")"), paste0(c(dtt1,"dtWV"), dttL)),
       lty=1, lwd=c(1,1,3), col=c("black", cc), bty = "n")
par(opa) # reset

# }