Fix Binomial GLM

In glm module binomial familiy is actually Bernoulli. This PR fix the Binomial likelihood, with the flexibility of specifying the `n`. 
Default n set to 1 for backward compatibility.
Also added test.

Author: Junpeng Lao

RELEASE-NOTES.md

@@ -16,6 +16,7 @@
 - Plots of discrete distributions in the docstrings
 - Add logitnormal distribution
 - Densityplot: add support for discrete variables
+- Fix the Binomial likelihood in `.glm.families.Binomial`, with the flexibility of specifying the `n`. 
 
 ### Fixes
 

pymc3/glm/families.py

@@ -1,4 +1,5 @@
 import numbers
+import numpy as np
 from copy import copy
 
 import theano.tensor as tt
@@ -51,7 +52,7 @@ def _get_priors(self, model=None, name=''):
         model = modelcontext(model)
         priors = {}
         for key, val in self.priors.items():
-            if isinstance(val, numbers.Number):
+            if isinstance(val, (numbers.Number, np.ndarray, np.generic)):
                 priors[key] = val
             else:
                 priors[key] = model.Var('{}{}'.format(name, key), val)
@@ -99,8 +100,9 @@ class Normal(Family):
 
 class Binomial(Family):
     link = logit
-    likelihood = pm_dists.Bernoulli
+    likelihood = pm_dists.Binomial
     parent = 'p'
+    priors = {'n': 1}
 
 
 class Poisson(Family):

pymc3/tests/test_glm.py

@@ -26,9 +26,13 @@ def setup_class(cls):
 
         x_logistic, y_logistic = generate_data(cls.intercept, cls.slope, size=3000)
         y_logistic = 1 / (1 + np.exp(-y_logistic))
-        bern_trials = [np.random.binomial(1, i) for i in y_logistic]
+        bern_trials = np.random.binomial(1, y_logistic)
         cls.data_logistic = dict(x=x_logistic, y=bern_trials)
 
+        n_trials = np.random.randint(1, 20, size=y_logistic.shape)
+        binom_trials = np.random.binomial(n_trials, y_logistic)
+        cls.data_logistic2 = dict(x=x_logistic, y=binom_trials, n=n_trials)
+
     def test_linear_component(self):
         with Model() as model:
             lm = LinearComponent.from_formula('y ~ x', self.data_linear)
@@ -65,6 +69,16 @@ def test_glm_link_func(self):
             assert round(abs(np.mean(trace['Intercept'])-self.intercept), 1) == 0
             assert round(abs(np.mean(trace['x'])-self.slope), 1) == 0
 
+    def test_glm_link_func2(self):
+        with Model() as model:
+            GLM.from_formula('y ~ x', self.data_logistic2,
+                    family=families.Binomial(priors={'n': self.data_logistic2['n']}))
+            trace = sample(1000, progressbar=False,
+                           random_seed=self.random_seed)
+
+            assert round(abs(np.mean(trace['Intercept'])-self.intercept), 1) == 0
+            assert round(abs(np.mean(trace['x'])-self.slope), 1) == 0
+
     def test_more_than_one_glm_is_ok(self):
         with Model():
             GLM.from_formula('y ~ x', self.data_logistic,