pt_Witkovsky_Tab1: Viktor Witosky's Table_1 pt() Examples

In DPQ: Density, Probability, Quantile ('DPQ') Computations

Viktor Witosky's Table_1 pt() Examples

Description

A data frame with 17 pt() examples from Witosky (2013)'s ‘Table 1’. We provide the results for the FOSS Softwares, additionally including octave's, running the original 2013 matlab code, and the corrected one from 2022.

Usage

data(pt_Witkovsky_Tab1)

Format

A data frame with 17 observations on the following numeric variables.

x

the abscissa, called q in pt().

nu

the positive degrees of freedom, called df in pt().

delta

the noncentrality parameter, called ncp in pt().

true_pnt

“true” values (computed via higher precision, see Witkovsky(2013)).

NCTCDFVW

the pt() values computed with Witkovsky's matlab implementation. Confirmed by using octave (on Fedora 40 Linux). These correspond to our R (package DPQ) pntVW13() values.

Boost

computed via the Boost C++ library; reported by Witkovsky.

R_3.3.0

computed by R version 3.3.0; confirmed to be identical using R 4.4.1

NCT2013_octave_7.3.0

values computed using Witkovsky's original matlab code, by octave 7.3.0

NCT2022_octave_8.4.0

values computed using Witkovsky's 2022 corrected matlab code, by octave 8.4.0

Source

The table was extracted (by MM) from the result of pdftotext --layout <*>.pdf from the publication. The NCT2013_octave_7.3.0 column was computed from the 2013 code, using GPL octave 7.3.0 on Linux Fedora 38, whereas NCT2013_octave_8.4..0 from the 2022 code, using GPL octave 8.4.0 on Linux Fedora 40.

Note that the ‘arXiv’ pre-publication has very slightly differing numbers in its ⁠R⁠ column, e.g., first entry ending in 00200 instead of 00111.

References

Witkovský, Viktor (2013) A Note on Computing Extreme Tail Probabilities of the Noncentral T Distribution with Large Noncentrality Parameter, Acta Universitatis Palackianae Olomucensis, Facultas Rerum Naturalium, Mathematica 52(2), 131–143.

Examples

data(pt_Witkovsky_Tab1)
stopifnot(is.data.frame(d.W <- pt_Witkovsky_Tab1), # shorter
          nrow(d.W) >= 17)
mW <- as.matrix(d.W); row.names(mW) <- NULL # more efficient
colnames(mW)[1:3] #  "x" "nu" "delta"
## use 'R pt() - compatible' names:
(n3 <- names(formals(pt)[1:3]))# "q" "df" "ncp"
colnames(mW)[1:3] <- n3
ptR <- apply(mW[, 1:3], 1, \(a3) unname(do.call(pt, as.list(a3))))
cNm <- paste0("R_", with(R.version, paste(major, minor, sep=".")))
mW <- cbind(mW, `colnames<-`(cbind(ptR), cNm),
            relErr = sfsmisc::relErrV(mW[,"true_pnt"], ptR))
mW
## is current R better than R 3.3.0?  -- or even "the same"?
all.equal(ptR, mW[,"R_3.3.0"])                    # still true in R 4.4.1
all.equal(ptR, mW[,"R_3.3.0"], tolerance = 1e-14) # (ditto)
table(ptR == mW[,"R_3.3.0"]) # {see only 4 (out of 17) *exactly* equal ??}


## How close to published NCTCDFVW is octave's run of the 2022 code?
with(d.W, all.equal(NCTCDFVW, NCT2022_octave_8.4.0, tolerance = 0)) # 3.977e-16

pVW <- apply(unname(mW[, 1:3]), 1, \(a3) unname(do.call(pntVW13, as.list(a3))))
all.equal(pVW, d.W$NCT2013_oct, tolerance = 0)# 2013-based pntVW13() --> 5.6443e-16
all.equal(pVW, d.W$NCT2022_oct, tolerance = 0)

DPQ documentation built on Sept. 11, 2024, 8:37 p.m.