Merge pull request #81 from joemoorhouse/main

Expose coastal inundation for real estate model

src/physrisk/data/event_provider.py

@@ -64,6 +64,22 @@ def _wri_inundation_prefix():
     return "inundation/wri/v2"
 
 
+_percentiles_map = {"95": "0", "5": "0_perc_05", "50": "0_perc_50"}
+_subsidence_set = {"wtsub", "nosub"}
+
+
+def get_source_path_wri_coastal_inundation(*, model: str, scenario: str, year: int):
+    type = "coast"
+    # model is expected to be of the form subsidence/percentile, e.g. wtsub/95
+    # if percentile is omitted then 95th percentile is used
+    model_components = model.split("/")
+    sub = model_components[0]
+    if sub not in _subsidence_set:
+        raise ValueError("expected model input of the form {subsidence/percentile}, e.g. wtsub/95, nosub/5, wtsub/50")
+    perc = "95" if len(model_components) == 1 else model_components[1]
+    return os.path.join(_wri_inundation_prefix(), f"inun{type}_{scenario}_{sub}_{year}_{_percentiles_map[perc]}")
+
+
 def get_source_path_wri_riverine_inundation(*, model: str, scenario: str, year: int):
     type = "river"
     return os.path.join(_wri_inundation_prefix(), f"inun{type}_{scenario}_{model}_{year}")

src/physrisk/data/inventory.py

@@ -117,48 +117,41 @@ def __init__(self):
                 "id": "nosub",
                 "display_name": "Baseline no subsidence",
                 "description": "Baseline condition; no subsidence",
-                "filename": "",
-                "scenarios": [
-                    {"id": "historical", "years": [1980]},
-                    {"id": "rcp4p5", "years": [2030, 2050, 2080]},
-                    {"id": "rcp8p5", "years": [2030, 2050, 2080]},
-                ],
+                "filename": "inuncoast_{scenario}_nosub_{year}_0",
+                "scenarios": [{"id": "historical", "years": [1980]}],
             },
             {
                 "event_type": "CoastalInundation",
                 "path": "coastal_inundation/wri/v2",
-                "id": "nosub_95",
+                "id": "nosub/95",
                 "display_name": "95% no subsidence",
-                "description": "No subsidence; 95th percentile sea rise",
-                "filename": "",
+                "description": "No subsidence; 95th percentile sea level rise",
+                "filename": "inuncoast_{scenario}_nosub_{year}_0",
                 "scenarios": [
-                    {"id": "historical", "years": [1980]},
                     {"id": "rcp4p5", "years": [2030, 2050, 2080]},
                     {"id": "rcp8p5", "years": [2030, 2050, 2080]},
                 ],
             },
             {
                 "event_type": "CoastalInundation",
                 "path": "coastal_inundation/wri/v2",
-                "id": "nosub_05",
+                "id": "nosub/5",
                 "display_name": "5% no subsidence",
-                "description": "No subsidence; 5th percentile sea rise",
-                "filename": "",
+                "description": "No subsidence; 5th percentile sea level rise",
+                "filename": "inuncoast_{scenario}_nosub_{year}_0_perc_05",
                 "scenarios": [
-                    {"id": "historical", "years": [1980]},
                     {"id": "rcp4p5", "years": [2030, 2050, 2080]},
                     {"id": "rcp8p5", "years": [2030, 2050, 2080]},
                 ],
             },
             {
                 "event_type": "CoastalInundation",
                 "path": "coastal_inundation/wri/v2",
-                "id": "nosub_50",
+                "id": "nosub/50",
                 "display_name": "50% no subsidence",
-                "description": "No subsidence; 50th percentile sea rise",
-                "filename": "",
+                "description": "No subsidence; 50th percentile sea level rise",
+                "filename": "inuncoast_{scenario}_nosub_{year}_0_perc_50",
                 "scenarios": [
-                    {"id": "historical", "years": [1980]},
                     {"id": "rcp4p5", "years": [2030, 2050, 2080]},
                     {"id": "rcp8p5", "years": [2030, 2050, 2080]},
                 ],
@@ -169,48 +162,41 @@ def __init__(self):
                 "id": "wtsub",
                 "display_name": "Baseline with subsidence",
                 "description": "Baseline condition; with subsidence",
-                "filename": "",
-                "scenarios": [
-                    {"id": "historical", "years": [1980]},
-                    {"id": "rcp4p5", "years": [2030, 2050, 2080]},
-                    {"id": "rcp8p5", "years": [2030, 2050, 2080]},
-                ],
+                "filename": "inuncoast_{scenario}_wtsub_{year}_0",
+                "scenarios": [{"id": "historical", "years": [1980]}],
             },
             {
                 "event_type": "CoastalInundation",
                 "path": "coastal_inundation/wri/v2",
-                "id": "wtsub_95",
+                "id": "wtsub/95",
                 "display_name": "95% with subsidence",
-                "description": "With subsidence; 95th percentile sea rise",
-                "filename": "",
+                "description": "With subsidence; 95th percentile sea level rise",
+                "filename": "inuncoast_{scenario}_wtsub_{year}_0",
                 "scenarios": [
-                    {"id": "historical", "years": [1980]},
                     {"id": "rcp4p5", "years": [2030, 2050, 2080]},
                     {"id": "rcp8p5", "years": [2030, 2050, 2080]},
                 ],
             },
             {
                 "event_type": "CoastalInundation",
                 "path": "coastal_inundation/wri/v2",
-                "id": "wtsub_05",
+                "id": "wtsub/5",
                 "display_name": "5% with subsidence",
-                "description": "With subsidence; 5th percentile sea rise",
-                "filename": "",
+                "description": "With subsidence; 5th percentile sea level rise",
+                "filename": "inuncoast_{scenario}_wtsub_{year}_0_perc_05",
                 "scenarios": [
-                    {"id": "historical", "years": [1980]},
                     {"id": "rcp4p5", "years": [2030, 2050, 2080]},
                     {"id": "rcp8p5", "years": [2030, 2050, 2080]},
                 ],
             },
             {
                 "event_type": "CoastalInundation",
                 "path": "coastal_inundation/wri/v2",
-                "id": "wtsub_50",
+                "id": "wtsub/50",
                 "display_name": "50% with subsidence",
-                "description": "With subsidence; 50th percentile sea rise",
-                "filename": "",
+                "description": "With subsidence; 50th percentile sea level rise",
+                "filename": "inuncoast_{scenario}_wtsub_{year}_0_perc_50",
                 "scenarios": [
-                    {"id": "historical", "years": [1980]},
                     {"id": "rcp4p5", "years": [2030, 2050, 2080]},
                     {"id": "rcp8p5", "years": [2030, 2050, 2080]},
                 ],

src/physrisk/kernel/calculation.py

@@ -2,12 +2,12 @@
 from collections import defaultdict
 from typing import Any, Dict, List, Optional
 
-from ..data.event_provider import get_source_path_wri_riverine_inundation
+from ..data.event_provider import get_source_path_wri_coastal_inundation, get_source_path_wri_riverine_inundation
 from ..data.pregenerated_hazard_model import ZarrHazardModel
 from ..models import power_generating_asset_models as pgam
 from ..utils.helpers import get_iterable
 from .assets import Asset, PowerGeneratingAsset, TestAsset
-from .events import RiverineInundation
+from .events import CoastalInundation, RiverineInundation
 from .hazard_event_distrib import HazardEventDistrib
 from .hazard_model import HazardModel
 from .impact_distrib import ImpactDistrib
@@ -19,8 +19,8 @@ class AssetImpactResult:
     def __init__(
         self,
         impact: ImpactDistrib,
-        vulnerability: Optional[VulnerabilityDistrib] = None,
-        event: Optional[HazardEventDistrib] = None,
+        vulnerability: VulnerabilityDistrib,
+        event: HazardEventDistrib,
         hazard_data=None,
     ):
         self.impact = impact
@@ -31,7 +31,10 @@ def __init__(
 
 
 def get_default_zarr_source_paths():
-    return {RiverineInundation: get_source_path_wri_riverine_inundation}
+    return {
+        RiverineInundation: get_source_path_wri_riverine_inundation,
+        CoastalInundation: get_source_path_wri_coastal_inundation,
+    }
 
 
 def get_default_hazard_model():

src/physrisk/kernel/impact_curve.py

@@ -21,17 +21,15 @@ def __init__(
         impact_cdfs=None,
     ):
         """Create a new asset event distribution.
+
         Args:
             intensities: possible intensities of hazard event.
             impacts: fractional damage or fractional average disruption occurring as a result
                 of hazard event of given intensity.
             distributions: provides the pdf and optiononally cdf of the impact distribution
         """
 
-        # probabilities must be sorted and decreasing
-        # values must be sorted and non-decreasing (intens[i + 1] >= intens[i])
-
-        if not np.all(np.diff(intensities) > 0):
+        if not np.all(np.diff(intensities) >= 0):
             raise ValueError("intensities must be sorted and increasing")
 
         self.intensities = np.array(intensities)

src/physrisk/models/real_estate_models.py

@@ -9,18 +9,32 @@
 from physrisk.kernel.impact_curve import ImpactCurve
 from physrisk.kernel.vulnerability_model import VulnerabilityModel
 
-from ..kernel.events import RiverineInundation
+from ..kernel.events import CoastalInundation, RiverineInundation
 from ..kernel.vulnerability_model import applies_to_events, get_vulnerability_curves_from_resource
 
 
-@applies_to_events([RiverineInundation])
 class RealEstateInundationModel(VulnerabilityModel):
+    _default_impact_bin_edges = np.array([0, 0.01, 0.1, 0.2, 0.3, 0.4, 0.5, 0.6, 0.7, 0.8, 0.9, 1.0])
+    _default_resource = "EU JRC global flood depth-damage functions"
+
     def __init__(
         self,
         *,
-        resource: str = "EU JRC global flood depth-damage functions",
-        impact_bin_edges=np.array([0, 0.01, 0.1, 0.2, 0.3, 0.4, 0.5, 0.6, 0.7, 0.8, 0.9, 1.0])
+        event_type: type,
+        model: str,
+        resource: str = _default_resource,
+        impact_bin_edges=_default_impact_bin_edges
     ):
+        """
+        Inundation vulnerability model for real estates assets. Applies to both riverine and coastal inundation.
+
+        Args:
+            event_type: Event type.
+            model: optional identifier for hazard event model, passed to HazardModel.
+            resource: embedded resource identifier used to infer vulnerability matrix.
+            impact_bin_edges: specifies the impact (fractional damage/disruption bins).
+        """
+
         curve_set: VulnerabilityCurves = get_vulnerability_curves_from_resource(resource)
 
         # for this model, key for looking up curves is (location, asset_type), e.g. ('Asian', 'Building/Industrial')
@@ -30,7 +44,8 @@ def __init__(
         for item in curve_set.items:
             self.vuln_curves_by_type[item.asset_type].append(item)
 
-        super().__init__(model="MIROC-ESM-CHEM", event_type=RiverineInundation, impact_bin_edges=impact_bin_edges)
+        # global circulation parameter 'model' is a hint; can be overriden by hazard model
+        super().__init__(model=model, event_type=event_type, impact_bin_edges=impact_bin_edges)
 
     def get_impact_curve(self, intensities, asset: RealEstateAsset):
         # we interpolate the mean and standard deviation and use this to construct distributions
@@ -81,3 +96,33 @@ def beta_distrib(mean, std):
     a = (1 - mean) / (cv * cv) - mean
     b = a * (1 - mean) / mean
     return lambda x, a=a, b=b: stats.beta.cdf(x, a, b)
+
+
+@applies_to_events([CoastalInundation])
+class RealEstateCoastalInundationModel(RealEstateInundationModel):
+    def __init__(
+        self,
+        *,
+        model: str = "wtsub/95",
+        resource: str = RealEstateInundationModel._default_resource,
+        impact_bin_edges=RealEstateInundationModel._default_impact_bin_edges
+    ):
+        # by default include subsidence and 95% sea-level rise
+        super().__init__(
+            event_type=CoastalInundation, model=model, resource=resource, impact_bin_edges=impact_bin_edges
+        )
+
+
+@applies_to_events([RiverineInundation])
+class RealEstateRiverineInundationModel(RealEstateInundationModel):
+    def __init__(
+        self,
+        *,
+        model: str = "MIROC-ESM-CHEM",
+        resource: str = RealEstateInundationModel._default_resource,
+        impact_bin_edges=RealEstateInundationModel._default_impact_bin_edges
+    ):
+        # by default request HazardModel to use "MIROC-ESM-CHEM" GCM
+        super().__init__(
+            event_type=RiverineInundation, model=model, resource=resource, impact_bin_edges=impact_bin_edges
+        )

src/test/data/hazard_model_store.py

@@ -3,7 +3,7 @@
 import zarr.storage
 from affine import Affine
 
-from physrisk.data.event_provider import get_source_path_wri_riverine_inundation
+from physrisk.data.event_provider import get_source_path_wri_coastal_inundation, get_source_path_wri_riverine_inundation
 
 
 class TestData:
@@ -43,6 +43,9 @@ class TestData:
         -39.2145,
     ]
 
+    coastal_longitudes = [12.2, 50.5919, 90.3473, 90.4295, 90.4804, 90.3429, 90.5153, 90.6007]
+    coastal_latitudes = [-5.55, 26.1981, 23.6473, 23.6783, 23.5699, 23.9904, 23.59, 23.6112]
+
 
 def get_mock_hazard_model_store_single_curve():
     return_periods = [5.0, 10.0, 25.0, 50.0, 100.0, 250.0, 500.0, 1000.0]
@@ -73,28 +76,39 @@ def get_mock_hazard_model_store_single_curve():
 
 
 def get_mock_hazard_model_store(longitudes, latitudes, curve):
+    """Returns a Zarr store, to be used with unit tests,
+    with the specified longitudes and latitudes set to the curve supplied."""
+
     return_periods = [2.0, 5.0, 10.0, 25.0, 50.0, 100.0, 250.0, 500.0, 1000.0]
     t = [0.008333333333333333, 0.0, -180.0, 0.0, -0.008333333333333333, 90.0, 0.0, 0.0, 1.0]
     shape = (len(return_periods), 21600, 43200)
     store = zarr.storage.MemoryStore(root="hazard.zarr")
     root = zarr.open(store=store, mode="w")
     for model, scenario, year in [("MIROC-ESM-CHEM", "rcp8p5", 2080), ("000000000WATCH", "historical", 1980)]:
-        array_path = get_source_path_wri_riverine_inundation(model=model, scenario=scenario, year=year)
-        z = root.create_dataset(  # type: ignore
-            array_path, shape=(shape[0], shape[1], shape[2]), chunks=(shape[0], 1000, 1000), dtype="f4"
-        )
-        z.attrs["transform_mat3x3"] = t
-        z.attrs["index_values"] = return_periods
-
-        t = z.attrs["transform_mat3x3"]
-        transform = Affine(t[0], t[1], t[2], t[3], t[4], t[5])
-
-        coords = np.vstack((longitudes, latitudes, np.ones(len(longitudes))))
-        inv_trans = ~transform
-        mat = np.array(inv_trans).reshape(3, 3)
-        frac_image_coords = mat @ coords
-        image_coords = np.floor(frac_image_coords).astype(int)
-        for j in range(len(longitudes)):
-            z[:, image_coords[1, j], image_coords[0, j]] = curve
+        array_path_riverine = get_source_path_wri_riverine_inundation(model=model, scenario=scenario, year=year)
+        _add_curves(root, longitudes, latitudes, array_path_riverine, shape, curve, return_periods, t)
+
+    for model, scenario, year in [("wtsub/95", "rcp8p5", 2080), ("wtsub", "historical", 1980)]:
+        array_path_coastal = get_source_path_wri_coastal_inundation(model=model, scenario=scenario, year=year)
+        _add_curves(root, longitudes, latitudes, array_path_coastal, shape, curve, return_periods, t)
 
     return store
+
+
+def _add_curves(root, longitudes, latitudes, array_path, shape, curve, return_periods, t):
+    z = root.create_dataset(  # type: ignore
+        array_path, shape=(shape[0], shape[1], shape[2]), chunks=(shape[0], 1000, 1000), dtype="f4"
+    )
+    z.attrs["transform_mat3x3"] = t
+    z.attrs["index_values"] = return_periods
+
+    t = z.attrs["transform_mat3x3"]
+    transform = Affine(t[0], t[1], t[2], t[3], t[4], t[5])
+
+    coords = np.vstack((longitudes, latitudes, np.ones(len(longitudes))))
+    inv_trans = ~transform
+    mat = np.array(inv_trans).reshape(3, 3)
+    frac_image_coords = mat @ coords
+    image_coords = np.floor(frac_image_coords).astype(int)
+    for j in range(len(longitudes)):
+        z[:, image_coords[1, j], image_coords[0, j]] = curve

src/test/data/test_events_retrieval.py

@@ -38,6 +38,7 @@ def test_hazard_data_availability_summary(self):
         self.assertEqual(summary["RiverineInundation"].years, [1980, 2030, 2050, 2080])
 
     def test_zarr_reading(self):
+
         request_dict = {
             "items": [
                 {
@@ -72,7 +73,7 @@ def test_zarr_reading(self):
     def test_zarr_reading_live(self):
         # needs valid OSC_S3_BUCKET, OSC_S3_ACCESS_KEY, OSC_S3_SECRET_KEY
 
-        request = {
+        request1 = {
             "items": [
                 {
                     "request_item_id": "test_inundation",
@@ -85,7 +86,23 @@ def test_zarr_reading_live(self):
                 }
             ],
         }
-        response = requests.get(request_id="get_hazard_data", request_dict=request)
+        response = requests.get(request_id="get_hazard_data", request_dict=request1)
+        print(response)
+
+        request2 = {
+            "items": [
+                {
+                    "request_item_id": "test_inundation",
+                    "event_type": "CoastalInundation",
+                    "longitudes": TestData.coastal_longitudes,
+                    "latitudes": TestData.coastal_latitudes,
+                    "year": 2080,
+                    "scenario": "rcp8p5",
+                    "model": "wtsub/95",
+                }
+            ],
+        }
+        response = requests.get(request_id="get_hazard_data", request_dict=request2)
         print(response)
 
     @unittest.skip("includes download of large files; deprecated")