Pn function

Utilities for projected-ecdf statistics of spherical uniformity

Computation of the kernels [REMOVE_ME]$\psi_p^W(\theta) := \int_{-1}^1 A_x(\theta)\,\mathrm{d}W(F_p(x)),\psi_p^W(\theta) := \int_{-1}^1 A_x(\theta) dW(F_p(x)), [REMOVE_ME_2]$

where $A_x(\theta)$ is the proportion of area surface of $S^{p - 1}$ covered by the intersection of two hyperspherical caps with common solid angle $\pi - \cos^{-1}(x)$ and centers separated by an angle $\theta \in [0, \pi]$, $F_p$ is the distribution function of the projected spherical uniform distribution , and $W$ is a measure on $[0, 1]$.

Also, computation of the Gegenbauer coefficients of $\psi_p^W$: [REMOVE_ME]$b_{k, p}^W := \frac{1}{c_{k, p}}\int_0^\pi \psi_p^W(\theta)C_k^{p / 2 - 1}(\cos\theta)\,\mathrm{d}\theta.b_{k, p}^W := \frac{1}{c_{k, p}} \int_0^\pi \psi_p^W(\theta)C_k^(p / 2 - 1)(\cos\theta) d\theta. [REMOVE_ME_2]$

These coefficients can also be computed via [REMOVE_ME]$b_{k, p}^W = \int_{-1}^1 a_{k, p}^x\,\mathrm{d}W(F_p(x))b_{k, p}^W = \int_{-1}^1 a_{k, p}^x dW(F_p(x)) [REMOVE_ME_2]$

for a certain function $x\rightarrow a_{k, p}^x$. They serve to define projected alternatives to uniformity .

psi_Pn(theta, q, type, Rothman_t = 1/3, tilde = FALSE, psi_Gauss = TRUE,
  psi_N = 320, tol = 1e-06)

Gegen_coefs_Pn(k, p, type, Rothman_t = 1/3, Gauss = TRUE, N = 320,
  tol = 1e-06, verbose = FALSE)

akx(x, p, k, sqr = FALSE)

f_locdev_Pn(p, type, K = 1000, N = 320, K_max = 10000, thre = 0.001,
  Rothman_t = 1/3, verbose = FALSE)

Arguments

• theta: vector with values in $[0, \pi]$.

• q: integer giving the dimension of the sphere $S^q$.

• type: type of projected-ecdf test statistic. Must be either "PCvM" (Cramér--von Mises), "PAD" (Anderson--Darling), or "PRt" (Rothman).

• Rothman_t: $t$ parameter for the Rothman test, a real in $(0, 1)$. Defaults to 1 / 3.

• tilde: include the constant and bias term? Defaults to FALSE.

• psi_Gauss: use a Gauss-Legendre quadrature

  rule with psi_N nodes in the computation of the kernel function? Defaults to TRUE.

• psi_N: number of points used in the Gauss--Legendre quadrature for computing the kernel function. Defaults to 320.

• tol: tolerance passed to integrate's rel.tol and abs.tol if Gauss = FALSE. Defaults to 1e-6.

• k: vector with the index of coefficients.

• p: integer giving the dimension of the ambient space $R^p$ that contains $S^{p-1}$.

• Gauss: use a Gauss--Legendre quadrature rule of N nodes in the computation of the Gegenbauer coefficients? Otherwise, call integrate. Defaults to TRUE.

• N: number of points used in the Gauss-Legendre quadrature for computing the Gegenbauer coefficients. Defaults to 320.

• verbose: flag to print informative messages. Defaults to FALSE.

• x: evaluation points for $a_{k, p}^x$, a vector with values in $[-1, 1]$.

• sqr: return the signed square root of $a_{k, p}^x$? Defaults to FALSE.

• K: number of equispaced points on $[-1, 1]$ used for evaluating $f$ and then interpolating. Defaults to 1e3.

• K_max: integer giving the truncation of the series. Defaults to 1e4.

• thre: proportion of norm not explained by the first terms of the truncated series. Defaults to 1e-3.

Returns

• psi_Pn: a vector of size length(theta) with the evaluation of $\psi$.

• Gegen_coefs_Pn: a vector of size length(k) containing the coefficients $b_{k, p}^W$.

• akx: a matrix of size c(length(x), length(k))

  containing the coefficients $a_{k, p}^x$.

• f_locdev_Pn: the projected alternative $f$ as a function ready to be evaluated.

Description

Computation of the kernels

where $A_x(\theta)$ is the proportion of area surface of $S^{p - 1}$ covered by the intersection of two hyperspherical caps with common solid angle $\pi - \cos^{-1}(x)$ and centers separated by an angle $\theta \in [0, \pi]$, $F_p$ is the distribution function of the projected spherical uniform distribution , and $W$ is a measure on $[0, 1]$.

Also, computation of the Gegenbauer coefficients of $\psi_p^W$:

These coefficients can also be computed via

for a certain function $x\rightarrow a_{k, p}^x$. They serve to define projected alternatives to uniformity .

Details

The evaluation of $\psi_p^W$ and $b_{k, p}^W$ depends on the type of projected-ecdf statistic:

• PCvM: closed-form expressions for $\psi_p^W$ and $b_{k, p}^W$

  with $p = 2, 3, 4$, numerical integration required for $p \ge 5$.

• PAD: closed-form expressions for $\psi_2^W$ and $b_{k, 3}^W$, numerical integration required for $\psi_p^W$ with $p \ge 3$ and $b_{k, p}^W$ with $p = 2$ and $p \ge 4$.

• PRt: closed-form expressions for $\psi_p^W$ and $b_{k, p}^W$

  for any $p \ge 2$.

See García-Portugués et al. (2023) for more details.

Examples

# Kernels in the projected-ecdf test statistics
k <- 0:10
coefs <- list()
(coefs$PCvM <- t(sapply(2:5, function(p)
  Gegen_coefs_Pn(k = k, p = p, type = "PCvM"))))
(coefs$PAD <- t(sapply(2:5, function(p)
  Gegen_coefs_Pn(k = k, p = p, type = "PAD"))))
(coefs$PRt <- t(sapply(2:5, function(p)
  Gegen_coefs_Pn(k = k, p = p, type = "PRt"))))

# Gegenbauer expansion
th <- seq(0, pi, length.out = 501)[-501]
old_par <- par(mfrow = c(3, 4))
for (type in c("PCvM", "PAD", "PRt")) {

  for (p in 2:5) {

    plot(th, psi_Pn(theta = th, q = p - 1, type = type), type = "l",
         main = paste0(type, ", p = ", p), xlab = expression(theta),
         ylab = expression(psi(theta)), axes = FALSE, ylim = c(-1.5, 1))
    axis(1, at = c(0, pi / 4, pi / 2, 3 * pi / 4, pi),
         labels = expression(0, pi / 4, pi / 2, 3 * pi / 4, pi))
    axis(2); box()
    lines(th, Gegen_series(theta = th, coefs = coefs[[type]][p - 1, ],
                           k = k, p = p), col = 2)

  }

}
par(old_par)

# Analytical coefficients vs. numerical integration
test_coef <- function(type, p, k = 0:20) {

  plot(k, log1p(abs(Gegen_coefs_Pn(k = k, p = p, type = type))),
       ylab = "Coefficients", main = paste0(type, ", p = ", p))
  points(k, log1p(abs(Gegen_coefs(k = k, p = p, psi = psi_Pn, type = type,
                                  q = p - 1))), col = 2)
  legend("topright", legend = c("log(1 + Gegen_coefs_Pn))",
                                "log(1 + Gegen_coefs(psi_Pn))"),
         lwd = 2, col = 1:2)

}

# PCvM statistic
old_par <- par(mfrow = c(2, 2))
for (p in 2:5) test_coef(type = "PCvM", p = p)
par(old_par)

# PAD statistic
old_par <- par(mfrow = c(2, 2))
for (p in 2:5) test_coef(type = "PAD", p = p)
par(old_par)

# PRt statistic
old_par <- par(mfrow = c(2, 2))
for (p in 2:5) test_coef(type = "PRt", p = p)
par(old_par)

# akx
akx(x = seq(-1, 1, l = 5), k = 1:4, p = 2)
akx(x = 0, k = 1:4, p = 3)

# PRt alternative to uniformity
z <- seq(-1, 1, l = 1e3)
p <- c(2:5, 10, 15, 17)
col <- viridisLite::viridis(length(p))
plot(z, f_locdev_Pn(p = p[1], type = "PRt")(z), type = "s",
     col = col[1], ylim = c(0, 0.6), ylab = expression(f[Rt](z)))
for (k in 2:length(p)) {
  lines(z, f_locdev_Pn(p = p[k], type = "PRt")(z), type = "s", col = col[k])
}
legend("topleft", legend = paste("p =", p), col = col, lwd = 2)

References

García-Portugués, E., Navarro-Esteban, P., Cuesta-Albertos, J. A. (2023) On a projection-based class of uniformity tests on the hypersphere. Bernoulli, 29(1):181--204. tools:::Rd_expr_doi("10.3150/21-BEJ1454") .

Author(s)

Eduardo García-Portugués and Paula Navarro-Esteban.