nu_OPP_estimator {fitHeavyTail}

Estimate the degrees of freedom of a heavy-tailed t distribution based on the OPP estimator

Description

This function estimates the degrees of freedom of a heavy-tailed t distribution based on the OPP estimator from paper [Ollila-Palomar-Pascal, TSP2021, Alg. 1]. Traditional nonparametric methods or likelihood methods provide erratic estimations of the degrees of freedom unless the number of observations is very large. The POP estimator provides a stable estimator based on random matrix theory. A number of different versions are provided, but the default POP method will most likely be the desired choice.

Usage

nu_OPP_estimator(var_X, trace_scatter, r2, method = c("OPP", "OPP-harmonic"))

Arguments

var_X

Vector with the sample variance of the columns of the data matrix.

trace_scatter

Trace of the scatter matrix.

r2

Vector containing the values of diag( Xc %*% inv(scatter) %*% t(Xc) ), where Xc is the centered data matrix.

method

String indicating the version of the OPP estimator (default is just "OPP"). Other option is the variation: "OPP-harmonic".

Value

Estimated value of the degrees of freedom nu of a heavy-tailed t distribution.

Author(s)

Esa Ollila, Frédéric Pascal, and Daniel P. Palomar

References

Esa Ollila, Daniel P. Palomar, and Frédéric Pascal, "Shrinking the Eigenvalues of M-estimators of Covariance Matrix," IEEE Trans. on Signal Processing, vol. 69, pp. 256-269, Jan. 2021. <https://doi.org/10.1109/TSP.2020.3043952>

Examples

library(mvtnorm)       # to generate heavy-tailed data
library(fitHeavyTail)

# parameters
N <- 5
T <- 100
nu_true <- 4           # degrees of freedom
mu_true <- rep(0, N)   # mean vector
Sigma_true <- diag(N)  # scatter matrix

# generate data
X <- rmvt(n = T, sigma = Sigma_true, delta = mu_true, df = nu_true)  # generate Student's t data
mu <- colMeans(X)
Xc <- X - matrix(mu, T, N, byrow = TRUE)    # center data

# usage #1
nu_OPP_estimator(var_X = 1/(T-1)*colSums(Xc^2), trace_scatter = sum(diag(Sigma_true)))

# usage #2
r2 <- rowSums(Xc * (Xc %*% solve(Sigma_true)))
nu_OPP_estimator(var_X = 1/(T-1)*colSums(Xc^2), trace_scatter = sum(diag(Sigma_true)),
                 method = "OPP-harmonic", r2 = r2)