Edgeworth Expansions for Triangular Arrays

Why triangular arrays need their own bookkeeping

In a triangular array, the observations in a row are independent but need not share a distribution, and the distributions may change from one row to the next. Consequently, a correction to the normal approximation cannot generally be obtained from one common skewness or kurtosis. The relevant quantities are the cumulants of the whole row sum, standardized by that row’s total variance.

This tutorial develops the first two correction orders for such sums. The formulas have the same Hermite-polynomial form as in the identically distributed case, but their coefficients are assembled from the row’s component cumulants. A separate regularity condition is needed to turn that algebra into a distributional approximation.

Definition: A triangular array is a collection $\{X_{nj}:1\leq j\leq k_n,\ n\geq1\}$, where for each fixed $n$, the variables $X_{n1},\ldots,X_{nk_n}$ are independent. Suppose they have finite means and variances, and write $$ S_n=\sum_{j=1}^{k_n}X_{nj},\qquad \mu_n=E[S_n],\qquad s_n^2=\operatorname{Var}(S_n)=\sum_{j=1}^{k_n}\operatorname{Var}(X_{nj}). $$ When $s_n^2>0$, the standardized sum is $Z_n=(S_n-\mu_n)/s_n$.

Aggregate standardized cumulants

For a centered variable $W$ with finite fourth moment, use the cumulant conventions $\kappa_2(W)=E[W^2]$, $\kappa_3(W)=E[W^3]$, and $\kappa_4(W)=E[W^4]-3(E[W^2])^2$. More generally, a finite-order cumulant is defined from the corresponding derivative of the log characteristic function at zero. Equivalently, with moments through order $r$, its log characteristic function has the Taylor expansion through order $r$, with remainder $o(t^r)$. Moments alone do not justify treating an infinite Taylor series as convergent.

For each row, let $\kappa_{r,nj}$ be the order-$r$ cumulant of $X_{nj}-E[X_{nj}]$. Independence makes cumulants additive. Thus, for $r\geq3$, the standardized cumulant of $Z_n$ is

\lambda_{r,n} =\frac{\sum_{j=1}^{k_n}\kappa_{r,nj}}{s_n^r}.

In particular, $\lambda_{3,n}$ is the aggregate standardized skewness and $\lambda_{4,n}$ is the aggregate standardized fourth cumulant, or excess kurtosis. The denominator is the total row standard deviation raised to the appropriate power; it is not a sum of component standardized cumulants.

Result: If the variables in a row are independent and their cumulants through order $r$ exist, then the order-$r$ cumulant of their centered sum is the sum of their order-$r$ cumulants. If the sum is divided by $s_n$, its order-$r$ cumulant is divided by $s_n^r$.

To verify this, let $\psi_j(t)$ be the characteristic function of the centered $j$th variable. Independence gives the characteristic function of the centered sum as $\prod_j\psi_j(t)$. Near zero each factor is nonzero, so the local logarithm of the product is $\sum_j\log\psi_j(t)$. Taking the order-$r$ derivative at zero adds the cumulants. Replacing the sum by its quotient by $s_n$ replaces $t$ by $t/s_n$, multiplying that derivative by $s_n^{-r}$. This proves both claims. $\square$

Worked Example: Unequal Exponential Components

In row $n$, take $m$ independent exponential variables of rate 1 and $m$ independent exponential variables of rate 2. Center each variable, and let $S$ be their sum. Find the aggregate standardized cumulants through order four and the corresponding formal correction terms.

For an exponential variable of rate $r$, the cumulant of order $q$ is $(q-1)!/r^q$. Centering changes the first cumulant but not those of orders $q\geq2$. Therefore $$ s^2=m\left(1+\frac14\right)=\frac{5m}{4},\qquad \kappa_3(S)=2m\left(1+\frac18\right)=\frac{9m}{4}, $$ and $$ \kappa_4(S)=6m\left(1+\frac1{16}\right)=\frac{51m}{8}. $$

Standardizing gives $$ \lambda_3=\frac{9m/4}{(5m/4)^{3/2}} =\frac{18}{5\sqrt5}\,m^{-1/2},\qquad \lambda_4=\frac{51m/8}{(5m/4)^2} =\frac{102}{25}\,m^{-1}. $$ Thus the density correction coefficients are $\lambda_3/6=3/(5\sqrt5)m^{-1/2}$, $\lambda_4/24=17/100m^{-1}$, and $\lambda_3^2/72=9/250m^{-1}$. Their orders show how the row’s heterogeneous components combine; they are not obtained by averaging the two component skewnesses.

The first two correction orders

Let $\phi(x)=(2\pi)^{-1/2}e^{-x^2/2}$ and let $\Phi(x)$ be its distribution function. The probabilists’ Hermite polynomials needed here are

\begin{aligned} H_2(x)&=x^2-1,& H_3(x)&=x^3-3x,\\ H_4(x)&=x^4-6x^2+3,& H_5(x)&=x^5-10x^3+15x,\\ H_6(x)&=x^6-15x^4+45x^2-15. \end{aligned}

The first density correction is proportional to $H_3$; at the next order, there is a fourth-cumulant term proportional to $H_4$ and a squared-skewness term proportional to $H_6$. The squared term appears when the exponential of the log characteristic function is expanded.

Theorem: Let $a_n\to\infty$, and suppose $s_n>0$, the row variables are independent, and their centered cumulants through order four exist. Suppose also that $\lambda_{3,n}=O(a_n^{-1/2})$, $\lambda_{4,n}=O(a_n^{-1})$, and the following Fourier remainder condition holds, where $\psi_n(t)=E[e^{itZ_n}]$: $$ \int_{\mathbb R} \frac{\left|\psi_n(t)-e^{-t^2/2}\left\{1+\frac{\lambda_{3,n}(it)^3}{6} +\frac{\lambda_{4,n}(it)^4}{24} +\frac{\lambda_{3,n}^2(it)^6}{72}\right\}\right|}{|t|}\,dt =o(a_n^{-1}). $$ Then the CDF of $Z_n$ has the uniform approximation $$ P(Z_n\leq x)=\Phi(x)-\phi(x)\left\{ \frac{\lambda_{3,n}}{6}H_2(x) +\frac{\lambda_{4,n}}{24}H_3(x) +\frac{\lambda_{3,n}^2}{72}H_5(x) \right\}+o(a_n^{-1}). $$

The displayed Fourier condition is a strong, explicit sufficient condition: it controls the integrated difference between the true characteristic function and the proposed truncated expansion. In applications it is typically verified by combining a finite-order cumulant expansion near zero with suitable control of the characteristic function away from zero. It makes clear that cumulant algebra alone does not guarantee a distributional approximation.

To see why the condition gives the stated CDF, the Fourier transform of $\phi(x)H_q(x)$, under the convention $E[e^{itX}]$, is $(it)^q e^{-t^2/2}$. Consequently the expression inside braces in the Fourier condition is the transform of the signed density $$ g_n(x)=\phi(x)\left\{1+\frac{\lambda_{3,n}}6H_3(x) +\frac{\lambda_{4,n}}{24}H_4(x) +\frac{\lambda_{3,n}^2}{72}H_6(x)\right\}. $$ Its integral is one because the nonconstant Hermite terms integrate to zero against $\phi$. Fourier inversion for the difference of the true distribution and this signed approximation bounds the absolute CDF difference uniformly in $x$ by $1/(2\pi)$ times the displayed integral. That bound is $o(a_n^{-1})$.

Finally, $(\phi H_{q-1})'=-\phi H_q$, and $\phi(x)H_{q-1}(x)$ tends to zero as $x$ tends to $-\infty$. Integrating each density term from $-\infty$ to $x$ therefore changes $\phi H_q$ into $-\phi H_{q-1}$. Applying this for $q=3,4,6$ gives the stated CDF formula. $\square$

Worked Example: A Symmetric Bernoulli Triangular Array

In row $m$, take $m$ independent Bernoulli variables with success probability $1/4$ and $m$ with success probability $3/4$. Compute the aggregate standardized cumulants. This is a lattice example, so the calculation gives the smooth Edgeworth polynomial’s cumulant terms, not by itself an ordinary smooth CDF approximation.

For a Bernoulli variable with success probability $p$, the variance, third cumulant, and fourth cumulant are $pq$, $pq(1-2p)$, and $pq(1-6pq)$, where $q=1-p$. At $p=1/4$, these are $3/16,\ 3/32,\ -3/128$. At $p=3/4$, they are $3/16,\ -3/32,\ -3/128$. Adding within the row gives $$ s_m^2=\frac{3m}{8},\qquad \kappa_3(S_m)=0,\qquad \kappa_4(S_m)=-\frac{3m}{64}. $$

Hence $\lambda_{3,m}=0$, while $$ \lambda_{4,m} =\frac{-3m/64}{(3m/8)^2} =-\frac{1}{3m}. $$ The smooth density polynomial would have the fourth-cumulant term $-H_4(x)/(72m)$, and the corresponding formal CDF term is $+\phi(x)H_3(x)/(72m)$. Because the summands are lattice-valued, an actual CDF approximation must account for the lattice correction as well; the cumulant calculation does not remove that requirement.

Worked Example: Heterogeneous Uniform Components

In row $m$, take $m$ independent variables uniform on $[-1,1]$ and $m$ independent variables uniform on $[-2,2]$. Each is centered. Compute the leading nonzero aggregate correction.

A uniform variable on $[-a,a]$ has variance $a^2/3$, third cumulant zero, and fourth cumulant $$ \frac{a^4}{5}-3\left(\frac{a^2}{3}\right)^2=-\frac{2a^4}{15}. $$ Thus the total variance is $m(1/3+4/3)=5m/3$, the aggregate third cumulant is zero, and the aggregate fourth cumulant is $-2m(1+16)/15=-34m/15$. It follows that $$ \lambda_3=0,\qquad \lambda_4=\frac{-34m/15}{(5m/3)^2}=-\frac{102}{125m}. $$

The first nonzero density correction is therefore $-17H_4(x)\phi(x)/(500m)$, and the CDF correction is $17\phi(x)H_3(x)/(500m)$. Unlike the Bernoulli example, this sum has continuous components; under the stated Fourier regularity condition, the smooth CDF expansion applies. The example also shows that zero aggregate skewness removes both the first-order term and the squared-skewness term, but does not remove the kurtosis correction.

Assembling and interpreting the terms

For a practical calculation, first center each component and add its variance to obtain $s_n^2$. Then add the component third and fourth cumulants separately. Divide those totals by $s_n^3$ and $s_n^4$ to obtain $\lambda_{3,n}$ and $\lambda_{4,n}$. Only after that standardization should the coefficients $\lambda_{3,n}/6$, $\lambda_{4,n}/24$, and $\lambda_{3,n}^2/72$ be inserted into the density or CDF polynomial.

The scale assumptions $\lambda_{3,n}=O(a_n^{-1/2})$ and $\lambda_{4,n}=O(a_n^{-1})$ explain the names “first” and “second” correction orders. They describe the size of the aggregate corrections, not a guarantee that every array behaves like an identically distributed sample of size $n$. Unequal variances, skewnesses, and tail behavior can change the aggregate coefficients substantially.

Common pitfall

Pooling component skewnesses or kurtoses directly. The correction uses sums of unstandardized component cumulants divided by the total standard deviation raised to the appropriate power. Averaging component standardized cumulants generally gives the wrong result.
Assuming cumulant algebra guarantees an Edgeworth approximation. Additivity identifies the polynomial coefficients, but a characteristic-function regularity condition is needed to control the distributional error.
Using the smooth CDF formula for a lattice sum without adjustment. A lattice sum can require a lattice correction even when its aggregate cumulants have been calculated correctly.

Key takeaway: For an independent triangular array, build Edgeworth coefficients from the cumulants of the entire row sum, standardize by the row’s total variance, and check regularity separately. The first two orders are determined by aggregate skewness, aggregate fourth cumulant, and the square of aggregate skewness.

Check Your Understanding

Use the definitions and formulas above to answer these questions.

Why is the order-four standardized cumulant divided by $s_n^4$, rather than by $s_n^2$?
If the aggregate third cumulant is zero, which terms in the displayed density expansion disappear?
For independent centered components, what quantities must be added before standardization to obtain $\lambda_{3,n}$ and $\lambda_{4,n}$?
Why does a correct cumulant calculation for a Bernoulli array not by itself justify the smooth CDF expansion?
What role does the Fourier remainder condition play in the theorem?