*! xtflucbreak 1.0.0 07aug2026

Title

xtflucbreak methods -- equation-by-equation derivation of the fluctuation test
and the exact correspondence with Li, Xiao and Chen (2024)

The model and the hypotheses

The heterogeneous panel regression (LXC eq. 1, p.1185) is

yit = xit'{&beta}i + eit,
eit = {&gamma}i'ft + {&epsilon}it,
xit = {&Gamma}i'ft + vit,

with xit of dimension K{&times}1, ft an unobserved m{&times}1 factor, and
{&beta}i = {&beta} + v{&beta},i, v{&beta},i ~ IID(0, {&Sigma}{&beta})
(Assumption 3.3 / 4.7: a random-coefficient model). If the slopes change at an unknown
k0,

yit = xit'({&beta}i + {&delta}i·1{t > k0}) + eit.

H0: {&delta}i = 0 for all i.
HA: {&delta}i {&ne} 0 for i {&isin} {&Pi}, with |{&Pi}|/N {&rarr} c {&isin} (0,1].
The last condition is the panel-unit-root convention of Choi (2001) and
Im-Pesaran-Shin (2003): the break need not be in every panel, but it must be in a
non-vanishing fraction of them.

Timing. The paper's displayed model writes 1{t {&ge} k0}, but
HA states {&delta}i = 0 for t = 1,...,k0 and {&delta}i {&ne} 0 for
t = k0+1,...,T, and the post-break design matrix on p.1189 is
X1i(T) = (0,...,0, xi,k0+1,...,xiT)'. The operative convention is
therefore 1{t > k0}: k0 is the last pre-break period.
r(breakdate) reports that period and r(breakpost) the first post-break one.
The paper's own application is consistent with this reading: on 1996-2020 data it dates
the break at 2008 and attributes it to policies enacted after the 2008 crisis.

Section 3: no common correlated effects

Set {&gamma}i = 0, {&Gamma}i = 0. Per-panel OLS on the full sample and on the
first k observations (LXC eq. 2 and the display below it, p.1187):

{&beta}hati = (Xi'Xi)-1Xi'Yi,
{&beta}hati(k) = (Xi(k)'Xi(k))-1Xi(k)'Yi(k).

With Qhati = Xi'Xi/T and
sigmahati2 = (1/T){&Sigma}t(ehatit - ebari)2
computed from the full-sample residuals, the fluctuation process is

S(k) = N-1/2 {&Sigma}i (1/sigmahati)(k/{&radic}T)·Qhati1/2({&beta}hati(k) - {&beta}hati).

Why this is a bridge. Substituting the OLS formulas,
Qhati1/2(Xi(k)'Xi(k))-1 {&asymp} (1/k)Qhati-1/2,
so

(k/{&radic}T)·Qhati1/2({&beta}hati(k) - {&beta}hati)
{&asymp} T-1/2Qhati-1/2{&Sigma}t{&le}kxiteit
- (k/T)·T-1/2Qhati-1/2{&Sigma}t{&le}Txiteit.

By the functional CLT of Assumptions 3.1-3.2 the first term is Wi(s) with
s = k/T and the second is sWi(1), so the difference is a Brownian bridge
Bi(s). Dividing by sigmahati standardises the variance
{&Omega} = E(xitxit'{&epsilon}it2) = {&sigma}i2Qi
to the identity, so Bi is a standard K-dimensional bridge. Averaging over i
with the N-1/2 normalisation and a second CLT across panels (LXC Appendix, p.1205)
gives Theorem 3.4:

max1<k<T S(k) {&rArr} sup0{&le}s{&le}1 B(s).

Note the algebraic point: the paper writes Qhati1/2, not
Qhati-1/2. That is correct, because it multiplies a matrix that already
contains the inverse (Xi(k)'Xi(k))-1; the net standardiser is
Qhati-1/2, as the display above shows. The code applies the paper's form
literally.

Section 4: common correlated effects

When {&gamma}i {&ne} 0 and {&Gamma}i {&ne} 0, per-panel OLS is inconsistent.
Stack wit = (yit, xit')' and take cross-section averages. Because
wit = Ci'ft + zit with Ci = ({&gamma}i, {&Gamma}i)
· [1, 0; {&beta}i, IK], Baltagi-Feng-Kao (2016) show that under
Assumption 4.4 (rank Cbar = m {&le} K+1)

ft - (CbarCbar')-1Cbar wbart {&rarr}p 0.

Hence wbart spans the factor space asymptotically and

Mw = IT - Wbar(Wbar'Wbar)-1Wbar',
Wbar = (wbar1,...,wbarT)',

annihilates it up to Op(N-1/2). The whole section-3 construction is then
repeated on Ytili = MwYi, Xtili = MwXi (LXC eq. 4), with
Qchki = Xtili'Xtili/T. Theorem 4.8 gives the same Brownian-bridge limit
provided {&radic}T/N {&rarr} 0 -- the rate at which the factor-estimation error must vanish
relative to the T1/2 scaling of the fluctuation process. xtflucbreak prints
{&radic}T/N in the header and warns when it exceeds 0.15.

The intercept. Under Mw a constant column of Xi would be annihilated,
producing an exactly singular Xtili'Xtili. So in the CCE branch the intercept
is absorbed rather than tested, and K equals the number of listed regressors. This matches
the paper's own Model 2, whose DGP has {&alpha}i in the equation but
xit = {&Gamma}i'ft + vit with no constant column and K = 2.

Critical values

Remark 3.5 rejects H0 when, for some component j,

max1<k<T|S(k)(j)| {&ge} C1({&alpha}*),
{&alpha}* = 1 - (1-{&alpha})1/K,

where C1 inverts

P(sup0{&le}u{&le}1|B(u)| {&ge} x) = {&Sigma}k{&ne}0(-1)k+1exp(-2k2x2)
= 2{&Sigma}k{&ge}1(-1)k+1exp(-2k2x2).

This is the Kolmogorov distribution -- the same law as the one-sample
Kolmogorov-Smirnov statistic. {&alpha}* is the Sidak (not Bonferroni) correction: if the K
component tests were independent, controlling each at {&alpha}* controls the family at
exactly {&alpha}. The paper's worked value is reproduced: at {&alpha} = 0.05 and K = 2,
{&alpha}* = 0.025321 and C1 = 1.4781 (the paper prints 1.4782).

xtflucbreak inverts the series by bisection. For x < 1 the alternating series loses
precision, so the theta-function form

P(sup|B| {&le} x) = ({&radic}(2{&pi})/x){&Sigma}k{&ge}1exp(-(2k-1)2{&pi}2/(8x2))

is used there instead. The two agree to machine precision near x = 1.

The reported overall p-value is 1 - (1 - minjpj)K. This is exactly
equivalent to the rule above: minjpj {&le} {&alpha}* if and only if
1 - (1-minjpj)K {&le} {&alpha}.

The change-point estimator

Remark 3.7 defines V(k) = {&Sigma}i(k/sigmahati)Qhati1/2({&beta}hati(k)-{&beta}hati)
and khat = argmax1<k<T||V(k)||, with ||·|| the Euclidean norm of the K-vector
taken after summing over panels. The motivation is the expected-value calculation on
p.1188-1189: |E(V(k)(j))| increases in k up to k0 and decreases thereafter, so
its maximum locates the break.

Note that V(k) = {&radic}(NT)·S(k) identically, so the argmax is invariant to the
scaling and xtflucbreak computes it from S(k) directly. Theorem 3.8 (and 4.12 in the
CCE case) establishes khat {&rarr}p k0 following Proposition 1 of Feng, Kao and
Lazarova (2009).

No confidence interval is reported. The paper's own conclusion states this
explicitly: "this work only discusses the consistency of the change-point estimator,
derive its asymptotic distribution is the next research topic". Anything the command
printed as an interval would be invented. Use xtbfkbreak, which implements the
Baltagi-Feng-Kao estimator and its interval, if you need one.

Step-to-equation map

Every block of the engine against the paper it comes from.

Code Paper

F1 {&beta}hati = (Xi'Xi)-1Xi'Yi LXC eq.(2), p.1187

F2 {&beta}hati(k), recursive OLS LXC p.1187

F3 Qhati = Xi'Xi/T; sigmahati2 = (1/T){&Sigma}(ehat-ebar)2 LXC p.1187

F4 S(k) LXC p.1187

F5 Sidak rule, {&alpha}* = 1-(1-{&alpha})1/K LXC Remark 3.5

F6 Kolmogorov inversion for C1({&alpha}*) LXC Remark 3.5

F7 khat = argmax ||V(k)|| LXC Remark 3.7

F8 Mw from (ybar, xbar) LXC p.1190; BFK 2016 eq.20

F9 {&beta}tili, {&beta}tili(k), Qchki LXC eq.(4), p.1191

F10 Stil(k), rate condition {&radic}T/N {&rarr} 0 LXC Theorem 4.8

A1-A9 benchmark tests, see below Antoch et al. (2018)

Limit theory used, for the record: Theorem 3.4 (null, no CCE), Theorem 3.6 (consistency
under HA), Theorem 3.8 (khat consistent), Theorem 4.8 (null, CCE), Theorem 4.10
(consistency, CCE), Theorem 4.12 (khat consistent, CCE).

Departures from the printed paper

Two places in the article are internally inconsistent. In both, xtflucbreak follows
the theorems, because the theorems are what the critical value is derived from. Both
are switchable so the literal text can be reproduced.

(D1) Where the absolute value goes. Remarks 3.5 and 4.9 typeset

maxk N-1/2{&Sigma}i{(1/sigmahati)(k/{&radic}T)|(·)(j)|} {&ge} C1({&alpha}*),

with |·| inside the sum over i. Taken literally this cannot work: under
H0 each summand converges to |Bi(s)|, a non-negative random variable with
positive mean, so N-1/2{&Sigma}i|Bi(s)| {&rarr} {&infin} and the test
would reject with probability one for any fixed critical value. Theorem 3.4 states the
limit for S(k) without absolute values, and the quoted distribution function is that
of sup|B|, the two-sided bridge. The bars therefore belong outside the sum:

maxk |N-1/2{&Sigma}i{(1/sigmahati)(k/{&radic}T)(·)(j)}|.

That is what the command computes. There is no option to reproduce the literal form,
because it has no valid critical value.

(D2) The missing sigma in section 4. The section-4 statistic is printed as

Stil(k) = N-1/2{&Sigma}i{(k/{&radic}T)Qchki1/2({&beta}tili(k) - {&beta}tili)},

with no 1/sigmatili. But Theorem 4.8 asserts convergence to a standard
Brownian bridge and Remark 4.9 uses the same C1({&alpha}*) as section 3. Both
require the variance standardisation; without it the limit is a bridge with variance
E({&sigma}i2) per component. The omission is nearly invisible in the paper's
own Table 6 because its DGP draws
{&sigma}i2 ~ U(0.5,1.5), so E({&sigma}i) {&asymp} 0.98 and the scaling is
almost exactly one -- which is very likely why it survived. The default therefore applies
1/sigmatili in both branches; nosigmascale reproduces the printed form.

(D3) The finite-sample size problem, and the trimming the paper does not mention.
This is the one that matters most in practice.

S(k) is asymptotically a standard Brownian bridge, so Var(S(k)(j)) {&rarr} s(1-s)
with s = k/T. Its exact variance, conditional on X and under
Var(ei) = {&sigma}i2I, is

Var[(k/{&radic}T)Qhati1/2({&beta}hati(k)-{&beta}hati)]
= (k2/T)·Qhati1/2(Aik-1 - AiT-1)Qhati1/2·{&sigma}i2,
Aik = Xi(k)'Xi(k).

Replacing Aik-1 by (kQi)-1 gives the bridge variance s(1-s), but
that substitution is poor for small k. For Gaussian regressors Aik is Wishart with k
degrees of freedom and E[Aik-1] = {&Sigma}-1/(k-K-1), not
{&Sigma}-1/k. The variance is therefore inflated by a factor

[1/(k-K-1) - 1/(T-K-1)] / [1/k - 1/T],

which at k = T/2 = 25, T = 50, K = 2 equals 1.208. Simulation confirms this exactly (ratio
1.20). Near the start of the grid the inflation is catastrophic: at k = K the subsample fit
is exact and the measured ratio is 2.7{&times}106; at k = K+1 it is 18.5, at K+2 it
is 3.0.

The consequence is that max1<k<T, as printed, is dominated by the first few
points of the grid. Measured empirical size at the nominal 5% level, LXC Model 1, iid
errors, N = T = 50, 500 replications:

trimming(0) 0.674

trimming(0.05) 0.140

trimming(0.10) 0.050

trimming(0.15) 0.050

trimming(0.10) reproduces the paper's Table 1 across all three error processes
(0.050 / 0.070 / 0.062 measured, against 0.050 / 0.081 / 0.069 published) and its Table 2
power (0.988 against 1.000). That agreement is strong evidence that the published results
were produced with roughly 10% trimming, and it is the default here.

Trimming a fixed fraction of T is not enough on short panels, though, because the
blow-up is governed by k relative to K, not to T. At T = 25 -- the length of the
paper's own application -- trimming(0.10) still starts the grid at k = 3, and the
measured size is 0.255 (Model 1) and 0.388 (Model 2). The command therefore standardises
each component by its exact conditional variance,

S*(k)(j) = S(k)(j)·{&radic}( s(1-s) / Vhat(k)(j) ),
Vhat(k) = (1/N){&Sigma}i(k2/T)·diag(Qhati1/2(Aik-1-AiT-1)Qhati1/2),

which requires no distributional assumption beyond the homoskedasticity already imposed by
Assumption 3.1, and satisfies Vhat(k) {&rarr} s(1-s), so S*(k) and S(k) have the same limit
and the same Kolmogorov critical value. Measured size at the nominal 5% level:

literal standardised

Model 1, N = T = 50 0.050 0.038

Model 1, N = 50, T = 25 0.255 0.028

Model 2 (CCE), N = T = 50 0.100 0.062

Model 2 (CCE), N = 200, T = 50 0.098 0.045

In the CCE branch the standardised form is also the closer match to the paper's own
tables: size 0.062 against LXC Table 6's 0.059, and power 0.870 against Table 7's 0.870 for
{&delta} = (0.2,0)' (the literal form gives 0.100 and 0.912). The default is therefore the
standardised statistic; asymptotic restores the literal one, and the literal value
and p-value are printed and returned in every case so any published number stays
recoverable.

khat is left alone. The change-point estimator continues to maximise the raw
||V(k)|| of Remark 3.7, because that is the object Theorems 3.8 and 4.12 prove consistency
for. Measured break-date accuracy is indistinguishable between the two at T {&ge} 50
(fraction within {&plusmn}3 of k0 at N = T = 50: 0.925 raw, 0.930 standardised); only
at T = 25 does the standardised argmax do noticeably better (0.843 against 0.682). That
gap is documented rather than acted on, since departing from a proved result is not
warranted by a difference this size at usable sample sizes.

A fourth, minor one. The Appendix proof of Theorem 3.4 writes
T-1/2{&Sigma}iBi(s) {&rArr} B(s) as N {&rarr} {&infin} (p.1205). The
normalisation must be N-1/2, as it is in the statement of the theorem and everywhere
else. This is a typographical slip with no consequence for the code.

Implementation choices the paper leaves open

The lower end of the search grid. The paper searches 1 < k < T but
{&beta}hati(k) does not exist until Xi(k)'Xi(k) is invertible, which
needs at least K observations and can need more if the early regressor values are
collinear. xtflucbreak finds, for each panel, the first k at which the accumulated
cross-product is of full rank, takes the maximum over panels, and starts the grid
there. The bound is reported in the header. No trimming beyond that is applied by default:
the k/{&radic}T weight already down-weights small k, which is precisely what makes a
fluctuation test usable without a trimming parameter. trimming() is offered for
robustness checks, not as a default.

Which square root. Qhat1/2 is not unique. The default is the symmetric
positive semi-definite root computed from the eigendecomposition, which is the standard
reading of the notation. It has a useful property: permuting the regressors permutes the
components of S(k), because the symmetric root of PQP' is PQ1/2P'. So
maxj -- and therefore the test decision -- does not depend on the order in which
regressors are listed. The Cholesky alternative (cholesky) does not have this
property, but its lower-triangular structure means component 1 loads only on the first
coefficient, component 2 on the first two, and so on.

Components are rotated. Because Qhati1/2 differs across panels, component
j of S(k) is a panel-varying linear combination of coefficient deviations. A rejection on
component j does not identify a regressor. That is a property of the statistic, not of
the implementation. The per-panel shift table reports
{&delta}hati = {&beta}hati(khat+1..T) - {&beta}hati(1..khat) in the
original coefficient space, which is the quantity to quote when discussing which
coefficient moved.

The CCE constant. The paper's Mw is built from Wbar alone. Pesaran's (2006)
CCE augmentation, and Baltagi-Feng-Kao's implementation, include a constant column
alongside the averages; under a slope break this also spans the regime-split factor space
{1, f, f·1{t>k0}}, so both the factor and its
regime-split copy are removed. Omitting the constant leaves the regime-split factor in the
errors. The default therefore includes it, which also makes xtflucbreak and
xtbfkbreak numerically consistent on the same data. nocceconstant reproduces
the literal display.

Balanced panels. S(k) is a sum over i evaluated at a common k. With unequal T
the panels would be evaluated at different points of their own bridges and the aggregation
would be meaningless, so unbalanced data is refused rather than silently truncated.

The benchmark tests

compare implements the three statistics Li, Xiao and Chen benchmark against, from
Antoch, Hanousek, Horvath, Huskova and Wang (2018). Let
Zit = {&Sigma}v{&le}txivxiv',
ehativ = yiv - xiv'{&beta}hatiT, and
sit = {&Sigma}v{&le}txivehativ.

The general Wald-type process is
UN(t) = {&Sigma}i({&beta}hatit-{&beta}hatiT)'Cit({&beta}hatit-{&beta}hatiT).
Since Zit({&beta}hatit - {&beta}hatiT) = sit exactly, the two
weighting matrices used in the paper collapse to score forms:

Wald 1 Cit = ZitZit {&rArr} UN(t) = {&Sigma}i sit'sit (Antoch S1, d = 2)

Wald 2 Cit = ZitZiT-1Zit {&rArr} UN(t) = {&Sigma}i sit'ZiT-1sit (Antoch S2, d = 5)

CUSUM VN(t) = {&Sigma}i{&Sigma}s{&le}tehatis2 (Antoch eq. 2.5)

Each is centred by its mean under H0 (Antoch eq. 3.1-3.4), with
sigchki2 = (1/(T-d)){&Sigma}tehatit2 (eq. 3.10) -- note this
is a different variance estimator from LXC's, and each paper's own is used:

Ahat(1)N(t) = {&Sigma}isigchki2·tr(Cit(Zit-1 - ZiT-1)),

Ahat(2)N(t) = {&Sigma}isigchki2·(t - tr(ZitZiT-1)).

The reported statistic is
maxt|N-1/2(UN(t) - AhatN(t))| over d {&le} t {&le} T-d.

Wild bootstrap (Antoch sec. 4.1-4.2). With qit the panel-i contribution to the
process, put
{&phi}it = qit - N-1{&Sigma}jqjt,
draw {&zeta}i ~ N(0,1) independently across panels, and form
U*N(t) = N-1/2{&Sigma}i{&zeta}i{&phi}it. The bootstrap
statistic is maxt|U*N(t)|; the critical value is its (1-{&alpha}) quantile over
B replications. Because the centring is by the cross-sectional mean, the bootstrap
reproduces the centred observed process automatically -- AhatN(t) does not enter it.
The whole replication is a single matrix product {&Phi}'{&zeta}, so B = 1000 is cheap.

Raw data. The benchmarks are always computed on the untransformed data with an
intercept, even when cce is specified for the fluctuation test. That is deliberate:
it is the configuration in the paper's Tables 6-9, where the Wald and CUSUM tests are
applied to factor-contaminated data because no CCE variant of them exists. Their collapse
in power there is the comparison's whole point.

Author

Dr Merwan Roudane

merwanroudane920@gmail.com

github.com/merwanroudane

Also see

Online: xtflucbreak,
xtflucbreak postestimation,
xtbfkbreak