[2] Tests for Multivariate Mean#

The functions in pysht.mean cover classical low-dimensional inference, high-dimensional trace and diagonal tests, unequal-covariance procedures, sparse alternatives, randomized dimension reduction, and multi-group designs. Observations are rows and variables are columns.

Choosing a procedure#

Design and working assumptions	Function	R routine
One sample; Gaussian; invertible covariance	`hotelling_1samp`	`mean1.1931Hotelling`
Two samples; Gaussian; common invertible covariance	`hotelling_2samp`	`mean2.1931Hotelling`
One or two samples; Gaussian; covariance may be singular	`dempster_1samp`, `dempster_2samp`	`mean1.1958Dempster`, `mean2.1958Dempster`
One or two samples; dense high-dimensional shift; common covariance for two samples	`bs_1samp`, `bs_2samp`	`mean1.1996BS`, `mean2.1996BS`
One or two samples; diagonal standardization; common covariance for two samples	`sd_1samp`, `sd_2samp`	`mean1.2008SD`, `mean2.2008SD`
Two samples; unequal covariance; more rows than variables	`yao_2samp`, `johansen_2samp`, `nvm_2samp`, `ky_2samp`	`mean2.1965Yao`, `mean2.1980Johansen`, `mean2.1986NVM`, `mean2.2004KY`
Two samples; high dimension; one Gaussian projection	`ljw_2samp`	`mean2.2011LJW`
Two samples; high dimension; averaged random subspaces	`thulin_2samp`	`mean2.2014Thulin`
Two samples; sparse shift; Bayesian evidence	`lyl_2samp`	`mean2.mxPBF`
Several groups; common covariance	`schott_ksamp`	`meank.2007Schott`
Several groups; unequal covariance; Gaussian Scheffé transformation	`zx_ksamp`	`meank.2009ZX`
Several groups; unequal covariance; high-dimensional factor model	`cph_ksamp`	`meank.2019CPH`

An asymptotic p-value is not a guarantee that a test is suitable for a small dataset. Check the dimensional regime and covariance assumptions in the method ledger before interpreting the result. The reported alternative describes the scientific direction: small p-values are evidence that at least one mean vector differs. lyl_2samp instead returns maximum log Bayes-factor evidence and deliberately has no p-value.

For numerical stability, null vectors are removed before scaling in one-sample routines. Translation-invariant two- and multi-sample routines remove one deterministic feature-wise anchor, shared by every group, in the original input coordinates before choosing a global or per-feature scale. For procedures whose definitions do not split or pair rows, this also makes the deterministic arithmetic insensitive to row or group order and avoids letting a huge common location erase representable variation. The documented CPH split estimator and Zhang–Xu Scheffé transformation are exceptions because row splitting or pairing is part of their definitions. No numerical method can recover information already lost when values were rounded to float64.

Classical and unequal-covariance tests#

pysht.mean.hotelling_1samp(x, *, popmean=None)[source]#

Perform the exact one-sample Hotelling \(T^2\) test.

Exact F calibration requires multivariate normal observations, n > p, and a positive-definite sample covariance matrix.

Parameters:

x (ArrayLike)
popmean (ArrayLike | None)

Return type:

HypothesisTestResult

pysht.mean.hotelling_2samp(x, y, *, paired=False)[source]#

Perform the exact equal-covariance two-sample Hotelling test.

Independent samples are assumed to share a positive-definite covariance matrix. With paired=True, the exact one-sample test is applied to the row-wise differences.

Parameters:

x (ArrayLike)
y (ArrayLike)
paired (bool)

Return type:

HypothesisTestResult

pysht.mean.dempster_1samp(x, *, popmean=None)[source]#

Perform Dempster’s one-sample non-exact mean-vector test.

The test uses a Euclidean mean-square ratio and an approximate F law with trace-estimated effective degrees of freedom. It assumes multivariate normal observations but does not require an invertible covariance matrix.

Parameters:

x (ArrayLike)
popmean (ArrayLike | None)

Return type:

HypothesisTestResult

pysht.mean.dempster_2samp(x, y)[source]#

Perform Dempster’s equal-covariance two-sample non-exact test.

Unlike the legacy SHT routine, the positive mean-square ratio is calibrated with Dempster’s approximate F distribution, not a standard-normal tail.

Parameters:

x (ArrayLike)
y (ArrayLike)

Return type:

HypothesisTestResult

pysht.mean.yao_2samp(x, y)[source]#

Perform Yao’s approximate multivariate Behrens-Fisher test.

Both group covariance matrices may differ. The covariance of the mean difference must be positive definite.

Parameters:

x (ArrayLike)
y (ArrayLike)

Return type:

HypothesisTestResult

pysht.mean.johansen_2samp(x, y)[source]#

Perform Johansen’s Welch-James multivariate mean test.

Parameters:

x (ArrayLike)
y (ArrayLike)

Return type:

HypothesisTestResult

pysht.mean.nvm_2samp(x, y)[source]#

Perform the Nel-Van der Merwe multivariate Behrens-Fisher test.

Parameters:

x (ArrayLike)
y (ArrayLike)

Return type:

HypothesisTestResult

pysht.mean.ky_2samp(x, y)[source]#

Perform the Krishnamoorthy-Yu modified Nel-Van der Merwe test.

Parameters:

x (ArrayLike)
y (ArrayLike)

Return type:

HypothesisTestResult

Validation: classical mean, Dempster, and multivariate Behrens–Fisher.

High-dimensional trace and diagonal tests#

pysht.mean.bs_1samp(x, *, popmean=None)[source]#

Perform the one-sample Bai-Saranadasa trace test.

The denominator uses the unbiased estimator of tr(Sigma @ Sigma); this corrects the missing n + 2 factor in legacy SHT’s one-sample code.

Parameters:

x (ArrayLike)
popmean (ArrayLike | None)

Return type:

HypothesisTestResult

pysht.mean.bs_2samp(x, y)[source]#

Perform the equal-covariance Bai-Saranadasa two-sample test.

Parameters:

x (ArrayLike)
y (ArrayLike)

Return type:

HypothesisTestResult

pysht.mean.sd_1samp(x, *, popmean=None)[source]#

Perform the Srivastava-Du one-sample diagonal mean test.

Parameters:

x (ArrayLike)
popmean (ArrayLike | None)

Return type:

HypothesisTestResult

pysht.mean.sd_2samp(x, y)[source]#

Perform the equal-covariance Srivastava-Du two-sample test.

Parameters:

x (ArrayLike)
y (ArrayLike)

Return type:

HypothesisTestResult

Validation: Bai–Saranadasa and Srivastava–Du.

Randomized and sparse tests#

pysht.mean.ljw_2samp(x, y, *, calibration='asymptotic', n_resamples=999, rng=None)[source]#

Perform the Lopes-Jacob-Wainwright random-projection mean test.

A single Gaussian projection is drawn independently of the data, as in the defining procedure. Permutation calibration, when requested, holds that projection fixed for the observed and every permuted statistic. The pooled observations are canonicalized before seeded randomization, so row reordering or swapping sample labels cannot change a fixed-seed result.

The conditional F calibration assumes Gaussian sampling with a common covariance matrix. Unrestricted permutation calibration is exact under the stronger null that the pooled observations are exchangeable; equality of means alone is not sufficient when the distributions differ.

Parameters:

x (ArrayLike)
y (ArrayLike)
calibration (Literal['asymptotic', 'monte-carlo'])
n_resamples (int)
rng (Generator | int | integer | None)

Return type:

HypothesisTestResult | ResamplingTestResult

calibration="asymptotic" uses the conditional F law for the single projection. calibration="monte-carlo" holds that same projection fixed for the observed statistic and every permutation. The F law assumes Gaussian samples with a common covariance. The permutation option requires pooled exchangeability; equal means alone do not make unrestricted relabeling exact. With an integer seed, canonical pooling makes the full Monte Carlo result invariant to row reordering and sample exchange.

pysht.mean.thulin_2samp(x, y, *, n_subspaces=100, n_resamples=999, rng=None)[source]#

Perform Thulin’s random-subspaces two-sample mean test.

The sampled feature subsets are drawn once and kept fixed throughout the permutation test. This both follows the conditional randomization design and removes the auxiliary-randomness bug in legacy SHT. The pooled rows and sample labels are canonicalized before seeded randomization.

Unrestricted permutation calibration is exact only when observations are exchangeable under the null. In particular, equality of means without a common distribution is not enough for exact finite-sample inference.

Parameters:

x (ArrayLike)
y (ArrayLike)
n_subspaces (int)
n_resamples (int)
rng (Generator | int | integer | None)

Return type:

ResamplingTestResult

The subspace dimension is the paper’s floor((n_x + n_y - 2) / 2) and is reported as a diagnostic. The sampled subspaces are fixed across all permutations. The same pooled-exchangeability requirement applies. Canonical pooling makes fixed-integer-seed results invariant to row reordering and sample exchange.

pysht.mean.lyl_2samp(x, y, *, a0=0.0, b0=0.0, alpha=2.01, gamma=None)[source]#

Compute the Lee-You-Lin maximum pairwise Bayes factor for means.

The returned statistic is the maximum log Bayes factor. Component log Bayes factors are retained, and no frequentist p-value or automatic evidence threshold is manufactured. Residual sums of squares use the ordinary centered maximum-likelihood variances in Equation (4) of the primary paper. By default, gamma=max(n_x+n_y, p)**(-alpha) with the paper’s numerical choice alpha=2.01. A supplied gamma overrides that rate. Nonzero a0 or b0 selects pySHT’s optional shared inverse-gamma extension rather than the published Equation (4).

Parameters:

x (ArrayLike)
y (ArrayLike)
a0 (float)
b0 (float)
alpha (float)
gamma (float | None)

Return type:

BayesFactorTestResult

The primary statistic is the maximum component log Bayes factor. Positive values favor a coordinate-wise alternative over the null, but pySHT does not invent a universal evidence threshold. The published default is alpha=2.01 with gamma=max(n_x+n_y, p)**(-alpha); an explicit gamma overrides it. The paper’s Equation (4) corresponds to a0=b0=0. Nonzero a0 or b0 selects a clearly labeled pySHT inverse-gamma extension. The paper’s two-sample model uses a common covariance matrix.

Validation: Lopes–Jacob–Wainwright, Thulin, Lee–You–Lin.

Validation-blocked SHT identity#

SHT’s mean2.2014CLX (the Cai–Liu–Xia maximum mean test) has no public pySHT callable. Its practical-size Gumbel calibration failed the release gate, and neither estimated-precision branch has a passing advertised sparse-covariance scenario. The private implementation is retained only for scientific audit work; see the Cai–Liu–Xia validation ledger.

Multi-group tests#

pysht.mean.schott_ksamp(*samples)[source]#

Perform Schott’s high-dimensional one-way MANOVA mean test.

The error sum-of-products matrix uses (n_i - 1) S_i for every group. This is the defining MANOVA error matrix; legacy SHT instead used n_i S_i and thereby mismatched its claimed N - k error degrees of freedom.

Parameters:: samples (ArrayLike)
Return type:: HypothesisTestResult

pysht.mean.zx_ksamp(*samples, base_test='bai-saranadasa')[source]#

Perform the Zhang-Xu k-sample Behrens-Fisher mean test.

Scheffe’s transformation reduces the unequal-covariance k-sample problem to one sample. base_test="bai-saranadasa" applies the corrected Bai-Saranadasa calibration; base_test="hotelling" applies the exact Gaussian Hotelling calibration when the transformed dimension permits it.

The paper’s transformation uses the first n_min observations from larger groups. Their row order must therefore be arbitrary with respect to the measurements, as it is for an ordinary i.i.d. sample. If groups tie for the smallest size, their input order fixes the reference group and is part of the realized Scheffe transformation.

Parameters:

samples (ArrayLike)
base_test (Literal['bai-saranadasa', 'hotelling'])

Return type:

HypothesisTestResult

base_test="bai-saranadasa" is the high-dimensional default; base_test="hotelling" requires more transformed observations than transformed variables. The defining Scheffé construction uses the first n_min rows from larger groups. Consequently, changing their row pairing can change the finite-sample statistic even though its null law remains valid for i.i.d. rows. If groups tie for the smallest size, input order selects the reference group; changing the order of tied groups can also change the realized statistic.

pysht.mean.cph_ksamp(*samples, variance_estimator='original')[source]#

Perform the Cao-Park-He high-dimensional k-sample mean test.

variance_estimator selects the paper’s split-sample estimator or the Hu estimator. Both implementations use the current group in every loop; legacy SHT accidentally reused the last group in the split-sample branch.

Parameters:

samples (ArrayLike)
variance_estimator (Literal['original', 'hu'])

Return type:

HypothesisTestResult

variance_estimator="original" uses independent covariance estimates from the paper’s split; "hu" uses the corrected full-sample trace estimator. The split branch requires five rows per group and is intentionally sensitive to row order. Do not sort rows by the outcome before applying it.

Validation: Schott, Zhang–Xu, and Cao–Park–He.