fix: prevent overflow in gaussian log-likelihood (#59)

probcomp · Dec 12, 2023 · 0b015ab · 0b015ab
1 parent eefa9af
commit 0b015ab
Show file tree

Hide file tree

Showing 2 changed files with 7 additions and 7 deletions.
diff --git a/src/gen/distribution/math/log_likelihood.cljc b/src/gen/distribution/math/log_likelihood.cljc
@@ -160,12 +160,11 @@
   \\end{aligned}
   $$"
   [mu sigma v]
-  (let [v-mu     (- v mu)
-        v-mu-sq  (* v-mu v-mu)
-        variance (* sigma sigma)]
+  (let [v-mu:sigma (/ (- v mu) sigma)]
     (* -0.5 (+ log-2pi
-               (Math/log variance)
-               (/ v-mu-sq variance)))))
+               (* 2 (Math/log sigma))
+               (* v-mu:sigma
+                  v-mu:sigma)))))
 
 (defn uniform
   "Returns the log-likelihood of the continuous [uniform

diff --git a/test/gen/distribution_test.cljc b/test/gen/distribution_test.cljc
@@ -174,15 +174,16 @@
                       (dist/logpdf (->normal 0.0 sigma) (- v)))
                 "Normal is symmetric about the mean")
 
-            (with-comparator (within 1e-12)
+            (with-comparator (within 1e-10)
               (is (ish? (dist/logpdf (->normal mu sigma) v)
                         (dist/logpdf (->normal (+ mu shift) sigma) (+ v shift)))
                   "shifting by the mean is a symmetry")))
 
   (testing "spot checks"
     (is (= -1.0439385332046727 (dist/logpdf (->normal 0 1) 0.5)))
     (is (= -1.643335713764618 (dist/logpdf (->normal 0 2) 0.5)))
-    (is (= -1.612085713764618 (dist/logpdf (->normal 0 2) 0)))))
+    (is (= -1.612085713764618 (dist/logpdf (->normal 0 2) 0)))
+    (is (= -8.996964098844152E20 (dist/logpdf (->normal 4241959036 0.1) 33978)))))
 
 (defn uniform-tests [->uniform]
   (checking "(log of the) Beta function is symmetrical"