Merge pull request #3164 from josef-pkt/bug_glm_residworking

BUG: fix GLM resid_working
statsmodels · Aug 17, 2016 · 1a09acf · 1a09acf
2 parents 87fd2d3 + 75b45d6
commit 1a09acf
Show file tree

Hide file tree

Showing 3 changed files with 19 additions and 9 deletions.
diff --git a/statsmodels/genmod/generalized_linear_model.py b/statsmodels/genmod/generalized_linear_model.py
@@ -1309,7 +1309,7 @@ def resid_pearson(self):
     @cache_readonly
     def resid_working(self):
         # Isn't self.resid_response is already adjusted by _n_trials?
-        val = (self.resid_response / self.family.link.deriv(self.mu))
+        val = (self.resid_response * self.family.link.deriv(self.mu))
         val *= self._n_trials
         return val
 

diff --git a/statsmodels/genmod/tests/test_glm.py b/statsmodels/genmod/tests/test_glm.py
@@ -66,10 +66,19 @@ def test_standard_errors(self):
 
     decimal_resids = DECIMAL_4
     def test_residuals(self):
+        # fix incorrect numbers in resid_working results
+        # residuals for Poisson are also tested in test_glm_weights.py
+        import copy
+        # new numpy would have copy method
+        resid2 = copy.copy(self.res2.resids)
+        resid2[:, 2] *= self.res1.family.link.deriv(self.res1.mu)**2
+
+        atol = 10**(-self.decimal_resids)
         resids = np.column_stack((self.res1.resid_pearson,
                 self.res1.resid_deviance, self.res1.resid_working,
                 self.res1.resid_anscombe, self.res1.resid_response))
-        assert_almost_equal(resids, self.res2.resids, self.decimal_resids)
+        assert_allclose(resids, resid2, rtol=1e-6, atol=atol)
+
 
     decimal_aic_R = DECIMAL_4
     def test_aic_R(self):
@@ -1423,13 +1432,13 @@ def test_resid(self):
                         np.concatenate((self.res2.resid_deviance[:17],
                                         [self.res2.resid_deviance[l2]])),
                         rtol=1e-5, atol=1e-5)
-        # Not working either...
-        if not isinstance(self, (TestTweedieLog1, TestTweedieLog15Fair)):
-            assert_allclose(np.concatenate((self.res1.resid_working[:17],
-                                            [self.res1.resid_working[l]])),
-                            np.concatenate((self.res2.resid_working[:17],
-                                            [self.res2.resid_working[l2]])),
-                            rtol=1e-5, atol=1e-5)
+
+        assert_allclose(np.concatenate((self.res1.resid_working[:17],
+                                        [self.res1.resid_working[l]])),
+                        np.concatenate((self.res2.resid_working[:17],
+                                        [self.res2.resid_working[l2]])),
+                        rtol=1e-5, atol=1e-5)
+
 
     def test_bse(self):
         assert_allclose(self.res1.bse, self.res2.bse, atol=1e-6, rtol=1e6)

diff --git a/statsmodels/genmod/tests/test_glm_weights.py b/statsmodels/genmod/tests/test_glm_weights.py
@@ -50,6 +50,7 @@ def test_residuals(self):
         assert_allclose(res1.resid_pearson, resid_all['resid_pearson'], atol= 1e-6, rtol=2e-6)
         assert_allclose(res1.resid_deviance, resid_all['resid_deviance'], atol= 1e-6, rtol=2e-6)
         assert_allclose(res1.resid_anscombe, resid_all['resid_anscombe'], atol= 1e-6, rtol=2e-6)
+        assert_allclose(res1.resid_working, resid_all['resid_working'], atol= 1e-6, rtol=2e-6)
 
 
 class TestGlmPoissonPlain(CheckWeight):