tests: finishing b-value tests and removing assumption of binning 0.1

seismostats/catalogs/catalog.py

@@ -328,6 +328,7 @@ def estimate_mc(
                                  self.magnitude.max(),
                                  delta_m)
 
+        # TODO change once we have a global estimate_mc
         mc_est = mc_ks(self.magnitude,
                        mcs_test,
                        delta_m,
@@ -348,7 +349,7 @@ def estimate_b(
         weights: list | None = None,
         b_parameter: str = "b_value",
         return_std: bool = False,
-        method: str = "tinti",
+        method: str = "classic",
         return_n: bool = False,
     ) -> float | tuple[float, float] | tuple[float, float, float]:
         """
@@ -394,14 +395,17 @@ def estimate_b(
         if delta_m is None:
             delta_m = self.delta_m
 
+        # filter magnitudes above mc without changing the original dataframe
+        df = self[self.magnitude >= mc]
+
         if method == "positive":
             # dataframe needs 'time' column to be sorted
-            if 'time' not in self.columns:
+            if 'time' not in df.columns:
                 raise ValueError('"time" column needs to be set in order to \
                                  use the b-positive method')
-            mags = self.sort_values("time").magnitude
+            mags = df.sort_values("time").magnitude
         else:
-            mags = self.magnitude
+            mags = df.magnitude
 
         b_estimate = estimate_b(mags,
                                 mc,

seismostats/catalogs/tests/test_catalog.py

@@ -9,6 +9,7 @@
 from seismostats.catalogs.catalog import (REQUIRED_COLS_CATALOG, Catalog,
                                           ForecastCatalog)
 from seismostats.utils.binning import bin_to_precision
+from seismostats.analysis.estimate_beta import estimate_b
 
 RAW_DATA = {'name': ['Object 1', 'Object 2', 'Object 3'],
             'magnitude': [10.0, 12.5, 8.2],
@@ -84,8 +85,6 @@ def test_forecast_catalog_strip():
     [
         (np.array([0.235, -0.235, 4.499, 4.5, 6, 0.1, 1.6]),
          0.1),
-        (np.array([0.235, -0.235, 4.499, 4.5, 6, 0.1, 1.6]),
-         None),
         (np.array([0.235, -0.235, 4.499, 5.5, 6, 0.1, 1.6]),
          0.2),
         ([0.235, -0.235, 4.499, 5.5, 6, 0.1, 1.6],
@@ -98,43 +97,49 @@ def test_catalog_bin(mag_values: np.ndarray, delta_m: float):
     assert (catalog.bin_magnitudes(
         delta_m)['magnitude'].tolist()
         == bin_to_precision(mag_values, delta_m)).all()
+    assert (catalog.bin_magnitudes()['magnitude'].tolist()
+        == bin_to_precision(mag_values, 0.1)).all()
 
     return_value = catalog.bin_magnitudes(delta_m, inplace=True)
     assert (catalog['magnitude'].tolist()
             == bin_to_precision(mag_values, delta_m)).all()
     assert return_value is None
-
     assert catalog.delta_m == delta_m
 
 
-@pytest.mark.parametrize(
-    "mag_values, delta_m",
-    [
-        (np.array([0.235, -0.235, 4.499, 4.5, 6, 0.1, 1.6]),
-         0),
-    ]
-)
-def test_catalog_bin_none(mag_values: np.ndarray, delta_m: float):
-    catalog = Catalog({'magnitude': mag_values})
-
-    with pytest.raises(ValueError):
-        catalog.bin_magnitudes(delta_m=delta_m)
-
-
 def test_catalog_estimate_mc():
+    # TODO once global_mc method is implemented
     catalog = Catalog({'magnitude': [0.235, -0.235, 4.499, 4.5, 6, 0.1, 1.6]})
 
     with pytest.raises(ValueError):
         catalog.estimate_mc()
 
-
-def test_catalog_estimate_b():
-    catalog = Catalog({'magnitude': [0.235, -0.235, 4.499, 4.5, 6, 0.1, 1.6]})
+@pytest.mark.parametrize(
+    "mag_values, delta_m, mc",
+    [
+        (np.array([0.0, 0.235, 0.238, 4.499, 4.5, 6, 0.1, 1.6]),
+        0.001, 0.0)
+    ]
+)
+def test_catalog_estimate_b(mag_values, delta_m, mc):
+    catalog = Catalog({'magnitude': mag_values})
 
     with pytest.raises(ValueError):
         catalog.estimate_b(mc=None, delta_m=None)
         catalog.estimate_b(mc=1.0, delta_m=None)
         catalog.estimate_b(mc=None, delta_m=0.1)
+        catalog.estimate_b(mc=1.0, delta_m=0.1, method='positive')
+
+    b_value = estimate_b(catalog['magnitude'], mc=mc, delta_m=delta_m, method="classic")
+    return_value = catalog.estimate_b(mc=mc, delta_m=delta_m, method="classic")
+    assert catalog.b_value == b_value
+    assert return_value == b_value
+
+    catalog.mc = mc
+    catalog.delta_m = delta_m
+    return_value = catalog.estimate_b(method="classic")
+    assert catalog.b_value == b_value
+    assert return_value == b_value
 
 
 def test_to_quakeml():

seismostats/utils/binning.py

@@ -57,8 +57,7 @@ def bin_to_precision(x: np.ndarray | list, delta_x: float = 0.1) -> np.ndarray:
 
     if isinstance(x, list):
         x = np.array(x)
-    if delta_x is None:
-        delta_x = 0.1
+
     d = decimal.Decimal(str(delta_x))
     decimal_places = abs(d.as_tuple().exponent)
     return np.round(normal_round_to_int(x / delta_x) * delta_x, decimal_places)

seismostats/utils/tests/test_binning.py

@@ -30,9 +30,6 @@ def test_normal_round(x: float, n: int, rounded_value: float):
         (np.array([0.235, -0.235, 4.499, 4.5, 6, 0.1, 1.6]),
          0.1,
          np.array([0.2, -0.2, 4.5, 4.5, 6, 0.1, 1.6])),
-        (np.array([0.235, -0.235, 4.499, 4.5, 6, 0.1, 1.6]),
-         None,
-         np.array([0.2, -0.2, 4.5, 4.5, 6, 0.1, 1.6])),
         (np.array([0.235, -0.235, 4.499, 5.5, 6, 0.1, 1.6]),
          0.2,
          np.array([0.2, -0.2, 4.4, 5.6, 6, 0.2, 1.6])),
@@ -46,6 +43,10 @@ def test_bin_to_precision(x: np.ndarray, delta_x: float,
     y = bin_to_precision(x, delta_x)
     assert (y == rounded_value).all()
 
+    y = bin_to_precision(x)
+    z = bin_to_precision(x, 0.1)
+    assert (y == z).all()
+
 
 def test_bin_to_precision_none():
     with pytest.raises(ValueError):