Fixed a sampling bug. Made sampling not reach the maximum recursion d…

…epth. Updated notebooks

src/probabilistic_model/bayesian_network/distributions.py

@@ -94,12 +94,9 @@ def forward_pass(self, event: ComplexEvent):
             # calculate the probability of said state
             parent_state_probability = self.parent.forward_message.likelihood(parent_state)
 
-            event_for_parent = ComplexEvent([Event({self.parent.variable: parent_state})])
             # construct the conditional distribution
-            # conditional, current_probability = (self.conditional_probability_distributions[parent_state]
-            #                                     ._conditional(event))
-
-            conditional, current_probability = self.conditional_probability_distributions[parent_state].conditional(event_for_parent)
+            conditional, current_probability = (self.conditional_probability_distributions[parent_state]
+                                                ._conditional(event))
 
             # if the conditional is None, skip
             if conditional is None:

src/probabilistic_model/distributions/distributions.py

@@ -453,8 +453,8 @@ def __init__(self, variable: Continuous, location: float, density_cap: float = f
         self.density_cap = density_cap
 
     @property
-    def domain(self) -> Event:
-        return Event({self.variable: portion.singleton(self.location)})
+    def domain(self) -> ComplexEvent:
+        return ComplexEvent([Event({self.variable: portion.singleton(self.location)})])
 
     def _pdf(self, value: float) -> float:
         if value == self.location:
@@ -474,8 +474,8 @@ def _probability(self, event: EncodedEvent) -> float:
         else:
             return 0
 
-    def _mode(self):
-        return ComplexEvent([self.domain.encode()]), self.density_cap
+    def _mode(self) -> Tuple[ComplexEvent, float]:
+        return self.domain.encode(), self.density_cap
 
     def sample(self, amount: int) -> List[List[float]]:
         return [[self.location] for _ in range(amount)]

src/probabilistic_model/distributions/gaussian.py

@@ -330,11 +330,13 @@ def robert_rejection_sample(self, amount: int) -> np.ndarray:
 
         # enforce an upper bound if it is infinite
         if standard_distribution.interval.upper >= float("inf"):
-            standard_distribution.interval = standard_distribution.interval.replace(upper=10)
+            standard_distribution.interval = (standard_distribution.interval.
+                                              replace(upper=standard_distribution.interval.lower + 10))
 
         # enforce a lower bound if it is infinite
         if standard_distribution.interval.lower <= -float("inf"):
-            standard_distribution.interval = standard_distribution.interval.replace(lower=-10)
+            standard_distribution.interval = (standard_distribution.interval.
+                                              replace(lower=standard_distribution.interval.upper - 10))
 
         # sample from double truncated standard normal instead
         samples = standard_distribution.robert_rejection_sample_from_standard_normal_with_double_truncation(amount)
@@ -348,13 +350,28 @@ def robert_rejection_sample(self, amount: int) -> np.ndarray:
     def robert_rejection_sample_from_standard_normal_with_double_truncation(self, amount: int) -> np.ndarray:
         """
         Use robert rejection sampling to sample from the truncated standard normal distribution.
+        Resamples as long as the amount of samples is not reached.
 
         :param amount: The amount of samples to generate
         :return: The samples
         """
         assert self.scale == 1 and self.mean == 0
         # sample from uniform distribution over this distribution's interval
+        accepted_samples = np.array([])
+        while len(accepted_samples) < amount:
+            accepted_samples = np.append(
+                accepted_samples,
+                self.robert_rejection_sample_from_standard_normal_with_double_truncation_helper(amount - len(accepted_samples)))
+        return accepted_samples
+
+    def robert_rejection_sample_from_standard_normal_with_double_truncation_helper(self, amount: int) -> np.ndarray:
+        """
+        Use robert rejection sampling to sample from the truncated standard normal distribution.
 
+        :param amount: The maximum number of samples to generate. The actual number of samples can be lower due to
+            rejection sampling.
+        :return: The samples
+        """
         uniform_samples = np.random.uniform(self.interval.lower, self.interval.upper, amount)
 
         # if the mean in the interval
@@ -371,18 +388,11 @@ def robert_rejection_sample_from_standard_normal_with_double_truncation(self, am
         else:
             raise ValueError("This should never happen")
 
-        # generate standard uniform samples
-        different_uniform_samples = np.random.uniform(0, 1, amount)
+        # generate standard uniform samples as acceptance probabilities
+        acceptance_probabilities = np.random.uniform(0, 1, amount)
 
         # accept samples that are below the limiting function
-        accepted_samples = uniform_samples[different_uniform_samples <= limiting_function]
-
-        # if any got rejected
-        if len(accepted_samples) < amount:
-            # resample the rejected samples
-            accepted_samples = np.concatenate(
-                [accepted_samples, self.robert_rejection_sample(amount - len(accepted_samples))])
-
+        accepted_samples = uniform_samples[acceptance_probabilities <= limiting_function]
         return accepted_samples
 
     def sample(self, amount: int) -> List[List[float]]:

src/probabilistic_model/probabilistic_circuit/probabilistic_circuit.py

@@ -1,11 +1,9 @@
 from __future__ import annotations
 
-import itertools
 import random
 from typing import Tuple, Iterable, TYPE_CHECKING
 
 import networkx as nx
-import numpy as np
 import portion
 from random_events.events import EncodedEvent, VariableMap, Event, ComplexEvent
 from random_events.variables import Variable, Symbolic, Continuous
@@ -195,7 +193,6 @@ def filter_variable_map_by_self(self, variable_map: VariableMap):
         return variable_map.__class__(
             {variable: value for variable, value in variable_map.items() if variable in variables})
 
-    # @cache_inference_result
     def _conditional(self, event: ComplexEvent) -> Tuple[Optional[Self], float]:
 
         # skip trivial case
@@ -882,6 +879,7 @@ def merge_modes_if_one_dimensional(self, modes: List[EncodedEvent]) -> List[Enco
 
     @cache_inference_result
     def _mode(self) -> Tuple[ComplexEvent, float]:
+
         modes = []
         likelihoods = []
 
@@ -893,7 +891,6 @@ def _mode(self) -> Tuple[ComplexEvent, float]:
 
         # get the most likely result
         maximum_likelihood = max(likelihoods)
-
         mode_events = []
 
         # gather all results that are maximum likely
@@ -979,6 +976,14 @@ def _probability(self, event: EncodedEvent) -> float:
     def is_deterministic(self) -> bool:
         return True
 
+    def is_decomposable(self):
+        for index, subcircuit in enumerate(self.subcircuits):
+            variables = subcircuit.variables
+            for subcircuit_ in self.subcircuits[index+1:]:
+                if len(set(subcircuit_.variables).intersection(set(variables))) > 0:
+                    return False
+        return True
+
     @cache_inference_result
     def _mode(self) -> Tuple[ComplexEvent, float]:
 
@@ -987,8 +992,10 @@ def _mode(self) -> Tuple[ComplexEvent, float]:
 
         # gather all modes from the children
         for subcircuit in self.subcircuits[1:]:
+
             subcircuit_mode, subcircuit_likelihood = subcircuit._mode()
             mode = mode & subcircuit_mode
+
             likelihood *= subcircuit_likelihood
 
         return mode, likelihood
@@ -1097,27 +1104,6 @@ def __copy__(self):
             result.probabilistic_circuit.add_edge(result, copied_subcircuit)
         return result
 
-    def is_decomposable(self) -> bool:
-        """
-        Check if only this product unit is decomposable.
-
-        A product mode is decomposable iff all children have disjoint scopes.
-
-        :return: if this product unit is decomposable
-        """
-        # for every child pair
-        for subcircuits_a, subcircuits_b in itertools.combinations(self.subcircuits, 2):
-
-            # form the intersection of the scopes
-            scope_intersection = set(subcircuits_a.variables) & set(subcircuits_b.variables)
-
-            # if this not empty, the product unit is not decomposable
-            if len(scope_intersection) > 0:
-                return False
-
-        # if every pairwise intersection is empty, the product unit is decomposable
-        return True
-
     @cache_inference_result
     def simplify(self) -> Self:
 
@@ -1234,7 +1220,6 @@ def marginal(self, variables: Iterable[Variable]) -> Optional[Self]:
         root.reset_result_of_current_query()
         return result.probabilistic_circuit
 
-    # @graph_inference_caching_wrapper
     def sample(self, amount: int) -> Iterable:
         return self.root.sample(amount)
 

test/test_bayesian_network/test_bayesian_network.py

@@ -77,16 +77,6 @@ def test_likelihood(self):
         likelihood = self.model.likelihood(event)
         self.assertEqual(likelihood, 0.3 * 0.3)
 
-    def test_probability(self):
-        event = Event({self.x: [1], self.y: [0, 1]})
-        probability = self.model.probability(event)
-        self.assertEqual(probability, self.bf_distribution.probability(event))
-
-    def test_probability_empty_x(self):
-        event = Event({self.y: [1]})
-        probability = self.model.probability(event)
-        self.assertEqual(probability, self.bf_distribution.probability(event))
-
     def test_as_probabilistic_circuit(self):
         circuit = self.model.as_probabilistic_circuit()
         for event in itertools.product(self.x.domain, self.y.domain):
@@ -141,11 +131,6 @@ def test_likelihood(self):
         for event in itertools.product(*[variable.domain for variable in self.model.variables]):
             self.assertEqual(self.model.likelihood(event), self.bf_distribution.likelihood(event))
 
-    def test_probability(self):
-        event = Event({self.x: [1], self.y: [0, 1], self.z: [0]})
-        probability = self.model.probability(event)
-        self.assertEqual(probability, self.bf_distribution.probability(event))
-
     def test_as_probabilistic_circuit(self):
         circuit = self.model.as_probabilistic_circuit().simplify()
         self.assertLess(len(circuit.weighted_edges), math.prod([len(v.domain) for v in circuit.variables]))
@@ -193,16 +178,11 @@ def test_likelihood(self):
         likelihood = self.bayesian_network.likelihood(event)
         self.assertEqual(likelihood, 0.3 * 0.5 * 1 / 3)
 
-    def test_probability(self):
-        event = Event({self.x: [0, 1], self.y: portion.closed(1.5, 2)})
-        probability = self.bayesian_network.probability(event)
-        self.assertEqual(probability, 0.3 * 0.25)
-
     def test_as_probabilistic_circuit(self):
         circuit = self.bayesian_network.as_probabilistic_circuit().simplify()
         self.assertEqual(circuit.probability(Event()), 1.)
         event = Event({self.x: [0, 1], self.y: portion.closed(1.5, 2)})
-        self.assertAlmostEqual(self.bayesian_network.probability(event), circuit.probability(event))
+        self.assertAlmostEqual(0.075, circuit.probability(event))
 
 
 class BayesianNetworkWrongOrderTestCase(unittest.TestCase):

test/test_bayesian_network/test_distributions.py

@@ -141,18 +141,13 @@ def test_forward_pass(self):
         self.assertEqual(self.p_x.forward_probability, 1)
         self.assertEqual(self.p_yzx.forward_probability, 1)
 
-    def test_probability(self):
-        event = Event({self.x: 0, self.y: portion.closed(0, 0.5)})
-        probability = self.bayesian_network.probability(event)
-        self.assertEqual(probability, 0.7 * 0.5)
-
     def test_joint_distribution_with_parent(self):
         event = self.bayesian_network.preprocess_event(Event())
         self.bayesian_network.forward_pass(event)
 
         joint_distribution = self.p_yzx.joint_distribution_with_parent()
         event = Event({self.x: 0, self.y: portion.closed(0, 0.5)})
-        self.assertEqual(joint_distribution.probability(event), self.bayesian_network.probability(event))
+        self.assertEqual(joint_distribution.probability(event),0.35)
 
 
 if __name__ == '__main__':

test/test_distributions/test_gaussian.py

@@ -326,6 +326,11 @@ def test_with_zero_in_bound(self):
         expectation = model.expectation(model.variables)[model.variable]
         self.assertAlmostEqual(mean, expectation, delta=0.1)
 
+    def test_with_far_right_interval(self):
+        model = TruncatedGaussianDistribution(self.x, portion.closed(11, float("inf")), 0, 1)
+        samples = model.sample(1000)
+        self.assertEqual(len(samples), 1000)
+
 
 if __name__ == '__main__':
     unittest.main()

test/test_jpt/test_jpt.py

@@ -121,21 +121,6 @@ def setUp(self):
         self.real, self.integer, self.symbol = infer_variables_from_dataframe(self.data, scale_continuous_types=False)
         self.model = JPT([self.real, self.integer, self.symbol])
 
-    @unittest.skip("Scaled variables are kinda buggy.")
-    def test_preprocess_data(self):
-        preprocessed_data = self.model.preprocess_data(self.data)
-        mean = preprocessed_data[:, 1].mean()
-        std = preprocessed_data[:, 1].std(ddof=1)
-        self.assertEqual(self.real.mean, mean)
-        self.assertEqual(self.real.std, std)
-
-        # assert that this does not throw exceptions
-        for variable, column in zip(self.model.variables, preprocessed_data.T):
-            if isinstance(variable, random_events.variables.Discrete):
-                variable.decode_many(column.astype(int))
-            else:
-                variable.decode_many(column)
-
     def test_construct_impurity(self):
         impurity = self.model.construct_impurity()
         self.assertIsInstance(impurity, Impurity)
@@ -291,7 +276,7 @@ def setUpClass(cls):
 
         variables = infer_variables_from_dataframe(cls.data, scale_continuous_types=False)
 
-        cls.model = JPT(variables, min_samples_leaf=0.1)
+        cls.model = JPT(variables, min_samples_leaf=0.4)
         cls.model.fit(cls.data)
 
     def test_serialization(self):
@@ -344,6 +329,10 @@ def test_marginal_conditional_chain(self):
         x, y = self.model.variables[:2]
         conditional, probability = model.conditional(Event({x: portion.closed(0, 10)}))
 
+    def test_mode(self):
+        mode, likelihood = self.model.mode()
+        self.assertGreater(len(mode.events), 0)
+
 
 class MNISTTestCase(unittest.TestCase):
     model: JPT
@@ -525,21 +514,13 @@ class MaxProblemTestCase(unittest.TestCase):
     model: ProbabilisticCircuit
     relative_x: Continuous
     relative_y: Continuous
-    event: ComplexEvent
 
     @classmethod
     def setUpClass(cls):
-        with open(os.path.join(os.path.expanduser("~"), "Documents", "move_and_pick_up.pm"),  "r") as file:
+        with open(os.path.join(os.path.expanduser("~"), "Documents", "boob_cancer_jpt.json"),  "r") as file:
             loaded_json = json.load(file)
             cls.model = ProbabilisticCircuit.from_json(loaded_json)
 
-        with open(os.path.join(os.path.expanduser("~"), "Documents", "complex.event"),  "r") as file:
-            loaded_json = json.load(file)
-            cls.event = ComplexEvent.from_json(loaded_json)
-        cls.relative_x = [var for var in cls.model.variables if var.name == "relative_x"][0]
-        cls.relative_y = [var for var in cls.model.variables if var.name == "relative_y"][0]
-
     def test_inference(self):
-        p_xy = self.model.marginal([self.relative_x, self.relative_y])
-        conditional, probability = self.model.conditional(self.event)
-        print(conditional)
+        mode, _ = self.model.mode()
+        print(mode.events)