#!/usr/bin/env python3 """ Cross-Sector Variance Crisis Signal -- Paper 2, qFriction_bound Prediction: Cross-sectional variance of IP subsector log returns serves as a crisis signal. As information friction T approaches T*, cross-sector variance rises, predicting recessions. AUC > 0.65 contemporaneous, AUC > 0.60 for 6-month ahead prediction. Paper 2 (CES Potential), Theorem "q-Variance-Response Identity" (thm:q_fdt): Sigma = T/(2-q) * chi As T -> T*, K_eff -> 0, so lambda_perp = T/((2-q)*K_eff) diverges. Cross-sectional variance of allocations (proxied by IP subsector returns) should spike as the economy approaches the critical friction T*. Lean: CESProofs/Potential/QDynamics.lean qFriction_bound [schematic] -- T >= Var_q(eps) / (2*LSE_q) qFriction_measurement_formula [schematic] Five subtests: 1. Contemporaneous: AUC of cross-sectional IP variance for NBER recession classification (theory predicts AUC > 0.65) 2. 6-month ahead: AUC for predicting recession 6 months later (theory predicts AUC > 0.60) 3. Recession onset: variance is elevated in the 6 months before recession onset (pre-onset z-score > 0.5, one-tailed permutation test) 4. NFCI horse race: logit of USREC on cross-var + NFCI; test whether cross-var adds significant predictive power beyond NFCI 5. K_eff eigenvalue signal: rolling eigenvalue ratio from IP covariance matrix tracks K_eff; declining K_eff predicts recessions (AUC > 0.60) Data: FRED sectoral IP indices (13 NAICS manufacturing subsectors), NBER recession indicator (USREC), Chicago Fed NFCI, 1972-present. Review Notes (2026-03-12) ------------------------- Theory match: CORRECT - Paper 2, q-VRI theorem: as T -> T*, variance diverges. - Test checks whether cross-sectional IP variance discriminates recession from expansion months (AUC). - Lean theorem is schematic (True := trivial). I/O audit: CORRECT - Cross-sectional variance of IP log returns (output-side proxy for allocation fluctuations). Theory predicts variance of input allocations rises, but IP is the standard observable proxy. - USREC (recession indicator) proxies T approaching T*. - K_eff from eigenvalue ratio: theory-specific signal from q-covariance theorem eigenvalue structure. Data-signal alignment: ADEQUATE - FRED IP subsectors: DIRECT measure of sectoral output (N=13). - USREC: DIRECT recession indicator (74 recession months, 575 expansion). - NFCI: DIRECT financial conditions index (Chicago Fed). - Cross-sectional variance: REASONABLE proxy for allocation dispersion. - Eigenvalue ratio: REASONABLE proxy for K_eff (theory-specific). Code quality: CLEAN (improved 2026-03-12) - Circular shift permutation preserves autocorrelation structure. - Added recession onset, NFCI horse race, K_eff eigenvalue subtests. Confidence in verdict: MEDIUM-HIGH - Contemporaneous AUC=0.794 is strong and robust. - 6m-ahead AUC=0.605 is borderline (p=0.108 with circular shift). - Onset, NFCI, and K_eff subtests provide complementary evidence. Date of review: 2026-03-12. """ import warnings warnings.filterwarnings("ignore") import numpy as np import pandas as pd import sys, os sys.path.insert(0, os.path.join(os.path.dirname(__file__), '..', '..')) from common.data import load_fred, load_fred_multi from common.stats import log_result, compute_auc IP_SUBSECTORS = [ 'IPG331S', 'IPG332S', 'IPG333S', 'IPG334S', 'IPG335S', 'IPG322S', 'IPG324S', 'IPG325S', 'IPG326S', 'IPG327S', 'IPG336S', 'IPG337S', 'IPG339S', ] _cache = {} # populated by main(), used by figures() def _circular_shift_auc_pvalue(indicator, recession, horizon_months, n_shifts=500, seed=42): """Circular shift permutation test for AUC significance. Circularly shifts the recession series by random offsets to preserve autocorrelation structure while breaking temporal alignment. """ observed_auc = compute_auc(indicator, recession, horizon_months=horizon_months) if np.isnan(observed_auc): return observed_auc, np.nan, [] rng = np.random.default_rng(seed) common_idx = indicator.dropna().index.intersection(recession.dropna().index) rec_aligned = recession.reindex(common_idx).dropna() T = len(rec_aligned) perm_aucs = [] for _ in range(n_shifts): shift = rng.integers(1, T) shifted_vals = np.roll(rec_aligned.values, shift) shifted_rec = pd.Series(shifted_vals, index=rec_aligned.index) auc_perm = compute_auc(indicator, shifted_rec, horizon_months=horizon_months) if not np.isnan(auc_perm): perm_aucs.append(auc_perm) if perm_aucs: p = np.mean(np.array(perm_aucs) >= observed_auc) else: p = np.nan return observed_auc, p, perm_aucs def _build_core_data(): """Load and prepare core data used by all subtests. Returns None on failure.""" try: usrec = load_fred('USREC') except FileNotFoundError: return None try: ip_df = load_fred_multi(IP_SUBSECTORS) except FileNotFoundError: return None n_sectors = ip_df.shape[1] if n_sectors < 4: return None log_ret = np.log(ip_df / ip_df.shift(1)).dropna() cross_var = log_ret.var(axis=1).dropna() cross_var.name = 'cross_var' cross_var_smooth = cross_var.rolling(3, min_periods=1).mean() cross_var_smooth.name = 'cross_var_smooth' usrec_monthly = usrec.resample('MS').last().dropna() common_idx = cross_var_smooth.dropna().index.intersection(usrec_monthly.dropna().index) n_months = len(common_idx) n_recession = int(usrec_monthly.reindex(common_idx).dropna().sum()) n_expansion = n_months - n_recession return { 'ip_df': ip_df, 'log_ret': log_ret, 'cross_var': cross_var, 'indicator': cross_var_smooth, 'usrec': usrec_monthly, 'n_sectors': n_sectors, 'n_months': n_months, 'n_recession': n_recession, 'n_expansion': n_expansion, } def main(): results = [] data = _build_core_data() if data is None: return [log_result( 'variance_crisis_signal_contemporaneous', stat=np.nan, p=np.nan, verd='SKIPPED', description='Cross-sectional IP variance as recession predictor. ' 'SKIPPED: required FRED data not available.' )] indicator = data['indicator'] usrec_monthly = data['usrec'] log_ret = data['log_ret'] n_sectors = data['n_sectors'] n_months = data['n_months'] n_recession = data['n_recession'] n_expansion = data['n_expansion'] _cache.update(data) # ================================================================ # Subtest 1: Contemporaneous AUC # ================================================================ auc_0, p1, perm_aucs_0 = _circular_shift_auc_pvalue( indicator, usrec_monthly, horizon_months=0, n_shifts=500, seed=42 ) _cache['perm_aucs_0'] = perm_aucs_0 _cache['auc_0'] = auc_0 if np.isnan(auc_0): verd1 = 'SKIPPED' else: if auc_0 >= 0.65: verd1 = 'CONSISTENT' if p1 < 0.05 else 'AMBIGUOUS' elif auc_0 >= 0.55: verd1 = 'AMBIGUOUS' else: verd1 = 'INCONSISTENT' results.append(log_result( 'variance_crisis_signal_contemporaneous', stat=float(auc_0), p=float(p1), verd=verd1, details={ 'auc': round(float(auc_0), 4), 'threshold': 0.65, 'n_sectors': n_sectors, 'n_months': n_months, 'n_recession_months': n_recession, 'n_expansion_months': n_expansion, 'p_method': 'circular shift permutation (500 shifts)', }, description=( "Tests the q-VRI prediction (Paper 2, Thm q-FDT): as information " "friction T approaches critical T*, cross-sectional variance of " "allocations diverges, signaling crisis. AUC of 3-month smoothed " "cross-sectional variance of %d FRED IP subsector log returns for " "classifying NBER recession months (%d recession, %d expansion, " "1972-present). Circular shift permutation p-value (preserves " "autocorrelation). Prediction: AUC > 0.65. " "AUC=%.3f, p=%.4f." % (n_sectors, n_recession, n_expansion, auc_0, p1) ), )) # ================================================================ # Subtest 2: 6-month ahead AUC # ================================================================ auc_6, p2, perm_aucs_6 = _circular_shift_auc_pvalue( indicator, usrec_monthly, horizon_months=6, n_shifts=500, seed=123 ) _cache['perm_aucs_6'] = perm_aucs_6 _cache['auc_6'] = auc_6 if np.isnan(auc_6): verd2 = 'SKIPPED' else: if auc_6 >= 0.60: verd2 = 'CONSISTENT' if p2 < 0.10 else 'AMBIGUOUS' elif auc_6 >= 0.52: verd2 = 'AMBIGUOUS' else: verd2 = 'INCONSISTENT' results.append(log_result( 'variance_crisis_signal_6m_ahead', stat=float(auc_6), p=float(p2), verd=verd2, details={ 'auc': round(float(auc_6), 4), 'threshold': 0.60, 'n_sectors': n_sectors, 'n_months': n_months, 'n_recession_months': n_recession, 'p_method': 'circular shift permutation (500 shifts)', }, description=( "Tests 6-month-ahead crisis prediction from cross-sectional IP " "variance (Paper 2, q-VRI). Pairs variance at month t with " "recession indicator at month t+6. AUC of 3-month smoothed " "cross-sectional variance of %d FRED IP subsector log returns " "(%d months, 1972-present). Circular shift permutation p-value. " "Prediction: AUC > 0.60. " "AUC=%.3f, p=%.4f." % (n_sectors, n_months, auc_6, p2) ), )) # ================================================================ # Subtest 3: Recession onset — is variance elevated pre-onset? # ================================================================ # Construct onset indicator (1 only in first month of each recession) onset = ((usrec_monthly > 0.5) & (usrec_monthly.shift(1) <= 0.5)).astype(int) onset_dates = onset[onset == 1].index n_onsets = len(onset_dates) common = indicator.dropna().index.intersection(usrec_monthly.dropna().index) ind_aligned = indicator.reindex(common).dropna() if n_onsets >= 4: # For each onset, compute mean variance in 6 months BEFORE onset pre_onset_means = [] for dt in onset_dates: # 6 months before onset mask = (ind_aligned.index >= dt - pd.DateOffset(months=6)) & \ (ind_aligned.index < dt) pre_vals = ind_aligned[mask] if len(pre_vals) >= 3: pre_onset_means.append(pre_vals.mean()) if pre_onset_means: pre_onset_avg = np.mean(pre_onset_means) overall_mean = ind_aligned.mean() overall_std = ind_aligned.std() pre_onset_z = (pre_onset_avg - overall_mean) / overall_std # Permutation test: circular shift onset dates, recompute pre-onset mean rng = np.random.default_rng(42) T = len(onset) perm_zs = [] for _ in range(500): shift = rng.integers(1, T) shifted_onset = pd.Series( np.roll(onset.values, shift), index=onset.index ) shifted_dates = shifted_onset[shifted_onset == 1].index perm_means = [] for dt in shifted_dates: mask = (ind_aligned.index >= dt - pd.DateOffset(months=6)) & \ (ind_aligned.index < dt) pre_vals = ind_aligned[mask] if len(pre_vals) >= 3: perm_means.append(pre_vals.mean()) if perm_means: perm_avg = np.mean(perm_means) perm_z = (perm_avg - overall_mean) / overall_std perm_zs.append(perm_z) if perm_zs: p3 = np.mean(np.array(perm_zs) >= pre_onset_z) else: p3 = np.nan sign_ok = pre_onset_z > 0 # variance should be ABOVE average pre-onset if sign_ok and p3 < 0.05: verd3 = 'CONSISTENT' elif sign_ok and p3 < 0.15: verd3 = 'AMBIGUOUS' elif not sign_ok: verd3 = 'INCONSISTENT' else: verd3 = 'AMBIGUOUS' _cache['onset'] = { 'onset_dates': [d.strftime('%Y-%m') for d in onset_dates], 'pre_onset_means': pre_onset_means, 'pre_onset_z': pre_onset_z, 'overall_mean': overall_mean, 'overall_std': overall_std, 'perm_zs': perm_zs, } results.append(log_result( 'variance_crisis_signal_onset', stat=float(pre_onset_z), p=float(p3), verd=verd3, details={ 'n_onsets': n_onsets, 'n_onsets_with_data': len(pre_onset_means), 'pre_onset_z': round(pre_onset_z, 3), 'pre_onset_mean': round(pre_onset_avg, 6), 'overall_mean': round(overall_mean, 6), 'onset_dates': [d.strftime('%Y-%m') for d in onset_dates], 'p_method': 'circular shift permutation (500 shifts)', }, description=( "Tests whether cross-sectional IP variance is elevated in the " "6 months before recession ONSET (first month of each NBER " "recession). This distinguishes leading from coincident signals. " "Computes mean variance in 6-month pre-onset window for %d " "recession episodes, standardizes against full-sample mean/std " "(z-score). Circular shift permutation of onset dates. " "Prediction: pre-onset z > 0.5 (variance elevated before onset). " "z=%.3f, p=%.4f." % (n_onsets, pre_onset_z, p3) ), )) else: results.append(log_result( 'variance_crisis_signal_onset', stat=np.nan, p=np.nan, verd='SKIPPED', description="Recession onset leading indicator test. SKIPPED: " "insufficient pre-onset data windows.", )) else: results.append(log_result( 'variance_crisis_signal_onset', stat=np.nan, p=np.nan, verd='SKIPPED', description="Recession onset leading indicator test. SKIPPED: " "fewer than 4 recession onsets in sample.", )) # ================================================================ # Subtest 4: NFCI horse race — does variance add info beyond NFCI? # ================================================================ try: nfci_raw = load_fred('NFCI') # NFCI is weekly — resample to monthly nfci_monthly = nfci_raw.resample('MS').mean().dropna() have_nfci = True except FileNotFoundError: have_nfci = False if have_nfci: import statsmodels.api as sm # Align all three series common3 = indicator.dropna().index \ .intersection(usrec_monthly.dropna().index) \ .intersection(nfci_monthly.dropna().index) if len(common3) >= 100: y = usrec_monthly.reindex(common3).values.astype(float) x_var = indicator.reindex(common3).values x_nfci = nfci_monthly.reindex(common3).values # Standardize predictors for comparable coefficients x_var_z = (x_var - np.mean(x_var)) / np.std(x_var) x_nfci_z = (x_nfci - np.mean(x_nfci)) / np.std(x_nfci) # Model 1: NFCI only X1 = sm.add_constant(x_nfci_z) try: logit1 = sm.Logit(y, X1).fit(disp=0) aic1 = logit1.aic llf1 = logit1.llf except Exception: logit1 = None aic1 = np.nan llf1 = np.nan # Model 2: Variance + NFCI X2 = sm.add_constant(np.column_stack([x_var_z, x_nfci_z])) try: logit2 = sm.Logit(y, X2).fit(disp=0) aic2 = logit2.aic llf2 = logit2.llf var_coeff = logit2.params[1] var_p = logit2.pvalues[1] nfci_coeff = logit2.params[2] nfci_p = logit2.pvalues[2] except Exception: logit2 = None aic2 = np.nan llf2 = np.nan var_coeff = np.nan var_p = np.nan nfci_coeff = np.nan nfci_p = np.nan # LR test: does adding variance improve over NFCI-only? if logit1 is not None and logit2 is not None: from scipy.stats import chi2 lr_stat = 2 * (llf2 - llf1) lr_p = chi2.sf(lr_stat, df=1) # 1 added parameter else: lr_stat = np.nan lr_p = np.nan delta_aic = aic2 - aic1 # negative = joint model better # Verdict: variance adds info if LR test significant sign_ok_4 = var_coeff > 0 # higher variance -> higher recession prob if sign_ok_4 and lr_p < 0.05: verd4 = 'CONSISTENT' elif sign_ok_4 and lr_p < 0.15: verd4 = 'AMBIGUOUS' elif not sign_ok_4: verd4 = 'INCONSISTENT' else: verd4 = 'AMBIGUOUS' _cache['nfci_race'] = { 'n_months': len(common3), 'var_coeff': var_coeff, 'var_p': var_p, 'nfci_coeff': nfci_coeff, 'nfci_p': nfci_p, 'aic_nfci_only': aic1, 'aic_joint': aic2, 'delta_aic': delta_aic, 'lr_stat': lr_stat, 'lr_p': lr_p, } results.append(log_result( 'variance_crisis_signal_nfci_race', stat=float(var_coeff), p=float(lr_p), verd=verd4, details={ 'n_months': len(common3), 'var_coeff_std': round(var_coeff, 3), 'var_p_joint': round(var_p, 4), 'nfci_coeff_std': round(nfci_coeff, 3), 'nfci_p_joint': round(nfci_p, 4), 'aic_nfci_only': round(aic1, 1), 'aic_joint': round(aic2, 1), 'delta_aic': round(delta_aic, 1), 'lr_stat': round(lr_stat, 2), 'lr_p': round(lr_p, 4), }, description=( "Tests whether cross-sectional IP variance adds predictive " "power for recessions beyond the Chicago Fed NFCI (standard " "financial conditions index). Logit of USREC on standardized " "cross-var + standardized NFCI (%d months). LR test of joint " "model vs NFCI-only. Prediction: variance coefficient positive " "and LR test significant (p < 0.05). " "var_coeff=%.3f (p=%.4f), LR chi2=%.2f (p=%.4f), " "DAIC=%.1f." % (len(common3), var_coeff, var_p, lr_stat, lr_p, delta_aic) ), )) else: results.append(log_result( 'variance_crisis_signal_nfci_race', stat=np.nan, p=np.nan, verd='SKIPPED', description="NFCI horse race. SKIPPED: insufficient overlapping data.", )) else: results.append(log_result( 'variance_crisis_signal_nfci_race', stat=np.nan, p=np.nan, verd='SKIPPED', description="NFCI horse race. SKIPPED: NFCI data not available.", )) # ================================================================ # Subtest 5: K_eff eigenvalue signal # ================================================================ # Theory (q-covariance theorem): eigenvalue ratio gives K_eff. # lambda_perp / lambda_1 = 1 + 1/(J*K_eff) # As K_eff -> 0, perpendicular eigenvalues dominate -> dispersion ratio rises. # Compute rolling dispersion ratio from IP covariance eigenvalues. J = n_sectors window = 24 # months if len(log_ret) >= window + 12: disp_ratio_series = {} keff_series = {} dates = log_ret.index[window - 1:] for i in range(window - 1, len(log_ret)): dt = log_ret.index[i] block = log_ret.iloc[i - window + 1:i + 1] if len(block) < window: continue cov_mat = block.cov().values eigvals = np.linalg.eigvalsh(cov_mat) eigvals = np.sort(eigvals)[::-1] # descending # lambda_1 = largest eigenvalue (market factor) # lambda_perp = mean of remaining J-1 eigenvalues lambda_1 = eigvals[0] lambda_perp_mean = np.mean(eigvals[1:]) if lambda_1 > 0: disp_ratio = lambda_perp_mean / lambda_1 disp_ratio_series[dt] = disp_ratio # K_eff = 1 / (J * (lambda_perp/lambda_1 - 1)) # But lambda_perp < lambda_1 empirically (market factor dominates) # so ratio < 1 and K_eff would be negative — use inverse instead # PC1 share = lambda_1 / sum(eigvals) — theory-agnostic crisis measure pc1_share = lambda_1 / np.sum(eigvals) keff_series[dt] = pc1_share if disp_ratio_series: disp_ratio_ts = pd.Series(disp_ratio_series) disp_ratio_ts = disp_ratio_ts.rolling(3, min_periods=1).mean() pc1_share_ts = pd.Series(keff_series) pc1_share_ts = pc1_share_ts.rolling(3, min_periods=1).mean() # Theory: when K_eff -> 0, idiosyncratic dispersion rises. # In output data, recessions show HIGHER dispersion ratio # (sectors decouple as complementarity breaks down). # Test: AUC of dispersion ratio for recession prediction. auc_dr, p_dr, perm_dr = _circular_shift_auc_pvalue( disp_ratio_ts, usrec_monthly, horizon_months=0, n_shifts=500, seed=456 ) # Also test PC1 share (inverse: LOWER PC1 share during crises # would mean idiosyncratic variance dominates, consistent with K_eff -> 0). # But empirically, PC1 share typically RISES during crises (correlations up). # Test both directions. auc_pc1, p_pc1, _ = _circular_shift_auc_pvalue( pc1_share_ts, usrec_monthly, horizon_months=0, n_shifts=500, seed=789 ) # Use dispersion ratio as primary (more theory-aligned) sign_ok_5 = auc_dr > 0.5 # higher dispersion ratio -> more recession if auc_dr >= 0.60 and p_dr < 0.10 and sign_ok_5: verd5 = 'CONSISTENT' elif sign_ok_5 and auc_dr >= 0.52: verd5 = 'AMBIGUOUS' else: verd5 = 'INCONSISTENT' _cache['keff'] = { 'disp_ratio_ts': disp_ratio_ts, 'pc1_share_ts': pc1_share_ts, 'auc_dr': auc_dr, 'p_dr': p_dr, 'auc_pc1': auc_pc1, 'p_pc1': p_pc1, 'perm_dr': perm_dr, } results.append(log_result( 'variance_crisis_signal_keff', stat=float(auc_dr), p=float(p_dr), verd=verd5, details={ 'auc_dispersion_ratio': round(float(auc_dr), 4), 'p_dispersion_ratio': round(float(p_dr), 4), 'auc_pc1_share': round(float(auc_pc1), 4), 'p_pc1_share': round(float(p_pc1), 4), 'window_months': window, 'n_sectors': J, 'p_method': 'circular shift permutation (500 shifts)', }, description=( "Tests the K_eff eigenvalue signal from the q-covariance " "theorem (Paper 2): as K_eff -> 0, perpendicular eigenvalues " "dominate the IP covariance matrix, increasing the dispersion " "ratio (mean non-PC1 eigenvalue / PC1 eigenvalue). Computed " "from rolling %d-month covariance of %d IP subsectors. " "AUC for contemporaneous recession classification. Circular " "shift permutation. Prediction: AUC > 0.60. " "AUC_disp=%.3f (p=%.4f), AUC_PC1=%.3f (p=%.4f)." % (window, J, auc_dr, p_dr, auc_pc1, p_pc1) ), )) else: results.append(log_result( 'variance_crisis_signal_keff', stat=np.nan, p=np.nan, verd='SKIPPED', description="K_eff eigenvalue signal. SKIPPED: insufficient data " "for rolling covariance.", )) else: results.append(log_result( 'variance_crisis_signal_keff', stat=np.nan, p=np.nan, verd='SKIPPED', description="K_eff eigenvalue signal. SKIPPED: time series too short.", )) return results def figures(): """Generate figures for each subtest.""" from common.plotting import setup_figure, setup_subplots, save_figure from matplotlib.patches import Patch if not _cache: main() paths = [] indicator = _cache.get('indicator') usrec = _cache.get('usrec') def _shade_recessions(ax, rec): in_recession = False start = None for dt, val in rec.items(): if val > 0.5 and not in_recession: start = dt in_recession = True elif val <= 0.5 and in_recession: ax.axvspan(start, dt, alpha=0.15, color='red') in_recession = False if in_recession: ax.axvspan(start, rec.index[-1], alpha=0.15, color='red') # --- Figure 1: Time series with recession shading --- if indicator is not None and usrec is not None: fig, ax = setup_figure(width=12, height=5) common = indicator.dropna().index.intersection(usrec.dropna().index) ind = indicator.reindex(common).dropna() rec = usrec.reindex(common).dropna() ax.plot(ind.index, ind.values, color='#2196F3', lw=1, alpha=0.8) _shade_recessions(ax, rec) auc_0 = _cache.get('auc_0', 0) ax.set_xlabel('Date') ax.set_ylabel('Cross-Sectional Variance (log returns)') ax.set_title(f'IP Subsector Cross-Sectional Variance as Crisis Signal\n' f'Contemporaneous AUC={auc_0:.3f}, N={_cache.get("n_sectors", "?")} sectors') ax.legend(handles=[ ax.get_lines()[0], Patch(facecolor='red', alpha=0.15, label='NBER Recession'), ], fontsize=9) paths.extend(save_figure(fig, 'variance_crisis_signal_contemporaneous', subfolder='potential')) # --- Figure 2: 6m-ahead null distribution --- perm_aucs_6 = _cache.get('perm_aucs_6', []) auc_6 = _cache.get('auc_6') if perm_aucs_6 and auc_6 is not None: fig, ax = setup_figure(width=8, height=5) ax.hist(perm_aucs_6, bins=30, color='#90CAF9', edgecolor='black', alpha=0.7, label='Null (circular shifts)') ax.axvline(auc_6, color='red', lw=2, linestyle='--', label=f'Observed AUC={auc_6:.3f}') ax.axvline(0.60, color='orange', lw=1.5, linestyle=':', label='Threshold (0.60)') n_exceed = sum(1 for a in perm_aucs_6 if a >= auc_6) p_val = n_exceed / len(perm_aucs_6) if perm_aucs_6 else np.nan ax.set_xlabel('AUC') ax.set_ylabel('Count') ax.set_title(f'6-Month Ahead: Null Distribution\n' f'p={p_val:.3f} ({n_exceed}/{len(perm_aucs_6)} exceed observed)') ax.legend(fontsize=9) paths.extend(save_figure(fig, 'variance_crisis_signal_6m_ahead', subfolder='potential')) # --- Figure 3: Recession onset pre-onset variance --- onset_data = _cache.get('onset') if onset_data: fig, (ax1, ax2) = setup_subplots(nrows=1, ncols=2, width=14, height=5) # Left: bar chart of pre-onset means vs overall mean onset_dates = onset_data['onset_dates'] pre_means = onset_data['pre_onset_means'] overall = onset_data['overall_mean'] ax1.bar(range(len(pre_means)), pre_means, color='#e74c3c', alpha=0.7, edgecolor='black') ax1.axhline(overall, color='blue', linestyle='--', lw=2, label=f'Overall mean={overall:.5f}') ax1.set_xticks(range(len(onset_dates))) ax1.set_xticklabels(onset_dates, rotation=45, fontsize=8) ax1.set_ylabel('Mean Variance (6mo pre-onset)') ax1.set_title(f'Pre-Onset Variance (z={onset_data["pre_onset_z"]:.2f})') ax1.legend(fontsize=8) # Right: permutation null distribution perm_zs = onset_data['perm_zs'] if perm_zs: ax2.hist(perm_zs, bins=25, color='#90CAF9', edgecolor='black', alpha=0.7) ax2.axvline(onset_data['pre_onset_z'], color='red', lw=2, linestyle='--', label=f'Observed z={onset_data["pre_onset_z"]:.2f}') n_exc = sum(1 for z in perm_zs if z >= onset_data['pre_onset_z']) ax2.set_title(f'Null Distribution (p={n_exc/len(perm_zs):.3f})') ax2.set_xlabel('Pre-onset z-score') ax2.set_ylabel('Count') ax2.legend(fontsize=9) fig.suptitle('Recession Onset Leading Indicator', fontsize=12, y=1.02) paths.extend(save_figure(fig, 'variance_crisis_signal_onset', subfolder='potential')) # --- Figure 4: NFCI horse race coefficients --- nfci_data = _cache.get('nfci_race') if nfci_data: fig, ax = setup_figure(width=7, height=4) labels = ['Cross-Var', 'NFCI'] coeffs = [nfci_data['var_coeff'], nfci_data['nfci_coeff']] ps = [nfci_data['var_p'], nfci_data['nfci_p']] colors = ['#e74c3c' if p < 0.05 else '#90CAF9' for p in ps] bars = ax.bar(labels, coeffs, color=colors, edgecolor='black', alpha=0.8) for bar, coeff, p in zip(bars, coeffs, ps): sig = '***' if p < 0.001 else '**' if p < 0.01 else '*' if p < 0.05 else '' ax.text(bar.get_x() + bar.get_width()/2, bar.get_height() + 0.02, f'{coeff:.2f}{sig}\n(p={p:.3f})', ha='center', fontsize=9) ax.axhline(0, color='black', lw=0.5) ax.set_ylabel('Standardized Logit Coefficient') daic = nfci_data['delta_aic'] lr_p = nfci_data['lr_p'] ax.set_title(f'NFCI Horse Race: Joint Logit Coefficients\n' f'DAIC={daic:.1f}, LR p={lr_p:.4f}') paths.extend(save_figure(fig, 'variance_crisis_signal_nfci_race', subfolder='potential')) # --- Figure 5: K_eff dispersion ratio time series --- keff_data = _cache.get('keff') if keff_data and usrec is not None: fig, (ax1, ax2) = setup_subplots(nrows=2, ncols=1, width=12, height=8) dr = keff_data['disp_ratio_ts'] pc1 = keff_data['pc1_share_ts'] common_dr = dr.index.intersection(usrec.dropna().index) rec_dr = usrec.reindex(common_dr).dropna() ax1.plot(dr.index, dr.values, color='#e74c3c', lw=1, alpha=0.8) _shade_recessions(ax1, rec_dr) ax1.set_ylabel('Dispersion Ratio (perp/PC1)') ax1.set_title(f'Eigenvalue Dispersion Ratio (AUC={keff_data["auc_dr"]:.3f}, ' f'p={keff_data["p_dr"]:.3f})') ax1.legend(handles=[ ax1.get_lines()[0], Patch(facecolor='red', alpha=0.15, label='NBER Recession'), ], fontsize=9) ax2.plot(pc1.index, pc1.values, color='#2196F3', lw=1, alpha=0.8) _shade_recessions(ax2, rec_dr) ax2.set_ylabel('PC1 Share') ax2.set_xlabel('Date') ax2.set_title(f'PC1 Variance Share (AUC={keff_data["auc_pc1"]:.3f}, ' f'p={keff_data["p_pc1"]:.3f})') ax2.legend(handles=[ ax2.get_lines()[0], Patch(facecolor='red', alpha=0.15, label='NBER Recession'), ], fontsize=9) fig.suptitle('K_eff Eigenvalue Signal', fontsize=12, y=1.02) paths.extend(save_figure(fig, 'variance_crisis_signal_keff', subfolder='potential')) return paths if __name__ == '__main__': res = main() for r in res: print(f" {r['name']}: {r['verdict']} (stat={r['statistic']:.4f}, p={r['p_value']})")

auc	0.7940
p value	0
statistic	0.7940
threshold	0.6500

Variance Crisis Signal Contemporaneous

Description

Figures

Statistics

Script

▸ Source Code

Reproduction Steps