Merge pull request #95 from dgergel/add_analogdownscaling_qdm_options

Add analogdownscaling qdm options

Author: brews

HISTORY.rst

@@ -5,6 +5,7 @@ History
 
 0.X.X (XXXX-XX-XX)
 ------------------
+* Add ``include-quantiles`` flag to ``apply_qdm`` to allow for including quantile information in bias corrected output. (PR #95, @dgergel)
 * Add ``astype`` argument to ``regrid``. (PR #92, @brews)
 * Make ``dodola`` container's default CMD. (PR #90, @brews)
 * Add ``dask-kubernetes``, ``distributed`` to container environment. (PR #90, @brews)

dodola/cli.py

@@ -49,10 +49,20 @@ def dodola_cli(debug):
     required=True,
     help="URL to write NetCDF4 with adjusted simulation year to",
 )
-def apply_qdm(simulation, qdm, year, variable, out):
+@click.option(
+    "--include-quantiles",
+    is_flag=True,
+    help="Include simulation quantiles in output",
+)
+def apply_qdm(simulation, qdm, year, variable, out, include_quantiles):
     """Adjust simulation year with QDM bias correction method, outputting to local NetCDF4 file"""
     services.apply_qdm(
-        simulation=simulation, qdm=qdm, year=year, variable=variable, out=out
+        simulation=simulation,
+        qdm=qdm,
+        year=year,
+        variable=variable,
+        out=out,
+        include_quantiles=include_quantiles,
     )
 
 

dodola/core.py

@@ -8,7 +8,7 @@
 import logging
 from skdownscale.spatial_models import SpatialDisaggregator
 import xarray as xr
-from xclim import sdba
+from xclim import sdba, set_options
 from xclim.sdba.utils import equally_spaced_nodes
 from xclim.core.calendar import convert_calendar
 import xesmf as xe
@@ -53,7 +53,7 @@ def train_quantiledeltamapping(
 
 
 def adjust_quantiledeltamapping_year(
-    simulation, qdm, year, variable, halfyearwindow_n=10
+    simulation, qdm, year, variable, halfyearwindow_n=10, include_quantiles=False
 ):
     """Apply QDM to adjust a year within a simulation.
 
@@ -74,6 +74,9 @@ def adjust_quantiledeltamapping_year(
     halfyearwindow_n : int, optional
         Half-length of the annual rolling window to extract along either
         side of `year`.
+    include_quantiles : bool, optional
+        Whether or not to output quantiles (sim_q) as a coordinate on
+        the bias corrected data variable in output.
 
     Returns
     -------
@@ -96,7 +99,14 @@ def adjust_quantiledeltamapping_year(
     simulation = simulation[variable].sel(
         time=timeslice
     )  # TODO: Need a check to ensure we have all the data in this slice!
-    out = qdm.adjust(simulation, interp="nearest").sel(time=str(year))
+    if include_quantiles:
+        # include quantile information in output
+        with set_options(sdba_extra_output=True):
+            out = qdm.adjust(simulation, interp="nearest").sel(time=str(year))
+            # make quantiles a coordinate of bias corrected output variable
+            out = out["scen"].assign_coords(sim_q=out.sim_q)
+    else:
+        out = qdm.adjust(simulation, interp="nearest").sel(time=str(year))
 
     return out.to_dataset(name=variable)
 

dodola/services.py

@@ -64,7 +64,7 @@ def train_qdm(historical, reference, out, variable, kind):
 
 
 @log_service
-def apply_qdm(simulation, qdm, year, variable, out):
+def apply_qdm(simulation, qdm, year, variable, out, include_quantiles=False):
     """Apply trained QDM to adjust a year within a simulation, dump to NetCDF.
 
     Dumping to NetCDF is a feature likely to change in the near future.
@@ -84,6 +84,9 @@ def apply_qdm(simulation, qdm, year, variable, out):
     out : str
         fsspec-compatible path or URL pointing to NetCDF4 file where the
         QDM-adjusted simulation data will be written.
+    include_quantiles : bool
+        Flag to indicate whether bias-corrected quantiles should be
+        included in the QDM-adjusted output.
     """
     sim_ds = storage.read(simulation)
     qdm_ds = storage.read(qdm)
@@ -92,7 +95,11 @@ def apply_qdm(simulation, qdm, year, variable, out):
     variable = str(variable)
 
     adjusted_ds = adjust_quantiledeltamapping_year(
-        simulation=sim_ds, qdm=qdm_ds, year=year, variable=variable
+        simulation=sim_ds,
+        qdm=qdm_ds,
+        year=year,
+        variable=variable,
+        include_quantiles=include_quantiles,
     )
 
     # Write to NetCDF, usually on local disk, pooling and "fanning-in" NetCDFs is

dodola/tests/test_core.py

@@ -115,6 +115,36 @@ def test_adjust_quantiledeltamapping_year_kind(variable_kind, expected):
     assert all(adjusted_ds[target_variable] == expected)
 
 
+def test_adjust_quantiledeltamapping_include_quantiles():
+    """Test that include-quantiles flag results in bias corrected quantiles
+    included in output"""
+    target_variable = "fakevariable"
+    n_simdays = 5 * 365  # 100 years of daily simulation.
+
+    model_bias = 2.0
+    ts_sim = np.ones(n_simdays, dtype=np.float64)
+    sim = _timeseriesfactory(
+        ts_sim * model_bias, start_dt="2015-01-01", variable_name=target_variable
+    )
+
+    target_year = 2017
+
+    # Yes, I'm intentionally training the QDM to a different bias. This is to
+    # spurn a difference between "kind" adjustments...
+    qdm = _train_simple_qdm(
+        target_variable="fakevariable", kind="+", additive_bias=model_bias + 1
+    )
+    adjusted_ds = adjust_quantiledeltamapping_year(
+        simulation=sim,
+        qdm=qdm,
+        year=target_year,
+        variable=target_variable,
+        include_quantiles=True,
+    )
+    # check that quantiles are contained in output
+    assert "sim_q" in adjusted_ds[target_variable].coords
+
+
 def test_adjust_quantiledeltamapping_year_output_time():
     """Check 'time' year and edges of QDM adjusted output
 

environment.yaml

@@ -11,6 +11,7 @@ dependencies:
 - fsspec=2021.5.0  # Prevent azure blob errors, not hard req.
 - gcsfs=2021.5.0
 - numpy=1.20.3
+- pandas=1.2.5
 - pip=21.1.2
 - pytest=6.2.4
 - python=3.9
@@ -21,4 +22,4 @@ dependencies:
 - zarr=2.8.3
 - pip:
   - git+https://github.com/dgergel/xsd@458f292dab6e8a6584659e97a66c37e028da2b7a
-  - git+https://github.com/Ouranosinc/xclim@d1e654cf59a5d5761f5586bca1102a1590ae6e32
+  - xclim