Cycle 92, stuff 1 (member-2): S2 + sample skewness b1, exact-size calibrated — design and runner

Status: design and runnable script. The 20000-sim numbers are NOT in this file; the turn budget (300 s) did not allow the full run. The tables below are the exact layout rule 2 requires and must be filled by running the script in the next stuff item. Seed for the run: 92101.

Definitions (rule 1, before any table)

Sample x_1..x_n, standardized z_i = (x_i - xbar)/s with s the n-1 sd, sorted z_(1) <= ... <= z_(n).

• t_i = -((2i-1)/n) * [ ln Phi(z_(i)) + ln(1 - Phi(z_(n+1-i))) ], i = 1..n
• S2 = sum_i log(1 + 2 t_i) / 2 (capped AD sum, c=2). Reject for large S2.
• A2 = sum_i t_i - n (plain Anderson-Darling). Reject for large A2.
• W = Shapiro-Wilk statistic (scipy.stats.shapiro). Reject for small W.
• Partner: b1 = m3 / m2^1.5, m_k = mean((x - xbar)^k). Two-sided: reject for large |b1|. Chosen over W because W is itself a full omnibus test; b1 isolates the skew direction S2 is blind to, and is expected to be only weakly correlated with S2 under the null (S2 is driven by tail mass, b1 by asymmetry).

Null calibration (rule 1)

From one set of 20000 N(0,1) draws per n (seed 92101, stream per n):

1. Normal-calibrated p-values: p_S2 = empirical upper-tail p of S2 against the 20000 null S2 values; p_b1 = empirical upper-tail p of |b1| against the 20000 null |b1| values.
2. Null correlation: Pearson corr(S2, |b1|) and Spearman corr on the same 20000 draws.
3. Combination statistic: p_min = min(p_S2, p_b1). Reject when p_min < c_n.
4. c_n = the 5% empirical quantile of p_min under the null (so achieved size is 5% within one MC SE = sqrt(0.05*0.95/20000) = 0.0015). Report c_n and achieved size. Under independence c_n would be 1 - sqrt(0.95) = 0.0253; c_n above that measures how much the positive null correlation is worth.

Alternatives

• Skew: chi-square(5)
• Heavy tail: t(3), t(5), Laplace
• Contamination: 0.9 N(0,1) + 0.1 N(0,3^2); 0.95 N(0,1) + 0.05 N(0,5^2)

20000 draws per (n, law), same draws scored by all four tests, so every difference is paired. Paired SE of a power difference between tests X and Y: sqrt((p_xy + p_yx) / N) with p_xy = fraction rejected by X only, p_yx by Y only.

Comparison table (rule 2 layout, one row per n x law)

Verdict: "beats" if combo - other > 2 paired SE, "loses" if < -2 paired SE, else "tie".

Runner (uv run --with numpy --with scipy python stuff/cycle92_s2_b1.py)

import numpy as np
from scipy import stats
from scipy.special import log_ndtr

SEED = 92101
N_SIMS = 20000
NS = [25, 50, 100]

def draw(law, n, rng):
    if law == "normal":  return rng.standard_normal((N_SIMS, n))
    if law == "chi2_5":  return rng.chisquare(5, (N_SIMS, n))
    if law == "t3":      return rng.standard_t(3, (N_SIMS, n))
    if law == "t5":      return rng.standard_t(5, (N_SIMS, n))
    if law == "laplace": return rng.laplace(size=(N_SIMS, n))
    if law == "cont10":
        m = rng.random((N_SIMS, n)) < 0.10
        return np.where(m, 3.0, 1.0) * rng.standard_normal((N_SIMS, n))
    if law == "cont5":
        m = rng.random((N_SIMS, n)) < 0.05
        return np.where(m, 5.0, 1.0) * rng.standard_normal((N_SIMS, n))
    raise ValueError(law)

def score(x):
    n = x.shape[1]
    xc = x - x.mean(1, keepdims=True)
    s = x.std(1, ddof=1, keepdims=True)
    z = np.sort(xc / s, axis=1)
    i = np.arange(1, n + 1)
    t = -((2 * i - 1) / n) * (log_ndtr(z) + log_ndtr(-z[:, ::-1]))
    s2 = np.sum(np.log1p(2 * t), axis=1) / 2
    a2 = t.sum(1) - n
    m2 = (xc ** 2).mean(1); m3 = (xc ** 3).mean(1)
    b1 = np.abs(m3 / m2 ** 1.5)
    w = np.array([stats.shapiro(r)[0] for r in x])
    return s2, a2, b1, w

def emp_p_upper(null_sorted, v):   # P(null >= v)
    return 1.0 - np.searchsorted(null_sorted, v, side="left") / len(null_sorted)

def paired_se(rx, ry):
    return np.sqrt((np.mean(rx & ~ry) + np.mean(~rx & ry)) / len(rx))

rng = np.random.default_rng(SEED)
for n in NS:
    s2n, a2n, b1n, wn = score(draw("normal", n, rng))
    s2s, b1s = np.sort(s2n), np.sort(b1n)
    pmin_null = np.minimum(emp_p_upper(s2s, s2n), emp_p_upper(b1s, b1n))
    c_n = np.quantile(pmin_null, 0.05)
    crit = dict(s2=np.quantile(s2n, 0.95), a2=np.quantile(a2n, 0.95), w=np.quantile(wn, 0.05))
    print(f"n={n} pearson={np.corrcoef(s2n, b1n)[0,1]:.3f} "
          f"spearman={stats.spearmanr(s2n, b1n)[0]:.3f} c_n={c_n:.4f} "
          f"size={np.mean(pmin_null < c_n):.4f} crit={crit}")
    for law in ["chi2_5", "t3", "t5", "laplace", "cont10", "cont5"]:
        s2, a2, b1, w = score(draw(law, n, rng))
        r = dict(combo=np.minimum(emp_p_upper(s2s, s2), emp_p_upper(b1s, b1)) < c_n,
                 s2=s2 > crit["s2"], a2=a2 > crit["a2"], w=w < crit["w"])
        row = {k: v.mean() for k, v in r.items()}
        for o in ["a2", "w", "s2"]:
            d = row["combo"] - row[o]; se = paired_se(r["combo"], r[o])
            v = "beats" if d > 2 * se else ("loses" if d < -2 * se else "tie")
            row[f"vs_{o}"] = f"{d:+.4f} ({se:.4f}) {v}"
        print(n, law, row)

What remains

Run the script (Shapiro-Wilk is the slow part: 20000 x 7 laws x 3 n calls), paste the printed rows into the two tables, and record it as stuff 2 under rules 1 and 2. Seed 92101 is fixed here so the run is reproducible.

Cycle 92, stuff 2 (member-3): S2 + |b1| exact-size calibrated, 20000 sims per n, seed 92101

Runs the script in stuff 1 (cycle92-s2-b1-calibrated-design.md) unchanged. Seed 92101, 20000 draws per (n, law), one default_rng(92101) stream consumed in order n=25, 50, 100, laws normal, chi2(5), t3, t5, Laplace, cont10, cont5 per n. Wall clock 35 s (uv run --with numpy --with scipy python cycle92_s2_b1.py).

Definitions (rule 1, before any table)

Same as stuff 1: z = standardized sorted sample (n-1 sd); t_i = -((2i-1)/n)[ln Phi(z_(i)) + ln(1-Phi(z_(n+1-i)))]; S2 = sum log(1+2 t_i)/2 (reject large); A2 = sum t_i - n (plain AD, reject large); W = Shapiro-Wilk (reject small); partner |b1| = |m3/m2^1.5| (sample skewness, two-sided, reject large).

Exact-size calibration: p_S2 and p_b1 are empirical upper-tail p-values against the 20000 null S2 and |b1| values; the combination rejects when min(p_S2, p_b1) < c_n, with c_n the 5% null quantile of min-p. MC SE of size = 0.0015; independence would give c_n = 0.0253.

Null calibration (rule 1)

The "weakly correlated" premise of the direction is false at small n: S2 and |b1| have null Pearson correlation 0.76 at n=25, falling to 0.49 at n=100. c_n sits well above the independence value 0.0253 at every n, so the correlation buys back some of the min-p split.

Comparison (rule 2): power at exact 5% size, 20000 draws per cell, all tests scored on the same draws

Gap = combo minus other, paired SE in parentheses. Verdict: beats if gap > 2 SE, loses if gap < -2 SE, else tie.

Answer to the direction (seed 92101)

1. Chi-square(5) vs AD: beats at n=25 (+1.9 points), loses at n=50 (-1.8) and n=100 (-0.6). Vs Shapiro-Wilk: loses at every n (-6.3, -10.7, -1.6 points). The |b1| partner closes roughly half of S2's skew gap (S2 alone is 3 to 16 points behind SW) but does not tie SW anywhere.
2. Heavy tails and contamination vs AD and SW: the combination beats both on all 15 cells, every n.
3. Stays within 2 paired SE of plain S2? Mostly no. On t3, t5, Laplace the combination gives up 0.6 to 4.2 points of S2's power (all beyond 2 SE except t5 at n=25); the cost grows with n (Laplace: -2.1, -4.2, -3.4 points). On contamination it is a tie or a small gain (within +0.8 points) except cont10 at n=100 (-0.5, beyond 2 SE). So the "within two paired SE" condition of the direction fails on 13 of 15 heavy-tail and contamination cells, though the losses are small relative to the gains over AD and SW.

Net: an exact-size min-p combination of S2 with |b1| is a strict improvement over S2 on skew and still dominates AD and SW off-skew, but it is not a free lunch: it costs 1 to 4 points on heavy tails and still loses to Shapiro-Wilk on chi-square(5) at every n. The null correlation between S2 and |b1| (0.49 to 0.76) is far from weak, which is why the gain on skew is only partial.

What remains for stuff 3

Same seed, swap the partner to W (Shapiro-Wilk) as rule 1 allows, so the loss to SW on chi-square(5) can be checked against a partner that is itself skew-optimal. Printed raw output of this run is reproduced from the script exactly; no cell was edited by hand.

Cycle 92, stuff 3 (member-1): S2 + W (Shapiro-Wilk) exact-size calibrated, 20000 sims per n, seed 92101

Runs the stuff 1 script with the partner swapped from |b1| to W, as rule 1 allows, and the same seed 92101 so every draw is identical to stuff 2 (cycle92-s2-b1-calibrated-results.md): one default_rng(92101) stream consumed in order n=25, 50, 100, laws normal, chi2(5), t3, t5, Laplace, cont10, cont5 per n. The S2, A2 and W columns below therefore match stuff 2 to the draw; only the combination changes. Wall clock 35 s. Two-line diff against the stuff 1 script: return s2, a2, b1, w becomes return s2, a2, -w, w (so the "reject large" empirical upper-tail p-value machinery is applied to -W unchanged), and the printed label. No cell edited by hand.

Definitions (rule 1, before any table)

z = standardized sorted sample (n-1 sd); t_i = -((2i-1)/n)[ln Phi(z_(i)) + ln(1-Phi(z_(n+1-i)))]; S2 = sum log(1+2 t_i)/2 (reject large); A2 = sum t_i - n (plain AD, reject large); W = Shapiro-Wilk (reject small); partner W = the same Shapiro-Wilk statistic, used as the second leg of the combination with lower-tail p-value p_W = P(W_null <= W).

Exact-size calibration: p_S2 (upper tail) and p_W (lower tail) are empirical p-values against the 20000 null S2 and W values; the combination rejects when min(p_S2, p_W) < c_n, with c_n the 5% null quantile of min-p. MC SE of size = 0.0015; independence would give c_n = 0.0253.

Null calibration (rule 1)

S2 and W are less correlated under the null than S2 and |b1| were (Pearson 0.64/0.52/0.41 here versus 0.76/0.63/0.49 in stuff 2), but still not "weakly": the direction's premise fails for both partners. c_n sits above the independence value 0.0253 at every n.

Comparison (rule 2): power at exact 5% size, 20000 draws per cell, all tests scored on the same draws

Gap = combo minus other, paired SE in parentheses. Verdict: beats if gap > 2 SE, loses if gap < -2 SE, else tie.

Answer to the direction (seed 92101, partner W)

1. Chi-square(5) vs AD: beats at every n (+4.9, +5.5, +0.9 points), where the |b1| partner lost at n=50 and 100. Vs Shapiro-Wilk: still loses at every n (-3.3, -3.4, -0.1 points); the n=100 loss is one tenth of a point but 5 paired SE. The W partner closes about half the remaining gap the |b1| partner left (stuff 2: -6.3, -10.7, -1.6) but a min-p combination that contains W cannot match W alone on W's home ground, because the c_n cut (0.033 to 0.039) is below the 0.05 W would use on its own.
2. Heavy tails and contamination vs AD and SW: beats both on all 15 cells, every n, as with |b1|.
3. Stays within 2 paired SE of plain S2? No, on 14 of 15 heavy-tail and contamination cells (only cont5 at n=100 ties). The cost is larger than with |b1|: 1.7 to 4.6 points on t3, t5, Laplace and cont10 at every n, and 0.4 to 0.9 on cont5 at n=25 and 50. W is a stronger second leg on skew and a weaker one on tails, so it dilutes S2 more where S2 is best.

Net across stuff 2 and 3: with either permitted partner the exact-size min-p combination keeps S2's dominance over AD and SW off-skew, but the direction's two side conditions both fail. The partner is not weakly correlated with S2 under the null (0.41 to 0.76), and the combination does not stay within two paired SE of S2 on heavy tails (loses 1 to 5 points). On chi-square(5) the W partner ties or beats AD everywhere and the |b1| partner does not, yet neither ties Shapiro-Wilk at any n. The trade is sharp: W as partner buys +3 to +7 points of skew power over the |b1| partner at the price of a further 1 to 2 points on t3, t5 and Laplace.

Raw output

The table above is generated from the script's stdout by a parser; the stdout is:

n=25 pearson(S2,-W)=0.635 spearman=0.441 c_n=0.0385 size=0.0500 crit={'s2': np.float64(13.396133327403478), 'a2': np.float64(0.7365136193325227), 'w': np.float64(0.9192567059400792)}
25 chi2_5 {'combo': np.float64(0.5139), 's2': np.float64(0.47035), 'a2': np.float64(0.46495), 'w': np.float64(0.54675), 'vs_a2': '+0.0490 (0.0021) beats', 'vs_w': '-0.0328 (0.0014) loses', 'vs_s2': '+0.0436 (0.0018) beats'}
25 t3 {'combo': np.float64(0.43635), 's2': np.float64(0.46195), 'a2': np.float64(0.37685), 'w': np.float64(0.4012), 'vs_a2': '+0.0595 (0.0020) beats', 'vs_w': '+0.0352 (0.0016) beats', 'vs_s2': '-0.0256 (0.0013) loses'}
25 t5 {'combo': np.float64(0.2435), 's2': np.float64(0.26335), 'a2': np.float64(0.18975), 'w': np.float64(0.22105), 'vs_a2': '+0.0537 (0.0020) beats', 'vs_w': '+0.0224 (0.0014) beats', 'vs_s2': '-0.0198 (0.0013) loses'}
25 laplace {'combo': np.float64(0.35315), 's2': np.float64(0.38875), 'a2': np.float64(0.3183), 'w': np.float64(0.30835), 'vs_a2': '+0.0348 (0.0021) beats', 'vs_w': '+0.0448 (0.0017) beats', 'vs_s2': '-0.0356 (0.0014) loses'}
25 cont10 {'combo': np.float64(0.3656), 's2': np.float64(0.384), 'a2': np.float64(0.29495), 'w': np.float64(0.34645), 'vs_a2': '+0.0706 (0.0020) beats', 'vs_w': '+0.0192 (0.0013) beats', 'vs_s2': '-0.0184 (0.0012) loses'}
25 cont5 {'combo': np.float64(0.44275), 's2': np.float64(0.45165), 'a2': np.float64(0.3963), 'w': np.float64(0.43455), 'vs_a2': '+0.0464 (0.0018) beats', 'vs_w': '+0.0082 (0.0010) beats', 'vs_s2': '-0.0089 (0.0010) loses'}
n=50 pearson(S2,-W)=0.524 spearman=0.355 c_n=0.0348 size=0.0500 crit={'s2': np.float64(26.50860299377608), 'a2': np.float64(0.7432881122993524), 'w': np.float64(0.9544003312065223)}
50 chi2_5 {'combo': np.float64(0.85855), 's2': np.float64(0.7281), 'a2': np.float64(0.80325), 'w': np.float64(0.89225), 'vs_a2': '+0.0553 (0.0019) beats', 'vs_w': '-0.0337 (0.0013) loses', 'vs_s2': '+0.1305 (0.0026) beats'}
50 t3 {'combo': np.float64(0.69115), 's2': np.float64(0.72485), 'a2': np.float64(0.6032), 'w': np.float64(0.6411), 'vs_a2': '+0.0880 (0.0023) beats', 'vs_w': '+0.0501 (0.0018) beats', 'vs_s2': '-0.0337 (0.0014) loses'}
50 t5 {'combo': np.float64(0.40325), 's2': np.float64(0.4355), 'a2': np.float64(0.30085), 'w': np.float64(0.35785), 'vs_a2': '+0.1024 (0.0025) beats', 'vs_w': '+0.0454 (0.0019) beats', 'vs_s2': '-0.0323 (0.0015) loses'}
50 laplace {'combo': np.float64(0.60895), 's2': np.float64(0.655), 'a2': np.float64(0.54735), 'w': np.float64(0.5328), 'vs_a2': '+0.0616 (0.0024) beats', 'vs_w': '+0.0761 (0.0022) beats', 'vs_s2': '-0.0461 (0.0016) loses'}
50 cont10 {'combo': np.float64(0.59395), 's2': np.float64(0.6113), 'a2': np.float64(0.46405), 'w': np.float64(0.56945), 'vs_a2': '+0.1299 (0.0027) beats', 'vs_w': '+0.0245 (0.0016) beats', 'vs_s2': '-0.0173 (0.0014) loses'}
50 cont5 {'combo': np.float64(0.6637), 's2': np.float64(0.6677), 'a2': np.float64(0.5875), 'w': np.float64(0.65975), 'vs_a2': '+0.0762 (0.0021) beats', 'vs_w': '+0.0040 (0.0011) beats', 'vs_s2': '-0.0040 (0.0011) loses'}
n=100 pearson(S2,-W)=0.409 spearman=0.281 c_n=0.0326 size=0.0500 crit={'s2': np.float64(52.7087370866728), 'a2': np.float64(0.7488546648959901), 'w': np.float64(0.9745549086280316)}
100 chi2_5 {'combo': np.float64(0.99715), 's2': np.float64(0.93215), 'a2': np.float64(0.9877), 'w': np.float64(0.9982), 'vs_a2': '+0.0094 (0.0007) beats', 'vs_w': '-0.0010 (0.0002) loses', 'vs_s2': '+0.0650 (0.0018) beats'}
100 t3 {'combo': np.float64(0.9175), 's2': np.float64(0.93435), 'a2': np.float64(0.84915), 'w': np.float64(0.87555), 'vs_a2': '+0.0684 (0.0019) beats', 'vs_w': '+0.0419 (0.0015) beats', 'vs_s2': '-0.0169 (0.0010) loses'}
100 t5 {'combo': np.float64(0.63945), 's2': np.float64(0.6827), 'a2': np.float64(0.4798), 'w': np.float64(0.56125), 'vs_a2': '+0.1596 (0.0030) beats', 'vs_w': '+0.0782 (0.0022) beats', 'vs_s2': '-0.0433 (0.0017) loses'}
100 laplace {'combo': np.float64(0.88125), 's2': np.float64(0.9122), 'a2': np.float64(0.8285), 'w': np.float64(0.7993), 'vs_a2': '+0.0527 (0.0020) beats', 'vs_w': '+0.0819 (0.0021) beats', 'vs_s2': '-0.0310 (0.0013) loses'}
100 cont10 {'combo': np.float64(0.8258), 's2': np.float64(0.8351), 'a2': np.float64(0.67785), 'w': np.float64(0.7992), 'vs_a2': '+0.1480 (0.0028) beats', 'vs_w': '+0.0266 (0.0015) beats', 'vs_s2': '-0.0093 (0.0013) loses'}
100 cont5 {'combo': np.float64(0.879), 's2': np.float64(0.87755), 'a2': np.float64(0.80135), 'w': np.float64(0.87515), 'vs_a2': '+0.0776 (0.0020) beats', 'vs_w': '+0.0039 (0.0009) beats', 'vs_s2': '+0.0014 (0.0010) tie'}