Repository: statsmodels/statsmodels (Python)

Issue:
PACF Fails When Matrix for Yule-Walker Equations is Singular
When we try to use the `pacf` function and the matrix for the Yule-Walker equations is singular, we get a [`LinAlgError`](https://numpy.org/doc/stable/reference/generated/numpy.linalg.LinAlgError.html) and the code fails.  This also leads to downstream failures.  We could resolve this by instead using the psuedo-inverse function ([`numpy.linalg.pinv`](https://numpy.org/doc/stable/reference/generated/numpy.linalg.pinv.html)) when solving the equations [here](https://github.com/statsmodels/statsmodels/blob/main/statsmodels/regression/linear_model.py#L1471).  

Here is an example of code that will currently fail with a `LinAlgError`.

```python
import numpy as np 
from statsmodels.tsa.stattools import acf, kpss, pacf

pacf(np.ones(30), nlags=6)
```


Hints:

Patch:
diff --git a/statsmodels/regression/linear_model.py b/statsmodels/regression/linear_model.py
--- a/statsmodels/regression/linear_model.py
+++ b/statsmodels/regression/linear_model.py
@@ -49,7 +49,10 @@
 # need import in module instead of lazily to copy `__doc__`
 from statsmodels.regression._prediction import PredictionResults
 from statsmodels.tools.decorators import cache_readonly, cache_writable
-from statsmodels.tools.sm_exceptions import InvalidTestWarning, ValueWarning
+from statsmodels.tools.sm_exceptions import (
+    InvalidTestWarning,
+    ValueWarning,
+    )
 from statsmodels.tools.tools import pinv_extended
 from statsmodels.tools.typing import Float64Array
 from statsmodels.tools.validation import bool_like, float_like, string_like
@@ -1475,7 +1478,15 @@ def yule_walker(x, order=1, method="adjusted", df=None, inv=False,
         r[k] = (x[0:-k] * x[k:]).sum() / (n - k * adj_needed)
     R = toeplitz(r[:-1])
 
-    rho = np.linalg.solve(R, r[1:])
+    try:
+        rho = np.linalg.solve(R, r[1:])
+    except np.linalg.LinAlgError as err:
+        if 'Singular matrix' in str(err):
+            warnings.warn("Matrix is singular. Using pinv.", ValueWarning)
+            rho = np.linalg.pinv(R) @ r[1:]
+        else:
+            raise
+
     sigmasq = r[0] - (r[1:]*rho).sum()
     sigma = np.sqrt(sigmasq) if not np.isnan(sigmasq) and sigmasq > 0 else np.nan
     if inv:
diff --git a/statsmodels/stats/gof.py b/statsmodels/stats/gof.py
--- a/statsmodels/stats/gof.py
+++ b/statsmodels/stats/gof.py
@@ -144,9 +144,9 @@ def powerdiscrepancy(observed, expected, lambd=0.0, axis=0, ddof=0):
             if axis == 0:
                 e = e.T
 
-    if np.all(np.sum(e, axis=axis) == n):
+    if np.allclose(np.sum(e, axis=axis), n, rtol=1e-8, atol=0):
         p = e/(1.0*nt)
-    elif np.all(np.sum(e, axis=axis) == 1):
+    elif np.allclose(np.sum(e, axis=axis), 1, rtol=1e-8, atol=0):
         p = e
         e = nt * e
     else:
diff --git a/statsmodels/tsa/stattools.py b/statsmodels/tsa/stattools.py
--- a/statsmodels/tsa/stattools.py
+++ b/statsmodels/tsa/stattools.py
@@ -24,6 +24,7 @@
     InfeasibleTestError,
     InterpolationWarning,
     MissingDataError,
+    ValueWarning,
 )
 from statsmodels.tools.tools import Bunch, add_constant
 from statsmodels.tools.validation import (
@@ -750,8 +751,10 @@ def pacf_yw(x, nlags=None, method="adjusted"):
 
     method = string_like(method, "method", options=("adjusted", "mle"))
     pacf = [1.0]
-    for k in range(1, nlags + 1):
-        pacf.append(yule_walker(x, k, method=method)[0][-1])
+    with warnings.catch_warnings():
+        warnings.simplefilter("once", ValueWarning)
+        for k in range(1, nlags + 1):
+            pacf.append(yule_walker(x, k, method=method)[0][-1])
     return np.array(pacf)
 
 


Test patch:
diff --git a/statsmodels/tsa/tests/test_stattools.py b/statsmodels/tsa/tests/test_stattools.py
--- a/statsmodels/tsa/tests/test_stattools.py
+++ b/statsmodels/tsa/tests/test_stattools.py
@@ -25,6 +25,7 @@
     InfeasibleTestError,
     InterpolationWarning,
     MissingDataError,
+    ValueWarning,
 )
 # Remove imports when range unit root test gets an R implementation
 from statsmodels.tools.validation import array_like, bool_like
@@ -346,6 +347,10 @@ def test_yw(self):
         pacfyw = pacf_yw(self.x, nlags=40, method="mle")
         assert_almost_equal(pacfyw[1:], self.pacfyw, DECIMAL_8)
 
+    def test_yw_singular(self):
+        with pytest.warns(ValueWarning):
+            pacf(np.ones(30), nlags=6)
+
     def test_ld(self):
         pacfyw = pacf_yw(self.x, nlags=40, method="mle")
         pacfld = pacf(self.x, nlags=40, method="ldb")