Binomial testing (#208)

* added to unit test for binary_joint_LL

Author: khawajasim

csep/core/binomial_evaluations.py

@@ -88,13 +88,13 @@ def binary_joint_log_likelihood_ndarray(forecast, catalog):
                     It has to be a either zero or positive integer only (No Floating Point)
     """
     #First, we mask the forecast in cells where we could find log=0.0 singularities:
-    forecast_masked = np.ma.masked_where(forecast.ravel() <= 0.0, forecast.ravel()) 
+    forecast_masked = numpy.ma.masked_where(forecast.ravel() <= 0.0, forecast.ravel())
 
     #Then, we compute the log-likelihood of observing one or more events given a Poisson distribution, i.e., 1 - Pr(0) 
     target_idx = numpy.nonzero(catalog.ravel())
     y = numpy.zeros(forecast_masked.ravel().shape)
     y[target_idx[0]] = 1
-    first_term = y * (np.log(1.0 - np.exp(-forecast_masked.ravel())))
+    first_term = y * (numpy.log(1.0 - numpy.exp(-forecast_masked.ravel())))
 
     #Also, we estimate the log-likelihood in cells no events are observed:
     second_term = (1-y) * (-forecast_masked.ravel().data)

tests/test_evaluations.py

@@ -3,6 +3,8 @@
 import unittest
 
 from csep.core.poisson_evaluations import _simulate_catalog, _poisson_likelihood_test
+from csep.core.binomial_evaluations import binary_joint_log_likelihood_ndarray
+
 
 def get_datadir():
     root_dir = os.path.dirname(os.path.abspath(__file__))
@@ -71,3 +73,17 @@ def test_likelihood(self):
         # calculated by hand given the expected data, see explanation in zechar et al., 2010.
         numpy.testing.assert_allclose(simulated_ll[0], -7.178053830347945)
 
+
+class TestBinomialLikelihood(unittest.TestCase):
+    def __init__(self, *args, **kwargs):
+        super().__init__(*args, **kwargs)
+        self.forecast_data = numpy.array([[0.1, 0.3, 0.4], [0.2, 0.1, 0.1]])
+        self.observed_data = numpy.array([[0, 1, 2], [1, 1, 0]])
+
+    def test_likelihood(self):
+        bill = binary_joint_log_likelihood_ndarray(self.forecast_data, self.observed_data)
+
+        numpy.testing.assert_allclose(bill, -6.7197988064)
+
+if __name__ == '__main__':
+    unittest.main()