diff --git a/climada/engine/impact_forecast.py b/climada/engine/impact_forecast.py index 2160b18da7..9e8d2af0ba 100644 --- a/climada/engine/impact_forecast.py +++ b/climada/engine/impact_forecast.py @@ -22,6 +22,7 @@ import logging import numpy as np +import scipy.sparse as sparse from ..util import log_level from ..util.checker import size @@ -185,6 +186,123 @@ def _check_sizes(self): size(exp_len=num_entries, var=self.member, var_name="Forecast.member") size(exp_len=num_entries, var=self.lead_time, var_name="Forecast.lead_time") + def _reduce_attrs(self, event_name: str): + """ + Reduce the attributes of an ImpactForecast to a single value. + + Attributes are modified as follows: + - lead_time: set to NaT + - member: set to -1 + - event_id: set to 0 + - event_name: set to the name of the reduction method (default) + - date: set to 0 + - frequency: set to 1 + + Parameters + ---------- + event_name : str + The event name given to the reduced data. + """ + reduced_attrs = { + "lead_time": np.array([np.timedelta64("NaT")]), + "member": np.array([-1]), + "event_id": np.array([0]), + "event_name": np.array([event_name]), + "date": np.array([0]), + "frequency": np.array([1]), + } + + return reduced_attrs + + def min(self): + """ + Reduce the impact matrix and at_event of an ImpactForecast to the minimum + value. + + Parameters + ---------- + None + + Returns + ------- + ImpactForecast + An ImpactForecast object with the min impact matrix and at_event. + """ + red_imp_mat = self.imp_mat.min(axis=0).tocsr() + red_at_event = np.array([red_imp_mat.sum()]) + return ImpactForecast( + frequency_unit=self.frequency_unit, + coord_exp=self.coord_exp, + crs=self.crs, + eai_exp=self.eai_exp, + at_event=red_at_event, + tot_value=self.tot_value, + aai_agg=self.aai_agg, + unit=self.unit, + imp_mat=red_imp_mat, + haz_type=self.haz_type, + **self._reduce_attrs("min"), + ) + + def max(self): + """ + Reduce the impact matrix and at_event of an ImpactForecast to the maximum + value. + + Parameters + ---------- + None + + Returns + ------- + ImpactForecast + An ImpactForecast object with the max impact matrix and at_event. + """ + red_imp_mat = self.imp_mat.max(axis=0).tocsr() + red_at_event = np.array([red_imp_mat.sum()]) + return ImpactForecast( + frequency_unit=self.frequency_unit, + coord_exp=self.coord_exp, + crs=self.crs, + eai_exp=self.eai_exp, + at_event=red_at_event, + tot_value=self.tot_value, + aai_agg=self.aai_agg, + unit=self.unit, + imp_mat=red_imp_mat, + haz_type=self.haz_type, + **self._reduce_attrs("max"), + ) + + def mean(self): + """ + Reduce the impact matrix and at_event of an ImpactForecast to the mean value. + + Parameters + ---------- + None + + Returns + ------- + ImpactForecast + An ImpactForecast object with the mean impact matrix and at_event. + """ + red_imp_mat = sparse.csr_matrix(self.imp_mat.mean(axis=0)) + red_at_event = np.array([red_imp_mat.sum()]) + return ImpactForecast( + frequency_unit=self.frequency_unit, + coord_exp=self.coord_exp, + crs=self.crs, + eai_exp=self.eai_exp, + at_event=red_at_event, + tot_value=self.tot_value, + aai_agg=self.aai_agg, + unit=self.unit, + imp_mat=red_imp_mat, + haz_type=self.haz_type, + **self._reduce_attrs("mean"), + ) + def select( self, event_ids=None, @@ -205,10 +323,6 @@ def select( lead_time : Sequence of numpy.timedelta64 Lead times to select - Returns - ------- - ImpactForecast - See Also -------- :py:meth:`~climada.engine.impact.Impact.select` diff --git a/climada/engine/test/test_impact_forecast.py b/climada/engine/test/test_impact_forecast.py index 94655c8b17..5742422f96 100644 --- a/climada/engine/test/test_impact_forecast.py +++ b/climada/engine/test/test_impact_forecast.py @@ -212,7 +212,7 @@ def test_impact_forecast_concat(impact_forecast, member): def test_impact_forecast_blocked_methods(impact_forecast): - """Check if blocked methods raise NotImplementedError""" + """Check if ImpactForecast.exceedance_freq_curve raises NotImplementedError""" with pytest.raises(NotImplementedError): impact_forecast.local_exceedance_impact(np.array([10, 50, 100])) @@ -221,3 +221,40 @@ def test_impact_forecast_blocked_methods(impact_forecast): with pytest.raises(NotImplementedError): impact_forecast.calc_freq_curve(np.array([10, 50, 100])) + + +@pytest.fixture +def impact_forecast_stats(impact_kwargs, lead_time, member): + max_index = 4 + for key, val in impact_kwargs.items(): + if isinstance(val, (np.ndarray, list)): + impact_kwargs[key] = val[:max_index] + elif isinstance(val, csr_matrix): + impact_kwargs[key] = val[:max_index, :] + impact_kwargs["imp_mat"] = csr_matrix([[1, 0], [0, 1], [3, 2], [2, 3]]) + impact_kwargs["at_event"] = np.array([1, 1, 5, 5]) + return ImpactForecast( + lead_time=lead_time[:max_index], member=member[:max_index], **impact_kwargs + ) + + +@pytest.mark.parametrize("attr", ["min", "mean", "max"]) +def test_impact_forecast_min_mean_max(impact_forecast_stats, attr): + """Check mean, min, and max methods for ImpactForecast""" + imp_fc_reduced = getattr(impact_forecast_stats, attr)() + + # assert imp_mat + npt.assert_array_equal( + imp_fc_reduced.imp_mat.todense(), + getattr(impact_forecast_stats.imp_mat.todense(), attr)(axis=0), + ) + at_event_expected = {"min": [0], "mean": [3], "max": [6]} + npt.assert_array_equal(imp_fc_reduced.at_event, at_event_expected[attr]) + + # check that attributes where reduced correctly + npt.assert_array_equal(np.isnat(imp_fc_reduced.lead_time), [True]) + npt.assert_array_equal(imp_fc_reduced.member, [-1]) + npt.assert_array_equal(imp_fc_reduced.event_name, [attr]) + npt.assert_array_equal(imp_fc_reduced.event_id, [0]) + npt.assert_array_equal(imp_fc_reduced.frequency, [1]) + npt.assert_array_equal(imp_fc_reduced.date, [0]) diff --git a/climada/hazard/forecast.py b/climada/hazard/forecast.py index a8c5cdc543..abab0d0787 100644 --- a/climada/hazard/forecast.py +++ b/climada/hazard/forecast.py @@ -22,6 +22,7 @@ import logging import numpy as np +import scipy.sparse as sparse from ..util.checker import size from ..util.forecast import Forecast @@ -105,6 +106,115 @@ def _check_sizes(self): size(exp_len=num_entries, var=self.member, var_name="Forecast.member") size(exp_len=num_entries, var=self.lead_time, var_name="Forecast.lead_time") + def _reduce_attrs(self, event_name: str): + """ + Reduce the attributes of a HazardForecast to a single value. + + Attributes are modified as follows: + - lead_time: set to NaT + - member: set to -1 + - event_id: set to 0 + - event_name: set to the name of the reduction method (default) + - date: set to 0 + - frequency: set to 1 + + Parameters + ---------- + event_name : str + The event_name given to the reduced data. + """ + reduced_attrs = { + "lead_time": np.array([np.timedelta64("NaT")]), + "member": np.array([-1]), + "event_id": np.array([0]), + "event_name": np.array([event_name]), + "date": np.array([0]), + "frequency": np.array([1]), + "orig": np.array([True]), + } + + return reduced_attrs + + def min(self): + """ + Reduce the intensity and fraction of a HazardForecast to the minimum + value. + + Parameters + ---------- + None + + Returns + ------- + HazardForecast + A HazardForecast object with the min intensity and fraction. + """ + red_intensity = self.intensity.min(axis=0).tocsr() + red_fraction = self.fraction.min(axis=0).tocsr() + return HazardForecast( + haz_type=self.haz_type, + pool=self.pool, + units=self.units, + centroids=self.centroids, + frequency_unit=self.frequency_unit, + intensity=red_intensity, + fraction=red_fraction, + **self._reduce_attrs("min"), + ) + + def max(self): + """ + Reduce the intensity and fraction of a HazardForecast to the maximum + value. + + Parameters + ---------- + None + + Returns + ------- + HazardForecast + A HazardForecast object with the min intensity and fraction. + """ + red_intensity = self.intensity.max(axis=0).tocsr() + red_fraction = self.fraction.max(axis=0).tocsr() + return HazardForecast( + haz_type=self.haz_type, + pool=self.pool, + units=self.units, + centroids=self.centroids, + frequency_unit=self.frequency_unit, + intensity=red_intensity, + fraction=red_fraction, + **self._reduce_attrs("max"), + ) + + def mean(self): + """ + Reduce the intensity and fraction of a HazardForecast to the mean value. + + Parameters + ---------- + None + + Returns + ------- + HazardForecast + A HazardForecast object with the min intensity and fraction. + """ + red_intensity = sparse.csr_matrix(self.intensity.mean(axis=0)) + red_fraction = sparse.csr_matrix(self.fraction.mean(axis=0)) + return HazardForecast( + haz_type=self.haz_type, + pool=self.pool, + units=self.units, + centroids=self.centroids, + frequency_unit=self.frequency_unit, + intensity=red_intensity, + fraction=red_fraction, + **self._reduce_attrs("mean"), + ) + @classmethod def concat(cls, haz_list: list): """Concatenate multiple HazardForecast instances and return a new object""" @@ -138,10 +248,6 @@ def select( lead_time : Sequence of numpy.timedelta64 Lead times to select - Returns - ------- - HazardForecast - See Also -------- :py:meth:`~climada.hazard.base.Hazard.select` diff --git a/climada/hazard/test/test_forecast.py b/climada/hazard/test/test_forecast.py index cb94767241..99bc1d6d73 100644 --- a/climada/hazard/test/test_forecast.py +++ b/climada/hazard/test/test_forecast.py @@ -233,3 +233,28 @@ def test_write_read_hazard_forecast(haz_fc, tmp_path): else: # npt.assert_array_equal also works for comparing int, float or list npt.assert_array_equal(haz_fc.__dict__[key], haz_fc_read.__dict__[key]) + + +@pytest.mark.parametrize("attr", ["min", "mean", "max"]) +def test_hazard_forecast_mean_min_max(haz_fc, attr): + """Check mean, min, and max methods for ImpactForecast""" + haz_fcst_reduced = getattr(haz_fc, attr)() + + # Assert sparse matrices + npt.assert_array_equal( + haz_fcst_reduced.intensity.todense(), + getattr(haz_fc.intensity.todense(), attr)(axis=0), + ) + npt.assert_array_equal( + haz_fcst_reduced.fraction.todense(), + getattr(haz_fc.fraction.todense(), attr)(axis=0), + ) + + # Check that attributes where reduced correctly + npt.assert_array_equal(np.isnat(haz_fcst_reduced.lead_time), [True]) + npt.assert_array_equal(haz_fcst_reduced.member, [-1]) + npt.assert_array_equal(haz_fcst_reduced.event_name, [attr]) + npt.assert_array_equal(haz_fcst_reduced.event_id, [0]) + npt.assert_array_equal(haz_fcst_reduced.frequency, [1]) + npt.assert_array_equal(haz_fcst_reduced.date, [0]) + npt.assert_array_equal(haz_fcst_reduced.orig, [True])