Learn R Programming

sphunif (version 1.4.3)

Sobolev: Asymptotic distributions of Sobolev statistics of spherical uniformity

Description

Approximated density, distribution, and quantile functions for the asymptotic null distributions of Sobolev statistics of uniformity on \(S^{p-1}:=\{{\bf x}\in R^p:||{\bf x}||=1\}\). These asymptotic distributions are infinite weighted sums of (central) chi squared random variables: $$\sum_{k = 1}^\infty v_k^2 \chi^2_{d_{p, k}},$$ where $$d_{p, k} := {{p + k - 3}\choose{p - 2}} + {{p + k - 2}\choose{p - 2}}$$ is the dimension of the space of eigenfunctions of the Laplacian on \(S^{p-1}\), \(p\ge 2\), associated to the \(k\)-th eigenvalue, \(k\ge 1\).

Usage

d_p_k(p, k, log = FALSE)

weights_dfs_Sobolev(p, K_max = 1000, thre = 0.001, type, Poisson_rho = 0.5, Pycke_q = 0.5, Riesz_s = 1, Rothman_t = 1/3, Sobolev_vk2 = c(0, 0, 1), Softmax_kappa = 1, Stein_cf = FALSE, Stereo_a = 0, log = FALSE, verbose = TRUE, Gauss = TRUE, N = 320, tol = 1e-06, force_positive = TRUE, x_tail = NULL)

d_Sobolev(x, p, type, method = c("I", "SW", "HBE")[1], K_max = 1000, thre = 0.001, Rothman_t = 1/3, Pycke_q = 0.5, Riesz_s = 1, Poisson_rho = 0.5, Softmax_kappa = 1, Stereo_a = 0, Sobolev_vk2 = c(0, 0, 1), ncps = 0, verbose = TRUE, N = 320, x_tail = NULL, ...)

p_Sobolev(x, p, type, method = c("I", "SW", "HBE", "MC")[1], K_max = 1000, thre = 0.001, Rothman_t = 1/3, Pycke_q = 0.5, Riesz_s = 1, Poisson_rho = 0.5, Softmax_kappa = 1, Stereo_a = 0, Sobolev_vk2 = c(0, 0, 1), ncps = 0, verbose = TRUE, N = 320, x_tail = NULL, ...)

q_Sobolev(u, p, type, method = c("I", "SW", "HBE", "MC")[1], K_max = 1000, thre = 0.001, Rothman_t = 1/3, Pycke_q = 0.5, Riesz_s = 1, Poisson_rho = 0.5, Softmax_kappa = 1, Stereo_a = 0, Sobolev_vk2 = c(0, 0, 1), ncps = 0, verbose = TRUE, N = 320, x_tail = NULL, ...)

Value

  • d_p_k: a vector of size length(k) with the evaluation of \(d_{p,k}\).

  • weights_dfs_Sobolev: a list with entries weights and dfs, automatically truncated according to K_max and thre (see details).

  • d_Sobolev: density function evaluated at x, a vector.

  • p_Sobolev: distribution function evaluated at x, a vector.

  • q_Sobolev: quantile function evaluated at u, a vector.

Arguments

p

integer giving the dimension of the ambient space \(R^p\) that contains \(S^{p-1}\).

k

sequence of integer indexes.

log

compute the logarithm of \(d_{p,k}\)? Defaults to FALSE.

K_max

integer giving the truncation of the series that compute the asymptotic p-value of a Sobolev test. Defaults to 1e3.

thre

error threshold for the tail probability given by the the first terms of the truncated series of a Sobolev test. Defaults to 1e-3.

type

name of the Sobolev statistic, using the naming from avail_cir_tests and avail_sph_tests.

Poisson_rho

\(\rho\) parameter for the Poisson test, a real in \([0, 1)\). Defaults to 0.5.

Pycke_q

\(q\) parameter for the Pycke "\(q\)-test", a real in \((0, 1)\). Defaults to 1 / 2.

Riesz_s

\(s\) parameter for the \(s\)-Riesz test, a real in \((0, 2)\). Defaults to 1.

Rothman_t

\(t\) parameter for the Rothman test, a real in \((0, 1)\). Defaults to 1 / 3.

Sobolev_vk2

weights for the finite Sobolev test. A non-negative vector or matrix. Defaults to c(0, 0, 1).

Softmax_kappa

\(\kappa\) parameter for the Softmax test, a non-negative real. Defaults to 1.

Stein_cf

logical indicating whether to use the characteristic function in the Stein test. Defaults to FALSE (moment generating function).

Stereo_a

\(a\) parameter for the Stereo test, a real in \([-1, 1]\). Defaults to 0.

verbose

output information about the truncation? Defaults to TRUE.

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.

tol

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

force_positive

set negative weights to zero? Defaults to TRUE.

x_tail

scalar evaluation point for determining the upper tail probability. If NULL, set to the 0.90 quantile of the whole series, computed by the "HBE" approximation.

x

vector of quantiles.

method

method for approximating the density, distribution, or quantile function of the weighted sum of chi squared random variables. Must be "I" (Imhof), "SW" (Satterthwaite--Welch), "HBE" (Hall--Buckley--Eagleson), or "MC" (Monte Carlo; only for distribution or quantile functions). Defaults to "I".

ncps

non-centrality parameters. Either 0 (default) or a vector with the same length as weights.

...

further parameters passed to *_wschisq.

u

vector of probabilities.

Details

The truncation of \(\sum_{k = 1}^\infty v_k^2 \chi^2_{d_{p, k}}\) is done to the first K_max terms and then up to the index such that the first terms explain the tail probability at the x_tail with an absolute error smaller than thre (see details in cutoff_wschisq). This automatic truncation takes place when calling *_Sobolev. Setting thre = 0 truncates to K_max terms exactly. If the series only contains odd or even non-zero terms, then only K_max / 2 addends are effectively taken into account in the first truncation.

Examples

Run this code
# Circular-specific statistics
curve(p_Sobolev(x = x, p = 2, type = "Watson", method = "HBE"),
      n = 2e2, ylab = "Distribution", main = "Watson")
curve(p_Sobolev(x = x, p = 2, type = "Rothman", method = "HBE"),
      n = 2e2, ylab = "Distribution", main = "Rothman")
curve(p_Sobolev(x = x, p = 2, type = "Pycke_q", method = "HBE"), to = 10,
      n = 2e2, ylab = "Distribution", main = "Pycke_q")
curve(p_Sobolev(x = x, p = 2, type = "Hermans_Rasson", method = "HBE"),
      to = 10, n = 2e2, ylab = "Distribution", main = "Hermans_Rasson")

# Statistics for arbitrary dimensions
test_statistic <- function(type, to = 1, pmax = 5, M = 1e3, ...) {

  col <- viridisLite::viridis(pmax - 1)
  curve(p_Sobolev(x = x, p = 2, type = type, method = "MC", M = M,
                  ...), to = to, n = 2e2, col = col[pmax - 1],
                  ylab = "Distribution", main = type, ylim = c(0, 1))
  for (p in 3:pmax) {
    curve(p_Sobolev(x = x, p = p, type = type, method = "MC", M = M,
                    ...), add = TRUE, n = 2e2, col = col[pmax - p + 1])
  }
  legend("bottomright", legend = paste("p =", 2:pmax), col = rev(col),
         lwd = 2)

}

# Ajne
test_statistic(type = "Ajne")
# \donttest{
# Gine_Gn
test_statistic(type = "Gine_Gn", to = 1.5)

# Gine_Fn
test_statistic(type = "Gine_Fn", to = 2)

# Bakshaev
test_statistic(type = "Bakshaev", to = 3)

# Riesz
test_statistic(type = "Riesz", Riesz_s = 0.5, to = 3)

# PCvM
test_statistic(type = "PCvM", to = 0.6)

# PAD
test_statistic(type = "PAD", to = 3)

# PRt
test_statistic(type = "PRt", Rothman_t = 0.5)

# Quantiles
p <- c(2, 3, 4, 11)
t(sapply(p, function(p) q_Sobolev(u = c(0.10, 0.05, 0.01), p = p,
                                  type = "PCvM")))
t(sapply(p, function(p) q_Sobolev(u = c(0.10, 0.05, 0.01), p = p,
                                  type = "PAD")))
t(sapply(p, function(p) q_Sobolev(u = c(0.10, 0.05, 0.01), p = p,
                                  type = "PRt")))

# Series truncation for thre = 1e-5
sapply(p, function(p) length(weights_dfs_Sobolev(p = p, type = "PCvM")$dfs))
sapply(p, function(p) length(weights_dfs_Sobolev(p = p, type = "PRt")$dfs))
sapply(p, function(p) length(weights_dfs_Sobolev(p = p, type = "PAD")$dfs))
# }

Run the code above in your browser using DataLab