Repository: scikit-learn/scikit-learn
Pull request: [MRG+1] FIX convert to CSR in decision_function
<!--
Thanks for contributing a pull request! Please ensure you have taken a look at
the contribution guidelines: https://github.com/scikit-learn/scikit-learn/blob/master/CONTRIBUTING.md#pull-request-checklist
-->

#### Reference Issues/PRs
<!--
Example: Fixes #1234. See also #3456.
Please use keywords (e.g., Fixes) to create link to the issues or pull requests
you resolved, so that they will automatically be closed when your pull request
is merged. See https://github.com/blog/1506-closing-issues-via-pull-requests
-->

closes #9991 

#### What does this implement/fix? Explain your changes.

Make an additional `check_input` in `_decision_function` to be sure that
`X` is converted to a CSR matrix when `warm_start=True`.

ping @jnothman @amueller 

NB: I made the checking in the `_decision_function` because each time `_init_decision_function` is actually doing the same and we also do the this in `decision_function`. So it looks like it would be the central place to make this checking. WDYT?

#### Any other comments?


<!--
Please be aware that we are a loose team of volunteers so patience is
necessary; assistance handling other issues is very welcome. We value
all user contributions, no matter how minor they are. If we are slow to
review, either the pull request needs some benchmarking, tinkering,
convincing, etc. or more likely the reviewers are simply busy. In either
case, we ask for your understanding during the review process.
For more information, see our FAQ on this topic:
http://scikit-learn.org/dev/faq.html#why-is-my-pull-request-not-getting-any-attention.

Thanks for contributing!
-->


Problem statement:
GradientBoostingClassifier: sparse data and warmstarting do not work together
#### Description
Using warmstarting and sparse data for gradient boosting do not work together.

#### Steps/Code to Reproduce
```
import numpy as np
import scipy.sparse
import sklearn.datasets
import sklearn.ensemble


def get_dataset(dataset='iris', make_sparse=False):
    iris = getattr(sklearn.datasets, "load_%s" % dataset)()
    X = iris.data.astype(np.float32)
    Y = iris.target
    rs = np.random.RandomState(42)
    indices = np.arange(X.shape[0])
    train_size = int(len(indices) / 3. * 2.)
    rs.shuffle(indices)
    X = X[indices]
    Y = Y[indices]
    X_train = X[:train_size]
    Y_train = Y[:train_size]
    X_test = X[train_size:]
    Y_test = Y[train_size:]

    if make_sparse:
        X_train[:, 0] = 0
        X_train[rs.random_sample(X_train.shape) > 0.5] = 0
        X_train = scipy.sparse.csc_matrix(X_train)
        X_train.eliminate_zeros()
        X_test[:, 0] = 0
        X_test[rs.random_sample(X_test.shape) > 0.5] = 0
        X_test = scipy.sparse.csc_matrix(X_test)
        X_test.eliminate_zeros()

    return X_train, Y_train, X_test, Y_test


X_train, Y_train, _, _ = get_dataset(dataset='iris')
print(type(X_train))

classifier = sklearn.ensemble.GradientBoostingClassifier(warm_start=True)
classifier.fit(X_train, Y_train)
classifier.n_estimators += 1
classifier.fit(X_train, Y_train)
print('Fitted dense data', flush=True)

X_train, Y_train, _, _ = get_dataset(dataset='iris', make_sparse=True)
print(type(X_train))

classifier = sklearn.ensemble.GradientBoostingClassifier(warm_start=True)
classifier.fit(X_train, Y_train)
classifier.n_estimators += 1
classifier.fit(X_train, Y_train)
print('Fitted sparse data', flush=True)
```

#### Expected Results
```
<class 'numpy.ndarray'>
Fitted dense data
<class 'scipy.sparse.csc.csc_matrix'>
Fitted Sparse data
```

#### Actual Results
```
<class 'numpy.ndarray'>
Fitted dense data
<class 'scipy.sparse.csc.csc_matrix'>
```

#### Versions
Linux-4.4.0-96-generic-x86_64-with-debian-stretch-sid
Python 3.6.2 |Continuum Analytics, Inc.| (default, Jul 20 2017, 13:51:32) 
[GCC 4.4.7 20120313 (Red Hat 4.4.7-1)]
NumPy 1.13.3
SciPy 0.19.1
Scikit-Learn 0.19.1

#### Update

The same issue happens for regression:
```
import numpy as np
import scipy.sparse
import sklearn.datasets
import sklearn.ensemble


def get_dataset(dataset='boston', make_sparse=False):
    iris = getattr(sklearn.datasets, "load_%s" % dataset)()
    X = iris.data.astype(np.float32)
    Y = iris.target
    rs = np.random.RandomState(42)
    indices = np.arange(X.shape[0])
    train_size = int(len(indices) / 3. * 2.)
    rs.shuffle(indices)
    X = X[indices]
    Y = Y[indices]
    X_train = X[:train_size]
    Y_train = Y[:train_size]
    X_test = X[train_size:]
    Y_test = Y[train_size:]

    if make_sparse:
        X_train[:, 0] = 0
        X_train[rs.random_sample(X_train.shape) > 0.5] = 0
        X_train = scipy.sparse.csc_matrix(X_train)
        X_train.eliminate_zeros()
        X_test[:, 0] = 0
        X_test[rs.random_sample(X_test.shape) > 0.5] = 0
        X_test = scipy.sparse.csc_matrix(X_test)
        X_test.eliminate_zeros()

    return X_train, Y_train, X_test, Y_test


X_train, Y_train, _, _ = get_dataset(dataset='boston')
print(type(X_train))

classifier = sklearn.ensemble.GradientBoostingRegressor(warm_start=True)
classifier.fit(X_train, Y_train)
classifier.n_estimators += 1
classifier.fit(X_train, Y_train)
print('Fitted dense data', flush=True)

X_train, Y_train, _, _ = get_dataset(dataset='boston', make_sparse=True)
print(type(X_train))

classifier = sklearn.ensemble.GradientBoostingRegressor(warm_start=True)
classifier.fit(X_train, Y_train)
classifier.n_estimators += 1
classifier.fit(X_train, Y_train)
print('Fitted sparse data', flush=True)
```

Hints from issue discussion:
What's the difference between your expected and actual results?
My bad, there was none in the snippet I pasted here. I fixed the snippet. That's what the script results in:
```
<class 'numpy.ndarray'>
Fitted dense data
<class 'scipy.sparse.csc.csc_matrix'>
```
vs what I expect:
```
<class 'numpy.ndarray'>
Fitted dense data
<class 'scipy.sparse.csc.csc_matrix'>
Fitted Sparse data
```
So there's a silent crash??
yes!? (could reproduce)
@jnothman please excuse my description missing **THE** relevant information. Here's the full output:
```
<class 'numpy.ndarray'>
Fitted dense data
<class 'scipy.sparse.csc.csc_matrix'>

Process finished with exit code 139 (interrupted by signal 11: SIGSEGV)
```
The bug here is that predict_stages
(_predict_regression_tree_stages_sparse) expects (but does not assert) a
CSR matrix, which predict, decision_function will give it, but `fit` is
ordinarily agnostic to sparse format (but not with warm_start)
.

Hi, 
I will try to solve this.
Thanks
thanks.

Commit message: address joels comments

Patch to review:
diff --git a/doc/whats_new/v0.20.rst b/doc/whats_new/v0.20.rst
index 9f3e755573ce6..20990f68add51 100644
--- a/doc/whats_new/v0.20.rst
+++ b/doc/whats_new/v0.20.rst
@@ -185,6 +185,12 @@ Decomposition, manifold learning and clustering
   points could be generated. :issue:`10037` by :user:`Christian Braune
   <christianbraune79>`_.
 
+- Fixed a bug when fitting :class:`ensemble.GradientBoostingClassifier` or
+  :class:`ensemble.GradientBoostingRegressor` with ``warm_start=True`` which
+  previously raised a segmentation fault due to a non-conversion of CSC
+  matrix into CSR format expected by ``decision_function``. :issue:`9991` by
+  :user:`Guillaume Lemaitre <glemaitre>`.
+
 Metrics
 
 - Fixed a bug due to floating point error in :func:`metrics.roc_auc_score` with
diff --git a/sklearn/ensemble/_gradient_boosting.pyx b/sklearn/ensemble/_gradient_boosting.pyx
index 71371f5c24a48..2edc73053dc42 100644
--- a/sklearn/ensemble/_gradient_boosting.pyx
+++ b/sklearn/ensemble/_gradient_boosting.pyx
@@ -203,6 +203,7 @@ def predict_stages(np.ndarray[object, ndim=2] estimators,
     cdef Tree tree
 
     if issparse(X):
+        assert X.format == 'csr'
         _predict_regression_tree_stages_sparse(estimators, X, scale, out)
     else:
         if not isinstance(X, np.ndarray):
diff --git a/sklearn/ensemble/gradient_boosting.py b/sklearn/ensemble/gradient_boosting.py
index 31a82b9ce2859..4048d33e17186 100644
--- a/sklearn/ensemble/gradient_boosting.py
+++ b/sklearn/ensemble/gradient_boosting.py
@@ -1034,6 +1034,9 @@ def fit(self, X, y, sample_weight=None, monitor=None):
                                  % (self.n_estimators,
                                     self.estimators_.shape[0]))
             begin_at_stage = self.estimators_.shape[0]
+            # The requirements of _decision_function are more constrained than
+            # fit. It accepts only CSR matrices.
+            X = check_array(X, dtype=DTYPE, order="C", accept_sparse='csr')
             y_pred = self._decision_function(X)
             self._resize_state()
 
@@ -1181,7 +1184,6 @@ def _decision_function(self, X):
         predict_stages(self.estimators_, X, self.learning_rate, score)
         return score
 
-
     def _staged_decision_function(self, X):
         """Compute decision function of ``X`` for each iteration.
 
@@ -1567,7 +1569,7 @@ def decision_function(self, X):
             Regression and binary classification produce an array of shape
             [n_samples].
         """
-        X = check_array(X, dtype=DTYPE, order="C",  accept_sparse='csr')
+        X = check_array(X, dtype=DTYPE, order="C", accept_sparse='csr')
         score = self._decision_function(X)
         if score.shape[1] == 1:
             return score.ravel()
@@ -2000,7 +2002,7 @@ def predict(self, X):
         y : array of shape = [n_samples]
             The predicted values.
         """
-        X = check_array(X, dtype=DTYPE, order="C",  accept_sparse='csr')
+        X = check_array(X, dtype=DTYPE, order="C", accept_sparse='csr')
         return self._decision_function(X).ravel()
 
     def staged_predict(self, X):
diff --git a/sklearn/ensemble/tests/test_gradient_boosting.py b/sklearn/ensemble/tests/test_gradient_boosting.py
index f4594529e034b..e8ff2d215ee66 100644
--- a/sklearn/ensemble/tests/test_gradient_boosting.py
+++ b/sklearn/ensemble/tests/test_gradient_boosting.py
@@ -849,6 +849,34 @@ def test_warm_start_oob():
                                   est.oob_improvement_[:100])
 
 
+def test_warm_start_sparse():
+    # Test that all sparse matrix types are supported
+    X, y = datasets.make_hastie_10_2(n_samples=100, random_state=1)
+    sparse_matrix_type = [csr_matrix, csc_matrix, coo_matrix]
+    for Cls in [GradientBoostingRegressor, GradientBoostingClassifier]:
+        est_dense = Cls(n_estimators=100, max_depth=1, subsample=0.5,
+                        random_state=1, warm_start=True)
+        est_dense.fit(X, y)
+        est_dense.predict(X)
+        est_dense.set_params(n_estimators=200)
+        est_dense.fit(X, y)
+        est_dense.predict(X)
+
+        for sparse_constructor in sparse_matrix_type:
+            X_sparse = sparse_constructor(X)
+
+            est_sparse = Cls(n_estimators=100, max_depth=1, subsample=0.5,
+                             random_state=1, warm_start=True)
+            est_sparse.fit(X_sparse, y)
+            est_sparse.predict(X)
+            est_sparse.set_params(n_estimators=200)
+            est_sparse.fit(X_sparse, y)
+            est_sparse.predict(X)
+
+            assert_array_almost_equal(est_dense.oob_improvement_[:100],
+                                      est_sparse.oob_improvement_[:100])
+
+
 def early_stopping_monitor(i, est, locals):
     """Returns True on the 10th iteration. """
     if i == 9:


Code review comments:
diff_hunk: @@ -203,6 +203,7 @@ def predict_stages(np.ndarray[object, ndim=2] estimators,
     cdef Tree tree
 
     if issparse(X):
+        assert X.format == 'csr'
path: sklearn/ensemble/_gradient_boosting.pyx
text: @user1:
Can you use a ValueError to have a better error message? There is one error message below that maybe you can reuse.

```py
raise ValueError("X should be in np.ndarray or csr_matrix format,"
                              "got %s" % type(X))
```

@user1:
Also I guess `predict_stages` can be imported, so you could unit-test the error message.

@user2:
Fwiw I suggested an assert because this is an internal problem, not a user value error

diff_hunk: @@ -185,6 +185,12 @@ Decomposition, manifold learning and clustering
   points could be generated. :issue:`10037` by :user:`Christian Braune
   <christianbraune79>`_.
 
+- Fixed a bug when fitting :class:`ensemble.GradientBoostingClassifier` or
path: doc/whats_new/v0.20.rst
text: @user1:
This does not seem the right section. At the moment it is in "Decomposition, manifold learning and clustering".