Data validation notebook and utilities - #45
Conversation
|
Check out this pull request on See visual diffs & provide feedback on Jupyter Notebooks. Powered by ReviewNB |
Codecov Report✅ All modified and coverable lines are covered by tests. Additional details and impacted files@@ Coverage Diff @@
## main #45 +/- ##
=======================================
Coverage 76.19% 76.19%
=======================================
Files 3 3
Lines 84 84
=======================================
Hits 64 64
Misses 20 20 ☔ View full report in Codecov by Harness. 🚀 New features to boost your workflow:
|
| # Weighted by cohort_population | ||
| pop_weight = xr.concat( | ||
| [ | ||
| socioeconomics["pop0to4"], |
There was a problem hiding this comment.
are these actual pop numbers or are they pop shares?
There was a problem hiding this comment.
They are actual pop numbers in the cohort for each region.
Total Deaths = da_cohort1 [death rate/100k] * pop_cohort1/100k + da_cohort2 [death rate/100k] * pop_cohort2/100k...
They become pop shares in the rate computation with age_weighted_total = age_weighted_total * 100000 / socioeconomics["pop"]
There was a problem hiding this comment.
if they are pop numbers, then the "weight" needs to be normalized by the total pop, I think?
There was a problem hiding this comment.
OH, I think I see now, this is working w/ mortality, not mortality rate
There was a problem hiding this comment.
Yeah, the general structure is to take everything to mortality for each age_cohort in age_weighted_total=da*pop_weight and then divide by total pop to take it back to an age-weighted rate if needed.
| # Mortality Rate = Total Age-Weighted Mortality/Total Pop (deaths/100k) | ||
| age_weighted_total = age_weighted_total * 100000 / socioeconomics["pop"] | ||
| # Sum across age cohorts for each Impact Region | ||
| regional_sum = age_weighted_total.sum(dim="age_cohort") |
There was a problem hiding this comment.
I can't tell, is this weighted by pop_share within the impact region if rate is True.
There was a problem hiding this comment.
looks like hotonly can be "hotonly", "coldonly", or else net
There was a problem hiding this comment.
Clunky change to add in "coldonly" functionality. Updated the docstring and added a validity check for the hotonly variable.
valid_hotonly = ("net", "hotonly", "coldonly") if hotonly not in valid_hotonly: raise ValueError(f"Invalid hotonly: {hotonly!r}. Must be one of {valid_hotonly}.")
| projected["/forecast_hotonly"]["effect"] | ||
| .chunk(chunks) | ||
| .groupby("time.month") | ||
| .mean() |
There was a problem hiding this comment.
what is this a mean over? Safer to define over what dim(s)
There was a problem hiding this comment.
oh right, this is a mean over "month" b/c of the groupby...I think you can ignore this comment.
But it seems these variables are all averaged over the months automatically?
| .chunk(chunks) | ||
| .mean(dim="number") | ||
| .groupby("time.month") | ||
| .mean() |
There was a problem hiding this comment.
same question as above, and also for "coldonly" and net below
| ) | ||
| return regional_sum | ||
|
|
||
| def compute_global_impact(impact, socioeconomics, rate=False): |
There was a problem hiding this comment.
does this work still if impact has age cohorts? this function assumes impact is a xr.DataArray ?
There was a problem hiding this comment.
Good call! I don't think I've ever actually preserved age-cohorts through all the processing steps. I've either done an age-weighting or selected a single cohort... This would be good to add!
| return regional_sum | ||
|
|
||
| def compute_global_impact(impact, socioeconomics, rate=False): | ||
| if rate: |
There was a problem hiding this comment.
this is a flag for whether impact contains a mortality rate or just mortality, right?
There was a problem hiding this comment.
Correct. This assumes that impact is already age_weighted, which it generally is for all of our current analysis but isn't guaranteed... I think a more robust way to handle this would be through attributes in the 'impact' xarray.
There was a problem hiding this comment.
best to set this to elif "net" and then have one more else that throws a ValueError for not having a correct option for hotonly.
|
|
||
| if "global_monthly" in output_scope: | ||
| global_monthly = dataset_to_dataframe(compute_stats(global_impact, dim=dims)) | ||
| stat_col_names = global_monthly.columns.difference(["region", "ISO"]) |
|
|
||
|
|
||
| ##### Plotting Functions ##### | ||
| from functools import lru_cache |
There was a problem hiding this comment.
might as well move these imports to the top of the script and clear these linting errors
analysis_utils.py: Functions to compute the mortality impact from a baseline year including age-weighting, population-weighted totals, and automatically computing and generating CSVs of statistics across ensemble members and gamma samplesc00_validate_data.ipynb: Simple notebook to plot impacts, print descriptive statistics for validation, and produce CSVs of quantile stats for impact regions and global scale.