support strength=0 PYPs, final notation clean-up

7 files changed, 131 insertions, 112 deletions

diff --git a/gi/pf/align-lexonly-pyp.cc b/gi/pf/align-lexonly-pyp.cc
index 4ce7cf62..87f7f6b5 100644
--- a/gi/pf/align-lexonly-pyp.cc
+++ b/gi/pf/align-lexonly-pyp.cc
@@ -104,7 +104,7 @@ struct HierarchicalWordBase {
   }
 
   void Summary() const {
-    cerr << "NUMBER OF CUSTOMERS: " << r.num_customers() << "  (d=" << r.discount() << ",\\alpha=" << r.alpha() << ')' << endl;
+    cerr << "NUMBER OF CUSTOMERS: " << r.num_customers() << "  (d=" << r.discount() << ",s=" << r.strength() << ')' << endl;
     for (MFCR<vector<WordID> >::const_iterator it = r.begin(); it != r.end(); ++it)
       cerr << "   " << it->second.total_dish_count_ << " (on " << it->second.table_counts_.size() << " tables)" << TD::GetString(it->first) << endl;
   }
diff --git a/gi/pf/conditional_pseg.h b/gi/pf/conditional_pseg.h
index f9841cbf..86403d8d 100644
--- a/gi/pf/conditional_pseg.h
+++ b/gi/pf/conditional_pseg.h
@@ -22,7 +22,7 @@ struct MConditionalTranslationModel {
   void Summary() const {
     std::cerr << "Number of conditioning contexts: " << r.size() << std::endl;
     for (RuleModelHash::const_iterator it = r.begin(); it != r.end(); ++it) {
-      std::cerr << TD::GetString(it->first) << "   \t(d=" << it->second.discount() << ",\\alpha = " << it->second.alpha() << ") --------------------------" << std::endl;
+      std::cerr << TD::GetString(it->first) << "   \t(d=" << it->second.discount() << ",s=" << it->second.strength() << ") --------------------------" << std::endl;
       for (MFCR<TRule>::const_iterator i2 = it->second.begin(); i2 != it->second.end(); ++i2)
         std::cerr << "   " << -1 << '\t' << i2->first << std::endl;
     }
diff --git a/gi/pf/learn_cfg.cc b/gi/pf/learn_cfg.cc
index 5b748311..bf157828 100644
--- a/gi/pf/learn_cfg.cc
+++ b/gi/pf/learn_cfg.cc
@@ -183,9 +183,9 @@ struct HieroLMModel {
       nts[i].resample_hyperparameters(rng);
     if (kHIERARCHICAL_PRIOR) {
       q0.resample_hyperparameters(rng);
-      cerr << "[base d=" << q0.discount() << ", alpha=" << q0.alpha() << "]";
+      cerr << "[base d=" << q0.discount() << ", s=" << q0.strength() << "]";
     }
-    cerr << " d=" << nts[0].discount() << ", alpha=" << nts[0].alpha() << endl;
+    cerr << " d=" << nts[0].discount() << ", s=" << nts[0].strength() << endl;
   }
 
   const BaseRuleModel base;
diff --git a/gi/pf/pyp_lm.cc b/gi/pf/pyp_lm.cc
index e5c44c8b..7ebada13 100644
--- a/gi/pf/pyp_lm.cc
+++ b/gi/pf/pyp_lm.cc
@@ -78,14 +78,14 @@ template <unsigned N> struct PYPLM {
       backoff(vs, da, db, ss, sr),
       discount_a(da), discount_b(db),
       strength_s(ss), strength_r(sr),
-      d(0.8), alpha(1.0), lookup(N-1) {}
+      d(0.8), strength(1.0), lookup(N-1) {}
   void increment(WordID w, const vector<WordID>& context, MT19937* rng) {
     const double bo = backoff.prob(w, context);
     for (unsigned i = 0; i < N-1; ++i)
       lookup[i] = context[context.size() - 1 - i];
     typename unordered_map<vector<WordID>, CCRP<WordID>, boost::hash<vector<WordID> > >::iterator it = p.find(lookup);
     if (it == p.end())
-      it = p.insert(make_pair(lookup, CCRP<WordID>(d,alpha))).first;
+      it = p.insert(make_pair(lookup, CCRP<WordID>(d,strength))).first;
     if (it->second.increment(w, bo, rng))
       backoff.increment(w, context, rng);
   }
@@ -107,7 +107,7 @@ template <unsigned N> struct PYPLM {
   }
 
   double log_likelihood() const {
-    return log_likelihood(d, alpha) + backoff.log_likelihood();
+    return log_likelihood(d, strength) + backoff.log_likelihood();
   }
 
   double log_likelihood(const double& dd, const double& aa) const {
@@ -125,15 +125,15 @@ template <unsigned N> struct PYPLM {
     DiscountResampler(const PYPLM& m) : m_(m) {}
     const PYPLM& m_;
     double operator()(const double& proposed_discount) const {
-      return m_.log_likelihood(proposed_discount, m_.alpha);
+      return m_.log_likelihood(proposed_discount, m_.strength);
     }
   };
 
   struct AlphaResampler {
     AlphaResampler(const PYPLM& m) : m_(m) {}
     const PYPLM& m_;
-    double operator()(const double& proposed_alpha) const {
-      return m_.log_likelihood(m_.d, proposed_alpha);
+    double operator()(const double& proposed_strength) const {
+      return m_.log_likelihood(m_.d, proposed_strength);
     }
   };
 
@@ -141,25 +141,25 @@ template <unsigned N> struct PYPLM {
     DiscountResampler dr(*this);
     AlphaResampler ar(*this);
     for (int iter = 0; iter < nloop; ++iter) {
-      alpha = slice_sampler1d(ar, alpha, *rng, 0.0,
+      strength = slice_sampler1d(ar, strength, *rng, 0.0,
                               std::numeric_limits<double>::infinity(), 0.0, niterations, 100*niterations);
       d = slice_sampler1d(dr, d, *rng, std::numeric_limits<double>::min(),
                           1.0, 0.0, niterations, 100*niterations);
     }
-    alpha = slice_sampler1d(ar, alpha, *rng, 0.0,
+    strength = slice_sampler1d(ar, strength, *rng, 0.0,
                             std::numeric_limits<double>::infinity(), 0.0, niterations, 100*niterations);
     typename unordered_map<vector<WordID>, CCRP<WordID>, boost::hash<vector<WordID> > >::iterator it;
-    cerr << "PYPLM<" << N << ">(d=" << d << ",a=" << alpha << ") = " << log_likelihood(d, alpha) << endl;
+    cerr << "PYPLM<" << N << ">(d=" << d << ",a=" << strength << ") = " << log_likelihood(d, strength) << endl;
     for (it = p.begin(); it != p.end(); ++it) {
       it->second.set_discount(d);
-      it->second.set_alpha(alpha);
+      it->second.set_strength(strength);
     }
     backoff.resample_hyperparameters(rng, nloop, niterations);
   }
 
   PYPLM<N-1> backoff;
   double discount_a, discount_b, strength_s, strength_r;
-  double d, alpha;
+  double d, strength;
   mutable vector<WordID> lookup;  // thread-local
   unordered_map<vector<WordID>, CCRP<WordID>, boost::hash<vector<WordID> > > p;
 };
diff --git a/phrasinator/gibbs_train_plm.cc b/phrasinator/gibbs_train_plm.cc
index 54861dcb..3b99e1b6 100644
--- a/phrasinator/gibbs_train_plm.cc
+++ b/phrasinator/gibbs_train_plm.cc
@@ -252,7 +252,7 @@ struct UniphraseLM {
   void ResampleHyperparameters(MT19937* rng) {
     phrases_.resample_hyperparameters(rng);
     gen_.resample_hyperparameters(rng);
-    cerr << " d=" << phrases_.discount() << ",a=" << phrases_.alpha();
+    cerr << " d=" << phrases_.discount() << ",s=" << phrases_.strength();
   }
 
   CCRP<vector<int> > phrases_;
diff --git a/utils/ccrp.h b/utils/ccrp.h
index c883c027..5f9db7a6 100644
--- a/utils/ccrp.h
+++ b/utils/ccrp.h
@@ -18,27 +18,27 @@
 template <typename Dish, typename DishHash = boost::hash<Dish> >
 class CCRP {
  public:
-  CCRP(double disc, double alpha) :
+  CCRP(double disc, double strength) :
       num_tables_(),
       num_customers_(),
       discount_(disc),
-      alpha_(alpha),
-      discount_prior_alpha_(std::numeric_limits<double>::quiet_NaN()),
+      strength_(strength),
+      discount_prior_strength_(std::numeric_limits<double>::quiet_NaN()),
       discount_prior_beta_(std::numeric_limits<double>::quiet_NaN()),
-      alpha_prior_shape_(std::numeric_limits<double>::quiet_NaN()),
-      alpha_prior_rate_(std::numeric_limits<double>::quiet_NaN()) {
+      strength_prior_shape_(std::numeric_limits<double>::quiet_NaN()),
+      strength_prior_rate_(std::numeric_limits<double>::quiet_NaN()) {
     check_hyperparameters();
   }
 
-  CCRP(double d_alpha, double d_beta, double c_shape, double c_rate, double d = 0.9, double c = 1.0) :
+  CCRP(double d_strength, double d_beta, double c_shape, double c_rate, double d = 0.9, double c = 1.0) :
       num_tables_(),
       num_customers_(),
       discount_(d),
-      alpha_(c),
-      discount_prior_alpha_(d_alpha),
+      strength_(c),
+      discount_prior_strength_(d_strength),
       discount_prior_beta_(d_beta),
-      alpha_prior_shape_(c_shape),
-      alpha_prior_rate_(c_rate) {
+      strength_prior_shape_(c_shape),
+      strength_prior_rate_(c_rate) {
     check_hyperparameters();
   }
 
@@ -47,23 +47,23 @@ class CCRP {
       std::cerr << "Bad discount: " << discount_ << std::endl;
       abort();
     }
-    if (alpha_ <= -discount_) {
-      std::cerr << "Bad strength: " << alpha_ << " (discount=" << discount_ << ")" << std::endl;
+    if (strength_ <= -discount_) {
+      std::cerr << "Bad strength: " << strength_ << " (discount=" << discount_ << ")" << std::endl;
       abort();
     }
   }
 
   double discount() const { return discount_; }
-  double alpha() const { return alpha_; }
+  double strength() const { return strength_; }
   void set_discount(double d) { discount_ = d; check_hyperparameters(); }
-  void set_alpha(double a) { alpha_ = a; check_hyperparameters(); }
+  void set_strength(double a) { strength_ = a; check_hyperparameters(); }
 
   bool has_discount_prior() const {
-    return !std::isnan(discount_prior_alpha_);
+    return !std::isnan(discount_prior_strength_);
   }
 
-  bool has_alpha_prior() const {
-    return !std::isnan(alpha_prior_shape_);
+  bool has_strength_prior() const {
+    return !std::isnan(strength_prior_shape_);
   }
 
   void clear() {
@@ -97,7 +97,7 @@ class CCRP {
     DishLocations& loc = dish_locs_[dish];
     bool share_table = false;
     if (loc.total_dish_count_) {
-      const double p_empty = (alpha_ + num_tables_ * discount_) * p0;
+      const double p_empty = (strength_ + num_tables_ * discount_) * p0;
       const double p_share = (loc.total_dish_count_ - loc.table_counts_.size() * discount_);
       share_table = rng->SelectSample(p_empty, p_share);
     }
@@ -131,7 +131,7 @@ class CCRP {
     DishLocations& loc = dish_locs_[dish];
     bool share_table = false;
     if (loc.total_dish_count_) {
-      const T p_empty = T(alpha_ + num_tables_ * discount_) * p0;
+      const T p_empty = T(strength_ + num_tables_ * discount_) * p0;
       const T p_share = T(loc.total_dish_count_ - loc.table_counts_.size() * discount_);
       share_table = rng->SelectSample(p_empty, p_share);
     }
@@ -198,47 +198,47 @@ class CCRP {
 
   double prob(const Dish& dish, const double& p0) const {
     const typename std::tr1::unordered_map<Dish, DishLocations, DishHash>::const_iterator it = dish_locs_.find(dish);
-    const double r = num_tables_ * discount_ + alpha_;
+    const double r = num_tables_ * discount_ + strength_;
     if (it == dish_locs_.end()) {
-      return r * p0 / (num_customers_ + alpha_);
+      return r * p0 / (num_customers_ + strength_);
     } else {
       return (it->second.total_dish_count_ - discount_ * it->second.table_counts_.size() + r * p0) /
-               (num_customers_ + alpha_);
+               (num_customers_ + strength_);
     }
   }
 
   template <typename T>
   T probT(const Dish& dish, const T& p0) const {
     const typename std::tr1::unordered_map<Dish, DishLocations, DishHash>::const_iterator it = dish_locs_.find(dish);
-    const T r = T(num_tables_ * discount_ + alpha_);
+    const T r = T(num_tables_ * discount_ + strength_);
     if (it == dish_locs_.end()) {
-      return r * p0 / T(num_customers_ + alpha_);
+      return r * p0 / T(num_customers_ + strength_);
     } else {
       return (T(it->second.total_dish_count_ - discount_ * it->second.table_counts_.size()) + r * p0) /
-               T(num_customers_ + alpha_);
+               T(num_customers_ + strength_);
     }
   }
 
   double log_crp_prob() const {
-    return log_crp_prob(discount_, alpha_);
+    return log_crp_prob(discount_, strength_);
   }
 
   // taken from http://en.wikipedia.org/wiki/Chinese_restaurant_process
   // does not include P_0's
-  double log_crp_prob(const double& discount, const double& alpha) const {
+  double log_crp_prob(const double& discount, const double& strength) const {
     double lp = 0.0;
     if (has_discount_prior())
-      lp = Md::log_beta_density(discount, discount_prior_alpha_, discount_prior_beta_);
-    if (has_alpha_prior())
-      lp += Md::log_gamma_density(alpha + discount, alpha_prior_shape_, alpha_prior_rate_);
+      lp = Md::log_beta_density(discount, discount_prior_strength_, discount_prior_beta_);
+    if (has_strength_prior())
+      lp += Md::log_gamma_density(strength + discount, strength_prior_shape_, strength_prior_rate_);
     assert(lp <= 0.0);
     if (num_customers_) {
       if (discount > 0.0) {
         const double r = lgamma(1.0 - discount);
-        if (alpha)
-          lp += lgamma(alpha) - lgamma(alpha / discount);
-        lp += - lgamma(alpha + num_customers_)
-             + num_tables_ * log(discount) + lgamma(alpha / discount + num_tables_);
+        if (strength)
+          lp += lgamma(strength) - lgamma(strength / discount);
+        lp += - lgamma(strength + num_customers_)
+             + num_tables_ * log(discount) + lgamma(strength / discount + num_tables_);
         assert(std::isfinite(lp));
         for (typename std::tr1::unordered_map<Dish, DishLocations, DishHash>::const_iterator it = dish_locs_.begin();
              it != dish_locs_.end(); ++it) {
@@ -247,8 +247,16 @@ class CCRP {
             lp += lgamma(*ti - discount) - r;
           }
         }
+      } else if (!discount) { // discount == 0.0
+        lp += lgamma(strength) + num_tables_ * log(strength) - lgamma(strength + num_tables_);
+        assert(std::isfinite(lp));
+        for (typename std::tr1::unordered_map<Dish, DishLocations, DishHash>::const_iterator it = dish_locs_.begin();
+             it != dish_locs_.end(); ++it) {
+          const DishLocations& cur = it->second;
+          lp += lgamma(cur.table_counts_.size());
+        }
       } else {
-        assert(!"not implemented yet");
+        assert(!"discount less than 0 detected!");
       }
     }
     assert(std::isfinite(lp));
@@ -256,22 +264,22 @@ class CCRP {
   }
 
   void resample_hyperparameters(MT19937* rng, const unsigned nloop = 5, const unsigned niterations = 10) {
-    assert(has_discount_prior() || has_alpha_prior());
+    assert(has_discount_prior() || has_strength_prior());
     DiscountResampler dr(*this);
     StrengthResampler sr(*this);
     for (int iter = 0; iter < nloop; ++iter) {
-      if (has_alpha_prior()) {
-        alpha_ = slice_sampler1d(sr, alpha_, *rng, -discount_,
+      if (has_strength_prior()) {
+        strength_ = slice_sampler1d(sr, strength_, *rng, -discount_,
                                std::numeric_limits<double>::infinity(), 0.0, niterations, 100*niterations);
       }
       if (has_discount_prior()) {
         double min_discount = std::numeric_limits<double>::min();
-        if (alpha_ < 0.0) min_discount = -alpha_;
+        if (strength_ < 0.0) min_discount = -strength_;
         discount_ = slice_sampler1d(dr, discount_, *rng, min_discount,
                                1.0, 0.0, niterations, 100*niterations);
       }
     }
-    alpha_ = slice_sampler1d(sr, alpha_, *rng, -discount_,
+    strength_ = slice_sampler1d(sr, strength_, *rng, -discount_,
                              std::numeric_limits<double>::infinity(), 0.0, niterations, 100*niterations);
   }
 
@@ -279,15 +287,15 @@ class CCRP {
     DiscountResampler(const CCRP& crp) : crp_(crp) {}
     const CCRP& crp_;
     double operator()(const double& proposed_discount) const {
-      return crp_.log_crp_prob(proposed_discount, crp_.alpha_);
+      return crp_.log_crp_prob(proposed_discount, crp_.strength_);
     }
   };
 
   struct StrengthResampler {
     StrengthResampler(const CCRP& crp) : crp_(crp) {}
     const CCRP& crp_;
-    double operator()(const double& proposed_alpha) const {
-      return crp_.log_crp_prob(crp_.discount_, proposed_alpha);
+    double operator()(const double& proposed_strength) const {
+      return crp_.log_crp_prob(crp_.discount_, proposed_strength);
     }
   };
 
@@ -299,7 +307,7 @@ class CCRP {
   };
 
   void Print(std::ostream* out) const {
-    std::cerr << "PYP(d=" << discount_ << ",c=" << alpha_ << ") customers=" << num_customers_ << std::endl;
+    std::cerr << "PYP(d=" << discount_ << ",c=" << strength_ << ") customers=" << num_customers_ << std::endl;
     for (typename std::tr1::unordered_map<Dish, DishLocations, DishHash>::const_iterator it = dish_locs_.begin();
          it != dish_locs_.end(); ++it) {
       (*out) << it->first << " (" << it->second.total_dish_count_ << " on " << it->second.table_counts_.size() << " tables): ";
@@ -324,15 +332,15 @@ class CCRP {
   std::tr1::unordered_map<Dish, DishLocations, DishHash> dish_locs_;
 
   double discount_;
-  double alpha_;
+  double strength_;
 
   // optional beta prior on discount_ (NaN if no prior)
-  double discount_prior_alpha_;
+  double discount_prior_strength_;
   double discount_prior_beta_;
 
-  // optional gamma prior on alpha_ (NaN if no prior)
-  double alpha_prior_shape_;
-  double alpha_prior_rate_;
+  // optional gamma prior on strength_ (NaN if no prior)
+  double strength_prior_shape_;
+  double strength_prior_rate_;
 };
 
 template <typename T,typename H>
diff --git a/utils/mfcr.h b/utils/mfcr.h
index df988f51..aeaf599d 100644
--- a/utils/mfcr.h
+++ b/utils/mfcr.h
@@ -39,37 +39,37 @@ template <typename Dish, typename DishHash = boost::hash<Dish> >
 class MFCR {
  public:
 
-  MFCR(unsigned num_floors, double d, double alpha) :
+  MFCR(unsigned num_floors, double d, double strength) :
     num_floors_(num_floors),
     num_tables_(),
     num_customers_(),
     discount_(d),
-    alpha_(alpha),
-    discount_prior_alpha_(std::numeric_limits<double>::quiet_NaN()),
+    strength_(strength),
+    discount_prior_strength_(std::numeric_limits<double>::quiet_NaN()),
     discount_prior_beta_(std::numeric_limits<double>::quiet_NaN()),
-    alpha_prior_shape_(std::numeric_limits<double>::quiet_NaN()),
-    alpha_prior_rate_(std::numeric_limits<double>::quiet_NaN()) {}
+    strength_prior_shape_(std::numeric_limits<double>::quiet_NaN()),
+    strength_prior_rate_(std::numeric_limits<double>::quiet_NaN()) {}
 
-  MFCR(unsigned num_floors, double discount_alpha, double discount_beta, double alpha_shape, double alpha_rate, double d = 0.9, double alpha = 10.0) :
+  MFCR(unsigned num_floors, double discount_strength, double discount_beta, double strength_shape, double strength_rate, double d = 0.9, double strength = 10.0) :
     num_floors_(num_floors),
     num_tables_(),
     num_customers_(),
     discount_(d),
-    alpha_(alpha),
-    discount_prior_alpha_(discount_alpha),
+    strength_(strength),
+    discount_prior_strength_(discount_strength),
     discount_prior_beta_(discount_beta),
-    alpha_prior_shape_(alpha_shape),
-    alpha_prior_rate_(alpha_rate) {}
+    strength_prior_shape_(strength_shape),
+    strength_prior_rate_(strength_rate) {}
 
   double discount() const { return discount_; }
-  double alpha() const { return alpha_; }
+  double strength() const { return strength_; }
 
   bool has_discount_prior() const {
-    return !std::isnan(discount_prior_alpha_);
+    return !std::isnan(discount_prior_strength_);
   }
 
-  bool has_alpha_prior() const {
-    return !std::isnan(alpha_prior_shape_);
+  bool has_strength_prior() const {
+    return !std::isnan(strength_prior_shape_);
   }
 
   void clear() {
@@ -122,7 +122,7 @@ class MFCR {
     int floor = -1;
     bool share_table = false;
     if (loc.total_dish_count_) {
-      const double p_empty = (alpha_ + num_tables_ * discount_) * marg_p0;
+      const double p_empty = (strength_ + num_tables_ * discount_) * marg_p0;
       const double p_share = (loc.total_dish_count_ - loc.table_counts_.size() * discount_);
       share_table = rng->SelectSample(p_empty, p_share);
     }
@@ -206,44 +206,53 @@ class MFCR {
     const double marg_p0 = std::inner_product(p0s.begin(), p0s.end(), lambdas.begin(), 0.0);
     assert(marg_p0 <= 1.0);
     const typename std::tr1::unordered_map<Dish, DishLocations, DishHash>::const_iterator it = dish_locs_.find(dish);
-    const double r = num_tables_ * discount_ + alpha_;
+    const double r = num_tables_ * discount_ + strength_;
     if (it == dish_locs_.end()) {
-      return r * marg_p0 / (num_customers_ + alpha_);
+      return r * marg_p0 / (num_customers_ + strength_);
     } else {
       return (it->second.total_dish_count_ - discount_ * it->second.table_counts_.size() + r * marg_p0) /
-               (num_customers_ + alpha_);
+               (num_customers_ + strength_);
     }
   }
 
   double log_crp_prob() const {
-    return log_crp_prob(discount_, alpha_);
+    return log_crp_prob(discount_, strength_);
   }
 
   // taken from http://en.wikipedia.org/wiki/Chinese_restaurant_process
   // does not include draws from G_w's
-  double log_crp_prob(const double& d, const double& alpha) const {
+  double log_crp_prob(const double& discount, const double& strength) const {
     double lp = 0.0;
     if (has_discount_prior())
-      lp = Md::log_beta_density(d, discount_prior_alpha_, discount_prior_beta_);
-    if (has_alpha_prior())
-      lp += Md::log_gamma_density(alpha, alpha_prior_shape_, alpha_prior_rate_);
+      lp = Md::log_beta_density(discount, discount_prior_strength_, discount_prior_beta_);
+    if (has_strength_prior())
+      lp += Md::log_gamma_density(strength + discount, strength_prior_shape_, strength_prior_rate_);
     assert(lp <= 0.0);
     if (num_customers_) {
-      if (d > 0.0) {
-        const double r = lgamma(1.0 - d);
-        lp += lgamma(alpha) - lgamma(alpha + num_customers_)
-             + num_tables_ * log(d) + lgamma(alpha / d + num_tables_)
-             - lgamma(alpha / d);
+      if (discount > 0.0) {
+        const double r = lgamma(1.0 - discount);
+        if (strength)
+          lp += lgamma(strength) - lgamma(strength / discount);
+        lp += - lgamma(strength + num_customers_)
+             + num_tables_ * log(discount) + lgamma(strength / discount + num_tables_);
         assert(std::isfinite(lp));
         for (typename std::tr1::unordered_map<Dish, DishLocations, DishHash>::const_iterator it = dish_locs_.begin();
              it != dish_locs_.end(); ++it) {
           const DishLocations& cur = it->second;
           for (std::list<TableCount>::const_iterator ti = cur.table_counts_.begin(); ti != cur.table_counts_.end(); ++ti) {
-            lp += lgamma(ti->count - d) - r;
+            lp += lgamma(ti->count - discount) - r;
           }
         }
+      } else if (!discount) { // discount == 0.0
+        lp += lgamma(strength) + num_tables_ * log(strength) - lgamma(strength + num_tables_);
+        assert(std::isfinite(lp));
+        for (typename std::tr1::unordered_map<Dish, DishLocations, DishHash>::const_iterator it = dish_locs_.begin();
+             it != dish_locs_.end(); ++it) {
+          const DishLocations& cur = it->second;
+          lp += lgamma(cur.table_counts_.size());
+        }
       } else {
-        assert(!"not implemented yet");
+        assert(!"discount less than 0 detected!");
       }
     }
     assert(std::isfinite(lp));
@@ -251,20 +260,22 @@ class MFCR {
   }
 
   void resample_hyperparameters(MT19937* rng, const unsigned nloop = 5, const unsigned niterations = 10) {
-    assert(has_discount_prior() || has_alpha_prior());
+    assert(has_discount_prior() || has_strength_prior());
     DiscountResampler dr(*this);
-    ConcentrationResampler cr(*this);
+    StrengthResampler sr(*this);
     for (int iter = 0; iter < nloop; ++iter) {
-      if (has_alpha_prior()) {
-        alpha_ = slice_sampler1d(cr, alpha_, *rng, 0.0,
+      if (has_strength_prior()) {
+        strength_ = slice_sampler1d(sr, strength_, *rng, -discount_,
                                std::numeric_limits<double>::infinity(), 0.0, niterations, 100*niterations);
       }
       if (has_discount_prior()) {
-        discount_ = slice_sampler1d(dr, discount_, *rng, std::numeric_limits<double>::min(),
+        double min_discount = std::numeric_limits<double>::min();
+        if (strength_ < 0.0) min_discount = -strength_;
+        discount_ = slice_sampler1d(dr, discount_, *rng, min_discount,
                                1.0, 0.0, niterations, 100*niterations);
       }
     }
-    alpha_ = slice_sampler1d(cr, alpha_, *rng, 0.0,
+    strength_ = slice_sampler1d(sr, strength_, *rng, -discount_,
                              std::numeric_limits<double>::infinity(), 0.0, niterations, 100*niterations);
   }
 
@@ -272,15 +283,15 @@ class MFCR {
     DiscountResampler(const MFCR& crp) : crp_(crp) {}
     const MFCR& crp_;
     double operator()(const double& proposed_d) const {
-      return crp_.log_crp_prob(proposed_d, crp_.alpha_);
+      return crp_.log_crp_prob(proposed_d, crp_.strength_);
     }
   };
 
-  struct ConcentrationResampler {
-    ConcentrationResampler(const MFCR& crp) : crp_(crp) {}
+  struct StrengthResampler {
+    StrengthResampler(const MFCR& crp) : crp_(crp) {}
     const MFCR& crp_;
-    double operator()(const double& proposediscount_alpha) const {
-      return crp_.log_crp_prob(crp_.discount_, proposediscount_alpha);
+    double operator()(const double& proposediscount_strength) const {
+      return crp_.log_crp_prob(crp_.discount_, proposediscount_strength);
     }
   };
 
@@ -292,7 +303,7 @@ class MFCR {
   };
 
   void Print(std::ostream* out) const {
-    (*out) << "MFCR(d=" << discount_ << ",alpha=" << alpha_ << ") customers=" << num_customers_ << std::endl;
+    (*out) << "MFCR(d=" << discount_ << ",strength=" << strength_ << ") customers=" << num_customers_ << std::endl;
     for (typename std::tr1::unordered_map<Dish, DishLocations, DishHash>::const_iterator it = dish_locs_.begin();
          it != dish_locs_.end(); ++it) {
       (*out) << it->first << " (" << it->second.total_dish_count_ << " on " << it->second.table_counts_.size() << " tables): ";
@@ -318,15 +329,15 @@ class MFCR {
   std::tr1::unordered_map<Dish, DishLocations, DishHash> dish_locs_;
 
   double discount_;
-  double alpha_;
+  double strength_;
 
   // optional beta prior on discount_ (NaN if no prior)
-  double discount_prior_alpha_;
+  double discount_prior_strength_;
   double discount_prior_beta_;
 
-  // optional gamma prior on alpha_ (NaN if no prior)
-  double alpha_prior_shape_;
-  double alpha_prior_rate_;
+  // optional gamma prior on strength_ (NaN if no prior)
+  double strength_prior_shape_;
+  double strength_prior_rate_;
 };
 
 template <typename T,typename H>