implemented changes suggested by Sarah H

climada/engine/impact.py

@@ -473,10 +473,9 @@ def local_exceedance_impact(
             self,
             return_periods=(25, 50, 100, 250),
             method = 'interpolate',
-            impact_cutoff=0,
+            min_impact=0,
             log_frequency=True,
-            log_impact=True,
-            extrapolation = True
+            log_impact=True
         ):
         """Compute local exceedance impact for given return periods. The default method
         is fitting the ordered impacts per centroid to the corresponding cummulated
@@ -488,24 +487,22 @@ def local_exceedance_impact(
             User-specified return periods for which the exceedance intensity should be calculated
             locally (at each centroid). Defaults to (25, 50, 100, 250).
         method : str
-            Method to interpolate to new return periods. Currently available are "interpolate" and 
-            "stepfunction". Defauls to "interpolate".
-        impact_cutoff : float, optional
-            Minimal threshold to filter the impact. Defaults to 0.
+            Method to interpolate to new return periods. Currently available are "interpolate",
+            "extrapolate" and "stepfunction". If set to "interpolate" or "stepfunction",
+            return periods larger than the Impact object's observed local return periods will be
+            assigned the largest local impact, and return periods smaller than the Impact object's
+            observed local return periods will be assigned 0. If set to "extrapolate", local
+            exceedance impacts will be extrapolated (and interpolated). Defauls to "interpolate".
+        min_impact : float, optional
+            Minimum threshold to filter the impact. Defaults to 0.
         log_frequency : bool, optional
-            This parameter is only used if method is set to "interpolate". If set to True,
-            (cummulative) frequency values are converted to log scale before inter- and
+            This parameter is only used if method is set to "extrapolate" or "interpolate". If set
+            to True, (cummulative) frequency values are converted to log scale before inter- and
             extrapolation. Defaults to True.
         log_impact : bool, optional
-            This parameter is only used if method is set to "interpolate". If set to True,
-            impact values are converted to log scale before inter- and extrapolation.
+            This parameter is only used if method is set to "extrapolate" or "interpolate". If set
+            to True, impact values are converted to log scale before inter- and extrapolation.
             Defaults to True.
-        extrapolation : bool, optional
-            This parameter is only used if method is set to "interpolate". If set to True, local
-            exceedance impacts will be extrapolated. If set to False, return periods larger than
-            the Impact object's observed local return periods will be assigned the largest
-            local impact, and return periods smaller than the Impact object's observed local
-            return periods will be assigned 0. Defaults to True.
 
         Returns
         -------
@@ -552,16 +549,17 @@ def local_exceedance_impact(
 
             # fit intensities to cummulative frequencies
             frequency = np.cumsum(frequency[::-1])[::-1]
-            if method == 'interpolate':
-                imp_stats[:,i] = u_interp.interpolate_ev(
-                    1/np.array(return_periods), frequency[::-1], impact[::-1], logx=log_frequency,
-                    logy=log_impact, y_threshold=impact_cutoff, extrapolation=extrapolation, y_asymptotic=0.
-                )
-            elif method == 'stepfunction':
+            if method == 'stepfunction':
                 imp_stats[:,i] = u_interp.stepfunction_ev(
-                    1/np.array(return_periods), frequency[::-1], impact[::-1], y_threshold=impact_cutoff,
+                    1/np.array(return_periods), frequency[::-1], impact[::-1], y_threshold=min_impact,
                     y_asymptotic=0.
                 )
+            elif method == 'extrapolate' or method == 'interpolate':
+                extrapolation = (method == 'extrapolate')
+                imp_stats[:,i] = u_interp.interpolate_ev(
+                    1/np.array(return_periods), frequency[::-1], impact[::-1], logx=log_frequency,
+                    logy=log_impact, y_threshold=min_impact, extrapolation=extrapolation, y_asymptotic=0.
+                )
             else:
                 raise ValueError(f"Unknown method: {method}")
 
@@ -584,8 +582,11 @@ def local_exceedance_imp(
             return_periods=(25, 50, 100, 250)
         ):
         """This function is deprecated, use Impact.local_exceedance_impact instead."""
-        LOGGER.warning("The use of Impact.local_exceedance_imp is deprecated."
-                       "Use Impact.local_exceedance_impact instead.")
+        LOGGER.warning(
+            "The use of Impact.local_exceedance_imp is deprecated. Use "
+            "Impact.local_exceedance_impact instead. Some errors in the previous calculation "
+            "in Impact.local_exceedance_imp have been corrected. To reproduce data with the "
+            "previous calculation, use CLIMADA v5.0.0 or less.")
 
         return self.local_exceedance_impact(return_periods)[0].values[:,1:].T.astype(float)
 

climada/engine/test/test_impact.py

@@ -517,7 +517,7 @@ def test_local_exceedance_impact(self):
         # fourth centroid has intensities 1,2,3,4 with cum frequencies 4,3,2,1
         # testing at frequencies 5, 2, 1, 0.5
         impact_stats, _, _ = impact.local_exceedance_impact(
-                return_periods=(.2, .5, 1, 2), log_frequency=False, log_impact=False)
+                return_periods=(.2, .5, 1, 2), method="extrapolate", log_frequency=False, log_impact=False)
         np.testing.assert_allclose(
             impact_stats.values[:,1:].astype(float),
             np.array([

climada/hazard/base.py

@@ -452,10 +452,9 @@ def local_exceedance_intensity(
             self,
             return_periods=(25, 50, 100, 250),
             method='interpolate',
-            intensity_cutoff=None,
+            min_intensity=None,
             log_frequeny=True,
-            log_intensity=True,
-            extrapolation=True
+            log_intensity=True
     ):
         """Compute local exceedance intensity for given return periods. The default method
         is fitting the ordered intensitites per centroid to the corresponding cummulated
@@ -467,10 +466,15 @@ def local_exceedance_intensity(
             User-specified return periods for which the exceedance intensity should be calculated
             locally (at each centroid). Defaults to (25, 50, 100, 250).
         method : str
-            Method to interpolate to new return periods. Currently available are "interpolate" and 
-            "stepfunction". Defauls to "interpolate".
-        intensity_cutoff : float, optional
-            Minimal threshold to filter the hazard intensity. If set to None, self.intensity_thres
+            Method to interpolate to new return periods. Currently available are "interpolate",
+            "extrapolate" and "stepfunction". If set to "interpolate" or "stepfunction", return
+            periods larger than the Hazard object's observed local return periods will be assigned
+            the largest local intensity, and return periods smaller than the Hazard object's
+            observed local return periods will be assigned 0. If set to "extrapolate", local
+            exceedance intensities will be extrapolated (and interpolated).
+            Defauls to "interpolate".
+        min_intensity : float, optional
+            Minimum threshold to filter the hazard intensity. If set to None, self.intensity_thres
             will be used. Defaults to None.
         log_frequency : bool, optional
             This parameter is only used if method is set to "interpolate". If set to True,
@@ -480,12 +484,6 @@ def local_exceedance_intensity(
             This parameter is only used if method is set to "interpolate". If set to True,
             intensity values are converted to log scale before inter- and extrapolation.
             Defaults to True.
-        extrapolation : bool, optional
-            This parameter is only used if method is set to "interpolate". If set to True, local
-            exceedance intensities will be extrapolated. If set to False, return periods larger than
-            the Hazard object's observed local return periods will be assigned the largest
-            local intensity, and return periods smaller than the Impact object's observed local
-            return periods will be assigned 0. Defaults to True.
 
         Returns
         -------
@@ -499,8 +497,8 @@ def local_exceedance_intensity(
         column_label : function
             Column-label-generating function, for reporting and plotting
         """
-        if not intensity_cutoff and intensity_cutoff != 0:
-            intensity_cutoff = self.intensity_thres
+        if not min_intensity and min_intensity != 0:
+            min_intensity = self.intensity_thres
         #check frequency unit
         if self.frequency_unit in ['1/year', 'annual', '1/y', '1/a']:
             return_period_unit = 'years'
@@ -528,16 +526,17 @@ def local_exceedance_intensity(
 
             # fit intensities to cummulative frequencies
             frequency = np.cumsum(frequency[::-1])[::-1]
-            if method == 'interpolate':
-                inten_stats[:,i] = u_interp.interpolate_ev(
-                    1/np.array(return_periods), frequency[::-1], intensity[::-1], logx=log_frequeny,
-                    logy=log_intensity, y_threshold=intensity_cutoff, y_asymptotic=0., extrapolation=extrapolation
-                )
-            elif method == 'stepfunction':
+            if method == 'stepfunction':
                 inten_stats[:,i] = u_interp.stepfunction_ev(
-                    1/np.array(return_periods), frequency[::-1], intensity[::-1], y_threshold=intensity_cutoff,
+                    1/np.array(return_periods), frequency[::-1], intensity[::-1], y_threshold=min_intensity,
                     y_asymptotic=0.
                 )
+            elif method == 'interpolate' or method == 'extrapolate':
+                extrapolation = (method == 'extrapolate')
+                inten_stats[:,i] = u_interp.interpolate_ev(
+                    1/np.array(return_periods), frequency[::-1], intensity[::-1], logx=log_frequeny,
+                    logy=log_intensity, y_threshold=min_intensity, y_asymptotic=0., extrapolation=extrapolation
+                )
             else:
                 raise ValueError(f"Unknown method: {method}")
 
@@ -556,8 +555,11 @@ def local_exceedance_intensity(
 #TODO: note different calculation in changelog 
     def local_exceedance_inten(self, return_period=(25, 50, 100, 250)):
         """This function is deprecated, use Hazard.local_exceedance_intensity instead."""
-        LOGGER.warning("The use of Hazard.local_exceedance_inten is deprecated."
-                       "Use Hazard.local_exceedance_intensity instead.")
+        LOGGER.warning(
+            "The use of Hazard.local_exceedance_inten is deprecated. Use "
+            "Hazard.local_exceedance_intensity instead. Some errors in the previous calculation "
+            "in Hazard.local_exceedance_inten have been corrected. To reproduce data with the "
+            "previous calculation, use CLIMADA v5.0.0 or less.")
         return self.local_exceedance_intensity(return_period)[0].values[:,1:].T.astype(float)
 
     def sanitize_event_ids(self):
@@ -570,8 +572,7 @@ def local_return_period(
             self,
             threshold_intensities=(10., 20.),
             method='interpolate',
-            intensity_cutoff = None,
-            extrapolation=False,
+            min_intensity = None,
             log_frequency=True,
             log_intensity=True
         ):
@@ -585,10 +586,15 @@ def local_return_period(
             User-specified hazard intensities for which the return period should be calculated
             locally (at each centroid). Defaults to (10, 20)
         method : str
-            Method to interpolate to new intensity values. Currently available are "interpolate" and 
-            "stepfunction". Defauls to "interpolate".
-        intensity_cutoff : float, optional
-            Minimal threshold to filter the hazard intensity. If set to None, self.intensity_thres
+            Method to interpolate to new threshold intensities. Currently available are
+            "interpolate", "extrapolate" and "stepfunction". If set to "interpolate" or
+            "stepfunction", threshold intensities larger than the Hazard object's local
+            intensities will be assigned NaN, and threshold intensities smaller than the Hazard
+            object's local intensities will be assigned the smallest observed local return period.
+            If set to "extrapolate", local return periods will be extrapolated (and interpolated).
+            Defaults to "interpolate".
+        min_intensity : float, optional
+            Minimum threshold to filter the hazard intensity. If set to None, self.intensity_thres
             will be used. Defaults to None.
         log_frequency : bool, optional
             This parameter is only used if method is set to "interpolate". If set to True,
@@ -598,12 +604,6 @@ def local_return_period(
             This parameter is only used if method is set to "interpolate". If set to True,
             intensity values are converted to log scale before inter- and extrapolation.
             Defaults to True.
-        extrapolation : bool, optional
-            This parameter is only used if method is set to "interpolate". If set to True, local
-            return periods will be extrapolated. If set to False, threshold intensities larger than
-            the Hazard object's local intensities will be assigned NaN, and threshold intensities
-            smaller than the Hazard object's local intensities will be assigned the smallest
-            observed local return period. Defaults to False.
 
         Returns
         -------
@@ -617,8 +617,8 @@ def local_return_period(
         column_label : function
             Column-label-generating function, for reporting and plotting
         """
-        if not intensity_cutoff and intensity_cutoff != 0:
-            intensity_cutoff = self.intensity_thres
+        if not min_intensity and min_intensity != 0:
+            min_intensity = self.intensity_thres
         #check frequency unit
         if self.frequency_unit in ['1/year', 'annual', '1/y', '1/a']:
             return_period_unit = 'Years'
@@ -649,14 +649,15 @@ def local_return_period(
 
             # fit intensities to cummulative frequencies
             frequency = np.cumsum(frequency[::-1])[::-1]
-            if method == 'interpolate':
+            if method == 'stepfunction':
+                return_periods[:,i] = u_interp.stepfunction_ev(
+                    threshold_intensities, intensity, frequency, x_threshold=min_intensity
+                )
+            elif method == 'interpolate' or method == "extrapolate":
+                extrapolation = (method == "extrapolate")
                 return_periods[:,i] = u_interp.interpolate_ev(
                     threshold_intensities, intensity, frequency, logx=log_intensity,
-                    logy=log_frequency, x_threshold=intensity_cutoff, extrapolation=extrapolation, y_asymptotic=np.nan
-                )
-            elif method == 'stepfunction':
-                return_periods[:,i] = u_interp.stepfunction_ev(
-                    threshold_intensities, intensity, frequency, x_threshold=intensity_cutoff
+                    logy=log_frequency, x_threshold=min_intensity, extrapolation=extrapolation, y_asymptotic=np.nan
                 )
             else:
                 raise ValueError(f"Unknown method: {method}")

climada/hazard/test/test_base.py

@@ -1060,7 +1060,7 @@ def test_local_exceedance_intensity(self):
         )
 
     def test_local_return_period(self):
-        """Test local return periods easy lin lin interpolation"""
+        """Test local return periods with lin lin interpolation"""
         haz = dummy_hazard()
         haz.intensity = sparse.csr_matrix([
             [1., 4., 1.],
@@ -1072,7 +1072,8 @@ def test_local_return_period(self):
         # second centroid has intensities 2, 4 with cum frequencies 2, 1
         # third centroid has intensities 1 with cum frequencies 1 (0 intensity is neglected)
         # testing at intensities 1, 2, 3
-        return_stats, _, _ = haz.local_return_period(threshold_intensities, log_frequency=False, log_intensity=False, intensity_cutoff=0)
+        return_stats, _, _ = haz.local_return_period(
+            threshold_intensities, log_frequency=False, log_intensity=False, min_intensity=0)
         np.testing.assert_allclose(
             return_stats[return_stats.columns[1:]].values,
             np.array([

climada/test/test_engine.py

@@ -291,8 +291,8 @@ def test_calc_sensitivity_pass(self):
 class TestImpactRPCals(unittest.TestCase):
     """Test the return periods and exceedance impact calculation of Impact objects"""
 
-    def test_local_exceedance_impact_options(self):
-        """Test local exceedance impacts per return period with different options"""
+    def test_local_exceedance_impact_methods(self):
+        """Test local exceedance impacts per return period with different methods"""
         impact = dummy_impact()
         impact.coord_exp = np.array([np.arange(4), np.arange(4)]).T
         impact.imp_mat = sp.sparse.csr_matrix(
@@ -324,7 +324,7 @@ def test_local_exceedance_impact_options(self):
 
         # test log log extrapolation
         impact_stats, _, _ = impact.local_exceedance_impact(
-                return_periods=(1000, 30, .1))
+                return_periods=(1000, 30, .1), method = "extrapolate")
         np.testing.assert_allclose(
             impact_stats.values[:,1:].astype(float),
             np.array([
@@ -337,7 +337,7 @@ def test_local_exceedance_impact_options(self):
 
         # test log log interpolation and no extrapolation
         impact_stats, _, _ = impact.local_exceedance_impact(
-                return_periods=(1000, 30, .1), extrapolation=False)
+                return_periods=(1000, 30, .1))
         np.testing.assert_allclose(
             impact_stats.values[:,1:].astype(float),
             np.array([
@@ -351,7 +351,7 @@ def test_local_exceedance_impact_options(self):
         # test lin lin interpolation with no extrapolation
         impact_stats, _, _ = impact.local_exceedance_impact(
                 return_periods=(1000, 30, .1),
-                log_frequency=False, log_impact=False, extrapolation=False)
+                log_frequency=False, log_impact=False)
         np.testing.assert_allclose(
             impact_stats.values[:,1:].astype(float),
             np.array([