Skip to content

Stabilize modeled-climatology change factors and fix CONUS CSV/JSON divergence #739

New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

Sign up for GitHub

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Merged
cstephen merged 3 commits into from
Aug 11, 2026
Merged

Stabilize modeled-climatology change factors and fix CONUS CSV/JSON divergence #739

Show file tree
Hide file tree
Changes from all commits
Commits
Show all changes
3 commits
Select commit Hold shift + click to select a range
ee73666
Apply Maurer adjustment before CSV export in CONUS modeled climatology
brucecrevensten Jul 31, 2026
d323457
Stabilize GCM change factors against near-zero historical baselines
brucecrevensten Jul 31, 2026
1d7167a
Enforce min <= mean <= max in adjusted climatology statistics
brucecrevensten Jul 31, 2026
File filter

Filter by extension

Filter by extension
Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
45 changes: 45 additions & 0 deletions postprocessing.py
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -267,3 +267,48 @@ def merge_dicts(dict1, dict2):
else:
merged[key] = value
return merged


# Stabilization constants for GCM change factors ("deltas") that are applied to
# historical baselines as ratios. The epsilon fraction sets the offset relative
# to the stream's historical mean flow; the floor keeps the offset nonzero for
# all-zero streams; the cap bounds how much any single statistic can be scaled.
RATIO_EPSILON_FRACTION = 0.01
RATIO_EPSILON_FLOOR = 0.0001
RATIO_CAP = 10.0
Comment on lines +277 to +278

@charparr charparr Aug 11, 2026

Copy link
Copy Markdown
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

One immediate question I have is how frequently are the floor / cap values triggered? If these values are returned frequently then they might need adjustment. It seems like they are intended to catch outliers which I think is the right approach but it could be worth an empirical test before all is said and done. I do see that Claude called these out already as "judgement calls" that Claude is "happy to tune".



def scale_aware_epsilon(values):
"""
Compute a stabilizing offset for ratio-based change factors, sized to the
magnitude of the historical data so that the same code works for large
rivers and small creeks alike.
Args:
values (list of float): Historical baseline values (e.g. all doy_mean
values for one model), used to establish the stream's flow scale
Returns:
float: RATIO_EPSILON_FRACTION times the mean of values, but never less
than RATIO_EPSILON_FLOOR
"""
if not values:
return RATIO_EPSILON_FLOOR
epsilon = RATIO_EPSILON_FRACTION * (sum(values) / len(values))

@charparr charparr Aug 11, 2026

Copy link
Copy Markdown
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

One choice I see embedded here is that the arithmetic mean is the value used to scale the epsilon. This makes sense to me because it is easy to explain and to conceptualize: big historical flows get bigger epsilons, and smaller historical flows get smaller epsilons.

But, are there any streamflow regimes for which this might have unanticipated consequences? Imagine a stream where the flow is: [small, small, small, BIG BIG GUSH O WATER, small, small, ...] : the mean may not scale the epsilon value in a desirable way because the big flow might get only get scaled a tiny bit.

Down the road it may be worthwhile exploring the use of the median the scale factor, or some other statistic.

And/or just adding some more text along the lines of
"Hey for a stream that is highly variable or seasonal, a spike in flow is going to raise the mean, which will raise the epsilon for the entire year's worth of flows.

return max(epsilon, RATIO_EPSILON_FLOOR)


def stabilized_ratio(projected, historical, epsilon):
"""
Ratio of projected to historical, stabilized for near-zero baselines.
Adding the same epsilon to numerator and denominator makes the ratio
approach 1 (no change) as both values approach zero, instead of exploding
when only the denominator is small. The result is clamped to
[1/RATIO_CAP, RATIO_CAP] as a backstop against noisy statistics.
Args:
projected (float): Future (projected) statistic value
historical (float): Historical statistic value for the same model
epsilon (float): Positive stabilizing offset from scale_aware_epsilon()
Returns:
float: Clamped change factor suitable for scaling a historical baseline
"""
ratio = (projected + epsilon) / (historical + epsilon)
return min(max(ratio, 1.0 / RATIO_CAP), RATIO_CAP)
34 changes: 29 additions & 5 deletions routes/arctic_hydrology.py
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -23,7 +23,11 @@
)
from fetch_data import fetch_data, fetch_layer_data, describe_via_wcps
from validate_request import get_axis_encodings
from postprocessing import prune_nulls_with_max_intensity
from postprocessing import (
prune_nulls_with_max_intensity,
scale_aware_epsilon,
stabilized_ratio,
)
from csv_functions import create_csv
from config import RAS_BASE_URL
from . import routes
Expand Down Expand Up @@ -490,6 +494,11 @@ def calculate_and_apply_gcm_diffs_to_blaskey_climatology(data_dict):
"""
Function to calculate GCM-projected changes in streamflow and apply them to the historical Blaskey climatology.
Models without a '1990-2021' era (e.g. PGW models with no historical baseline) are silently skipped.

Change factors are stabilized ratios (see postprocessing.stabilized_ratio):
a symmetric offset scaled to the model's historical mean flow keeps the
factor near 1 when both values are near zero (e.g. frozen winter minimums),
and factors are clamped to [1/RATIO_CAP, RATIO_CAP].
Args:
data_dict (dict): Climatology data dict keyed by model then era
Returns:
Expand All @@ -502,6 +511,13 @@ def calculate_and_apply_gcm_diffs_to_blaskey_climatology(data_dict):
continue
if "1990-2021" not in data_dict[model]:
continue
epsilon = scale_aware_epsilon(
[
row["doy_mean"]
for row in data_dict[model]["1990-2021"]
if "doy_mean" in row
]
)
adjusted_data_dict[model] = {}
for era in data_dict[model].keys():
if era == "1990-2021":
Expand All @@ -519,12 +535,20 @@ def calculate_and_apply_gcm_diffs_to_blaskey_climatology(data_dict):
blaskey_historical = data_dict["historical"]["1990-2021"][i][stat]
gcm_historical = data_dict[model]["1990-2021"][i][stat]
gcm_projected = entry[stat]
denominator = gcm_historical
if denominator == 0:
denominator = 0.0001
projected_quotient = gcm_projected / denominator
projected_quotient = stabilized_ratio(
gcm_projected, gcm_historical, epsilon
)
blaskey_adjusted = round(blaskey_historical * projected_quotient, 3)
doy_stats[stat] = blaskey_adjusted
# each stat is scaled by its own ratio, so the adjusted values
# can cross; clamp to preserve min <= mean <= max
if all(k in doy_stats for k in ("doy_min", "doy_mean", "doy_max")):

@charparr charparr Aug 11, 2026

Copy link
Copy Markdown
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Good comment on the "why" -- feels like we'll come back to this

doy_stats["doy_min"] = min(
doy_stats["doy_min"], doy_stats["doy_mean"]
)
doy_stats["doy_max"] = max(
doy_stats["doy_max"], doy_stats["doy_mean"]
)
adjusted_data_dict[model][era].append(doy_stats)
return adjusted_data_dict

Expand Down
56 changes: 41 additions & 15 deletions routes/conus_hydrology.py
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -25,7 +25,11 @@
)
from fetch_data import fetch_data, fetch_layer_data, describe_via_wcps
from validate_request import get_axis_encodings
from postprocessing import prune_nulls_with_max_intensity
from postprocessing import (
prune_nulls_with_max_intensity,
scale_aware_epsilon,
stabilized_ratio,
)
from csv_functions import create_csv
from config import RAS_BASE_URL
from . import routes
Expand Down Expand Up @@ -784,6 +788,11 @@ def calculate_and_apply_gcm_diffs_to_maurer_climatology(data_dict):
Function to calculate the GCM-projected changes in streamflow stats and apply those changes to the historical Maurer climatology stats.
This is done by first calculating the ratio between the GCM-projected future stat values and the GCM historical stat values,
then applying that scaling factor to the Maurer historical stat values.

Change factors are stabilized ratios (see postprocessing.stabilized_ratio):
a symmetric offset scaled to the model's historical mean flow keeps the
factor near 1 when both values are near zero (e.g. dry-season minimums),
and factors are clamped to [1/RATIO_CAP, RATIO_CAP].
Args:
data_dict (dict): Data dictionary with the hydrology data populated
Returns:
Expand All @@ -795,6 +804,13 @@ def calculate_and_apply_gcm_diffs_to_maurer_climatology(data_dict):
if model == "Maurer":
adjusted_data_dict[model] = data_dict[model]
continue
epsilon = scale_aware_epsilon(
[
row["doy_mean"]
for row in data_dict[model]["historical"]["1976-2005"]
if "doy_mean" in row
]
)
if model not in adjusted_data_dict:
adjusted_data_dict[model] = {}
for scenario in data_dict[model].keys():
Expand All @@ -819,14 +835,24 @@ def calculate_and_apply_gcm_diffs_to_maurer_climatology(data_dict):
stat
]
gcm_projected = entry[stat]
denominator = gcm_historical
if denominator == 0:
denominator = 0.0001
projected_quotient = gcm_projected / denominator
projected_quotient = stabilized_ratio(
gcm_projected, gcm_historical, epsilon
)
maurer_adjusted = round(
maurer_historical * projected_quotient, 3
)
doy_stats[stat] = maurer_adjusted
# each stat is scaled by its own ratio, so the adjusted
# values can cross; clamp to preserve min <= mean <= max
if all(
k in doy_stats for k in ("doy_min", "doy_mean", "doy_max")
):
doy_stats["doy_min"] = min(
doy_stats["doy_min"], doy_stats["doy_mean"]
)
doy_stats["doy_max"] = max(
doy_stats["doy_max"], doy_stats["doy_mean"]
)
Comment on lines +845 to +855

@charparr charparr Aug 11, 2026

Copy link
Copy Markdown
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

I think this is identical to R543 to R551 in routes/arctic_hydrology.py

And thematically it feels like postprocessing. Consider consolidating and moving to the post-processing module. If nothing else the route is already 1000+ LOC so it could be nice. However, this suggestion is not functional so just take it under consideration.

adjusted_data_dict[model][scenario][era].append(doy_stats)
return adjusted_data_dict

Expand Down Expand Up @@ -962,6 +988,16 @@ def run_get_conus_hydrology_modeled_climatology(stream_id):
data_dict = populate_feature_name_and_location_attributes(data_dict, gdf)
data_dict = prune_nulls_with_max_intensity(data_dict)

# apply GCM-projected changes to Maurer climatology stats if source is "gcm_diff_applied_to_maurer"
# otherwise, if source is "original_gcm", then we are just returning the original GCM stats with no adjustments
if source == "gcm_diff_applied_to_maurer":
for landcover in data_dict["data"]:
data_dict["data"][landcover] = (
calculate_and_apply_gcm_diffs_to_maurer_climatology(
data_dict["data"][landcover]
)
)

if request.args.get("format") == "csv":
try:
return create_csv(
Expand All @@ -976,16 +1012,6 @@ def run_get_conus_hydrology_modeled_climatology(stream_id):
except Exception as exc:
return render_template("500/server_error.html"), 500

# apply GCM-projected changes to Maurer climatology stats if source is "gcm_diff_applied_to_maurer"
# otherwise, if source is "original_gcm", then we are just returning the original GCM stats with no adjustments
if source == "gcm_diff_applied_to_maurer":
for landcover in data_dict["data"]:
data_dict["data"][landcover] = (
calculate_and_apply_gcm_diffs_to_maurer_climatology(
data_dict["data"][landcover]
)
)

return jsonify(data_dict)

except Exception as exc:
Expand Down
Loading
Loading