CH4 emission accounting from Material Reporting

[Tyler-lc]

Currently the Material reporting has a bug when calculating CH4 emissions.
The file run_reporting.py loads the
CH4-emi.yaml. The method pyam_df_from_rep requests a standard set of indexes for all reporting keys.
The genno reporter is not filtering properly, because the requested keys do not all exist (commodity and level are not part of the key emi), which means that Genno will return ALL indexes.
The .join() applied to DF, will correctly apply the reporting variable to each technology.
But then the .groupby is where things break. Beacuse we are grouping by ["nl", "ya", "iamc_variable"], we are actually summing up ALL available emissions. That means that if a technology has in the parameter emission_fator :

• TCE
• PM2p5
• CO2_industry
• CH4_IND_ENE_COAL
• N2O_IND_ENE_COAL
• TCE_non-CO2
• TCE_other

All of these entries will be summed up. I think the reason is that we are not filtering by emission type. And there is no information anywhere for the reporter that provides that emission type.
Most likely, the solution would be to have not just the

var: emi

but to spell out the whole key:

var: emi:t-ya-yv-m-e  #or similar

Plus, we would need to specify in the yaml file, the actualy emission type we are looking for, so that we can filter for it properly.
For example, in the case of Industrial Processes|Iron and Steel

  Industrial Processes|Iron and Steel:
    filter:
      { technology: [DUMMY_limestone_supply_steel, prod_charcoal_steel, sinter_steel, pellet_steel, bf_steel, bof_steel,  eaf_steel, finishing_steel], mode: [M1, M2, M3, M4], emission[ADD EMISSIONS TYPE HERE]}
    short:
      emi_ch4_ind_proc_steel

This way, technologies that do not actually contain any CH4 (such as pellet_steel) will be filtered out. While those that do contain any of the emissions in the specified emission will be included, and only that type of emission will be kept in the DF.
I'm not entirely sure if this is the best way to go about it. There are several types of CH4 emissions to account for different things (leakages on the pipelines, industrial use, non-energy use, etc., pure CH4 emissions). I suppose that each of these emission types should end up in a different iamc variable anyway.
@macflo8 and I briefly discussed this issue already. Currently this may not be a big issue. I have found this to happen only for the CH4_emi.yaml. But it may affect all emission reporting done through the genno reporter that follow this structure.
I know @khaeru has definitely more experience with the Genno reporting. I don't know if you have any suggestion on the matter.

[khaeru]

The description of the intended behaviour and the failure mode is not entirely clear to me, but I can comment on the code that is pointed to. We have code like:

message-ix-models/message_ix_models/model/material/report/run_reporting.py

Lines 141 to 145 in f83d830

	def run_ch4_reporting(rep, model_name: str, scen_name: str) -> pyam.IamDataFrame:
	"""Generate reporting for industry methane emissions."""
	var = "ch4_emi"
	config = load_config(var)
	df = pyam_df_from_rep(rep, config.var, config.mapping)

that invoke this:

message-ix-models/message_ix_models/model/material/report/run_reporting.py

Lines 31 to 36 in f83d830

	filters_dict = {
	col: list(mapping_df.index.get_level_values(col).unique())
	for col in mapping_df.index.names
	}
	rep.set_filters(**filters_dict)
	df_var = pd.DataFrame(rep.get(f"{reporter_var}:nl-t-ya-m-c-l"))

In the present example, config.var and reporter_var are "emi".

If there are multiple functions like run_ch4_reporting() (hypothetically run_nox_reporting()), then Reporter.set_filters() and Reporter.get() are called multiple times. Within any one call of Reporter.get() intermediate results are stored and reused but they do not persist across multiple calls. (This is behaviour not from genno but from the popular dask package whose core get() mechanism we re-use.) So what would happen is:

1. Call run_ch4_reporting().
   1. Call pyam_df_from_rep().
      1. Call Reporter.set_filters().
      2. Call Reporter.get()
         1. Reporter retrieves emission_factor with the given filters.
         2. Reporter retrieves ACT with the given filters.
         3. emi = ACT × emission_factor
         4. The desired partial sum "emi:nl-t-ya-m-c-l" is computed and returned.
         5. Call to .get() finishes, so intermediate results like (a), (b), (c) are discarded.
      3. (Other operations)
2. Call run_nox_reporting().
   1. Call pyam_df_from_rep().
      1. Call Reporter.set_filters().
      2. Call Reporter.get()
         1. Reporter retrieves emission_factor with the given filters. ← repeats 1(i)(b)(a)
         2. Reporter retrieves ACT with the given filters. ← repeats 1(i)(b)(b)
         3. emi = ACT × emission_factor
         4. The desired partial sum "emi:nl-t-ya-m-c-l" is computed and returned.
      3. (Other operations)

This is why the preferred way to use the Reporter is to "push into the graph" all operations such that they happen with a single (or fewer) Reporter.get() call and we make best use of those dask/genno features. We would want this to happen:

1. Retrieve emission_factor with some (broad) filters.
2. Retrieve ACT with some (board) filters.
3. Compute emi = ACT × emission_factor (full dimensionality).
4. From (3), select a (narrower) subset, e.g. only those with e=="CH4".
5. Compute a partial sum of (4) like "emi:nl-t-ya-m-c-l".
6. Convert (5) to pd.DataFrame or pyam.IamDataFrame.
7. Repeat (4–6) as necessary, e.g. with different filters at step (4).
8. Concatenate results from (6).

If we set it up this way, then the intermediate results (1–3) are done only once, and reused. In adding (4) (multiple times, per (7)) to the Reporter, one could provide different set of indexers for each call.

I need to change to a different task, but perhaps we can continue in Slack where I can advise on how to rearrange the existing code in this way.

[Tyler-lc]

So what happens is:

ch4_emi.yaml defines:

var: emi
iamc_prefix: Emissions|CH4|
unit: Mt CH4/yr

vars:
  Industrial Processes|Iron and Steel:
    filter:
      { technology: [DUMMY_limestone_supply_steel, prod_charcoal_steel, sinter_steel,
                     pellet_steel, bf_steel, bof_steel, eaf_steel, finishing_steel],
        mode: [M1, M2, M3, M4] }

Then pyam_df_from_rep (run_reporting.py:36) requests:

df_var = pd.DataFrame(rep.get(f"{reporter_var}:nl-t-ya-m-c-l"))   # in this case this becomes "emi:nl-t-ya-m-c-l"

But in emi there are no "c" and "l" so genno uses as a fallback the entire dimensionality of emi, and it silently returns emi:nl-t-yv-ya-m-e-h with the emission species dimension e intact. For example:

nl	t	ya	m	e	...	value
R12_CHN	bf_steel	2020	M1	SO2		...
R12_CHN	bf_steel	2020	M1	NOx		...
R12_CHN	bf_steel	2020	M1	CH4_IND_ENE_COAL		...
R12_CHN	bf_steel	2020	M1	TCE		...
R12_CHN	bf_steel	2020	M1	PM2p5		...

The emission species were never filtered because the yaml schema only specifies the technology, mode and the iamc_variable name.

After this we call join, and here we would get:

nl	ya	t	e	iamc_name
R12_CHN	2020	bf_steel	SO2	Industrial Processes|Iron and Steel
R12_CHN	2020	bf_steel	NOx	Industrial Processes|Iron and Steel
R12_CHN	2020	bf_steel	CH4_IND_ENE_COAL	Industrial Processes|Iron and Steel
R12_CHN	2020	bf_steel	TCE	Industrial Processes|Iron and Steel

Basically, this should be a CH4 emission. But so far we still have SO2, NOx, TCE as well.

Now we call the gruopby

.groupby(["nl", "ya", "iamc_name"]).sum(numeric_only=True)

and this will sum over yv, h, m, t, and e. The result is that now, what we have labeled as Mt CH4/yr is the Σ of all emission species, because again there has been no per-emission filtering.
So we end up with something that looks like this:

emission species	contribution
SO2	2267.93
NOx	762.89
CH4_IND_ENE_COAL	32.37
N2O_IND_ENE_COAL	15.78
CO2	10.26
TCE / TCE_non-CO2 / TCE_other	4.52
PM2p5	0.69
CH4_IND_ENE_GAS_TD	0.62
CH4_IND_ENE_GAS	0.01
total	3095.07

I hope this is clearer now.