diff --git a/bedrock/utils/validation/__init__.py b/bedrock/utils/validation/__init__.py index b9b550ec..eb17aec8 100644 --- a/bedrock/utils/validation/__init__.py +++ b/bedrock/utils/validation/__init__.py @@ -3,16 +3,22 @@ from bedrock.utils.validation.eeio_diagnostics import ( DiagnosticResult, + assert_eeio_year_alignment_precondition, compare_commodity_output_to_domestics_use_plus_exports, + compare_E_and_LCI_result, compare_output_vs_leontief_x_demand, + eeio_year_alignment_precondition_ok, format_diagnostic_result, run_all_diagnostics, ) __all__ = [ "DiagnosticResult", + "assert_eeio_year_alignment_precondition", + "compare_E_and_LCI_result", "compare_commodity_output_to_domestics_use_plus_exports", "compare_output_vs_leontief_x_demand", + "eeio_year_alignment_precondition_ok", "format_diagnostic_result", "run_all_diagnostics", ] diff --git a/bedrock/utils/validation/__tests__/test_eeio_diagnostics.py b/bedrock/utils/validation/__tests__/test_eeio_diagnostics.py index 13bf11ab..11143e4f 100644 --- a/bedrock/utils/validation/__tests__/test_eeio_diagnostics.py +++ b/bedrock/utils/validation/__tests__/test_eeio_diagnostics.py @@ -1,6 +1,11 @@ # ruff: noqa: PLC0415 """Unit tests for the EEIO diagnostics module.""" +from __future__ import annotations + +import dataclasses as dc + +import numpy as np import pandas as pd import pytest @@ -23,11 +28,19 @@ derive_cornerstone_y_nab, derive_cornerstone_Ytot_matrix_set, ) +from bedrock.utils.config.usa_config import ( + USAConfig, + reset_usa_config, + set_global_usa_config, +) from bedrock.utils.math.formulas import compute_y_imp from bedrock.utils.validation.eeio_diagnostics import ( DiagnosticResult, + assert_eeio_year_alignment_precondition, compare_commodity_output_to_domestics_use_plus_exports, + compare_E_and_LCI_result, compare_output_vs_leontief_x_demand, + eeio_year_alignment_precondition_ok, format_diagnostic_result, run_all_diagnostics, validate_result, @@ -454,3 +467,194 @@ def test_compare_output_and_L_y( output=output, L=L, y=y, tolerance=0.01, include_details=True ) assert len(r_output_L_y_validation.failing_sectors) == 0 + + +# --------------------------------------------------------------------------- +# useeior-parity LCI ≈ E + year-alignment precondition +# --------------------------------------------------------------------------- + + +@dc.dataclass(frozen=True) +class _ToyLciMatrices: + B: pd.DataFrame + L: pd.DataFrame + y: pd.Series[float] + E_ind: pd.DataFrame + V: pd.DataFrame + x: pd.Series[float] + q: pd.Series[float] + + +def _toy_lci_matrices() -> _ToyLciMatrices: + """Minimal commodity LCI≈E system with C_m = I and L = I.""" + industries = ['i1', 'i2'] + commodities = ['c1', 'c2'] + flows = ['CO2'] + + V = pd.DataFrame( + [[10.0, 0.0], [0.0, 20.0]], + index=industries, + columns=commodities, + ) + x = V.sum(axis=1) + q = pd.Series([10.0, 20.0], index=commodities) + E_ind = pd.DataFrame([[100.0, 200.0]], index=flows, columns=industries) + # With C_m = I, E_c = E_ind; B = E_c / q + B = pd.DataFrame([[10.0, 10.0]], index=flows, columns=commodities) + L = pd.DataFrame( + np.eye(2), + index=commodities, + columns=commodities, + ) + y = q.copy() + return _ToyLciMatrices(B=B, L=L, y=y, E_ind=E_ind, V=V, x=x, q=q) + + +def test_eeio_year_alignment_precondition_ok() -> None: + cfg = USAConfig.model_validate( + { + 'model_base_year': 2024, + 'usa_ghg_data_year': 2024, + 'use_E_data_year_for_x_in_B': True, + }, + strict=True, + ) + assert eeio_year_alignment_precondition_ok(cfg) + assert_eeio_year_alignment_precondition(cfg) + + +def test_eeio_year_alignment_precondition_fails_loud() -> None: + cfg = USAConfig.model_validate( + { + 'model_base_year': 2024, + 'usa_ghg_data_year': 2023, + 'use_E_data_year_for_x_in_B': True, + }, + strict=True, + ) + assert not eeio_year_alignment_precondition_ok(cfg) + with pytest.raises(ValueError, match='precondition failed'): + assert_eeio_year_alignment_precondition(cfg) + + +def test_toy_lci_passes_when_LCI_equals_E_c() -> None: + m = _toy_lci_matrices() + result = compare_E_and_LCI_result( + B=m.B, + L=m.L, + y=m.y, + E_ind=m.E_ind, + V=m.V, + x=m.x, + tolerance=0.01, + check_precondition=False, + ) + assert result.passed + assert result.failing_sectors == [] + + +def test_toy_lci_fails_when_E_perturbed() -> None: + m = _toy_lci_matrices() + E_ind = m.E_ind.copy() + E_ind.iloc[0, 0] = 1000.0 + result = compare_E_and_LCI_result( + B=m.B, + L=m.L, + y=m.y, + E_ind=E_ind, + V=m.V, + x=m.x, + tolerance=0.01, + check_precondition=False, + ) + assert not result.passed + assert len(result.failing_sectors) > 0 + + +@pytest.mark.eeio_integration +def test_v0_3_domestic_lci_equals_e() -> None: + """Live domestic LCI≈E under aligned v0.3 (Vnorm path, χ=1 precondition).""" + reset_usa_config(should_reset_env_var=True) + set_global_usa_config('2025_usa_cornerstone_v0_3.yaml') + try: + from bedrock.transform.allocation.derived import derive_E_usa + from bedrock.transform.eeio.derived_cornerstone import ( + derive_cornerstone_Aq_scaled, + derive_cornerstone_B_non_finetuned, + derive_cornerstone_Vnorm_scrap_corrected, + derive_cornerstone_x_after_redefinition, + derive_cornerstone_y_nab, + ) + from bedrock.utils.math.formulas import compute_L_matrix + + assert_eeio_year_alignment_precondition() + + Aq = derive_cornerstone_Aq_scaled() + L_d = compute_L_matrix(A=Aq.Adom) + y_d = derive_cornerstone_y_nab() + result = compare_E_and_LCI_result( + B=derive_cornerstone_B_non_finetuned(), + L=L_d, + y=y_d, + E_ind=derive_E_usa(), + x=derive_cornerstone_x_after_redefinition(), + Vnorm=derive_cornerstone_Vnorm_scrap_corrected(), + q=Aq.scaled_q, + tolerance=0.01, + include_details=True, + check_precondition=False, + ) + assert result.passed, ( + f'domestic LCI≈E failed ({len(result.failing_sectors)} cells): ' + f'{result.failing_sectors[:20]}' + ) + finally: + reset_usa_config(should_reset_env_var=True) + + +@pytest.mark.eeio_integration +@pytest.mark.xfail( + reason=( + 'Cornerstone total LCI≈E uses y_trade (ytot/trade), not IO-balanced ' + 'throughput; same imbalance as total Ly≈q (~1749 cells fail).' + ), +) +def test_v0_3_total_lci_equals_e() -> None: + """Live total LCI≈E under aligned v0.3 (Vnorm path; L_tot / y_trade pairing).""" + reset_usa_config(should_reset_env_var=True) + set_global_usa_config('2025_usa_cornerstone_v0_3.yaml') + try: + from bedrock.transform.allocation.derived import derive_E_usa + from bedrock.transform.eeio.derived_cornerstone import ( + derive_cornerstone_Aq_scaled, + derive_cornerstone_B_non_finetuned, + derive_cornerstone_Vnorm_scrap_corrected, + derive_cornerstone_x_after_redefinition, + derive_cornerstone_Y_and_trade_scaled, + ) + from bedrock.utils.math.formulas import compute_L_matrix + + assert_eeio_year_alignment_precondition() + + Aq = derive_cornerstone_Aq_scaled() + L_tot = compute_L_matrix(A=Aq.Adom + Aq.Aimp) + y_trade = derive_cornerstone_Y_and_trade_scaled() + y_tot = y_trade.ytot + y_trade.exports - y_trade.imports + result = compare_E_and_LCI_result( + B=derive_cornerstone_B_non_finetuned(), + L=L_tot, + y=y_tot, + E_ind=derive_E_usa(), + x=derive_cornerstone_x_after_redefinition(), + Vnorm=derive_cornerstone_Vnorm_scrap_corrected(), + q=Aq.scaled_q, + tolerance=0.01, + include_details=True, + check_precondition=False, + ) + assert result.passed, ( + f'total LCI≈E failed ({len(result.failing_sectors)} cells): ' + f'{result.failing_sectors[:20]}' + ) + finally: + reset_usa_config(should_reset_env_var=True) diff --git a/bedrock/utils/validation/eeio_diagnostics.py b/bedrock/utils/validation/eeio_diagnostics.py index a8dfe94d..a5a2a081 100644 --- a/bedrock/utils/validation/eeio_diagnostics.py +++ b/bedrock/utils/validation/eeio_diagnostics.py @@ -14,9 +14,11 @@ import numpy as np import pandas as pd +from bedrock.utils.config.usa_config import USAConfig, get_usa_config from bedrock.utils.math.formulas import ( backcompute_q_from_L_and_y, compute_commodity_mix_matrix, + compute_E_from_BLy, ) from bedrock.utils.schemas.single_region_types import SingleRegionYtotAndTradeVectorSet @@ -465,3 +467,202 @@ def compare_output_from_make_and_use( 'invalid output parameter requested for comparison between make and use, select commodity or industry' ) return d_result + + +def eeio_year_alignment_precondition_ok(cfg: USAConfig) -> bool: + """Return True when χ=1 is a valid surrogate for useeior Chi year alignment. + + Requires matching model/GHG years, ``use_E_data_year_for_x_in_B``, and no + scaled/deflated B paths that introduce intermediate dollar years. + """ + if cfg.model_base_year != cfg.usa_ghg_data_year: + return False + if not cfg.use_E_data_year_for_x_in_B: + return False + if cfg.use_scaled_x_and_scaled_Vnorm_for_B: + return False + if cfg.deflate_x_to_detail_io_year_for_B: + return False + return True + + +def assert_eeio_year_alignment_precondition( + cfg: ta.Optional[USAConfig] = None, +) -> None: + """Raise if the active (or given) config is not aligned for χ=1 LCI≈E checks. + + Replaces useeior ``generateChiMatrix``: when this passes, ``compare_E_and_LCI_result`` + may use χ=1. Failures are loud (``ValueError``), never silent skips. + + Parameters + ---------- + cfg + Config to check; defaults to ``get_usa_config()``. + """ + if cfg is None: + cfg = get_usa_config() + if eeio_year_alignment_precondition_ok(cfg): + return + + reasons: list[str] = [] + if cfg.model_base_year != cfg.usa_ghg_data_year: + reasons.append( + f'model_base_year ({cfg.model_base_year}) != ' + f'usa_ghg_data_year ({cfg.usa_ghg_data_year})' + ) + if not cfg.use_E_data_year_for_x_in_B: + reasons.append('use_E_data_year_for_x_in_B is False') + if cfg.use_scaled_x_and_scaled_Vnorm_for_B: + reasons.append( + 'use_scaled_x_and_scaled_Vnorm_for_B is True ' + '(intermediate dollar years break χ=1)' + ) + if cfg.deflate_x_to_detail_io_year_for_B: + reasons.append( + 'deflate_x_to_detail_io_year_for_B is True ' + '(intermediate dollar years break χ=1)' + ) + raise ValueError( + 'EEIO year-alignment precondition failed for χ=1 LCI≈E validation: ' + + '; '.join(reasons) + ) + + +def _flatten_matrix_for_validate(df: pd.DataFrame) -> pd.Series[float]: + """Stack a flow×sector matrix to a Series with MultiIndex labels ``flow|sector``.""" + stacked = df.stack() + stacked.index = stacked.index.map( + lambda idx: f'{idx[0]}|{idx[1]}' if isinstance(idx, tuple) else str(idx) + ) + return ta.cast('pd.Series[float]', stacked.astype(float)) + + +def compare_E_and_LCI_result( + *, + B: pd.DataFrame, + L: pd.DataFrame, + y: pd.Series[float], + E_ind: pd.DataFrame, + V: ta.Optional[pd.DataFrame] = None, + x: ta.Optional[pd.Series[float]] = None, + Vnorm: ta.Optional[pd.DataFrame] = None, + q: ta.Optional[pd.Series[float]] = None, + tolerance: float = 0.01, + include_details: bool = False, + check_precondition: bool = True, + cfg: ta.Optional[USAConfig] = None, +) -> DiagnosticResult: + """Compare direct-perspective LCI to commodity-transformed satellite totals. + + Port of useeior ``compareEandLCIResult`` with χ=1 (no Chi matrix). Requires + :func:`assert_eeio_year_alignment_precondition` unless + ``check_precondition=False`` (unit tests with synthetic matrices). + + Commodity ``E_c`` (pick one path):: + + # A) Vnorm path (Cornerstone B = (E/x) @ Vnorm_scrap): preferred when + # scrap-corrected Vnorm is used — C_m(V, x_Make) will not match. + E_c = (E_ind / x @ Vnorm) · diag(q) + + # B) useeior C_m path (no scrap adjustment on market shares): + E_c = (C_m @ E_ind.T).T # C_m from Make V, x + + c = L @ y + LCI = B · diag(c) # B already commodity — do not @ V_n again + compare LCI to E_c cell-wise + + Provide either (``Vnorm``, ``x``, ``q``) or (``V``, ``x``). Distinct from + NAB national ``sum(diag(D) @ L @ y) ≈ sum(E)``. + """ + if check_precondition: + assert_eeio_year_alignment_precondition(cfg) + + use_vnorm_path = Vnorm is not None and q is not None and x is not None + use_cm_path = V is not None and x is not None and not use_vnorm_path + if not use_vnorm_path and not use_cm_path: + raise ValueError( + 'compare_E_and_LCI_result requires either (Vnorm, x, q) or (V, x)' + ) + + flows = B.index.intersection(E_ind.index) + if len(flows) == 0: + return DiagnosticResult( + name='compare_E_and_LCI_result: empty flow intersection', + passed=False, + tolerance=tolerance, + max_rel_diff=float('inf'), + failing_sectors=[], + details=None, + ) + + if use_vnorm_path: + assert Vnorm is not None and q is not None and x is not None + industries = Vnorm.index.intersection(E_ind.columns).intersection(x.index) + commodities = ( + B.columns.intersection(Vnorm.columns) + .intersection(q.index) + .intersection(L.index) + .intersection(y.index) + ) + if len(industries) == 0 or len(commodities) == 0: + return DiagnosticResult( + name='compare_E_and_LCI_result: empty industry/commodity intersection', + passed=False, + tolerance=tolerance, + max_rel_diff=float('inf'), + failing_sectors=[], + details=None, + ) + Bi = ( + E_ind.loc[flows, industries] + .divide(x.reindex(industries).fillna(0.0), axis=1) + .fillna(0.0) + ) + E_c = Bi @ Vnorm.loc[industries, commodities] + E_c = E_c.multiply(q.reindex(commodities).fillna(0.0), axis=1) + else: + assert V is not None and x is not None + C_m = compute_commodity_mix_matrix(V=V, x=x) + industries = C_m.columns.intersection(E_ind.columns) + if len(industries) == 0: + return DiagnosticResult( + name='compare_E_and_LCI_result: empty industry intersection', + passed=False, + tolerance=tolerance, + max_rel_diff=float('inf'), + failing_sectors=[], + details=None, + ) + E_aligned = E_ind.loc[flows, industries] + E_c = (C_m.loc[:, industries] @ E_aligned.T).T + commodities = ( + B.columns.intersection(E_c.columns) + .intersection(L.index) + .intersection(y.index) + ) + if len(commodities) == 0: + return DiagnosticResult( + name='compare_E_and_LCI_result: empty commodity intersection', + passed=False, + tolerance=tolerance, + max_rel_diff=float('inf'), + failing_sectors=[], + details=None, + ) + E_c = E_c.loc[flows, commodities] + + B_c = B.loc[flows, commodities] + L_a = L.loc[commodities, commodities] + y_a = y.reindex(commodities).fillna(0.0) + + LCI = compute_E_from_BLy(B=B_c, L=L_a, y=y_a) + LCI = LCI.reindex(index=E_c.index, columns=E_c.columns).fillna(0.0) + + # useeior: (LCI - E) / E → value=E_c, value_check=LCI + return validate_result( + 'compare_E_and_LCI_result', + _flatten_matrix_for_validate(E_c), + _flatten_matrix_for_validate(LCI), + tolerance=tolerance, + include_details=include_details, + )