wstdiff (version 1.0.0)

ws_tdiff_multivariate_general: Welch-Satterthwaite Approximation for General Multivariate t-Differences

Description

Approximates the distribution of differences between two independent multivariate t-distributed random vectors with arbitrary covariance structure. This implements Theorem 3 from Yamaguchi et al. (2025).

Usage

ws_tdiff_multivariate_general(
  mu1,
  Sigma1,
  nu1,
  mu2,
  Sigma2,
  nu2,
  max_iter = 10,
  tol = 0.001
)

Value

An S3 object of class "ws_tdiff_multivariate_general" containing:

mu_diff

Location vector of difference

Sigma_star

Effective scale matrix

nu_star

Effective degrees of freedom (scalar)

converged

Logical indicating convergence

iterations

Number of iterations performed

method

Character string "multivariate_general"

Arguments

mu1

Location vector of first distribution (length p)

Sigma1

Scale matrix of first distribution (p x p, positive definite)

nu1

Degrees of freedom of first distribution (must be > 4)

mu2

Location vector of second distribution (length p)

Sigma2

Scale matrix of second distribution (p x p, positive definite)

nu2

Degrees of freedom of second distribution (must be > 4)

max_iter

Maximum iterations for convergence (default: 10)

tol

Convergence tolerance (default: 1e-6)

Details

This function handles the general case where components may be correlated within each multivariate t-distribution. The approximation uses a single scalar degrees of freedom parameter to capture the overall tail behavior.

The iterative algorithm (Section 4.3 of the paper):

1. Initialize with sum of covariance matrices

2. Compute effective degrees of freedom using trace formulas

3. Update scale matrix

4. Iterate until convergence

Note: For high dimensions with heterogeneous component behaviors, consider using ws_tdiff_multivariate_independent instead.

Examples

Sigma1 <- matrix(c(1, 0.3, 0.3, 1), 2, 2)
Sigma2 <- matrix(c(1.5, 0.5, 0.5, 1.2), 2, 2)
result <- ws_tdiff_multivariate_general(
  mu1 = c(0, 1), Sigma1 = Sigma1, nu1 = 10,
  mu2 = c(0, 0), Sigma2 = Sigma2, nu2 = 15
)
print(result)