Introduction to mmtdiff

library(mmtdiff)

Overview

The mmtdiff package implements the moment-matching approximation for differences of non-standardized t-distributed random variables.

Background

This package apply the moment-matching approximation to t-distribution differences, which arise in:

Univariate Example

Consider two independent t-distributed variables with different parameters:

# X1 ~ t(mu=0, sigma=1, nu=10)
# X2 ~ t(mu=0, sigma=1.5, nu=15)
result <- mm_tdiff_univariate(
  mu1 = 0, sigma1 = 1, nu1 = 10,
  mu2 = 0, sigma2 = 1.5, nu2 = 15
)

print(result)
#> Moment-Matching Approximation (Univariate)
#> ==============================================
#> Location (mu1 - mu2): 0.0000
#> Scale (sigma*): 1.8652
#> df (nu*): 20.9421

# The difference Z = X1 - X2 is approximately:
# Z ~ t(mu_diff, sigma_star^2, nu_star)

Using Distribution Functions

Once we have the approximation, we can compute probabilities and quantiles:

# Probability that difference is negative
p_negative <- ptdiff(0, result)
print(paste("P(X1 - X2 < 0) =", round(p_negative, 4)))
#> [1] "P(X1 - X2 < 0) = 0.5"

# 95% confidence interval for the difference
ci_95 <- qtdiff(c(0.025, 0.975), result)
print(paste("95% CI:", round(ci_95[1], 3), "to", round(ci_95[2], 3)))
#> [1] "95% CI: -3.879 to 3.879"

# Generate random samples
samples <- rtdiff(1000, result)
hist(samples, breaks = 30, main = "Simulated Difference Distribution",
     xlab = "X1 - X2", probability = TRUE)

# Overlay theoretical density
x_seq <- seq(min(samples), max(samples), length.out = 100)
lines(x_seq, dtdiff(x_seq, result), col = "red", lwd = 2)

Multivariate Example

For multiple independent components:

result_multi <- mm_tdiff_multivariate_independent(
  mu1 = c(0, 1, 2),
  sigma1 = c(1, 1.5, 2),
  nu1 = c(10, 12, 15),
  mu2 = c(0, 0, 0),
  sigma2 = c(1.2, 1, 1.8),
  nu2 = c(15, 20, 25)
)

print(result_multi)
#> Moment-Matching Approximation (Multivariate Independent)
#> =============================================================
#> Location difference:
#> [1] 0 1 2
#> 
#> Scale:
#> [1] 1.621278 1.844829 2.756493
#> 
#> df:
#> [1] 20.57649 18.69487 30.20081