#!/usr/bin/env python3 """ Business Cycle Predictions #9-18 ================================= Tests ten predictions from the Business Cycles / Landscape Dynamics papers: #9-11 rho-ordering of recession entry/recovery (Kendall tau pooled) #12 Expansion/contraction duration asymmetry ratio ~ 1/eps ~ 10 #13 Critical slowing down (rolling AR1 rise before recessions) #14 Information temperature T = sigma^2/chi as leading indicator #15 Regulation reduces cycle amplitude (Great Moderation) #16 Phillips curve slope proportional to K_bar (flattening) #17 rho-diversity predicts recession resilience #18 Power-law recession depths Data: 7 FRED manufacturing IP subsectors (monthly, 1972+), INDPRO, USREC, GDPC1, CPIAUCSL, DGS2, DGS10, VIXCLS Migrated from scripts/test_business_cycle_predictions.py Review Notes (2026-03-11) ------------------------- Theory match: CORRECT for #9-11 (rho_ordering), #14 (info_temperature). PARTIAL for #12 (asymmetry — not uniquely CES), #15 (regulation — many explanations). CORRECT for #16 (phillips_flattening — proved in Lean). N/A for #17 (n=2), #18 (schematic). I/O audit: CORRECT. All variables correctly oriented: - rho is input parameter (calibrated), lead_months is output (measured) - T = sigma^2/chi correctly computes fluctuation-response ratio - Phillips slope uses inflation ~ output_gap (correct direction) Data-signal alignment: ADEQUATE overall. - #9-11: 7 FRED manufacturing IP sectors, N=49 pooled obs. ADEQUATE. - #12: 11 NBER cycles. DIRECT (NBER dates are definitive). - #13: WEAK. 4 recessions × 7 sectors = 21 pairs. Great Moderation confound. - #14: ADEQUATE. Probit on ~400 monthly obs, pseudo_R²=0.035 (low). - #15: STRONG. Pre/post 1984 with >60 months each side. - #16: WEAK. Paper acknowledges K_bar declined only ~1pp over available data. - #17/#18: INSUFFICIENT (n=2 recessions, n=7 depths). Code quality: BUGS FOUND (v1 fixes applied 2026-03-11): - Bug 1: Test #16 had SIGN ERROR — checked tau > 0 but theory predicts flattening (tau < 0). Fixed: now checks tau < 0. - Bug 2: Test #12 could never return CONSISTENT — verdict logic only had AMBIGUOUS/INCONSISTENT paths even when ratio was in predicted range. Fixed: ratio in [3,15] → CONSISTENT. - Bug 3: Test #16 had no INCONSISTENT path (falsifiability failure). Fixed: tau > 0 and p < 0.10 → INCONSISTENT. Confidence: MEDIUM. - 4 CONSISTENT: rho_ordering (HIGH), expansion_asymmetry (MEDIUM — not uniquely CES), info_temperature (MEDIUM — low pseudo_R²), regulation_amplitude (MEDIUM — many explanations). - 1 INCONSISTENT: phillips_slope (tau=+0.084, wrong sign, p=0.033). Root cause: CONFOUNDED. Post-COVID supply-side inflation re-steepened the Phillips curve. K_bar channel (~1pp decline) is too weak to dominate. - 3 AMBIGUOUS: critical_slowing (Great Moderation confound), rho_diversity (n=2), power_law_depths (n=7). Test #13 (Critical Slowing Down): - Only 4 target recessions × 7 sectors = 21 pairs. - 120-month rolling window; results sensitive to window choice. - Great Moderation confound reduces autocorrelation mechanically. - Recommendation: test on FINANCIAL data (return correlations). - Current verdict AMBIGUOUS is appropriate. """ import warnings warnings.filterwarnings("ignore") import numpy as np import pandas as pd from scipy import stats import sys, os sys.path.insert(0, os.path.join(os.path.dirname(__file__), '..', '..')) from common.data import load_fred from common.stats import log_result from common.plotting import setup_figure, setup_subplots, save_figure # Manufacturing subsectors with calibrated rho (ordinal ranking). # Sources: Primary Metals (σ≈1.2, ρ≈-0.2) from Oberfield-Raval (2021) Table 2, # highly complementary raw inputs. Food Manufacturing (σ≈3.3, ρ≈0.7) from # Broda-Weinstein (2006) consumer goods median. Intermediate values interpolated # by input complementarity: heavy industry (low ρ) → electronics/food (high ρ). # Exact values matter less than ordinal ranking for Kendall tau test. MFG_SECTORS = { "IPG331S": ("Primary Metals", -0.2), # Strong complements (ore + energy + labor) "IPG333S": ("Machinery", 0.2), # Moderate complements "IPG336S": ("Transport Equip", 0.3), # Complex assembly, many inputs "IPG325S": ("Chemicals", 0.4), # Process industry, some substitution "IPG335S": ("Electrical Equip", 0.5), # Mixed complement/substitute "IPG334S": ("Computer/Electronics", 0.6), # More substitutable components "IPG311A2S": ("Food Manufacturing", 0.7), # Highly substitutable consumer goods } # NBER recession dates (post-war) RECESSIONS = [ ("1948-11-01", "1949-10-01"), ("1953-07-01", "1954-05-01"), ("1957-08-01", "1958-04-01"), ("1960-04-01", "1961-02-01"), ("1969-12-01", "1970-11-01"), ("1973-11-01", "1975-03-01"), ("1980-01-01", "1980-07-01"), ("1981-07-01", "1982-11-01"), ("1990-07-01", "1991-03-01"), ("2001-03-01", "2001-11-01"), ("2007-12-01", "2009-06-01"), ("2020-02-01", "2020-04-01"), ] def _load_mfg_data(): """Load manufacturing IP series and compute log-growth.""" raw, growth = {}, {} for sid in list(MFG_SECTORS.keys()) + ["INDPRO"]: try: s = load_fred(sid) if len(s) >= 60: raw[sid] = s growth[sid] = np.log(s).diff().dropna() except FileNotFoundError: pass return raw, growth def _find_peak_timing(growth, recession_peak, window_before=18, window_after=6): """Find cumulative growth peak for each sector relative to NBER peak.""" peak_dt = pd.Timestamp(recession_peak) results = {} for sid, (name, rho) in MFG_SECTORS.items(): if sid not in growth: continue g = growth[sid] start = peak_dt - pd.DateOffset(months=window_before) end = peak_dt + pd.DateOffset(months=window_after) window = g.loc[start:end] if len(window) < 6: continue cum = window.cumsum() peak_idx = cum.idxmax() lead_months = (peak_dt - peak_idx).days / 30.44 results[sid] = {"name": name, "rho": rho, "lead_months": lead_months} return results # ── Test #9-11: rho-ordering of recession entry ── def test_rho_ordering(growth): """Pooled Kendall tau of (rho, lead_months) across mfg recessions.""" mfg_recessions = [r for r in RECESSIONS if r[0] >= "1973-01-01"] all_leads = [] for peak, trough in mfg_recessions: peaks = _find_peak_timing(growth, peak) if len(peaks) >= 4: for sid, info in peaks.items(): all_leads.append({ 'rho': info['rho'], 'lead_months': info['lead_months'], 'recession': peak[:4], }) if len(all_leads) < 10: return log_result('business_cycle_rho_ordering', np.nan, np.nan, 'AMBIGUOUS', description="Tests whether low-rho manufacturing sectors lead high-rho sectors into recessions, using pooled OLS and Kendall tau on (rho, lead_months) across NBER recession peaks since 1973.") df = pd.DataFrame(all_leads) # Pooled OLS: lead_months = b0 + b1*rho from numpy.polynomial.polynomial import polyfit as _pf rho_arr = df['rho'].values lead_arr = df['lead_months'].values # OLS via numpy X = np.column_stack([np.ones(len(rho_arr)), rho_arr]) beta = np.linalg.lstsq(X, lead_arr, rcond=None)[0] beta1 = beta[1] resid = lead_arr - X @ beta n = len(rho_arr) se_beta1 = np.sqrt(np.sum(resid**2) / (n - 2) / np.sum((rho_arr - rho_arr.mean())**2)) t_stat = beta1 / se_beta1 from scipy.stats import t as t_dist p_ols = 2 * (1 - t_dist.cdf(abs(t_stat), n - 2)) tau, p_tau = stats.kendalltau(rho_arr, lead_arr) # Verdict: beta1 < 0 (low-rho leads) and p < 0.10 if beta1 < 0 and p_ols < 0.10: verd = 'CONSISTENT' elif beta1 < 0: verd = 'AMBIGUOUS' else: verd = 'INCONSISTENT' return log_result('business_cycle_rho_ordering', tau, p_tau, verd, details={'pooled_beta1': beta1, 'pooled_p': p_ols, 'kendall_tau': tau, 'n': n}, description="Tests whether low-rho manufacturing sectors lead high-rho sectors into recessions. Uses pooled OLS (lead_months = b0 + b1*rho) and Kendall tau across 7 FRED manufacturing IP subsectors over post-1973 NBER recessions. Prediction: beta1 < 0.") # ── Test #12: Expansion/contraction asymmetry ── def test_asymmetry(): """Duration ratio of expansion/contraction from NBER dates.""" ratios = [] for i in range(1, len(RECESSIONS)): peak, trough = RECESSIONS[i] peak_dt = pd.Timestamp(peak) trough_dt = pd.Timestamp(trough) contraction = (trough_dt - peak_dt).days / 30.44 prev_trough = pd.Timestamp(RECESSIONS[i - 1][1]) expansion = (peak_dt - prev_trough).days / 30.44 if contraction > 0: ratios.append(expansion / contraction) if not ratios: return log_result('business_cycle_expansion_asymmetry', np.nan, np.nan, 'AMBIGUOUS', description="Tests whether expansion/contraction duration ratio is approximately 1/epsilon ~ 10, as predicted by the slow-manifold asymmetry of the CES landscape dynamics.") mean_ratio = np.mean(ratios) median_ratio = np.median(ratios) # Prediction: ratio ~ 1/eps ~ 10; within order of magnitude [3, 15] if 3 <= mean_ratio <= 15: verd = 'CONSISTENT' elif 1.5 <= mean_ratio <= 30: verd = 'AMBIGUOUS' else: verd = 'INCONSISTENT' return log_result('business_cycle_expansion_asymmetry', mean_ratio, np.nan, verd, details={'mean_ratio': mean_ratio, 'median_ratio': median_ratio, 'n_cycles': len(ratios)}, description="Tests expansion/contraction duration asymmetry from NBER recession dates. Computes ratio of preceding expansion to contraction duration for each post-war cycle. Prediction: mean ratio in [3, 15], consistent with 1/epsilon ~ 10.") # ── Test #13: Critical slowing down ── def _rolling_ar1(series, window=120): """Compute rolling AR(1) via OLS linregress.""" n = len(series) ar1 = np.full(n, np.nan) vals = series.values for t in range(window, n): y = vals[t - window:t] if np.std(y) < 1e-10: continue slope, _, _, _, _ = stats.linregress(y[:-1], y[1:]) ar1[t] = slope return pd.Series(ar1, index=series.index) def test_critical_slowing(growth): """Rolling AR(1) rises in 36-month pre-crisis window. Fraction-based verdict.""" target_recessions = [ ("1973-11-01", "1975-03-01"), ("1990-07-01", "1991-03-01"), ("2001-03-01", "2001-11-01"), ("2007-12-01", "2009-06-01"), ] # Compute rolling AR(1) for each sector sector_ar1 = {} for sid in MFG_SECTORS: if sid not in growth: continue g = growth[sid] if len(g) < 180: continue sector_ar1[sid] = _rolling_ar1(g, window=120) if not sector_ar1: return log_result('business_cycle_critical_slowing', np.nan, np.nan, 'AMBIGUOUS', description="Tests whether rolling AR(1) coefficients rise in 36-month pre-crisis windows, indicating critical slowing down before regime shifts.") csd_results = [] for peak, trough in target_recessions: peak_dt = pd.Timestamp(peak) pre_start = peak_dt - pd.DateOffset(months=36) for sid in sector_ar1: ar1 = sector_ar1[sid] pre_window = ar1.loc[pre_start:peak_dt].dropna() if len(pre_window) < 12: continue # Kendall tau of AR(1) against time in pre-crisis window t_vals = np.arange(len(pre_window)) tau, p = stats.kendalltau(t_vals, pre_window.values) csd_results.append({'tau': tau, 'recession': peak[:4]}) if not csd_results: return log_result('business_cycle_critical_slowing', np.nan, np.nan, 'AMBIGUOUS', description="Tests whether rolling AR(1) coefficients rise in 36-month pre-crisis windows across manufacturing IP sectors.") n_rising = sum(1 for r in csd_results if r['tau'] > 0) n_total = len(csd_results) frac_rising = n_rising / n_total mean_tau = np.mean([r['tau'] for r in csd_results]) # Sign test try: p_sign = stats.binomtest(n_rising, n_total, 0.5).pvalue except AttributeError: p_sign = np.nan # Verdict: fraction > 0.6 = CONSISTENT, > 0.5 = AMBIGUOUS (was DIRECTIONAL) if frac_rising > 0.6 and (np.isnan(p_sign) or p_sign < 0.10): verd = 'CONSISTENT' elif frac_rising > 0.5: verd = 'AMBIGUOUS' else: verd = 'INCONSISTENT' return log_result('business_cycle_critical_slowing', frac_rising, p_sign, verd, details={'frac_rising': frac_rising, 'mean_tau': mean_tau, 'n_pairs': n_total}, description="Tests critical slowing down before recessions: rolling AR(1) coefficients (120-month window) should trend upward in 36-month pre-crisis windows. Uses Kendall tau of AR(1) vs time for each sector x recession pair, then binomial sign test on fraction with positive trend. Data: 7 FRED manufacturing IP subsectors over 4 major recessions.") # ── Test #14: Information temperature T = sigma^2/chi ── def test_temperature(growth): """Sector-level T = rolling_var / rolling_corr. Probit predicting recession 12mo ahead.""" try: usrec = load_fred('USREC') except FileNotFoundError: return log_result('business_cycle_info_temperature', np.nan, np.nan, 'AMBIGUOUS', description="Tests whether information temperature T = sigma^2/chi predicts recessions 12 months ahead via probit regression.") if 'INDPRO' not in growth: return log_result('business_cycle_info_temperature', np.nan, np.nan, 'AMBIGUOUS', description="Tests whether information temperature T = sigma^2/chi predicts recessions via probit. Insufficient INDPRO data.") g_indpro = growth['INDPRO'] # Compute sector-level T = sigma^2 / chi sector_T = {} for sid in MFG_SECTORS: if sid not in growth: continue g_s = growth[sid] common = g_s.dropna().index.intersection(g_indpro.dropna().index) if len(common) < 120: continue g_s_aligned = g_s.reindex(common) g_i_aligned = g_indpro.reindex(common) sigma2 = g_s_aligned.rolling(60, min_periods=30).var() chi_raw = g_s_aligned.rolling(60, min_periods=30).corr(g_i_aligned) chi = chi_raw.clip(lower=0.01) T_n = (sigma2 / chi).dropna() if len(T_n) > 60: sector_T[sid] = T_n if len(sector_T) < 3: return log_result('business_cycle_info_temperature', np.nan, np.nan, 'AMBIGUOUS', description="Tests whether information temperature T = sigma^2/chi predicts recessions. Insufficient sector-level T series.") # Aggregate T = mean across sectors T_df = pd.DataFrame(sector_T).dropna(how='all') T_agg = T_df.mean(axis=1).dropna() # Probit: P(recession_{t+12}) = Phi(alpha + beta*T) rec = usrec.reindex(T_agg.index, method='ffill').dropna() common = T_agg.index.intersection(rec.index) if len(common) < 60: return log_result('business_cycle_info_temperature', np.nan, np.nan, 'AMBIGUOUS', description="Tests whether information temperature T = sigma^2/chi predicts recessions. Insufficient aligned observations.") T_aligned = T_agg.reindex(common) rec_aligned = rec.reindex(common) # 12-month ahead recession indicator rec_ahead = rec_aligned.shift(-12).dropna() T_for_probit = T_aligned.reindex(rec_ahead.index) valid = T_for_probit.dropna().index.intersection(rec_ahead.dropna().index) if len(valid) < 60: return log_result('business_cycle_info_temperature', np.nan, np.nan, 'AMBIGUOUS', description="Tests whether information temperature T = sigma^2/chi predicts recessions. Insufficient valid probit observations.") y = (rec_ahead.reindex(valid).values > 0).astype(float) X_T = np.column_stack([np.ones(len(valid)), T_for_probit.reindex(valid).values]) try: import statsmodels.api as sm probit = sm.Probit(y, X_T).fit(disp=0) coef = probit.params[1] pval = probit.pvalues[1] pseudo_r2 = probit.prsquared except Exception: return log_result('business_cycle_info_temperature', np.nan, np.nan, 'AMBIGUOUS', description="Tests whether information temperature T = sigma^2/chi predicts recessions. Probit estimation failed.") if coef > 0 and pval < 0.10: verd = 'CONSISTENT' elif coef > 0: verd = 'AMBIGUOUS' else: verd = 'INCONSISTENT' return log_result('business_cycle_info_temperature', coef, pval, verd, details={'probit_coef': coef, 'probit_pval': pval, 'pseudo_r2': pseudo_r2}, description="Tests whether information temperature T = sigma^2/chi (rolling variance divided by rolling correlation with INDPRO) predicts recessions 12 months ahead. Uses probit regression P(recession_{t+12}) = Phi(alpha + beta*T) on aggregate T across 7 FRED manufacturing sectors. Prediction: beta > 0.") # ── Test #15: Regulation reduces amplitude (Great Moderation) ── def test_regulation(): """Pre-1984 vs 1984-2007 INDPRO volatility.""" try: indpro = load_fred('INDPRO') except FileNotFoundError: return log_result('business_cycle_regulation_amplitude', np.nan, np.nan, 'AMBIGUOUS', description="Tests whether post-Great-Moderation INDPRO volatility is lower than pre-1984 volatility, consistent with regulation reducing cycle amplitude.") g = np.log(indpro).diff().dropna() * 100 # percent pre_84 = g.loc[:'1983-12-31'] post_84 = g.loc['1984-01-01':'2007-11-30'] if len(pre_84) < 60 or len(post_84) < 60: return log_result('business_cycle_regulation_amplitude', np.nan, np.nan, 'AMBIGUOUS', description="Tests Great Moderation amplitude reduction. Insufficient INDPRO data in one or both sub-periods.") vol_pre = pre_84.std() vol_post = post_84.std() mean_pre = pre_84.mean() mean_post = post_84.mean() vol_ratio = vol_post / vol_pre # Verdict: vol ratio < 0.8 and mean growth roughly unchanged if vol_ratio < 0.8 and abs(mean_post - mean_pre) < 0.15: verd = 'CONSISTENT' elif vol_ratio < 0.8: verd = 'AMBIGUOUS' else: verd = 'INCONSISTENT' return log_result('business_cycle_regulation_amplitude', vol_ratio, np.nan, verd, details={'vol_pre': vol_pre, 'vol_post': vol_post, 'mean_pre': mean_pre, 'mean_post': mean_post}, description="Tests whether regulation reduces business cycle amplitude (Great Moderation). Compares INDPRO log-growth volatility pre-1984 vs 1984-2007, checking that the ratio < 0.8 while mean growth remains roughly unchanged.") # ── Test #16: Phillips slope proportional to K_bar ── def test_phillips_slope(): """Rolling Phillips slope vs time trend (secular flattening).""" try: cpi = load_fred('CPIAUCSL') gdp = load_fred('GDPC1') except FileNotFoundError: return log_result('business_cycle_phillips_slope', np.nan, np.nan, 'AMBIGUOUS', description="Tests whether the Phillips curve slope has secularly declined as predicted by falling K_bar. Requires CPI and GDP data.") # YoY inflation cpi_q = cpi.resample('QS').last() inflation = cpi_q.pct_change(4) * 100 # HP-filtered output gap from statsmodels.tsa.filters.hp_filter import hpfilter log_gdp = np.log(gdp.dropna()) cycle, _ = hpfilter(log_gdp, lamb=1600) gap = cycle * 100 common = inflation.dropna().index.intersection(gap.dropna().index) if len(common) < 40: return log_result('business_cycle_phillips_slope', np.nan, np.nan, 'AMBIGUOUS', description="Tests secular Phillips curve flattening. Insufficient common observations between inflation and output gap.") inf_a = inflation.reindex(common) gap_a = gap.reindex(common) # Rolling 5-year OLS slopes window_q = 20 rolling_slopes = [] for i in range(window_q, len(common)): sub_inf = inf_a.iloc[i - window_q:i] sub_gap = gap_a.iloc[i - window_q:i] if sub_gap.std() < 1e-6: continue slope, _, _, _, _ = stats.linregress(sub_gap.values, sub_inf.values) rolling_slopes.append(slope) if len(rolling_slopes) < 20: return log_result('business_cycle_phillips_slope', np.nan, np.nan, 'AMBIGUOUS', description="Tests secular Phillips curve flattening. Insufficient rolling windows.") # Secular trend: Kendall tau of slope vs time t_vals = np.arange(len(rolling_slopes)) tau, p = stats.kendalltau(t_vals, rolling_slopes) # Prediction: negative tau (K_bar declining → slope declining = flattening) # Phillips_flattening theorem: higher rho → lower K → lower |slope| # As economy shifts to services (higher rho), K_bar falls, slope flattens if tau < 0 and p < 0.10: verd = 'CONSISTENT' elif tau < 0: verd = 'AMBIGUOUS' else: verd = 'INCONSISTENT' if p < 0.10 else 'AMBIGUOUS' return log_result('business_cycle_phillips_slope', tau, p, verd, details={'kendall_tau': tau, 'n_windows': len(rolling_slopes)}, description="Tests whether the Phillips curve slope has secularly flattened, as predicted by declining aggregate K_bar (Corollary cor:flat-phillips). Computes rolling 5-year OLS slopes of YoY CPI inflation on HP-filtered output gap (GDPC1), then tests for a negative secular trend via Kendall tau of slope vs time. Prediction: tau < 0 (slope declining). Note: Paper's own remark acknowledges K_bar declined only ~1pp (0.555→0.545) over 1997-2024, likely too small for detection.") # ── Test #16a: Phillips slope pre-2020 (diagnostic) ── def test_phillips_pre2020(): """Same as test_phillips_slope but truncated before COVID to isolate secular trend.""" try: cpi = load_fred('CPIAUCSL') gdp = load_fred('GDPC1') except FileNotFoundError: return log_result('business_cycle_phillips_pre2020', np.nan, np.nan, 'AMBIGUOUS', description="Tests Phillips flattening pre-COVID (truncated at 2019Q4).") cpi_q = cpi.resample('QS').last() inflation = cpi_q.pct_change(4) * 100 from statsmodels.tsa.filters.hp_filter import hpfilter gdp_pre = gdp.loc[:'2019-12-31'].dropna() if len(gdp_pre) < 80: return log_result('business_cycle_phillips_pre2020', np.nan, np.nan, 'AMBIGUOUS', description="Tests Phillips flattening pre-COVID. Insufficient GDP data.") log_gdp = np.log(gdp_pre) cycle, _ = hpfilter(log_gdp, lamb=1600) gap = cycle * 100 inflation_pre = inflation.loc[:'2019-12-31'] common = inflation_pre.dropna().index.intersection(gap.dropna().index) if len(common) < 40: return log_result('business_cycle_phillips_pre2020', np.nan, np.nan, 'AMBIGUOUS', description="Tests Phillips flattening pre-COVID. Insufficient aligned data.") inf_a = inflation_pre.reindex(common) gap_a = gap.reindex(common) window_q = 20 rolling_slopes = [] for i in range(window_q, len(common)): sub_inf = inf_a.iloc[i - window_q:i] sub_gap = gap_a.iloc[i - window_q:i] if sub_gap.std() < 1e-6: continue slope, _, _, _, _ = stats.linregress(sub_gap.values, sub_inf.values) rolling_slopes.append(slope) if len(rolling_slopes) < 20: return log_result('business_cycle_phillips_pre2020', np.nan, np.nan, 'AMBIGUOUS', description="Tests Phillips flattening pre-COVID. Insufficient rolling windows.") t_vals = np.arange(len(rolling_slopes)) tau, p = stats.kendalltau(t_vals, rolling_slopes) if tau < 0 and p < 0.10: verd = 'CONSISTENT' elif tau < 0: verd = 'AMBIGUOUS' else: verd = 'INCONSISTENT' if p < 0.10 else 'AMBIGUOUS' _cache['phillips_pre2020'] = {'slopes': rolling_slopes, 'tau': tau, 'p': p} return log_result('business_cycle_phillips_pre2020', tau, p, verd, details={'kendall_tau': tau, 'n_windows': len(rolling_slopes)}, description="Diagnostic: Tests Phillips flattening truncated at 2019Q4 to remove post-COVID supply-side re-steepening. Same method as business_cycle_phillips_slope but on pre-pandemic sample only. Prediction: tau < 0.") # ── Test #16b: Cross-country Phillips slope vs manufacturing share ── def test_phillips_cross_country(): """Cross-country: countries with higher manufacturing share should have steeper Phillips curves.""" try: from common.data import load_oecd_stan, load_worldbank except ImportError: return log_result('business_cycle_phillips_cross_country', np.nan, np.nan, 'AMBIGUOUS', description="Tests cross-country Phillips: manufacturing share vs Phillips slope.") # Load manufacturing share from OECD STAN try: stan = load_oecd_stan() except FileNotFoundError: return log_result('business_cycle_phillips_cross_country', np.nan, np.nan, 'AMBIGUOUS', description="Tests cross-country Phillips. OECD STAN data not available.") # Load WB GDP growth (proxy for output gap) and CPI inflation try: gdp_g = load_worldbank('NY.GDP.MKTP.KD.ZG') cpi_inf = load_worldbank('FP.CPI.TOTL.ZG') except FileNotFoundError: return log_result('business_cycle_phillips_cross_country', np.nan, np.nan, 'AMBIGUOUS', description="Tests cross-country Phillips. World Bank data not available.") # OECD STAN uses ISO3 (USA, AUS, DEU); WB uses ISO2 (US, AU, DE) ISO3_TO_ISO2 = { 'AUS': 'AU', 'AUT': 'AT', 'BEL': 'BE', 'CAN': 'CA', 'CHE': 'CH', 'CHL': 'CL', 'COL': 'CO', 'CRI': 'CR', 'CZE': 'CZ', 'DEU': 'DE', 'DNK': 'DK', 'ESP': 'ES', 'EST': 'EE', 'FIN': 'FI', 'FRA': 'FR', 'GBR': 'GB', 'GRC': 'GR', 'HUN': 'HU', 'IRL': 'IE', 'ISL': 'IS', 'ISR': 'IL', 'ITA': 'IT', 'JPN': 'JP', 'KOR': 'KR', 'LTU': 'LT', 'LUX': 'LU', 'LVA': 'LV', 'MEX': 'MX', 'NLD': 'NL', 'NOR': 'NO', 'NZL': 'NZ', 'POL': 'PL', 'PRT': 'PT', 'SVK': 'SK', 'SVN': 'SI', 'SWE': 'SE', 'TUR': 'TR', 'USA': 'US', } # Compute manufacturing VA share per country (ISIC 'C' = manufacturing) # _T = total economy in OECD STAN (NOT BTE which is a narrow subset) mfg = stan[stan['activity'] == 'C'].copy() total = stan[stan['activity'] == '_T'].copy() if total.empty: total = stan[stan['activity'].isin(['TOTAL', 'TOT'])].copy() if mfg.empty or total.empty: return log_result('business_cycle_phillips_cross_country', np.nan, np.nan, 'AMBIGUOUS', description="Tests cross-country Phillips. Cannot identify manufacturing or total in STAN data.") # Use nominal VA (price_base='V') for consistent share computation if 'price_base' in mfg.columns: mfg = mfg[mfg['price_base'] == 'V'] total = total[total['price_base'] == 'V'] # Average manufacturing share 2000-2019 per country mfg_avg = mfg[(mfg['year'] >= 2000) & (mfg['year'] <= 2019)].groupby('country')['value'].mean() tot_avg = total[(total['year'] >= 2000) & (total['year'] <= 2019)].groupby('country')['value'].mean() common_countries = mfg_avg.index.intersection(tot_avg.index) mfg_share = (mfg_avg[common_countries] / tot_avg[common_countries]).dropna() mfg_share = mfg_share[(mfg_share > 0.01) & (mfg_share < 0.80)] # sanity filter # For each country with WB data, estimate Phillips slope # Phillips: inflation_t = a + b * gdp_growth_t (annual) # Use 2000-2019 to avoid COVID phillips_slopes = {} for country_iso3 in mfg_share.index: country_iso2 = ISO3_TO_ISO2.get(country_iso3) if country_iso2 is None: continue g = gdp_g[(gdp_g['country_code'] == country_iso2) & (gdp_g['year'] >= 2000) & (gdp_g['year'] <= 2019)].set_index('year')['value'].dropna() inf = cpi_inf[(cpi_inf['country_code'] == country_iso2) & (cpi_inf['year'] >= 2000) & (cpi_inf['year'] <= 2019)].set_index('year')['value'].dropna() common_yrs = g.index.intersection(inf.index) if len(common_yrs) < 10: continue g_a = g.reindex(common_yrs).values inf_a = inf.reindex(common_yrs).values if np.std(g_a) < 1e-6: continue slope, _, _, _, _ = stats.linregress(g_a, inf_a) phillips_slopes[country_iso3] = slope if len(phillips_slopes) < 8: return log_result('business_cycle_phillips_cross_country', np.nan, np.nan, 'AMBIGUOUS', details={'n_countries': len(phillips_slopes)}, description="Tests cross-country Phillips. Insufficient countries with both manufacturing share and Phillips slope data.") # ISO3 vs ISO3 alignment — OECD STAN may use 2-letter codes # Try matching slopes_s = pd.Series(phillips_slopes) common_final = mfg_share.index.intersection(slopes_s.index) if len(common_final) < 8: # OECD might use 2-letter codes while WB uses 3-letter — need concordance # Try: OECD 3-letter codes that match WB 3-letter codes return log_result('business_cycle_phillips_cross_country', np.nan, np.nan, 'AMBIGUOUS', details={'n_mfg': len(mfg_share), 'n_slopes': len(slopes_s), 'n_matched': len(common_final)}, description="Tests cross-country Phillips. Country code mismatch between OECD STAN and World Bank. N matched = %d." % len(common_final)) x_mfg = mfg_share[common_final].values y_slope = slopes_s[common_final].values # Spearman: manufacturing share vs Phillips slope # Theory predicts POSITIVE: more manufacturing → higher K_bar → steeper Phillips rho_sp, p_sp = stats.spearmanr(x_mfg, y_slope) if rho_sp > 0 and p_sp < 0.10: verd = 'CONSISTENT' elif rho_sp > 0: verd = 'AMBIGUOUS' else: verd = 'INCONSISTENT' if p_sp < 0.10 else 'AMBIGUOUS' _cache['phillips_cross_country'] = { 'mfg_share': x_mfg, 'slopes': y_slope, 'countries': list(common_final), 'rho': rho_sp, 'p': p_sp } return log_result('business_cycle_phillips_cross_country', rho_sp, p_sp, verd, details={'spearman_rho': rho_sp, 'n_countries': len(common_final)}, description="Tests cross-country Phillips curve prediction: countries with higher manufacturing share (proxy for higher K_bar) should have steeper Phillips curves. Computes manufacturing VA share from OECD STAN (2000-2019 avg) and Phillips slope (OLS of WB CPI inflation on GDP growth, 2000-2019) per country, then Spearman rank correlation. Prediction: rho > 0.") # ── Test #17: rho-diversity predicts resilience ── def test_rho_diversity(): """Higher rho-diversity -> shallower recessions. (n=2 post-2005, AMBIGUOUS.)""" try: indpro = load_fred('INDPRO') except FileNotFoundError: return log_result('business_cycle_rho_diversity', np.nan, np.nan, 'AMBIGUOUS', description="Tests whether higher cross-sector rho-diversity predicts shallower recessions. Requires INDPRO data.") # Only 2 recessions in post-2005 era with sector data # n=2 is too few for inference return log_result('business_cycle_rho_diversity', np.nan, np.nan, 'AMBIGUOUS', details={'reason': 'n=2 recessions with diversity data, too few'}, description="Tests whether higher cross-sector rho-diversity predicts shallower recessions. Only 2 post-2005 recessions have sector-level diversity data, which is insufficient for statistical inference.") # ── Test #18: Power-law recession depths ── def test_power_law_depths(): """Hill estimator of recession depth tail exponent.""" try: indpro = load_fred('INDPRO') except FileNotFoundError: return log_result('business_cycle_power_law_depths', np.nan, np.nan, 'AMBIGUOUS', description="Tests whether recession depths follow a power-law distribution using Hill estimation on INDPRO peak-to-trough declines.") depths = [] for peak, trough in RECESSIONS: peak_dt = pd.Timestamp(peak) trough_dt = pd.Timestamp(trough) window = indpro.loc[peak_dt - pd.DateOffset(months=3):trough_dt + pd.DateOffset(months=3)] if len(window) > 0: peak_val = window.max() trough_val = window.min() depths.append((peak_val - trough_val) / peak_val * 100) if len(depths) < 5: return log_result('business_cycle_power_law_depths', np.nan, np.nan, 'AMBIGUOUS', description="Tests power-law recession depths. Insufficient recession episodes for Hill estimation.") depth_vals = np.sort(depths)[::-1] n = len(depth_vals) # Hill estimator x_min = np.median(depth_vals) tail = depth_vals[depth_vals >= x_min] if len(tail) >= 3: alpha_hill = 1 + len(tail) / np.sum(np.log(tail / x_min)) else: alpha_hill = np.nan # Verdict: alpha in plausible range [1, 4] = AMBIGUOUS (suggestive, n too small) if not np.isnan(alpha_hill) and 1.0 < alpha_hill < 4.0: verd = 'AMBIGUOUS' # old: SUGGESTIVE (n=7 too few) else: verd = 'AMBIGUOUS' return log_result('business_cycle_power_law_depths', alpha_hill, np.nan, verd, details={'hill_alpha': alpha_hill, 'n_recessions': n}, description="Tests whether recession depths follow a power-law distribution. Applies Hill estimator to INDPRO peak-to-trough percentage declines across post-war NBER recessions. Prediction: tail exponent alpha in [1, 4], consistent with CES landscape Kramers escape statistics.") # ── Main ── _cache = {} def main(): """Run all business cycle tests. Returns list of result dicts.""" raw, growth = _load_mfg_data() if len(growth) < 5: return [log_result('business_cycle_predictions', np.nan, np.nan, 'AMBIGUOUS', details={'reason': f'Only {len(growth)} sectors loaded'}, description="Combined business cycle predictions test. Insufficient manufacturing IP sectors loaded for testing.")] results = [] results.append(test_rho_ordering(growth)) results.append(test_asymmetry()) results.append(test_critical_slowing(growth)) results.append(test_temperature(growth)) results.append(test_regulation()) results.append(test_phillips_slope()) results.append(test_phillips_pre2020()) results.append(test_phillips_cross_country()) results.append(test_rho_diversity()) results.append(test_power_law_depths()) # Combined verdict n_consistent = sum(1 for r in results if r['verdict'] == 'CONSISTENT') n_total = len(results) if n_consistent >= 3: combined = 'CONSISTENT' elif n_consistent >= 1: combined = 'AMBIGUOUS' else: combined = 'INCONSISTENT' results.append(log_result('business_cycle_combined', float(n_consistent), np.nan, combined, details={'n_consistent': n_consistent, 'n_total': n_total}, description=f"Combined verdict across all business cycle sub-tests. {n_consistent} of {n_total} CONSISTENT.")) _cache['results'] = results _cache['growth'] = growth _cache['raw'] = raw # Cache rho-ordering scatter data mfg_recessions = [r for r in RECESSIONS if r[0] >= "1973-01-01"] all_leads = [] for peak, trough in mfg_recessions: peaks = _find_peak_timing(growth, peak) if len(peaks) >= 4: for sid, info in peaks.items(): all_leads.append({'rho': info['rho'], 'lead_months': info['lead_months'], 'recession': peak[:4]}) if all_leads: _cache['rho_ordering'] = pd.DataFrame(all_leads) # Cache expansion/contraction ratios ratios = [] for i in range(1, len(RECESSIONS)): peak, trough = RECESSIONS[i] peak_dt = pd.Timestamp(peak) trough_dt = pd.Timestamp(trough) contraction = (trough_dt - peak_dt).days / 30.44 prev_trough = pd.Timestamp(RECESSIONS[i - 1][1]) expansion = (peak_dt - prev_trough).days / 30.44 if contraction > 0: ratios.append({'cycle': i, 'ratio': expansion / contraction, 'peak': peak[:4]}) if ratios: _cache['asymmetry'] = pd.DataFrame(ratios) return results def figures(): """Generate per-subtest figures for business cycle predictions.""" if not _cache: main() paths = [] results = _cache.get('results', []) growth = _cache.get('growth', {}) if not results: return paths import matplotlib.pyplot as plt # --- Figure 1: rho-ordering scatter --- rho_df = _cache.get('rho_ordering') if rho_df is not None and len(rho_df) > 0: fig, ax = setup_figure() rho_r = [r for r in results if r['name'] == 'business_cycle_rho_ordering'] ax.scatter(rho_df['rho'], rho_df['lead_months'], alpha=0.5, s=30, c='steelblue') # Add OLS line x = rho_df['rho'].values y = rho_df['lead_months'].values m, b = np.polyfit(x, y, 1) x_line = np.linspace(x.min(), x.max(), 50) ax.plot(x_line, m * x_line + b, 'r--', alpha=0.7) if rho_r: tau = rho_r[0].get('kendall_tau', np.nan) p = rho_r[0].get('p_value', np.nan) ax.text(0.05, 0.95, f'Kendall tau={tau:.3f}\np={p:.1e}\nN={len(rho_df)}', transform=ax.transAxes, va='top', fontsize=9, bbox=dict(boxstyle='round', facecolor='wheat', alpha=0.5)) ax.set_xlabel(r'Calibrated $\hat\rho$') ax.set_ylabel('Lead months before NBER peak') ax.set_title('Rho-Ordering of Recession Entry') paths.extend(save_figure(fig, 'business_cycle_rho_ordering', subfolder='dynamics')) # --- Figure 2: expansion/contraction asymmetry --- asym_df = _cache.get('asymmetry') if asym_df is not None and len(asym_df) > 0: fig, ax = setup_figure() ax.bar(range(len(asym_df)), asym_df['ratio'].values, color='steelblue', alpha=0.7) ax.axhline(10, color='red', linestyle='--', linewidth=2, label='Predicted ~10') ax.axhspan(3, 15, alpha=0.1, color='green', label='CONSISTENT range [3,15]') ax.set_xticks(range(len(asym_df))) ax.set_xticklabels(asym_df['peak'].values, rotation=45, fontsize=8) ax.set_ylabel('Expansion / Contraction duration ratio') ax.set_title('Duration Asymmetry by Cycle') ax.text(0.95, 0.95, f'Mean={asym_df["ratio"].mean():.1f}\nMedian={asym_df["ratio"].median():.1f}', transform=ax.transAxes, va='top', ha='right', fontsize=9, bbox=dict(boxstyle='round', facecolor='wheat', alpha=0.5)) ax.legend(loc='upper left') paths.extend(save_figure(fig, 'business_cycle_expansion_asymmetry', subfolder='dynamics')) # --- Figure 3: critical slowing down --- if growth: sector_ar1 = {} for sid in MFG_SECTORS: if sid not in growth: continue g = growth[sid] if len(g) < 180: continue sector_ar1[sid] = _rolling_ar1(g, window=120) if sector_ar1: fig, ax = setup_figure() for sid, ar1 in sector_ar1.items(): name = MFG_SECTORS[sid][0] ar1_clean = ar1.dropna() if len(ar1_clean) > 0: ax.plot(ar1_clean.index, ar1_clean.values, alpha=0.5, linewidth=0.8, label=name) # Shade recessions for peak, trough in RECESSIONS: ax.axvspan(pd.Timestamp(peak), pd.Timestamp(trough), alpha=0.15, color='red') ax.set_ylabel('Rolling AR(1) coefficient') ax.set_title('Critical Slowing Down: Rolling AR(1)') ax.legend(fontsize=7, loc='lower left') csd_r = [r for r in results if r['name'] == 'business_cycle_critical_slowing'] if csd_r: frac = csd_r[0].get('frac_rising', np.nan) ax.text(0.95, 0.95, f'Fraction rising\npre-crisis: {frac:.2f}', transform=ax.transAxes, va='top', ha='right', fontsize=9, bbox=dict(boxstyle='round', facecolor='wheat', alpha=0.5)) paths.extend(save_figure(fig, 'business_cycle_critical_slowing', subfolder='dynamics')) # --- Figure 4: information temperature --- temp_r = [r for r in results if r['name'] == 'business_cycle_info_temperature'] if temp_r and not np.isnan(temp_r[0].get('probit_coef', np.nan)): fig, ax = setup_figure() # Reconstruct T_agg for plotting try: usrec = load_fred('USREC') g_indpro = growth.get('INDPRO') if g_indpro is not None: sector_T = {} for sid in MFG_SECTORS: if sid not in growth: continue g_s = growth[sid] common = g_s.dropna().index.intersection(g_indpro.dropna().index) if len(common) < 120: continue g_s_a = g_s.reindex(common) g_i_a = g_indpro.reindex(common) sigma2 = g_s_a.rolling(60, min_periods=30).var() chi_raw = g_s_a.rolling(60, min_periods=30).corr(g_i_a) chi = chi_raw.clip(lower=0.01) T_n = (sigma2 / chi).dropna() if len(T_n) > 60: sector_T[sid] = T_n if sector_T: T_agg = pd.DataFrame(sector_T).dropna(how='all').mean(axis=1).dropna() ax.plot(T_agg.index, T_agg.values, 'b-', alpha=0.7, linewidth=1) for peak, trough in RECESSIONS: ax.axvspan(pd.Timestamp(peak), pd.Timestamp(trough), alpha=0.15, color='red') ax.set_ylabel(r'$T = \sigma^2 / \chi$') ax.set_title('Information Temperature as Leading Indicator') coef = temp_r[0]['probit_coef'] pval = temp_r[0].get('probit_pval', np.nan) ax.text(0.05, 0.95, f'Probit coef={coef:.1f}\np={pval:.4f}', transform=ax.transAxes, va='top', fontsize=9, bbox=dict(boxstyle='round', facecolor='wheat', alpha=0.5)) paths.extend(save_figure(fig, 'business_cycle_info_temperature', subfolder='dynamics')) else: plt.close(fig) else: plt.close(fig) except Exception: plt.close(fig) # --- Figure 5: regulation amplitude (Great Moderation) --- reg_r = [r for r in results if r['name'] == 'business_cycle_regulation_amplitude'] if reg_r and not np.isnan(reg_r[0].get('statistic', np.nan)): try: indpro = load_fred('INDPRO') g = np.log(indpro).diff().dropna() * 100 fig, ax = setup_figure() rolling_vol = g.rolling(60).std() ax.plot(rolling_vol.index, rolling_vol.values, 'b-', alpha=0.7, linewidth=1) ax.axvline(pd.Timestamp('1984-01-01'), color='red', linestyle='--', linewidth=2, label='Great Moderation (1984)') ax.set_ylabel('Rolling 5-year volatility (%)') ax.set_title('Great Moderation: INDPRO Volatility') vol_ratio = reg_r[0]['statistic'] ax.text(0.95, 0.95, f'Post/Pre ratio={vol_ratio:.3f}', transform=ax.transAxes, va='top', ha='right', fontsize=9, bbox=dict(boxstyle='round', facecolor='wheat', alpha=0.5)) ax.legend() paths.extend(save_figure(fig, 'business_cycle_regulation_amplitude', subfolder='dynamics')) except Exception: pass # --- Figure 6: Phillips slope trend --- phil_r = [r for r in results if r['name'] == 'business_cycle_phillips_slope'] if phil_r and not np.isnan(phil_r[0].get('kendall_tau', np.nan)): try: cpi = load_fred('CPIAUCSL') gdp = load_fred('GDPC1') cpi_q = cpi.resample('QS').last() inflation = cpi_q.pct_change(4) * 100 from statsmodels.tsa.filters.hp_filter import hpfilter log_gdp = np.log(gdp.dropna()) cycle, _ = hpfilter(log_gdp, lamb=1600) gap = cycle * 100 common = inflation.dropna().index.intersection(gap.dropna().index) inf_a = inflation.reindex(common) gap_a = gap.reindex(common) window_q = 20 slopes, dates = [], [] for i in range(window_q, len(common)): sub_inf = inf_a.iloc[i - window_q:i] sub_gap = gap_a.iloc[i - window_q:i] if sub_gap.std() < 1e-6: continue slope, _, _, _, _ = stats.linregress(sub_gap.values, sub_inf.values) slopes.append(slope) dates.append(common[i]) if slopes: fig, ax = setup_figure() ax.plot(dates, slopes, 'b-', alpha=0.7, linewidth=1) ax.axhline(0, color='gray', linestyle='-', alpha=0.3) for peak, trough in RECESSIONS: ax.axvspan(pd.Timestamp(peak), pd.Timestamp(trough), alpha=0.1, color='red') ax.set_ylabel('Rolling 5-year Phillips slope') ax.set_title('Phillips Curve Slope Over Time') tau = phil_r[0]['kendall_tau'] p = phil_r[0].get('p_value', np.nan) ax.text(0.05, 0.95, f'Kendall tau={tau:.3f}\np={p:.3f}\nPrediction: tau < 0', transform=ax.transAxes, va='top', fontsize=9, bbox=dict(boxstyle='round', facecolor='wheat', alpha=0.5)) paths.extend(save_figure(fig, 'business_cycle_phillips_slope', subfolder='dynamics')) except Exception: pass # --- Figure 7: power-law depths --- try: indpro = load_fred('INDPRO') depths = [] for peak, trough in RECESSIONS: peak_dt = pd.Timestamp(peak) trough_dt = pd.Timestamp(trough) window = indpro.loc[peak_dt - pd.DateOffset(months=3):trough_dt + pd.DateOffset(months=3)] if len(window) > 0: peak_val = window.max() trough_val = window.min() depths.append((peak_val - trough_val) / peak_val * 100) if depths: fig, ax = setup_figure() depths_sorted = np.sort(depths)[::-1] ranks = np.arange(1, len(depths_sorted) + 1) ax.scatter(np.log10(depths_sorted), np.log10(ranks / len(depths_sorted)), s=50, c='steelblue', zorder=5) ax.set_xlabel('log10(Recession depth %)') ax.set_ylabel('log10(Exceedance probability)') ax.set_title('Recession Depth Distribution') pl_r = [r for r in results if r['name'] == 'business_cycle_power_law_depths'] if pl_r: alpha = pl_r[0].get('hill_alpha', np.nan) ax.text(0.95, 0.95, f'Hill alpha={alpha:.2f}\nn={len(depths)}', transform=ax.transAxes, va='top', ha='right', fontsize=9, bbox=dict(boxstyle='round', facecolor='wheat', alpha=0.5)) paths.extend(save_figure(fig, 'business_cycle_power_law_depths', subfolder='dynamics')) except Exception: pass # --- Figure 7b: Phillips pre-2020 diagnostic --- pre2020 = _cache.get('phillips_pre2020') if pre2020 and pre2020.get('slopes'): fig, ax = setup_figure() n_s = len(pre2020['slopes']) ax.plot(range(n_s), pre2020['slopes'], 'b-', alpha=0.7, linewidth=1) ax.axhline(0, color='gray', linestyle='-', alpha=0.3) ax.set_ylabel('Rolling 5-year Phillips slope') ax.set_title('Phillips Slope Pre-2020 (COVID removed)') tau_p = pre2020['tau'] p_p = pre2020['p'] ax.text(0.05, 0.95, f'Kendall tau={tau_p:.3f}\np={p_p:.3f}\nPrediction: tau < 0', transform=ax.transAxes, va='top', fontsize=9, bbox=dict(boxstyle='round', facecolor='wheat', alpha=0.5)) paths.extend(save_figure(fig, 'business_cycle_phillips_pre2020', subfolder='dynamics')) # --- Figure 7c: Cross-country Phillips scatter --- xc = _cache.get('phillips_cross_country') if xc and len(xc.get('countries', [])) > 0: fig, ax = setup_figure() ax.scatter(xc['mfg_share'], xc['slopes'], s=40, c='steelblue', alpha=0.7, zorder=5) for i, c in enumerate(xc['countries']): ax.annotate(c, (xc['mfg_share'][i], xc['slopes'][i]), fontsize=6, alpha=0.7, xytext=(3, 3), textcoords='offset points') # Trend line m, b = np.polyfit(xc['mfg_share'], xc['slopes'], 1) x_line = np.linspace(min(xc['mfg_share']), max(xc['mfg_share']), 50) ax.plot(x_line, m * x_line + b, 'r--', alpha=0.5) ax.set_xlabel('Manufacturing VA share (2000-2019 avg)') ax.set_ylabel('Phillips slope (inflation ~ GDP growth)') ax.set_title('Cross-Country: Manufacturing Share vs Phillips Slope') ax.text(0.05, 0.95, f'Spearman rho={xc["rho"]:.3f}\np={xc["p"]:.3f}\nN={len(xc["countries"])}', transform=ax.transAxes, va='top', fontsize=9, bbox=dict(boxstyle='round', facecolor='wheat', alpha=0.5)) paths.extend(save_figure(fig, 'business_cycle_phillips_cross_country', subfolder='dynamics')) # --- Figure 8: combined verdict summary --- fig, ax = setup_figure(width=8, height=4) names = [r['name'].replace('business_cycle_', '') for r in results if 'combined' not in r['name']] verdicts = [r['verdict'] for r in results if 'combined' not in r['name']] colors = ['#2ecc71' if v == 'CONSISTENT' else '#f39c12' if v == 'AMBIGUOUS' else '#e74c3c' for v in verdicts] y_pos = range(len(names)) ax.barh(y_pos, [1]*len(names), color=colors, height=0.6) ax.set_yticks(y_pos) ax.set_yticklabels(names, fontsize=9) ax.set_xlim(0, 1.3) ax.set_xticks([]) for i, v in enumerate(verdicts): ax.text(1.05, i, v, va='center', fontsize=8) n_con = sum(1 for v in verdicts if v == 'CONSISTENT') ax.set_title(f'Business Cycle Predictions: {n_con}/{len(verdicts)} CONSISTENT') paths.extend(save_figure(fig, 'business_cycle_combined', subfolder='dynamics')) return paths if __name__ == "__main__": res = main() for r in res: print(f"{r['name']}: {r['verdict']} (stat={r['statistic']}, p={r['p_value']})")

statistic	3.5431
hill alpha	3.5431
n recessions	7

Business Cycle Power Law Depths

Description

Figures

Statistics

Script

▸ Source Code

Reproduction Steps