BUG: np.polynomial.Polynomail.fit gets failed Singular Value Decomposition for deg=0 when there is only 1 datapoint to fit #25982
Description
Describe the issue:
Hello,
I found that calling np.polynomial.Polynomial.fit fails with a LinalgError (see below) when x and y represent a single datapoint and deg=0 (average of the data). As soon as x and y have at least 2 entries, everything works again.
np.polyfit on the other hand, handles this very special case properly.
In the reproducible example, there are 3 combinations yielding:
- ❌ fails (
do_fail=True➡️ size 1,poly_degree=0)
- ✅ works (
do_fail=False➡️ size 2,poly_degree=0)
- ✅ works (
do_fail=False➡️ size 2,poly_degree=1)
Reproduce the code example:
### Imports ###
import traceback
import numpy as np
from matplotlib import pyplot as plt
### Test ###
# the xy-data to fit are set up
do_fail = True
poly_degree = 0 # 0 will fail if `do_fail=True`, 1 will work
if do_fail:
x_fit = np.array([9.0])
y_fit = np.array([1.0])
else:
x_fit = np.array([8.0, 9.0])
y_fit = np.array([1.0, 2.0])
# `np.polyfit` will not fail
poly_old = np.polyfit(x_fit, y_fit, deg=poly_degree)
# `np.polynomial.Polynomial.fit` will fail for size 1 and `deg`=0
try:
poly_new = np.polynomial.Polynomial.fit(x_fit, y_fit, deg=poly_degree)
except np.linalg.LinAlgError:
poly_new = None
traceback.print_exc()
### Visualisation ###
x_predict = np.linspace(start=0.0, stop=10.0, num=101)
plt.scatter(x_fit, y_fit, marker="o", c="blue", s=100.0, label="data to fit")
plt.plot(
x_predict,
np.polyval(poly_old, x_predict),
label="np.polyfit",
color="blue",
linewidth=5.0,
)
if poly_new is not None:
plt.plot(
x_predict,
poly_new(x_predict),
label="np.polynomial.Polynomial.fit",
color="cyan",
linestyle="--",
linewidth=3.0,
)
# end if
plt.legend()
plt.show()Error message:
SomeUserPath\.venv311\Lib\site-packages\numpy\polynomial\polyutils.py:303: RuntimeWarning: divide by zero encountered in scalar divide
off = (old[1]*new[0] - old[0]*new[1])/oldlen
SomeUserPath\.venv311\Lib\site-packages\numpy\polynomial\polyutils.py:304: RuntimeWarning: divide by zero encountered in scalar divide
scl = newlen/oldlen
SomeUserPath\.venv311\Lib\site-packages\numpy\polynomial\polyutils.py:372: RuntimeWarning: invalid value encountered in add
return off + scl*x
** On entry to DLASCLS parameter number 4 had an illegal value
Traceback (most recent call last):
File "SomeUserPath\numpy_poly.py", line 23, in <module>
poly_new = np.polynomial.Polynomial.fit(x_fit, y_fit, deg=poly_degree)
^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^
File "SomeUserPath\.venv311\Lib\site-packages\numpy\polynomial\_polybase.py", line 1037, in fit
res = cls._fit(xnew, y, deg, w=w, rcond=rcond, full=full)
^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^
File "SomeUserPath\.venv311\Lib\site-packages\numpy\polynomial\polynomial.py", line 1362, in polyfit
return pu._fit(polyvander, x, y, deg, rcond, full, w)
^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^
File "SomeUserPath\.venv311\Lib\site-packages\numpy\polynomial\polyutils.py", line 664, in _fit
c, resids, rank, s = np.linalg.lstsq(lhs.T/scl, rhs.T, rcond)
^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^
File "SomeUserPath\.venv311\Lib\site-packages\numpy\linalg\linalg.py", line 2326, in lstsq
x, resids, rank, s = gufunc(a, b, rcond, signature=signature, extobj=extobj)
^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^
File "SomeUserPath\.venv311\Lib\site-packages\numpy\linalg\linalg.py", line 124, in _raise_linalgerror_lstsq
raise LinAlgError("SVD did not converge in Linear Least Squares")
numpy.linalg.LinAlgError: SVD did not converge in Linear Least SquaresPython and NumPy Versions:
numpy: 1.26.4
python: 3.11.6 (tags/v3.11.6:8b6ee5b, Oct 2 2023, 14:57:12) [MSC v.1935 64 bit (AMD64)]
Runtime Environment:
[{'numpy_version': '1.26.4',
'python': '3.11.6 (tags/v3.11.6:8b6ee5b, Oct 2 2023, 14:57:12) [MSC v.1935 '
'64 bit (AMD64)]',
'uname': uname_result(system='Windows', node='SomeUsersPC', release='10', version='10.0.22631', machine='AMD64')},
{'simd_extensions': {'baseline': ['SSE', 'SSE2', 'SSE3'],
'found': ['SSSE3',
'SSE41',
'POPCNT',
'SSE42',
'AVX',
'F16C',
'FMA3',
'AVX2'],
'not_found': ['AVX512F',
'AVX512CD',
'AVX512_SKX',
'AVX512_CLX',
'AVX512_CNL',
'AVX512_ICL']}},
{'architecture': 'Haswell',
'filepath': SomeUserPath\\.venv311\\Lib\\site-packages\\numpy.libs\\libopenblas64__v0.3.23-293-gc2f4bdbb-gcc_10_3_0-2bde3a66a51006b2b53eb373ff767a3f.dll',
'internal_api': 'openblas',
'num_threads': 12,
'prefix': 'libopenblas',
'threading_layer': 'pthreads',
'user_api': 'blas',
'version': '0.3.23.dev'}]
NoneContext for the issue:
I ran into this issue when running a pytest on a project that caused this relatively unlikely edge case.
The priority for this is thus low I guess because this will probably not occur in real world applications.
However, np.polyfit handles it correctly, so it would be more a matter of consistency.
Thanks for your time! 🙃