Abstract tubes associated with perturbed polyhedra with applications to multidimensional normal probability computations

Let $K$ be a closed convex polyhedron defined by a finite number of linear inequalities. In this paper we refine the theory of abstract tubes (Naiman and Wynn, 1997) associated with $K$ when $K$ is perturbed. In particular, we focus on the perturbation that is lexicographic and in an outer direction. An algorithm for constructing the abstract tube by means of linear programming and its implementation are discussed. Using the abstract tube for perturbed $K$ combined with the recursive integration technique proposed by Miwa, Hayter and Kuriki (2003), we show that the multidimensional normal probability for a polyhedral region $K$ can be computed efficiently. In addition, abstract tubes and the distribution functions of studentized range statistics are exhibited as numerical examples.

💡 Research Summary

The paper addresses the problem of computing the probability that a standard n‑dimensional Gaussian vector falls inside a closed convex polyhedron K defined by a finite set of linear inequalities, P(K)=Pr(x∈K). While for simple cones (the case m=n and A nonsingular) recursive integration methods exist (Miwa, Hayter, and Kuriki, 2003), general polyhedra present two difficulties: (i) the complement of K is a union of half‑space complements whose intersections are not necessarily simple cones, and (ii) the half‑spaces may not be in general position, causing the inclusion‑exclusion expansion to involve many redundant or non‑simple terms.

The authors refine the abstract‑tube framework originally introduced by Naiman and Wynn (1997). An abstract tube is a pair (H,F) where H={H_i^c} are the complements of the defining half‑spaces and F is a simplicial complex of index sets J such that the indicator of the union of H_i^c can be expressed as a signed sum over intersections indexed by J∈F. In the original theory the complex F may be “weak” – it can contain index sets whose associated intersections are not simple cones, especially when the hyperplanes are not in general position.

To overcome this, the paper proposes a lexicographic perturbation of the right‑hand side vector: b(ε)=b+ε·(1,2,…,m)ᵀ with ε>0 infinitesimal. This pushes each hyperplane outward in a prescribed order, guaranteeing that for sufficiently small ε the perturbed hyperplanes are in strict general position. Lemma 2.1 proves that (i) the perturbed hyperplanes are indeed in general position and (ii) the resulting simplicial complex F(ε) stabilises for all sufficiently small ε; the limit is denoted F(0⁺). Consequently, every J∈F(0⁺) consists of linearly independent normal vectors, its cardinality never exceeds rank(A), and the intersection ∩_{i∈J}H_i^c is either an n‑dimensional simple cone or a direct sum of a simple cone with a linear subspace. This property makes each term amenable to the recursive integration algorithm of Miwa et al. (2003).

Theorem 2.1 establishes the strong abstract‑tube identity  1_{∪{i=1}^m H_i^c}(x)=∑{J∈F(0⁺)} (−1)^{|J|−1} 1_{∩_{i∈J} H_i^c}(x) for all x∈ℝⁿ, improving on the earlier weak version that held only almost everywhere. Because |F(0⁺)|≤|F|, the number of terms in the inclusion‑exclusion sum is reduced, and each term is computationally tractable.

The practical construction of F(0⁺) is reduced to a feasibility test for each candidate index set J. The authors formulate a linear programming problem (Problem 3.1) that asks whether there exists x satisfying a_iᵀx=b_i(ε) for i∈J and a_iᵀx≤b_i(ε) for i∉J. By building a tableau that incorporates the perturbation terms and applying Irie’s (1973) pivot algorithm (Algorithm 3.1), feasibility can be decided efficiently. Repeating this test over all subsets yields the complete complex F(0⁺). An R implementation of this algorithm, together with a call to an LP solver, is provided by the authors.

With F(0⁺) in hand, the probability of the polyhedron is obtained via  1−P(K)=∑{J∈F(0⁺)} (−1)^{|J|−1} P(∩{i∈J} H_i^c), where each term P(∩_{i∈J} H_i^c) is the Gaussian probability of a simple cone (or cone plus subspace) and can be evaluated by the recursive integration technique of Miwa et al. This yields an exact, efficient algorithm for P(K) even in dimensions up to about 20, far beyond what naïve Monte‑Carlo would allow.

The paper presents several numerical illustrations. In a three‑dimensional “pyramid” example, the unperturbed complex contains 15 index sets, including a non‑simple 4‑fold intersection; after perturbation, the complex shrinks to 11 sets, all of which correspond to simple cones, dramatically simplifying computation. A redundant‑inequality example shows that different orderings of the perturbation can lead to different complexes, highlighting the importance of the lexicographic scheme. Finally, the method is applied to compute distribution functions of Tukey’s studentized range statistics in multiple‑comparison settings, demonstrating both accuracy and speed compared with existing approximations.

In conclusion, the authors deliver a theoretically solid and computationally practical framework: lexicographic outer perturbation guarantees a strong abstract tube, linear programming provides an automated construction of the associated simplicial complex, and the recursive integration algorithm evaluates each cone probability efficiently. This combination enables exact multidimensional normal probability calculations for arbitrary polyhedral regions and opens avenues for extensions to non‑Gaussian distributions, optimal perturbation ordering, and large‑scale problems.