Repository: scikit-learn/scikit-learn
Pull request: [MRG] Changed VarianceThreshold behaviour when threshold is zero. See #13691
<!--
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
Fixes #13691

#### What does this implement/fix? Explain your changes.
Currently, VarianceThreshold can fail to remove features that only contain one repeated value when the threshold is 0 due to floating point errors introduced by division in `np.var`. As proposed in #13691, this changes the implementation of `VarianceThreshold.fit` to use `np.ptp` (difference between max and min) rather than `np.var` when the threshold is 0 (and the analogous `utils.sparsefuncs.min_max_axis` for sparse inputs). This avoids division and causes VarianceThreshold to behave as expected (verified by testing on the code in the issue).

#### 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:
VarianceThreshold doesn't remove feature with zero variance
#### Description
When calling VarianceThreshold().fit_transform() on certain inputs, it fails to remove a column that has only one unique value.

#### Steps/Code to Reproduce
```python
import numpy as np
from sklearn.feature_selection import VarianceThreshold

works_correctly = np.array([[-0.13725701,  7.        ],
                            [-0.13725701, -0.09853293],
                            [-0.13725701, -0.09853293],
                            [-0.13725701, -0.09853293],
                            [-0.13725701, -0.09853293],
                            [-0.13725701, -0.09853293],
                            [-0.13725701, -0.09853293],
                            [-0.13725701, -0.09853293],
                            [-0.13725701, -0.09853293]])

broken = np.array([[-0.13725701,  7.        ],
                   [-0.13725701, -0.09853293],
                   [-0.13725701, -0.09853293],
                   [-0.13725701, -0.09853293],
                   [-0.13725701, -0.09853293],
                   [-0.13725701, -0.09853293],
                   [-0.13725701, -0.09853293],
                   [-0.13725701, -0.09853293],
                   [-0.13725701, -0.09853293],
                   [-0.13725701, -0.09853293]])

selector = VarianceThreshold()
print(selector.fit_transform(works_correctly))

selector = VarianceThreshold()
print(selector.fit_transform(broken))
print(set(broken[:, 0]))
```

#### Expected Results
The Variance threshold should produce
```
[[ 7.        ]
 [-0.09853293]
 [-0.09853293]
 [-0.09853293]
 [-0.09853293]
 [-0.09853293]
 [-0.09853293]
 [-0.09853293]
 [-0.09853293]]
```
#### Actual Results
```
[[ 7.        ]
 [-0.09853293]
 [-0.09853293]
 [-0.09853293]
 [-0.09853293]
 [-0.09853293]
 [-0.09853293]
 [-0.09853293]
 [-0.09853293]]
[[-0.13725701  7.        ]
 [-0.13725701 -0.09853293]
 [-0.13725701 -0.09853293]
 [-0.13725701 -0.09853293]
 [-0.13725701 -0.09853293]
 [-0.13725701 -0.09853293]
 [-0.13725701 -0.09853293]
 [-0.13725701 -0.09853293]
 [-0.13725701 -0.09853293]
 [-0.13725701 -0.09853293]]
{-0.13725701}
```
This issue arose when I was using VarianceThreshold on a real dataset (of which this is a subset). It appears to work correctly in other situations (for instance I can't reproduce this behaviour if I replace the first column with 1's).

#### Versions
System
------
    python: 3.5.6 |Anaconda, Inc.| (default, Aug 26 2018, 16:30:03)  [GCC 4.2.1 Compatible Clang 4.0.1 (tags/RELEASE_401/final)]
   machine: Darwin-18.2.0-x86_64-i386-64bit
executable: /anaconda3/envs/tensorflow/bin/python3

BLAS
----
    macros: NO_ATLAS_INFO=3, HAVE_CBLAS=None
  lib_dirs:
cblas_libs: cblas

Python deps
-----------
setuptools: 40.2.0
     numpy: 1.15.4
   sklearn: 0.20.0
    Cython: None
     scipy: 1.1.0
    pandas: 0.24.0
       pip: 19.0.1

Hints from issue discussion:
On closer inspection this is just caused by floating point error in calculating the variance, and therefore not a bug with sklearn. It is resolvable by setting the variance threshold to e.g. 1e-33 rather than 0.
We should probably avoid 0 as a default. I'd be happy to deprecate the
current default and change it to np.finfo(X.dtype) or something smaller by
default.

Another option would be to use np.ptp rather than np.var when the threshold is 0.
Cool, here's a pull request that does the latter (and the analogous operation for sparse inputs).

Commit message: Update v0.22.rst

Patch to review:
diff --git a/doc/whats_new/v0.22.rst b/doc/whats_new/v0.22.rst
index 6cc93ecdfe358..dda9adbc2e79e 100644
--- a/doc/whats_new/v0.22.rst
+++ b/doc/whats_new/v0.22.rst
@@ -57,7 +57,6 @@ Changelog
   ``decision_function_shape='ovr'``, and the number of target classes > 2.
   :pr:`12557` by `Adrin Jalali`_.
 
-
 :mod:`sklearn.cluster`
 ..................
 
@@ -65,8 +64,15 @@ Changelog
   parameter. This parameter extends `SpectralClustering` class functionality to
   match `spectral_clustering`.
   :pr:`13726` by :user:`Shuzhe Xiao <fdas3213>`.
+
+:mod:`sklearn.feature_selection`
+................................
+- |Fix| :class:`feature_selection.VarianceThreshold` now uses range
+  rather than variance when :term:`threshold` is 0. This prevents constant
+  features from sometimes being incorrectly kept due to floating-point error
+  in variance calculation.
+  :pr:`13704` by `Roddy MacSween <rlms>`.
   
-	
 Miscellaneous
 .............
 
diff --git a/sklearn/feature_selection/tests/test_variance_threshold.py b/sklearn/feature_selection/tests/test_variance_threshold.py
index a40491302f350..fba4478a28e2f 100644
--- a/sklearn/feature_selection/tests/test_variance_threshold.py
+++ b/sklearn/feature_selection/tests/test_variance_threshold.py
@@ -1,3 +1,6 @@
+import numpy as np
+import pytest
+
 from sklearn.utils.testing import (assert_array_equal, assert_equal,
                                    assert_raises)
 
@@ -26,3 +29,17 @@ def test_variance_threshold():
     for X in [data, csr_matrix(data)]:
         X = VarianceThreshold(threshold=.4).fit_transform(X)
         assert_equal((len(data), 1), X.shape)
+
+
+def test_zero_variance_floating_point_error():
+    # Test that VarianceThreshold(0.0).fit eliminates features that have
+    # the same value in every sample, even when floating point errors
+    # cause np.var not to be 0 for the feature.
+    # See #13691
+
+    data = [[-0.13725701]] * 10
+    assert np.var(data) != 0
+    for X in [data, csr_matrix(data), csc_matrix(data), bsr_matrix(data)]:
+        msg = "No feature in X meets the variance threshold 0.00000"
+        with pytest.raises(ValueError, match=msg):
+            VarianceThreshold().fit(X)
diff --git a/sklearn/feature_selection/variance_threshold.py b/sklearn/feature_selection/variance_threshold.py
index c9e018d94a84e..7d98de82c9711 100644
--- a/sklearn/feature_selection/variance_threshold.py
+++ b/sklearn/feature_selection/variance_threshold.py
@@ -5,7 +5,7 @@
 from ..base import BaseEstimator
 from .base import SelectorMixin
 from ..utils import check_array
-from ..utils.sparsefuncs import mean_variance_axis
+from ..utils.sparsefuncs import mean_variance_axis, min_max_axis
 from ..utils.validation import check_is_fitted
 
 
@@ -65,8 +65,18 @@ def fit(self, X, y=None):
 
         if hasattr(X, "toarray"):   # sparse matrix
             _, self.variances_ = mean_variance_axis(X, axis=0)
+            if self.threshold == 0:
+                mins, maxes = min_max_axis(X, axis=0)
+                peak_to_peaks = maxes - mins
         else:
             self.variances_ = np.var(X, axis=0)
+            if self.threshold == 0:
+                peak_to_peaks = np.ptp(X, axis=0)
+
+        if self.threshold == 0:
+            # Use peak-to-peak to avoid numeric precision issues
+            # for constant features
+            self.variances_ = np.minimum(self.variances_, peak_to_peaks)
 
         if np.all(self.variances_ <= self.threshold):
             msg = "No feature in X meets the variance threshold {0:.5f}"


Code review comments:
diff_hunk: @@ -57,16 +57,22 @@ Changelog
   ``decision_function_shape='ovr'``, and the number of target classes > 2.
   :pr:`12557` by `Adrin Jalali`_.
 
-
 :mod:`sklearn.cluster`
 ..................
 
 - |Enhancement| :class:`cluster.SpectralClustering` now accepts a ``n_components`` 
   parameter. This parameter extends `SpectralClustering` class functionality to
   match `spectral_clustering`.
   :pr:`13726` by :user:`Shuzhe Xiao <fdas3213>`.
+
+:mod:`sklearn.feature_selection`
+................................
+- |Fix| :class:`feature_selection.VarianceThreshold` now uses range
+  rather than variance when :term:`threshold` is 0. This prevents constant
path: doc/whats_new/v0.22.rst
text: @user1:
I should have picked these things up before. `:term:` can only be used when there is an entry in the glossary.

it may also be more helpful to phrase this with respect to the user: "Fixed a bug where VarianceThreshold with threshold=0 did not remove constant features due to numerical instability."