tang-norm-decomp {rotasym}

Distributions based on the tangent-normal decomposition

Description

Density and simulation of a distribution on S^{p-1}:=\{\mathbf{x}\in R^p:||\mathbf{x}||=1\}, p\ge 2, obtained by the tangent-normal decomposition. The tangent-normal decomposition of the random vector \mathbf{X}\in S^{p-1} is

V\boldsymbol{\theta} + \sqrt{1 - V^2}\boldsymbol{\Gamma}_{\boldsymbol{\theta}}\mathbf{U}

where V := \mathbf{X}'\boldsymbol{\theta} is a random variable in [-1, 1] (the cosines of \mathbf{X}) and \mathbf{U} := \boldsymbol{\Gamma}_{\boldsymbol{\theta}}\mathbf{X}/ ||\boldsymbol{\Gamma}_{\boldsymbol{\theta}}\mathbf{X}|| is a random vector in S^{p-2} (the multivariate signs of \mathbf{X}) and \boldsymbol{\Gamma}_{\boldsymbol{\theta}} is the p\times(p-1) matrix computed by Gamma_theta.

The tangent-normal decomposition can be employed for constructing distributions for \mathbf{X} that arise for certain choices of V and \mathbf{U}. If V and \mathbf{U} are independent, then simulation from \mathbf{X} is straightforward using the tangent-normal decomposition. Also, the density of \mathbf{X} at \mathbf{x}\in S^{p-1}, f_\mathbf{X}(\mathbf{x}), is readily computed as

f_\mathbf{X}(\mathbf{x})= \omega_{p-1}c_g g(t)(1-t^2)^{(p-3)/2}f_\mathbf{U}(\mathbf{u})

where t:=\mathbf{x}'\boldsymbol{\theta}, \mathbf{u}:=\boldsymbol{\Gamma}_{\boldsymbol{\theta}}\mathbf{x}/ ||\boldsymbol{\Gamma}_{\boldsymbol{\theta}}\mathbf{x}||, f_\mathbf{U} is the density of \mathbf{U}, and f_V(v) := \omega_{p-1} c_g g(v) (1 - v^2)^{(p-3)/2} is the density of V for an angular function g with normalizing constant c_g. \omega_{p-1} is the surface area of S^{p-2}.

Usage

d_tang_norm(x, theta, g_scaled, d_V, d_U, log = FALSE)

r_tang_norm(n, theta, r_U, r_V)

Arguments

x

locations in S^{p-1} to evaluate the density. Either a matrix of size c(nx, p) or a vector of length p. Normalized internally if required (with a warning message).

theta

a unit norm vector of size p giving the axis of rotational symmetry.

g_scaled

the scaled angular density c_g g. In the form
g_scaled <- function(t, log = TRUE) {...}. See examples.

d_V

the density f_V. In the form d_V <- function(v, log = TRUE) {...}. See examples.

d_U

the density f_\mathbf{U}. In the form d_U <- function(u, log = TRUE) {...}. See examples.

log

flag to indicate if the logarithm of the density (or the normalizing constant) is to be computed.

n

sample size, a positive integer.

r_U

a function for simulating \mathbf{U}. Its first argument must be the sample size. See examples.

r_V

a function for simulating V. Its first argument must be the sample size. See examples.

Details

Either g_scaled or d_V can be supplied to d_tang_norm (the rest of the arguments are compulsory). One possible choice for g_scaled is g_vMF with scaled = TRUE. Another possible choice is the angular function g(t) = 1 - t^2, normalized by its normalizing constant c_g = (\Gamma(p/2) p) / (2\pi^{p/2} (p - 1)) (see examples). This angular function makes V^2 to be distributed as a \mathrm{Beta}(1/2,(p+1)/2).

The normalizing constants and densities are computed through log-scales for numerical accuracy.

Value

Depending on the function:

• d_tang_norm: a vector of length nx or 1 with the evaluated density at x.

• r_tang_norm: a matrix of size c(n, p) with the random sample.

Author(s)

Eduardo García-Portugués, Davy Paindaveine, and Thomas Verdebout.

References

García-Portugués, E., Paindaveine, D., Verdebout, T. (2020) On optimal tests for rotational symmetry against new classes of hyperspherical distributions. Journal of the American Statistical Association, 115(532):1873–1887. doi:10.1080/01621459.2019.1665527

See Also

Gamma_theta, signs, tangent-elliptical, tangent-vMF, vMF.

Examples

## Simulation and density evaluation for p = 2

# Parameters
n <- 1e3
p <- 2
theta <- c(rep(0, p - 1), 1)
mu <- c(rep(0, p - 2), 1)
kappa_V <- 2
kappa_U <- 0.1

# The vMF scaled angular function
g_scaled <- function(t, log) {
  g_vMF(t, p = p - 1, kappa = kappa_V, scaled = TRUE, log = log)
}

# Cosine density for the vMF distribution
d_V <- function(v, log) {
 log_dens <- g_scaled(v, log = log) + (p - 3)/2 * log(1 - v^2)
 switch(log + 1, exp(log_dens), log_dens)
}

# Multivariate signs density based on a vMF
d_U <- function(x, log) d_vMF(x = x, mu = mu, kappa = kappa_U, log = log)

# Simulation functions
r_V <- function(n) r_g_vMF(n = n, p = p, kappa = kappa_V)
r_U <- function(n) r_vMF(n = n, mu = mu, kappa = kappa_U)

# Sample and color according to density
x <- r_tang_norm(n = n, theta = theta, r_V = r_V, r_U = r_U)
r <- runif(n, 0.95, 1.05) # Radius perturbation to improve visualization
col <- viridisLite::viridis(n)
dens <- d_tang_norm(x = x, theta = theta, g_scaled = g_scaled, d_U = d_U)
# dens <- d_tang_norm(x = x, theta = theta, d_V = d_V, d_U = d_U) # The same
plot(r * x, pch = 16, col = col[rank(dens)])

## Simulation and density evaluation for p = 3

# Parameters
p <- 3
n <- 5e3
theta <- c(rep(0, p - 1), 1)
mu <- c(rep(0, p - 2), 1)
kappa_V <- 2
kappa_U <- 2

# Sample and color according to density
x <- r_tang_norm(n = n, theta = theta, r_V = r_V, r_U = r_U)
col <- viridisLite::viridis(n)
dens <- d_tang_norm(x = x, theta = theta, g_scaled = g_scaled, d_U = d_U)
if (requireNamespace("rgl")) {
  rgl::plot3d(x, col = col[rank(dens)], size = 5)
}

## A non-vMF angular function: g(t) = 1 - t^2. It is sssociated to the
## Beta(1/2, (p + 1)/2) distribution.

# Scaled angular function
g_scaled <- function(t, log) {
  log_c_g <- lgamma(0.5 * p) + log(0.5 * p / (p - 1)) - 0.5 * p * log(pi)
  log_g <- log_c_g + log(1 - t^2)
  switch(log + 1, exp(log_g), log_g)
}

# Cosine density
d_V <- function(v, log) {
  log_dens <- w_p(p = p - 1, log = TRUE) + g_scaled(t = v, log = TRUE) +
    (0.5 * (p - 3)) * log(1 - v^2)
  switch(log + 1, exp(log_dens), log_dens)
}

# Simulation
r_V <- function(n) {
  sample(x = c(-1, 1), size = n, replace = TRUE) *
    sqrt(rbeta(n = n, shape1 = 0.5, shape2 = 0.5 * (p + 1)))
}

# Sample and color according to density
r_U <- function(n) r_unif_sphere(n = n, p = p - 1)
x <- r_tang_norm(n = n, theta = theta, r_V = r_V, r_U = r_U)
col <- viridisLite::viridis(n)
dens <- d_tang_norm(x = x, theta = theta, d_V = d_V, d_U = d_unif_sphere)
# dens <- d_tang_norm(x = x, theta = theta, g_scaled = g_scaled,
#                     d_U = d_unif_sphere) # The same
if (requireNamespace("rgl")) {
  rgl::plot3d(x, col = col[rank(dens)], size = 5)
}