GLM tests probably pass (pymc-devs#1352)

Author: ColCarroll

pymc3/tests/helpers.py

@@ -3,7 +3,11 @@
 
 
 class SeededTest(unittest.TestCase):
-    random_seed = 20160907
+    random_seed = 20160911
+
+    @classmethod
+    def setUpClass(cls):
+        nr.seed(cls.random_seed)
 
     def setUp(self):
         nr.seed(self.random_seed)

pymc3/tests/test_glm.py

@@ -1,81 +1,61 @@
-import unittest
-from nose import SkipTest
 import numpy as np
-import sys
-try:
-    import statsmodels.api as sm
-except ImportError:
-    raise SkipTest("Test requires statsmodels.")
 
-from pymc3.examples import glm_linear, glm_robust
+from .helpers import SeededTest
+from pymc3 import glm, Model, Uniform, Normal, find_MAP, Slice, sample
 
 
-np.random.seed(1)
 # Generate data
-true_intercept = 0
-true_slope = 3
-
-
-def generate_data(size=700):
+def generate_data(intercept, slope, size=700):
     x = np.linspace(-1, 1, size)
-    y = true_intercept + x * true_slope
+    y = intercept + x * slope
     return x, y
 
-true_sd = .05
-x_linear, y_linear = generate_data(size=1000)
-y_linear += np.random.normal(size=1000, scale=true_sd)
-data_linear = dict(x=x_linear, y=y_linear)
 
-x_logistic, y_logistic = generate_data(size=3000)
-y_logistic = 1 / (1 + np.exp(-y_logistic))
-bern_trials = [np.random.binomial(1, i) for i in y_logistic]
-data_logistic = dict(x=x_logistic, y=bern_trials)
+class TestGLM(SeededTest):
+    @classmethod
+    def setUpClass(cls):
+        super(TestGLM, cls).setUpClass()
+        cls.intercept = 1
+        cls.slope = 3
+        cls.sd = .05
+        x_linear, cls.y_linear = generate_data(cls.intercept, cls.slope, size=1000)
+        cls.y_linear += np.random.normal(size=1000, scale=cls.sd)
+        cls.data_linear = dict(x=x_linear, y=cls.y_linear)
 
+        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]
+        cls.data_logistic = dict(x=x_logistic, y=bern_trials)
 
-class TestGLM(unittest.TestCase):
-
-    @unittest.skip("Fails only on travis. Investigate")
     def test_linear_component(self):
         with Model() as model:
-            y_est, coeffs = glm.linear_component('y ~ x', data_linear)
-            for coeff, true_val in zip(coeffs, [true_intercept, true_slope]):
-                self.assertAlmostEqual(coeff.tag.test_value, true_val, 1)
+            y_est, _ = glm.linear_component('y ~ x', self.data_linear)
             sigma = Uniform('sigma', 0, 20)
-            y_obs = Normal('y_obs', mu=y_est, sd=sigma, observed=y_linear)
+            Normal('y_obs', mu=y_est, sd=sigma, observed=self.y_linear)
             start = find_MAP(vars=[sigma])
             step = Slice(model.vars)
-            trace = sample(2000, step, start, progressbar=False)
+            trace = sample(500, step, start, progressbar=False, random_seed=self.random_seed)
 
-            self.assertAlmostEqual(
-                np.mean(trace['Intercept']), true_intercept, 1)
-            self.assertAlmostEqual(np.mean(trace['x']), true_slope, 1)
-            self.assertAlmostEqual(np.mean(trace['sigma']), true_sd, 1)
+            self.assertAlmostEqual(np.mean(trace['Intercept']), self.intercept, 1)
+            self.assertAlmostEqual(np.mean(trace['x']), self.slope, 1)
+            self.assertAlmostEqual(np.mean(trace['sigma']), self.sd, 1)
 
-    @unittest.skip("Fails only on travis. Investigate")
     def test_glm(self):
         with Model() as model:
-            vars = glm.glm('y ~ x', data_linear)
-            for coeff, true_val in zip(vars[1:], [true_intercept, true_slope, true_sd]):
-                self.assertAlmostEqual(coeff.tag.test_value, true_val, 1)
+            glm.glm('y ~ x', self.data_linear)
             step = Slice(model.vars)
-            trace = sample(2000, step, progressbar=False)
+            trace = sample(500, step, progressbar=False, random_seed=self.random_seed)
 
-            self.assertAlmostEqual(
-                np.mean(trace['Intercept']), true_intercept, 1)
-            self.assertAlmostEqual(np.mean(trace['x']), true_slope, 1)
-            self.assertAlmostEqual(np.mean(trace['sigma']), true_sd, 1)
+            self.assertAlmostEqual(np.mean(trace['Intercept']), self.intercept, 1)
+            self.assertAlmostEqual(np.mean(trace['x']), self.slope, 1)
+            self.assertAlmostEqual(np.mean(trace['sd']), self.sd, 1)
 
-    @unittest.skip("Was an error, then a fail, now a skip.")
     def test_glm_link_func(self):
         with Model() as model:
-            vars = glm.glm('y ~ x', data_logistic,
-                           family=glm.families.Binomial(link=glm.families.logit))
-
-            for coeff, true_val in zip(vars[1:], [true_intercept, true_slope]):
-                self.assertAlmostEqual(coeff.tag.test_value, true_val, 0)
+            glm.glm('y ~ x', self.data_logistic,
+                    family=glm.families.Binomial(link=glm.families.logit))
             step = Slice(model.vars)
-            trace = sample(2000, step, progressbar=False)
+            trace = sample(1000, step, progressbar=False, random_seed=self.random_seed)
 
-            self.assertAlmostEqual(
-                np.mean(trace['Intercept']), true_intercept, 1)
-            self.assertAlmostEqual(np.mean(trace['x']), true_slope, 0)
+            self.assertAlmostEqual(np.mean(trace['Intercept']), self.intercept, 1)
+            self.assertAlmostEqual(np.mean(trace['x']), self.slope, 1)