Skip to content

LogisticRegression fails in sklearn 1.1.1 with newton-cg solver when X only contains one predictor #23605

New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Closed
raphaelvallat opened this issue Jun 13, 2022 · 1 comment · Fixed by #23608
Labels

Comments

@raphaelvallat
Copy link

Describe the bug

I think this is related to #21808

As explained in the title, this bug appeared in version 1.1.1. LogisticRegression(solver="newton-cg") fails when X only contains a single predictor, unless using fit_intercept=False.

Steps/Code to Reproduce

import numpy as np
from sklearn.linear_model import LogisticRegression

y = np.array([1, 1, 0, 0, 1, 1, 0, 1])
X = np.array([[0.5, 0.65, 1.1, 1.25, 0.8, 0.54, 0.95, 0.7]]).T

# Newton-cg solver
LogisticRegression(solver="newton-cg").fit(X, y)  # FAILS in 1.1.1, works in previous versions
LogisticRegression(solver="newton-cg", fit_intercept=False).fit(X, y)  # Works in 1.1.1 and previous versions

# Default solver
LogisticRegression().fit(X, y)  # Works in 1.1.1 and previous versions

Expected Results

No error is thrown

Actual Results

---------------------------------------------------------------------------
ValueError                                Traceback (most recent call last)
<ipython-input-1-bc267bc834b6> in <module>
      4 y = np.array([1, 1, 0, 0, 1, 1, 0, 1])
      5 X = np.array([[0.5, 0.65, 1.1, 1.25, 0.8, 0.54, 0.95, 0.7]]).T
----> 6 LogisticRegression(solver="newton-cg").fit(X, y)  # Fails
      7 LogisticRegression().fit(X, y)

~/.pyenv/versions/3.8.3/lib/python3.8/site-packages/sklearn/linear_model/_logistic.py in fit(self, X, y, sample_weight)
   1231             n_threads = 1
   1232 
-> 1233         fold_coefs_ = Parallel(n_jobs=self.n_jobs, verbose=self.verbose, prefer=prefer)(
   1234             path_func(
   1235                 X,

~/.pyenv/versions/3.8.3/lib/python3.8/site-packages/joblib/parallel.py in __call__(self, iterable)
   1039             # remaining jobs.
   1040             self._iterating = False
-> 1041             if self.dispatch_one_batch(iterator):
   1042                 self._iterating = self._original_iterator is not None
   1043 

~/.pyenv/versions/3.8.3/lib/python3.8/site-packages/joblib/parallel.py in dispatch_one_batch(self, iterator)
    857                 return False
    858             else:
--> 859                 self._dispatch(tasks)
    860                 return True
    861 

~/.pyenv/versions/3.8.3/lib/python3.8/site-packages/joblib/parallel.py in _dispatch(self, batch)
    775         with self._lock:
    776             job_idx = len(self._jobs)
--> 777             job = self._backend.apply_async(batch, callback=cb)
    778             # A job can complete so quickly than its callback is
    779             # called before we get here, causing self._jobs to

~/.pyenv/versions/3.8.3/lib/python3.8/site-packages/joblib/_parallel_backends.py in apply_async(self, func, callback)
    206     def apply_async(self, func, callback=None):
    207         """Schedule a func to be run"""
--> 208         result = ImmediateResult(func)
    209         if callback:
    210             callback(result)

~/.pyenv/versions/3.8.3/lib/python3.8/site-packages/joblib/_parallel_backends.py in __init__(self, batch)
    570         # Don't delay the application, to avoid keeping the input
    571         # arguments in memory
--> 572         self.results = batch()
    573 
    574     def get(self):

~/.pyenv/versions/3.8.3/lib/python3.8/site-packages/joblib/parallel.py in __call__(self)
    260         # change the default number of processes to -1
    261         with parallel_backend(self._backend, n_jobs=self._n_jobs):
--> 262             return [func(*args, **kwargs)
    263                     for func, args, kwargs in self.items]
    264 

~/.pyenv/versions/3.8.3/lib/python3.8/site-packages/joblib/parallel.py in <listcomp>(.0)
    260         # change the default number of processes to -1
    261         with parallel_backend(self._backend, n_jobs=self._n_jobs):
--> 262             return [func(*args, **kwargs)
    263                     for func, args, kwargs in self.items]
    264 

~/.pyenv/versions/3.8.3/lib/python3.8/site-packages/sklearn/utils/fixes.py in __call__(self, *args, **kwargs)
    115     def __call__(self, *args, **kwargs):
    116         with config_context(**self.config):
--> 117             return self.function(*args, **kwargs)
    118 
    119 

~/.pyenv/versions/3.8.3/lib/python3.8/site-packages/sklearn/linear_model/_logistic.py in _logistic_regression_path(X, y, pos_class, Cs, fit_intercept, max_iter, tol, verbose, solver, coef, class_weight, dual, penalty, intercept_scaling, multi_class, random_state, check_input, max_squared_sum, sample_weight, l1_ratio, n_threads)
    452             l2_reg_strength = 1.0 / C
    453             args = (X, target, sample_weight, l2_reg_strength, n_threads)
--> 454             w0, n_iter_i = _newton_cg(
    455                 hess, func, grad, w0, args=args, maxiter=max_iter, tol=tol
    456             )

~/.pyenv/versions/3.8.3/lib/python3.8/site-packages/sklearn/utils/optimize.py in _newton_cg(grad_hess, func, grad, x0, args, tol, maxiter, maxinner, line_search, warn)
    191         # Inner loop: solve the Newton update by conjugate gradient, to
    192         # avoid inverting the Hessian
--> 193         xsupi = _cg(fhess_p, fgrad, maxiter=maxinner, tol=termcond)
    194 
    195         alphak = 1.0

~/.pyenv/versions/3.8.3/lib/python3.8/site-packages/sklearn/utils/optimize.py in _cg(fhess_p, fgrad, maxiter, tol)
     86             break
     87 
---> 88         Ap = fhess_p(psupi)
     89         # check curvature
     90         curv = np.dot(psupi, Ap)

~/.pyenv/versions/3.8.3/lib/python3.8/site-packages/sklearn/linear_model/_linear_loss.py in hessp(s)
    344                 if self.fit_intercept:
    345                     ret[:n_features] += s[-1] * hX_sum
--> 346                     ret[-1] = hX_sum @ s[:n_features] + hessian_sum * s[-1]
    347                 return ret
    348 

ValueError: matmul: Input operand 0 does not have enough dimensions (has 0, gufunc core with signature (n?,k),(k,m?)->(n?,m?) requires 1)

Versions

System:
    python: 3.8.3 (default, Jun 26 2020, 00:34:38)  [Clang 11.0.3 (clang-1103.0.32.62)]
executable: /Users/raphael/.pyenv/versions/3.8.3/bin/python3.8
   machine: macOS-10.16-x86_64-i386-64bit

Python dependencies:
      sklearn: 1.1.1
          pip: 22.1.2
   setuptools: 41.2.0
        numpy: 1.21.6
        scipy: 1.8.0
       Cython: 0.29.23
       pandas: 1.4.1
   matplotlib: 3.5.1
       joblib: 1.0.1
threadpoolctl: 2.1.0

Built with OpenMP: True
@TomDLT
Copy link
Member

TomDLT commented Jun 14, 2022

Thanks for the detailed report, it made the bug easy to fix.

Fix in #23608

@TomDLT TomDLT removed the Needs Triage Issue requires triage label Jun 14, 2022
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment
Labels
Projects
None yet
2 participants