diff --git a/bedrock/transform/eeio/__tests__/test_useeio_b_steps.py b/bedrock/transform/eeio/__tests__/test_useeio_b_steps.py index ba7d391da..13aa1d797 100644 --- a/bedrock/transform/eeio/__tests__/test_useeio_b_steps.py +++ b/bedrock/transform/eeio/__tests__/test_useeio_b_steps.py @@ -60,6 +60,12 @@ def test_useeio_b_adjust_divide_transform_steps( lambda original_year, target_year: ratio, ) monkeypatch.setattr(dc, "derive_Vnorm_scrap_corrected", lambda: vnorm) + # Gross output includes scrap, so this branch divides by 1 - s (#1293). + monkeypatch.setattr( + dc, + "scrap_share_of_output", + lambda: pd.Series([0.1, 0.0], index=["I1", "I2"]), + ) try: out = dc.derive_B_via_vnorm() @@ -69,10 +75,10 @@ def test_useeio_b_adjust_divide_transform_steps( # Step 1: adjusted x = x_nominal / ratio = [125, 160] # Step 2: Bi = E / adjusted x = [0.08, 0.125] - # Step 3: B = Bi @ I = Bi + # Step 3: B = Bi @ (I / (1 - s)), with s = [0.1, 0] expected = pd.DataFrame( [ - [0.08, 0.125], # CO2 + [0.08 / 0.9, 0.125], # CO2 [0.0, 0.0], # CH4 [0.0, 0.0], # N2O [0.0, 0.0], # HFCs diff --git a/bedrock/transform/eeio/usa_matrices.py b/bedrock/transform/eeio/usa_matrices.py index 61ed7a9c0..487c9a66b 100644 --- a/bedrock/transform/eeio/usa_matrices.py +++ b/bedrock/transform/eeio/usa_matrices.py @@ -308,12 +308,53 @@ def derive_q() -> pd.Series[float]: return q +def scrap_share_of_output() -> pd.Series[float]: + """Each industry's scrap (``S00401``) over its output **including** scrap. + + That is USEEIO's ``s``. The Cornerstone Make has no scrap, so + ``derive_Vnorm_scrap_corrected`` needs no adjustment where ``x`` is its row + sums. Where ``x`` is BEA gross output, which includes scrap, B divides the + market shares by ``1 - s`` (``derive_B_via_vnorm``, #1293). + + A BEA industry split into several Cornerstone industries (the waste + children of ``562000``, the electricity children) shares its scrap across + them **by output**, not by the correspondence's equal weights. So every + child gets its parent's scrap share, as in USEEIO's waste disaggregation. + """ + V = derive_V() + x = compute_x(V=V) + corresp = industry_corresp() + by_output = corresp.gt(0).mul(x.reindex(corresp.index).fillna(0.0), axis=0) + totals = by_output.sum(axis=0) + # A BEA industry whose children have no output keeps the correspondence. + weights = by_output.div(totals.where(totals > 0), axis=1).fillna(corresp) + scrap = weights @ load_detail_V_usa().loc[:, 'S00401'] + if get_usa_config().implement_electricity_disaggregation: + parent = float(scrap.get(ELECTRICITY_AGGREGATE_SECTOR, 0.0)) + scrap = scrap.drop(ELECTRICITY_AGGREGATE_SECTOR, errors='ignore') + children = [c for c in ELECTRICITY_DISAGG_SECTORS if c in x.index] + child_x = float(x.reindex(children).sum()) + for code in children: + scrap.loc[code] = parent * float(x[code]) / child_x if child_x else 0.0 + scrap = scrap.reindex(V.index, fill_value=0.0) + return (scrap / (x.reindex(V.index, fill_value=0.0) + scrap)).fillna(0.0) + + @functools.cache def derive_Vnorm_scrap_corrected( apply_inflation: bool | None = None, target_year: int = 0, ) -> pd.DataFrame: - """Scrap-corrected V norm. Inflation is applied via ``derive_V``. + """Market shares ``V q̂⁻¹``, net of scrap. Inflation is applied via ``derive_V``. + + USEEIO divides market shares by ``1 - s``, with ``s`` scrap's share of an + industry's output **including** scrap, so that the industry's inputs and + emissions land on its other products. The Cornerstone Make table has no + scrap: ``S00401`` maps with weight 0, so ``x`` is already output without + scrap, and ``U / x`` and ``E / x`` already land on the other products. + Dividing by ``1 - s`` again would count scrap twice (#1293). The plain + market shares give exactly USEEIO's A and B: + ``U x_incl⁻¹ V_n / (1 - scrap / x_incl) = U x_excl⁻¹ V_n``. When ``target_year`` is not positive and inflation is on, uses ``USAConfig.model_base_year``. @@ -327,27 +368,8 @@ def derive_Vnorm_scrap_corrected( V = derive_V(use_inflation, effective_year) q = compute_q(V=V) - x = compute_x(V=V) - - Vnorm = compute_Vnorm_matrix(V=V, q=q) - scrap_detail = load_detail_V_usa().loc[:, 'S00401'] - scrap_fraction = industry_corresp() @ scrap_detail - if get_usa_config().implement_electricity_disaggregation: - parent_scrap = float(scrap_fraction.get(ELECTRICITY_AGGREGATE_SECTOR, 0.0)) - scrap_fraction = scrap_fraction.drop( - ELECTRICITY_AGGREGATE_SECTOR, errors='ignore' - ) - for code in ELECTRICITY_DISAGG_SECTORS: - scrap_fraction.loc[code] = parent_scrap - # The delivery sector has no scrap: the reindex below gives it zero. - scrap_fraction = scrap_fraction.reindex(V.index, fill_value=0.0) - x_aligned = x.reindex(V.index, fill_value=0.0) - - V_scrap_corrected = Vnorm.divide( - (1.0 - (scrap_fraction / x_aligned).fillna(0.0)), axis=0 - ) - V_scrap_corrected = V_scrap_corrected.reindex( + V_scrap_corrected = compute_Vnorm_matrix(V=V, q=q).reindex( index=V.index, columns=V.columns, fill_value=0.0 ) validate_cornerstone(V_scrap_corrected, 'V') @@ -834,6 +856,12 @@ def derive_B_via_vnorm() -> pd.DataFrame: else: x = derive_x_after_redefinition() if cfg.use_ghg_year_x_in_B else derive_x() Vnorm = derive_Vnorm_scrap_corrected() + if cfg.usa_detail_io_source != 'nowcast' and ( + cfg.deflate_x_to_detail_io_year_for_B or cfg.use_ghg_year_x_in_B + ): + # BEA gross output includes scrap, so the scrap correction applies here, + # as in USEEIO: divide by 1 - s with s over output including scrap (#1293). + Vnorm = Vnorm.divide(1.0 - scrap_share_of_output(), axis=0) Bi = E.divide(x, axis=1).fillna(0.0) B = Bi @ Vnorm validate_cornerstone(B, 'B') diff --git a/bedrock/utils/economic/inflation_helpers_usa.py b/bedrock/utils/economic/inflation_helpers_usa.py index 81d8d096b..8b234ca69 100644 --- a/bedrock/utils/economic/inflation_helpers_usa.py +++ b/bedrock/utils/economic/inflation_helpers_usa.py @@ -423,12 +423,10 @@ def get_vnorm_adjusted_commodity_price_ratio( aligned = industry_ratio.reindex(Vnorm.index, fill_value=1.0) # Normalize V_norm columns to sum to 1 so the dot-product is a true weighted - # *average*, not a weighted *sum*. Scrap correction in - # `derive_Vnorm_scrap_corrected` is a row-axis scaling that - # leaves column sums slightly > 1 — using the un-renormalized form would - # produce a non-1 ratio at year=year (off by the column-sum drift, ~5–7%). - # Done locally here rather than upstream because A/B matrix consumers want - # the un-renormalized scrap-corrected V_norm. + # *average*, not a weighted *sum*. `derive_Vnorm_scrap_corrected` returns the + # plain market shares (#1293), whose columns already sum to 1, so this is a + # guard: it keeps the ratio at 1 for year = year if the shares ever change. + # B's gross-output path scales them by 1 / (1 - s) separately. column_sums = Vnorm.sum(axis=0) weights = Vnorm.divide(column_sums.where(column_sums > 1e-9, 1.0), axis=1) commodity_ratio = aligned @ weights diff --git a/bedrock/utils/validation/__tests__/test_eeio_diagnostics.py b/bedrock/utils/validation/__tests__/test_eeio_diagnostics.py index 0a45aebfd..7b52210a2 100644 --- a/bedrock/utils/validation/__tests__/test_eeio_diagnostics.py +++ b/bedrock/utils/validation/__tests__/test_eeio_diagnostics.py @@ -417,8 +417,8 @@ def test_compare_output_and_L_y( @pytest.mark.xfail( reason=( 'Make and Use disagree on commodity output. On MUT v1.0.0_8365fa4 ' - '(2024), with the F09000 sales in y: 334118 (+4.4%, imports above ' - 'use) and S00402 ($4M on q = 0).' + '(2024) one sector fails: 334118 (+4.4%), whose imports exceed its ' + 'use. 2017-2023 pass.' ), ) def test_domestic_L_y_from_use_final_demand_equals_q( @@ -569,6 +569,21 @@ def test_domestic_lci_equals_e(canonical_usa_config: None) -> None: ) +@pytest.mark.eeio_integration +def test_b_times_q_conserves_emissions(canonical_usa_config: None) -> None: + """``B @ q`` sums to E: commodity factors neither add nor lose emissions. + + It failed by 4.0 Mt when the scrap correction divided market shares by + ``1 - s`` on an ``x`` that already excluded scrap (#1293). The LCI≈E test + above can't see that, because it puts the same market shares on both sides. + """ + B = derive_B() + q = derive_Aq_scaled().scaled_q + assert float((B.reindex(columns=q.index, fill_value=0.0) @ q).sum()) == ( + pytest.approx(float(derive_E_usa().values.sum()), rel=1e-9) + ) + + @pytest.mark.eeio_integration @pytest.mark.xfail( reason=( @@ -600,3 +615,20 @@ def test_total_lci_equals_e(canonical_usa_config: None) -> None: f'total LCI≈E failed ({len(result.failing_sectors)} cells): ' f'{result.failing_sectors[:20]}' ) + + +def test_output_vs_leontief_skips_the_accounting_sectors() -> None: + """S00402 has no output, so q = L @ y isn't checked for it (#1270).""" + sectors = ['1111A0', 'S00402'] + L = pd.DataFrame(np.eye(2), index=sectors, columns=sectors) + q = pd.Series([100.0, 0.0], index=sectors) + + off_on_used_goods = compare_output_vs_leontief_x_demand( + output=q, L=L, y=pd.Series([100.0, -6.0e6], index=sectors) + ) + assert off_on_used_goods.passed + + off_on_oilseeds = compare_output_vs_leontief_x_demand( + output=q, L=L, y=pd.Series([103.0, 0.0], index=sectors) + ) + assert not off_on_oilseeds.passed diff --git a/bedrock/utils/validation/eeio_diagnostics.py b/bedrock/utils/validation/eeio_diagnostics.py index de24772d9..3d0c56639 100644 --- a/bedrock/utils/validation/eeio_diagnostics.py +++ b/bedrock/utils/validation/eeio_diagnostics.py @@ -28,6 +28,12 @@ # (e.g. S00402 used goods) compare |diff| against atol instead of (diff / q). _VALIDATE_RESULT_ATOL = 1e-4 +# BEA's accounting rows: scrap, used and secondhand goods, noncomparable +# imports and the rest-of-world adjustment. They have no output of their own +# (S00402's q is 0), so q = L @ y can't hold for them and isn't checked, as in +# USEEIO 2. Of the four, only S00402 is a Cornerstone commodity today (#1270). +ACCOUNTING_SECTORS: frozenset[str] = frozenset({'S00300', 'S00401', 'S00402', 'S00900'}) + @dc.dataclass class DiagnosticResult: @@ -347,7 +353,7 @@ def compare_output_vs_leontief_x_demand( Compares the total sector output (commodity or industry) against the model result calculation of L @ y. Pass/fail: ``|diff| <= tolerance * |value| + atol`` for all sectors (see - ``validate_result``). + ``validate_result``) except the ``ACCOUNTING_SECTORS``, which aren't checked. Parameters ---------- @@ -387,8 +393,14 @@ def compare_output_vs_leontief_x_demand( output_check = backcompute_q_from_L_and_y(L=L, y=y) name = "compare output and L * y" + # L @ y uses every sector; only the comparison leaves the accounting rows out. + checked = output.index.difference(list(ACCOUNTING_SECTORS), sort=False) d_result = validate_result( - name, output, output_check, tolerance=tolerance, include_details=include_details + name, + output.loc[checked], + output_check.loc[checked], + tolerance=tolerance, + include_details=include_details, ) return d_result