#4475 : Add a safe_pairwise_distances function, dealing with zero variance samples when using correlation metric.

Lowik CHANUSSOT · Lowik CHANUSSOT · commit 359fd1a21ec4 · 2015-04-02T14:49:40.000+02:00
The best fix would be to have the metric not returning NaN values, but as the correlation metric is actually computed by spicy, we can't modify it directly.
So, in case of metric=='correlation', we replace rows and cols corresponding to zero variance samples by the maximum distance (here 1.0).
diff --git a/sklearn/manifold/t_sne.py b/sklearn/manifold/t_sne.py
@@ -17,7 +17,7 @@
 from ..utils import check_random_state
 from ..utils.extmath import _ravel
 from ..decomposition import RandomizedPCA
-from ..metrics.pairwise import pairwise_distances
+from ..metrics.pairwise import safe_pairwise_distances
 from . import _utils
 
 
@@ -269,8 +269,8 @@ def trustworthiness(X, X_embedded, n_neighbors=5, precomputed=False):
     if precomputed:
         dist_X = X
     else:
-        dist_X = pairwise_distances(X, squared=True)
-    dist_X_embedded = pairwise_distances(X_embedded, squared=True)
+        dist_X = safe_pairwise_distances(X, squared=True)
+    dist_X_embedded = safe_pairwise_distances(X_embedded, squared=True)
     ind_X = np.argsort(dist_X, axis=1)
     ind_X_embedded = np.argsort(dist_X_embedded, axis=1)[:, 1:n_neighbors + 1]
 
@@ -438,9 +438,9 @@ def fit(self, X, y=None):
                 print("[t-SNE] Computing pairwise distances...")
 
             if self.metric == "euclidean":
-                distances = pairwise_distances(X, metric=self.metric, squared=True)
+                distances = safe_pairwise_distances(X, metric=self.metric, squared=True)
             else:
-                distances = pairwise_distances(X, metric=self.metric)
+                distances = safe_pairwise_distances(X, metric=self.metric)
 
         # Degrees of freedom of the Student's t-distribution. The suggestion
         # alpha = n_components - 1 comes from "Learning a Parametric Embedding
diff --git a/sklearn/manifold/tests/test_t_sne.py b/sklearn/manifold/tests/test_t_sne.py
@@ -271,3 +271,21 @@ def test_reduction_to_one_component():
     X = random_state.randn(5, 2)
     X_embedded = tsne.fit_transform(X)
     assert(np.all(np.isfinite(X_embedded)))
+
+
+def test_undefined_correlation():
+    #t_SNE throws an exception with undefined correlation (issue #4475)
+    from sklearn.manifold import TSNE
+    import numpy as np
+    np.random.seed(42)
+
+    data = np.random.rand(10, 3)
+    data[-1, :] = 0
+
+    try:
+        model = TSNE(metric="correlation")
+        model.fit_transform(data)
+        assert True
+    except ValueError as e:
+        assert False
+
diff --git a/sklearn/metrics/pairwise.py b/sklearn/metrics/pairwise.py
@@ -1113,6 +1113,32 @@ def pairwise_distances(X, Y=None, metric="euclidean", n_jobs=1, **kwds):
     return _parallel_pairwise(X, Y, func, n_jobs, **kwds)
 
 
+def safe_pairwise_distances(X, Y=None, metric="euclidean", n_jobs=1, **kwds):
+    """ Compute a 'safe' distance matrix from a vector array X and optional Y, especially when metric is 'correlation'.
+
+    When metric is 'correlation', scipy.distance.pdist returns NaN values when a sample has zero variance.
+    Those values will be replaced by 1.0.
+
+    See also
+    --------
+    pairwise_distances
+    """
+    distances = pairwise_distances(X, Y, metric, n_jobs, **kwds)
+    if metric == 'correlation':
+        invalid_peak_to_peak_value = 0.0
+        correlation_for_invalid_value = 1.0  # maximum value of scipy.spatial.distance.correlation
+        invalid_rows_X = np.where(np.ptp(X, axis=1) == invalid_peak_to_peak_value)
+        distances[invalid_rows_X, :] = correlation_for_invalid_value
+
+        if Y is not None:
+            invalid_rows_Y = np.where(np.ptp(Y, axis=1) == invalid_peak_to_peak_value)
+            distances[:, invalid_rows_Y] = correlation_for_invalid_value
+        else:
+            distances[:, invalid_rows_X] = correlation_for_invalid_value
+
+    return distances
+
+
 # Helper functions - distance
 PAIRWISE_KERNEL_FUNCTIONS = {
     # If updating this dictionary, update the doc in both distance_metrics()
