BUG/REF/DOC: more internal consistency fixes, update docstrings

Author: josef-pkt

statsmodels/robust/norms.py

@@ -406,16 +406,16 @@ def rho(self, z):
         Returns
         -------
         rho : ndarray
-            rho(z) = a*(1-cos(z/a))     for \|z\| <= a*pi
+            rho(z) = a**2 *(1-cos(z/a))     for \|z\| <= a*pi
 
-            rho(z) = 2*a                for \|z\| > a*pi
+            rho(z) = 2*a                    for \|z\| > a*pi
         """
 
         a = self.a
         z = np.asarray(z)
         test = self._subset(z)
-        return (test * a**2 * (0 - np.cos(z / a)) +
-                (1 - test) * 2 * a)
+        return (test * a**2 * (1 - np.cos(z / a)) +
+                (1 - test) * a**2 * 2)
 
     def psi(self, z):
         r"""
@@ -431,7 +431,7 @@ def psi(self, z):
         Returns
         -------
         psi : ndarray
-            psi(z) = sin(z/a)       for \|z\| <= a*pi
+            psi(z) = a * sin(z/a)   for \|z\| <= a*pi
 
             psi(z) = 0              for \|z\| > a*pi
         """
@@ -455,9 +455,9 @@ def weights(self, z):
         Returns
         -------
         weights : ndarray
-            weights(z) = sin(z/a)/(z/a)     for \|z\| <= a*pi
+            weights(z) = sin(z/a) / (z/a)     for \|z\| <= a*pi
 
-            weights(z) = 0                  for \|z\| > a*pi
+            weights(z) = 0                    for \|z\| > a*pi
         """
         a = self.a
         z = np.asarray(z)
@@ -527,12 +527,12 @@ def rho(self, z):
         rho : ndarray
             rho(z) = (1/2.)*z**2    for \|z\| <= c
 
-            rho(z) = 0              for \|z\| > c
+            rho(z) = (1/2.)*c**2              for \|z\| > c
         """
 
         z = np.asarray(z)
         test = self._subset(z)
-        return test * z**2 * 0.5
+        return test * z**2 * 0.5 + (1 - test) * self.c**2 * 0.5
 
     def psi(self, z):
         r"""
@@ -635,13 +635,13 @@ def rho(self, z):
         Returns
         -------
         rho : ndarray
-            rho(z) = (1/2.)*z**2                    for \|z\| <= a
+            rho(z) = z**2 / 2                     for \|z\| <= a
 
-            rho(z) = a*\|z\| - 1/2.*a**2              for a < \|z\| <= b
+            rho(z) = a*\|z\| - 1/2.*a**2               for a < \|z\| <= b
 
-            rho(z) = a*(c*\|z\|-(1/2.)*z**2)/(c-b)    for b < \|z\| <= c
+            rho(z) = a*(c - \|z\|)**2 / (c - b) / 2    for b < \|z\| <= c
 
-            rho(z) = a*(b + c - a)                  for \|z\| > c
+            rho(z) = a*(b + c - a) / 2                 for \|z\| > c
         """
         a, b, c = self.a, self.b, self.c
 
@@ -654,7 +654,8 @@ def rho(self, z):
         v = np.zeros(z.shape, dtype=dt)
         z = np.abs(z)
         v[t1] = z[t1]**2 * 0.5
-        v[t2] = (a * (z[t2] - a) + a**2 * 0.5)
+        # v[t2] = (a * (z[t2] - a) + a**2 * 0.5)
+        v[t2] = (a * z[t2] - a**2 * 0.5)
         v[t3] = a * (c - z[t3])**2 / (c - b) * (-0.5)
         v[t34] += a * (b + c - a) * 0.5
 
@@ -694,14 +695,11 @@ def psi(self, z):
         dt = np.promote_types(z.dtype, "float")
         v = np.zeros(z.shape, dtype=dt)
         s = np.sign(z)
-        z = np.abs(z)
+        za = np.abs(z)
 
-        v[t1] = z[t1] * s[t1]
+        v[t1] = z[t1]
         v[t2] = a * s[t2]
-        v[t3] = a * s[t3] * (c - z[t3]) / (c - b)
-        # v = (t1 * z*s +
-        #          t2 * a*s +
-        #          t3 * a*s * (c - z) / (c - b))
+        v[t3] = a * s[t3] * (c - za[t3]) / (c - b)
 
         if z_isscalar:
             v = v[0]
@@ -721,13 +719,13 @@ def weights(self, z):
         Returns
         -------
         weights : ndarray
-            weights(z) = 1                            for \|z\| <= a
+            weights(z) = 1                                for \|z\| <= a
 
-            weights(z) = a/\|z\|                        for a < \|z\| <= b
+            weights(z) = a/\|z\|                          for a < \|z\| <= b
 
             weights(z) = a*(c - \|z\|)/(\|z\|*(c-b))      for b < \|z\| <= c
 
-            weights(z) = 0                            for \|z\| > c
+            weights(z) = 0                                for \|z\| > c
         """
         a, b, c = self.a, self.b, self.c
 

statsmodels/robust/tests/test_norms.py

@@ -79,10 +79,15 @@ def test_norms_consistent(case):
     x = np.array([-9, -6, -2, -1, 0, 1, 2, 6, 9], dtype=float)
 
     weights = norm.weights(x)
-    # rho = norm.rho(x)  # not used
+    rho = norm.rho(x)  # not used
     psi = norm.psi(x)
     psi_deriv = norm.psi_deriv(x)
 
+    # check location and u-shape of rho
+    assert rho[4] == 0
+    assert np.all(np.diff(rho[4:]) >= 0)
+    assert np.all(np.diff(rho[:4]) <= 0)
+
     # avoid zero division nan:
     assert_allclose(weights, (psi + 1e-50) / (x + 1e-50), rtol=1e-6, atol=1e-8)
     psid = _approx_fprime_scalar(x, norm.rho)