From 89788800abf14b12a904d5ec28fb22ef944b828c Mon Sep 17 00:00:00 2001
From: martin_cerny <cerny.m@gmail.com>
Date: Thu, 16 Jan 2020 18:39:49 +0100
Subject: [PATCH] Quick fix for #1495 and the worst part of #1496

---
 .../prim/prob/neg_binomial_2_log_lpmf.hpp     | 31 ++++++++++-------
 stan/math/prim/prob/neg_binomial_2_lpmf.hpp   | 34 ++++++++++---------
 2 files changed, 37 insertions(+), 28 deletions(-)

diff --git a/stan/math/prim/prob/neg_binomial_2_log_lpmf.hpp b/stan/math/prim/prob/neg_binomial_2_log_lpmf.hpp
index b0a0eb5450a..91c069bb3c0 100644
--- a/stan/math/prim/prob/neg_binomial_2_log_lpmf.hpp
+++ b/stan/math/prim/prob/neg_binomial_2_log_lpmf.hpp
@@ -9,6 +9,7 @@
 #include <stan/math/prim/fun/lgamma.hpp>
 #include <stan/math/prim/fun/log_sum_exp.hpp>
 #include <stan/math/prim/fun/value_of.hpp>
+#include <stan/math/prim/prob/poisson_log_lpmf.hpp>
 #include <cmath>
 
 namespace stan {
@@ -84,19 +85,25 @@ return_type_t<T_log_location, T_precision> neg_binomial_2_log_lpmf(
   }
 
   for (size_t i = 0; i < max_size_seq_view; i++) {
-    if (include_summand<propto>::value) {
-      logp -= lgamma(n_vec[i] + 1.0);
+    if (phi__[i] > 1e10) {
+      // TODO(martinmodrak) This is wrong, but shouldn't brake most models.
+      // Will be adressed better in PR #1497
+      logp += poisson_log_lpmf(n_vec[i], eta__[i]);
+    } else {
+      if (include_summand<propto>::value) {
+        logp -= lgamma(n_vec[i] + 1.0);
+      }
+      if (include_summand<propto, T_precision>::value) {
+        logp += multiply_log(phi__[i], phi__[i]) - lgamma(phi__[i]);
+      }
+      if (include_summand<propto, T_log_location>::value) {
+        logp += n_vec[i] * eta__[i];
+      }
+      if (include_summand<propto, T_precision>::value) {
+        logp += lgamma(n_plus_phi[i]);
+      }
+      logp -= (n_plus_phi[i]) * logsumexp_eta_logphi[i];
     }
-    if (include_summand<propto, T_precision>::value) {
-      logp += multiply_log(phi__[i], phi__[i]) - lgamma(phi__[i]);
-    }
-    if (include_summand<propto, T_log_location>::value) {
-      logp += n_vec[i] * eta__[i];
-    }
-    if (include_summand<propto, T_precision>::value) {
-      logp += lgamma(n_plus_phi[i]);
-    }
-    logp -= (n_plus_phi[i]) * logsumexp_eta_logphi[i];
 
     if (!is_constant_all<T_log_location>::value) {
       ops_partials.edge1_.partials_[i]
diff --git a/stan/math/prim/prob/neg_binomial_2_lpmf.hpp b/stan/math/prim/prob/neg_binomial_2_lpmf.hpp
index 435ebeff2e3..e88d12423c8 100644
--- a/stan/math/prim/prob/neg_binomial_2_lpmf.hpp
+++ b/stan/math/prim/prob/neg_binomial_2_lpmf.hpp
@@ -79,23 +79,25 @@ return_type_t<T_location, T_precision> neg_binomial_2_lpmf(
   }
 
   for (size_t i = 0; i < max_size_seq_view; i++) {
-    if (include_summand<propto>::value) {
-      logp -= lgamma(n_vec[i] + 1.0);
-    }
-    if (include_summand<propto, T_precision>::value) {
-      logp += multiply_log(phi__[i], phi__[i]) - lgamma(phi__[i]);
-    }
-    if (include_summand<propto, T_location>::value) {
-      logp += multiply_log(n_vec[i], mu__[i]);
-    }
-    if (include_summand<propto, T_precision>::value) {
-      logp += lgamma(n_plus_phi[i]);
-    }
-    logp -= (n_plus_phi[i]) * log_mu_plus_phi[i];
-
     // if phi is large we probably overflow, defer to Poisson:
-    if (phi__[i] > 1e5) {
-      logp = poisson_lpmf(n_vec[i], mu__[i]);
+    if (phi__[i] > 1e10) {
+      // TODO(martinmodrak) This is wrong, but shouldn't brake most models.
+      // Will be adressed better in PR #1497
+      logp += poisson_lpmf(n_vec[i], mu__[i]);
+    } else {
+      if (include_summand<propto>::value) {
+        logp -= lgamma(n_vec[i] + 1.0);
+      }
+      if (include_summand<propto, T_precision>::value) {
+        logp += multiply_log(phi__[i], phi__[i]) - lgamma(phi__[i]);
+      }
+      if (include_summand<propto, T_location>::value) {
+        logp += multiply_log(n_vec[i], mu__[i]);
+      }
+      if (include_summand<propto, T_precision>::value) {
+        logp += lgamma(n_plus_phi[i]);
+      }
+      logp -= (n_plus_phi[i]) * log_mu_plus_phi[i];
     }
 
     if (!is_constant_all<T_location>::value) {