opened
reiinakano
username_0
<!--
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#Contributing-Pull-Requests
-->
#### Reference Issue
<!-- Example: Fixes #1234 -->
Fixes #9589 

#### What does this implement/fix? Explain your changes.
I figured that the purpose of this code (original)

```python
if method in ['decision_function', 'predict_proba', 'predict_log_proba']:
        n_classes = len(set(y))
        predictions_ = np.zeros((X_test.shape[0], n_classes))
        if method == 'decision_function' and len(estimator.classes_) == 2:
                predictions_[:, estimator.classes_[-1]] = predictions
        else:
            predictions_[:, estimator.classes_] = predictions
        predictions = predictions_
```

is to handle cases where the cross-validation done in `cross_val_predict` doesn't properly stratify the classes i.e. one fold in the CV may not have all classes present.

However, it misses an edge case for `decision_function` when `len(estimator.classes_) == 2` (train split has only 2 classes) and `n_classes == 2` (total data has only 2 classes). It wrongly assumes that the output will be of shape `(n_samples, n_classes)`, when in fact, this is only true for `decision_function` if `n_classes > 2`. For the correct behavior, we must first check if the total number of classes is greater than 2, and if it is not, we must stick to an output shape of `(n_samples,)`.

#### Any other comments?

While this PR fixes the case raised in #9589, it still doesn't take into account the fact that for `sklearn.svm.SVC`, `decision_function` does not always follow the shape `(n_samples, n_classes)` or `(n_samples,)`. If `decision_function_shape` is not set to `ovr`, the shape is `(n_samples, n_classes * (n_classes-1) / 2)`, and if different train splits have different numbers of classes (not stratified), this inconsistency can't be fixed by the approach currently taken by `cross_val_predict`.

I would also like to point out my personal opinion that `cross_val_predict` should probably just throw an error if the classes are not stratified properly (by this I mean some classes are entirely absent in some splits) instead of handling it automatically as it tries to do now. Metrics like log loss and ROC might be affected if for some splits, entire columns are zeroed out. Idk, it just feels inconsistent. Of course, you guys have the final say on this.

What I mean by inconsistent:

```python
from sklearn.model_selection import cross_val_predict, KFold
from sklearn.datasets import load_iris
from sklearn.linear_model import LogisticRegression

X, y = load_iris(True)
cross_val_predict(LogisticRegression(), X, y, cv=KFold().split(X), method='predict_proba')
```

returns (array cut down in size)

```
array([[  0.00000000e+00,   9.99987841e-01,   1.21585683e-05],
       [  0.00000000e+00,   9.99963190e-01,   3.68104260e-05],
       [  0.00000000e+00,   9.99970224e-01,   2.97757119e-05],
       [  0.00000000e+00,   9.99978481e-01,   2.15185386e-05],
       [  0.00000000e+00,   9.99993136e-01,   6.86428459e-06],
       [  1.20030529e-01,   0.00000000e+00,   8.79969471e-01],
       [  1.13183668e-01,   0.00000000e+00,   8.86816332e-01],
       [  6.75506225e-02,   0.00000000e+00,   9.32449377e-01],
       [  1.07273496e-01,   0.00000000e+00,   8.92726504e-01],
       [  1.26854747e-01,   0.00000000e+00,   8.73145253e-01],
       [  1.16881976e-01,   0.00000000e+00,   8.83118024e-01],
       [  2.45077442e-04,   9.99754923e-01,   0.00000000e+00],
       [  2.70475201e-04,   9.99729525e-01,   0.00000000e+00],
       [  6.59523734e-04,   9.99340476e-01,   0.00000000e+00],
       [  4.07574999e-04,   9.99592425e-01,   0.00000000e+00],
       [  2.42451670e-03,   9.97575483e-01,   0.00000000e+00],
       [  2.03959503e-03,   9.97960405e-01,   0.00000000e+00]])
```
Notice how there are "zeroed out" columns.

And even worse:

```python
from sklearn.model_selection import cross_val_predict, KFold
from sklearn.datasets import load_iris
from sklearn.linear_model import LogisticRegression

X, y = load_iris(True)
cross_val_predict(LogisticRegression(), X, y, cv=KFold().split(X), method='decision_function')
```
returns (array cut down in size)
```
array([[  0.        ,   0.        , -11.31746426],
       [  0.        ,   0.        , -10.20969263],
       [  0.        ,   0.        , -10.42178776],
       [  0.        ,   0.        ,  -9.60397656],
       [  0.        ,   0.        , -11.30005126],
       [  0.        ,   0.        , -11.1906238 ],
       [  0.        ,   0.        , -10.05510072],
       [  0.        ,   0.        , -10.74657422],
       [  0.        ,   0.        ,  -9.20295991],
       [  0.        ,   8.3136912 ,   0.        ],
       [  0.        ,   6.77281268,   0.        ],
       [  0.        ,   7.79874366,   0.        ],
       [  0.        ,   7.29615134,   0.        ],
       [  0.        ,   7.9199709 ,   0.        ],
       [  0.        ,   9.163118  ,   0.        ],
       [  0.        ,   5.8858718 ,   0.        ],
       [  0.        ,   6.37425957,   0.        ],
       [  0.        ,   6.66261924,   0.        ],
       [  0.        ,   6.19296233,   0.        ]])
```
which is hardly of use to anyone.

<!--
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!
-->
[WIP] Fix cross_val_predict behavior for binary classification in decision_function (Fixes #9589)
issue
created
jnothman
username_1
I think better not to error in the context of using this in ClassifierChain or a stacking meta-estimator.
comment
created
reiinakano
username_0
NaN makes it clearer that these columns weren't considered at all for this particular split, but wouldn't that break things? For instance, I frequently use `np.argmax(cross_val_predict(..., method='predict_proba))` to recover the true labels for calculating accuracy, etc. Of course, I always make sure my folds are stratified so I'm probably not going to face it, but still.. Also, wouldn't that break ClassifierChains and stacking meta-estimators?
comment
created
jnothman
username_1
yes, nan would probably break it all. a warning is probably a good idea.
But first we need to fix the regression.
comment
created
reiinakano
username_0
I've added the necessary unit tests.

Although at this point, I'd like to point out there's another regression for the case of SVC with `decision_function_shape=ovo`

Code to test:
```python
from sklearn.datasets import load_digits
from sklearn.svm import SVC
from sklearn.model_selection import cross_val_predict

X, y = load_digits(return_X_y=True)
cross_val_predict(SVC(kernel='linear', decision_function_shape='ovo'), X, y, method='decision_function').shape
```

**Expected Results**
In previous versions, this is `(1797, 45)`

**Actual Results**
```
ValueError: shape mismatch: value array of shape (602,45) could not be broadcast to indexing result of shape (10,602)
```

This error is due to the shape of an SVC decision function with `decision_function_shape` set to `'ovr'` which is `(n_samples, n_classes * (n_classes-1) / 2)`

Should I include the fix in this PR?
comment
created
jnothman
username_1
Thanks. Maybe fix that in a separate PR unless they're going to step on
each others' toes.
comment
created
reiinakano
username_0
Hi @username_1, I changed this back to WIP since there was a conflict with a recent merge. I also figured out a better fix that also fixes the SVC "ovo" regression.
comment
created
jnothman
username_1
Okay. Let me know when it's MRG, i.e. ready for review.
comment
created
reiinakano
username_0
@username_1 Ready for review
comment
created
reiinakano
username_0
Yup, sorry. Added unit test and addressed comments.
comment
created
amueller
username_2
A solution to some more of this will be in #9532 but that will probably still take a while.
comment
created
amueller
username_2
comment
created
amueller
username_2
This is the case when some classifiers have 1 class and some have 2 classes.
comment
created
reiinakano
username_0
Got it. Didn't know some classifiers could be fit on a single class. So that's *another* regression. How do you want to handle this case? It sounds similar to the case `if method == 'decision_function' and len(estimator.classes_) == 2`, which is currently handled by doing `predictions_[:, estimator.classes_[-1]] = predictions`.

If that's fine, I think the fix is to simply change the former line to `if method == 'decision_function' and len(estimator.classes_) <= 2`.

This would fix the errors, but frankly, I'd still prefer a stern warning that something went wrong during cross-validation using `cross_val_predict`.

For the ovo case, should I include that in this PR as well? I didn't include it previously since this PR primarily fixes regressions and I didn't want to include a new "feature".
comment
created
amueller
username_2
Not sure if the 1 class thing is a regression or was always broken ;)
And I think I'd like the error / test for OVO in this PR.

Sorry, I don't understand your proposal for the 1class / 2class case. The result in the 2 class case should be 1d, not 2d. That's the original bug this is trying to solve, right?
comment
created
reiinakano
username_0
The case I'm talking about is when the original dataset has >2 classes, and a particular fold is <=2 *and* the method is `decision_function`. The original bug is still fixed by the `if n_classes != len(estimator.classes_):` check.
comment
created
amueller
username_2
Ah. My breaking example is with the dataset having two classes, and some folds have one class.
comment
created
reiinakano
username_0
You're right. It's not a regression, but another edge case that needs to be fixed. Yes, we should add an explicit check that `n_classes > 2` but only if the method used is `decision_function`.
comment
created
reiinakano
username_0
@username_2 I've added both the error for OvO and the fix for a dataset having two classes, and some folds having one class. Added unit tests for both.
comment
created
reiinakano
username_0
@username_1 Addressed comments. Failure in Appveyor appears from an unrelated test.
comment
created
reiinakano
username_0
I'm not sure why an unrelated test in `sklearn.neighbors` is failing in Appveyor when the only difference in the commit that caused the failure is a few comments and additional unit tests.
comment
created
jnothman
username_1
That appears to be a very rare failure at random. I'll fix that test and restart appveyor.
comment
created
reiinakano
username_0
Finally green. Thanks! Ready for review.
comment
created
amueller
username_2
Ok, I think if one fold only contains one class, we should probably just error becaus it's a bit unclear to me what is supposed to happen. And this is such a weird edge-case that I don't think we should hold 0.19.1 for that.

For the if-else thing in the decision function case: maybe just add some comments? right now I find the code hard to read. But seems ok to me.
comment
created
reiinakano
username_0
@username_2 I've addressed your comments. So, ready for review again.

The test failures are due to something unrelated again..
comment
created
reiinakano
username_0
@username_2 addressed your comments. Ready for review again.

The test failures seem to be due to something unrelated again...
comment
created
jnothman
username_1
Add a blank commit and the failures should go away due to a patch in master
comment
created
reiinakano
username_0
So.. something failed. But again, I have no idea what...
comment
created
lesteve
username_3
You can see that this is lgtm failing (google it if you want to know more what lgtm does). For some reason it has done that recently, not sure why but you don't need to worry about this.
comment
created
jnothman
username_1
I'll add a what's new entry when throwing it into 0.19.1, to avoid having it incorrectly in 0.20
comment
created
lesteve
username_3
I looked at this for a bit, my main comment: with your change it seems like we are not very consistent  when we ave a training fold which does not have all the training labels present in the complete dataset:
* if prediction.ndim == 2, which I think for method='decision_function' means `len(est.classes_) >= 3`, we raise an error saying use a proper stratification technique
* if prediction.ndim < 2, which I think for method='decision_function' means `len(est.classes_) <= 2`, we do not raise an error.

Also it seems like the fix for cross_val_predict in binary classification as reported in #9589 is actually moving the from the logic from https://github.com/scikit-learn/scikit-learn/pull/7889 in a `if n_classes != len(estimator.classes_)`. I would like to take a closer look at this to understand better whether this is the right thing to do.
comment
created
reiinakano
username_0
I don't think there is ever a case when `decision_function` returns `predictions.ndim == 2`.
comment
created
lesteve
username_3
You are right, I read the code too quickly, which to my defence is a bit tricky to follow.

I am not convinced we are doing the right thing when we don't have all the classes in the training fold (for example setting decision function to 0 in the missing classes does not seem right).
One way to avoid the problem would be to raise an error when we don't have all the classes in the training fold? Is there a reason why we don't do this?
comment
created
reiinakano
username_0
Yup, I agree this does not make sense.

Perhaps what we could do is keep the special imbalanced handling code only for `predict_proba` and error out when it is imbalanced for `decision_function`.
comment
created
jnothman
username_1
The problem is that it does just happen in practice, when some classes are
infrequent. It might not make sense for decision_function to assume some
arbitrary value but it does make sense to zero out a column in
predict_proba.
comment
created
reiinakano
username_0
How about erroring out imbalanced classes in `decision_function` but zeroing out columns in `predict_proba` and `log_predict_proba`? Literally all the complication in this PR is because of the different weird ways `decision_function` behaves.
comment
created
jnothman
username_1
Essentially all we do in that case is make it harder for people to do the
same analysis on a highly imbalanced classification problem that they do
with a balanced one, which in turn makes it harder for people to provide
generic ML tools built on top of our library.

One example, again, is stacking.
comment
created
jnothman
username_1
I suppose that is, at least, more principled, and I can accept that. (And
by that I mean erroring for decision_function but not for
predict_{log_}proba).

Yes, these quirks and edge cases are annoying to deal with. But I feel like
sometimes it's our role as a library provider to handle them as best we
can, and raise an error only when the user has done something wrong or at
least very rare. I don't think this falls under very rare.
comment
created
lesteve
username_3
I read a bit the discussion in the related PRs, sorry it's taking me some time to catch up. Basically here is my summary:

* I would be a favour of a warning when the training fold does not contain all the classes. https://github.com/scikit-learn/scikit-learn/pull/8773#issuecomment-330362879 agrees with that as well.
* the fix for the wrong shape in cross_validation_predict with method=decision_function (which consists in only applying the logic from #7889 when `n_classes != n_classes_train`) is fine
* `method='decision_function'` has a number of edge cases. Raising an exception for `method='decision_function'` if `n_classes_train != n_classes` may be an easy way to avoid tackling all these special cases.
  1. SVC with ovo with uncommon number of columns for decision function `n_classes_train * (n_classes_train - 1) / 2` rather than `n_classes`. Tackled in this PR already.
  2. `n_classes_train <= 2` and `n_classes >= 3`, you get a 1d array as the output of decision_function and then it is not clear what to do with it. I would be in favour of raising an exception in this case.
  3. `n_classes_train == 1` and `n_classes == 2`. No idea what to do in this case. I don't know how common this is for classifier to be able to be fittable with a single class and have decision_function. I quickly looked using all_estimators and I found `['AdaBoostClassifier', 'LinearDiscriminantAnalysis',  'RidgeClassifier', 'RidgeClassifierCV']`.
* predictions should be initialized to 0 if `method='predict_proba'` but `-inf` (or maybe `np.finfo(X.dtype).min`) if `method='predict_log_proba'`.
comment
created
reiinakano
username_0
There is another case for `decision_function`

iv. `n_classes_train > 2` and `n_classes > 2` but `n_classes_train < n_classes`. Currently we set missing classes to 0.
comment
created
lesteve
username_3
For decision_function either we do the same as predict_log_proba or we decide that we always raise an exception, given that there are many edge cases for decision_function ...
comment
created
reiinakano
username_0
I think I've mostly sorted this out. For a non-special case of `decision_function` (`iv. n_classes_train > 2 and n_classes > 2 but n_classes_train != n_classes.`), I do the same as `predict_log_proba`. For the other special cases, there is an error message. 

Apologies for the confusing code previously, I should've done it this way in the first place, but I didn't expect `decision_function` to have so many edge cases.
comment
created
jnothman
username_1
I would like to see this merged and 0.19.1 released by the end of the week. Is that reasonable, @username_2?
comment