Skip to contents

Cumulative Distribution Function of the Difference of Two t-Distributed Variables by Monte Carlo Simulation

Source: R/pMC2tdiff.R
pMC2tdiff.Rd

This function calculates the cumulative distribution function (CDF) of the difference between two independent t-distributed random variables using Monte Carlo simulation. This method provides a flexible approach that can handle any combination of parameters, with accuracy improving as the number of simulations increases. Specifically, it computes P(T1 - T2 ≤ q) or P(T1 - T2 > q) where T1 and T2 follow t-distributions with potentially different location, scale, and degrees of freedom parameters.

Usage

pMC2tdiff(nMC, q, mu.t1, mu.t2, sd.t1, sd.t2, nu.t1, nu.t2, lower.tail = TRUE)

Arguments

nMC

A positive integer representing the number of Monte Carlo iterations for simulation. Typical values range from 10,000 (for quick estimates) to 100,000+ (for high precision). Larger values yield more accurate results but require more computational time.

q

A numeric value representing the quantile threshold.

mu.t1

A numeric value representing the location parameter (μ) of the first t-distribution.

mu.t2

A numeric value representing the location parameter (μ) of the second t-distribution.

sd.t1

A positive numeric value representing the scale parameter (σ) of the first t-distribution.

sd.t2

A positive numeric value representing the scale parameter (σ) of the second t-distribution.

nu.t1

A positive numeric value representing the degrees of freedom (ν) of the first t-distribution. Must be > 2 for finite variance.

nu.t2

A positive numeric value representing the degrees of freedom (ν) of the second t-distribution. Must be > 2 for finite variance.

lower.tail

A logical value; if TRUE (default), probabilities are P(T1 - T2 ≤ q), otherwise P(T1 - T2 > q).

Value

A numeric value in [0, 1] representing the estimated cumulative probability that the difference between the two t-distributed variables is below (if lower.tail = TRUE) or exceeds (if lower.tail = FALSE) the quantile q. The estimate is subject to Monte Carlo error that decreases as √(1/nMC).

Details

This function uses Monte Carlo simulation to approximate the distribution of the difference between two t-distributed variables. The algorithm consists of:

  • Step 1: Generate nMC random samples from the first t-distribution T1 ~ t(μ₁, σ₁², ν₁)

  • Step 2: Generate nMC random samples from the second t-distribution T2 ~ t(μ₂, σ₂², ν₂)

  • Step 3: Compute the difference D = T1 - T2 for each pair of samples

  • Step 4: Calculate the empirical probability as the proportion of differences that satisfy the condition (D > q or D ≤ q)

Monte Carlo error:

  • The standard error of the estimate is approximately √(p(1-p)/nMC), where p is the true probability

  • For a probability near 0.5, the standard error is roughly 0.5/√nMC

  • For nMC = 10,000: SE ≈ 0.005 (±0.01 with 95% confidence)

  • For nMC = 100,000: SE ≈ 0.0016 (±0.003 with 95% confidence)

Advantages:

  • Highly flexible - works with any parameter combination

  • Easy to understand and implement

  • Naturally handles complex scenarios

  • Provides approximate confidence intervals via bootstrap

Computational considerations:

  • Computational time scales linearly with nMC

  • Suitable for moderate precision requirements

  • For highest precision, consider numerical integration (pNIdifft)

  • Results are stochastic - different runs yield slightly different estimates

Examples

# Calculate P(t1 - t2 > 3) for equal parameters using 10,000 simulations
pMC2tdiff(nMC = 1e4, q = 3, mu.t1 = 2, mu.t2 = 0, sd.t1 = 1, sd.t2 = 1,
         nu.t1 = 17, nu.t2 = 17, lower.tail = FALSE)
#> [1] 0.2576

# Calculate P(t1 - t2 > 1) for unequal variances using 50,000 simulations
pMC2tdiff(nMC = 5e4, q = 1, mu.t1 = 5, mu.t2 = 3, sd.t1 = 2, sd.t2 = 1.5,
         nu.t1 = 10, nu.t2 = 15, lower.tail = FALSE)
#> [1] 0.65012

# Calculate P(t1 - t2 > 0) with high precision using 100,000 simulations
pMC2tdiff(nMC = 1e5, q = 0, mu.t1 = 1, mu.t2 = 0, sd.t1 = 1, sd.t2 = 1,
         nu.t1 = 5, nu.t2 = 20, lower.tail = FALSE)
#> [1] 0.74238

# Calculate lower tail probability P(t1 - t2 ≤ 2)
pMC2tdiff(nMC = 2e4, q = 2, mu.t1 = 3, mu.t2 = 0, sd.t1 = 1.5, sd.t2 = 1.2,
         nu.t1 = 12, nu.t2 = 15, lower.tail = TRUE)
#> [1] 0.31515