Aluode/PerceptionLabPortable
0
1import functools
2import warnings
3from typing import Any, List
4
5import numpy as np
6import pytest
7import scipy.sparse as sp
8
9from sklearn.exceptions import DataDimensionalityWarning, NotFittedError
10from sklearn.metrics import euclidean_distances
11from sklearn.random_projection import (
12 GaussianRandomProjection,
13 SparseRandomProjection,
14 _gaussian_random_matrix,
15 _sparse_random_matrix,
16 johnson_lindenstrauss_min_dim,
17)
18from sklearn.utils._testing import (
19 assert_allclose,
20 assert_allclose_dense_sparse,
21 assert_almost_equal,
22 assert_array_almost_equal,
23 assert_array_equal,
24)
25from sklearn.utils.fixes import COO_CONTAINERS
26
27all_sparse_random_matrix: List[Any] = [_sparse_random_matrix]
28all_dense_random_matrix: List[Any] = [_gaussian_random_matrix]
29all_random_matrix = all_sparse_random_matrix + all_dense_random_matrix
30
31all_SparseRandomProjection: List[Any] = [SparseRandomProjection]
32all_DenseRandomProjection: List[Any] = [GaussianRandomProjection]
33all_RandomProjection = all_SparseRandomProjection + all_DenseRandomProjection
34
35
36def make_sparse_random_data(
37 coo_container,
38 n_samples,
39 n_features,
40 n_nonzeros,
41 random_state=None,
42 sparse_format="csr",
43):
44 """Make some random data with uniformly located non zero entries with
45 Gaussian distributed values; `sparse_format` can be `"csr"` (default) or
46 `None` (in which case a dense array is returned).
47 """
48 rng = np.random.RandomState(random_state)
49 data_coo = coo_container(
50 (
51 rng.randn(n_nonzeros),
52 (
53 rng.randint(n_samples, size=n_nonzeros),
54 rng.randint(n_features, size=n_nonzeros),
55 ),
56 ),
57 shape=(n_samples, n_features),
58 )
59 if sparse_format is not None:
60 return data_coo.asformat(sparse_format)
61 else:
62 return data_coo.toarray()
63
64
65def densify(matrix):
66 if not sp.issparse(matrix):
67 return matrix
68 else:
69 return matrix.toarray()
70
71
72n_samples, n_features = (10, 1000)
73n_nonzeros = int(n_samples * n_features / 100.0)
74
75
76###############################################################################
77# test on JL lemma
78###############################################################################
79
80
81@pytest.mark.parametrize(
82 "n_samples, eps",
83 [
84 ([100, 110], [0.9, 1.1]),
85 ([90, 100], [0.1, 0.0]),
86 ([50, -40], [0.1, 0.2]),
87 ],
88)
89def test_invalid_jl_domain(n_samples, eps):
90 with pytest.raises(ValueError):
91 johnson_lindenstrauss_min_dim(n_samples, eps=eps)
92
93
94def test_input_size_jl_min_dim():
95 with pytest.raises(ValueError):
96 johnson_lindenstrauss_min_dim(3 * [100], eps=2 * [0.9])
97
98 johnson_lindenstrauss_min_dim(
99 np.random.randint(1, 10, size=(10, 10)), eps=np.full((10, 10), 0.5)
100 )
101
102
103###############################################################################
104# tests random matrix generation
105###############################################################################
106def check_input_size_random_matrix(random_matrix):
107 inputs = [(0, 0), (-1, 1), (1, -1), (1, 0), (-1, 0)]
108 for n_components, n_features in inputs:
109 with pytest.raises(ValueError):
110 random_matrix(n_components, n_features)
111
112
113def check_size_generated(random_matrix):
114 inputs = [(1, 5), (5, 1), (5, 5), (1, 1)]
115 for n_components, n_features in inputs:
116 assert random_matrix(n_components, n_features).shape == (
117 n_components,
118 n_features,
119 )
120
121
122def check_zero_mean_and_unit_norm(random_matrix):
123 # All random matrix should produce a transformation matrix
124 # with zero mean and unit norm for each columns
125
126 A = densify(random_matrix(10000, 1, random_state=0))
127
128 assert_array_almost_equal(0, np.mean(A), 3)
129 assert_array_almost_equal(1.0, np.linalg.norm(A), 1)
130
131
132def check_input_with_sparse_random_matrix(random_matrix):
133 n_components, n_features = 5, 10
134
135 for density in [-1.0, 0.0, 1.1]:
136 with pytest.raises(ValueError):
137 random_matrix(n_components, n_features, density=density)
138
139
140@pytest.mark.parametrize("random_matrix", all_random_matrix)
141def test_basic_property_of_random_matrix(random_matrix):
142 # Check basic properties of random matrix generation
143 check_input_size_random_matrix(random_matrix)
144 check_size_generated(random_matrix)
145 check_zero_mean_and_unit_norm(random_matrix)
146
147
148@pytest.mark.parametrize("random_matrix", all_sparse_random_matrix)
149def test_basic_property_of_sparse_random_matrix(random_matrix):
150 check_input_with_sparse_random_matrix(random_matrix)
151
152 random_matrix_dense = functools.partial(random_matrix, density=1.0)
153
154 check_zero_mean_and_unit_norm(random_matrix_dense)
155
156
157def test_gaussian_random_matrix():
158 # Check some statical properties of Gaussian random matrix
159 # Check that the random matrix follow the proper distribution.
160 # Let's say that each element of a_{ij} of A is taken from
161 # a_ij ~ N(0.0, 1 / n_components).
162 #
163 n_components = 100
164 n_features = 1000
165 A = _gaussian_random_matrix(n_components, n_features, random_state=0)
166
167 assert_array_almost_equal(0.0, np.mean(A), 2)
168 assert_array_almost_equal(np.var(A, ddof=1), 1 / n_components, 1)
169
170
171def test_sparse_random_matrix():
172 # Check some statical properties of sparse random matrix
173 n_components = 100
174 n_features = 500
175
176 for density in [0.3, 1.0]:
177 s = 1 / density
178
179 A = _sparse_random_matrix(
180 n_components, n_features, density=density, random_state=0
181 )
182 A = densify(A)
183
184 # Check possible values
185 values = np.unique(A)
186 assert np.sqrt(s) / np.sqrt(n_components) in values
187 assert -np.sqrt(s) / np.sqrt(n_components) in values
188
189 if density == 1.0:
190 assert np.size(values) == 2
191 else:
192 assert 0.0 in values
193 assert np.size(values) == 3
194
195 # Check that the random matrix follow the proper distribution.
196 # Let's say that each element of a_{ij} of A is taken from
197 #
198 # - -sqrt(s) / sqrt(n_components) with probability 1 / 2s
199 # - 0 with probability 1 - 1 / s
200 # - +sqrt(s) / sqrt(n_components) with probability 1 / 2s
201 #
202 assert_almost_equal(np.mean(A == 0.0), 1 - 1 / s, decimal=2)
203 assert_almost_equal(
204 np.mean(A == np.sqrt(s) / np.sqrt(n_components)), 1 / (2 * s), decimal=2
205 )
206 assert_almost_equal(
207 np.mean(A == -np.sqrt(s) / np.sqrt(n_components)), 1 / (2 * s), decimal=2
208 )
209
210 assert_almost_equal(np.var(A == 0.0, ddof=1), (1 - 1 / s) * 1 / s, decimal=2)
211 assert_almost_equal(
212 np.var(A == np.sqrt(s) / np.sqrt(n_components), ddof=1),
213 (1 - 1 / (2 * s)) * 1 / (2 * s),
214 decimal=2,
215 )
216 assert_almost_equal(
217 np.var(A == -np.sqrt(s) / np.sqrt(n_components), ddof=1),
218 (1 - 1 / (2 * s)) * 1 / (2 * s),
219 decimal=2,
220 )
221
222
223###############################################################################
224# tests on random projection transformer
225###############################################################################
226
227
228def test_random_projection_transformer_invalid_input():
229 n_components = "auto"
230 fit_data = [[0, 1, 2]]
231 for RandomProjection in all_RandomProjection:
232 with pytest.raises(ValueError):
233 RandomProjection(n_components=n_components).fit(fit_data)
234
235
236@pytest.mark.parametrize("coo_container", COO_CONTAINERS)
237def test_try_to_transform_before_fit(coo_container, global_random_seed):
238 data = make_sparse_random_data(
239 coo_container,
240 n_samples,
241 n_features,
242 n_nonzeros,
243 random_state=global_random_seed,
244 sparse_format=None,
245 )
246 for RandomProjection in all_RandomProjection:
247 with pytest.raises(NotFittedError):
248 RandomProjection(n_components="auto").transform(data)
249
250
251@pytest.mark.parametrize("coo_container", COO_CONTAINERS)
252def test_too_many_samples_to_find_a_safe_embedding(coo_container, global_random_seed):
253 data = make_sparse_random_data(
254 coo_container,
255 n_samples=1000,
256 n_features=100,
257 n_nonzeros=1000,
258 random_state=global_random_seed,
259 sparse_format=None,
260 )
261
262 for RandomProjection in all_RandomProjection:
263 rp = RandomProjection(n_components="auto", eps=0.1)
264 expected_msg = (
265 "eps=0.100000 and n_samples=1000 lead to a target dimension"
266 " of 5920 which is larger than the original space with"
267 " n_features=100"
268 )
269 with pytest.raises(ValueError, match=expected_msg):
270 rp.fit(data)
271
272
273@pytest.mark.parametrize("coo_container", COO_CONTAINERS)
274def test_random_projection_embedding_quality(coo_container):
275 data = make_sparse_random_data(
276 coo_container,
277 n_samples=8,
278 n_features=5000,
279 n_nonzeros=15000,
280 random_state=0,
281 sparse_format=None,
282 )
283 eps = 0.2
284
285 original_distances = euclidean_distances(data, squared=True)
286 original_distances = original_distances.ravel()
287 non_identical = original_distances != 0.0
288
289 # remove 0 distances to avoid division by 0
290 original_distances = original_distances[non_identical]
291
292 for RandomProjection in all_RandomProjection:
293 rp = RandomProjection(n_components="auto", eps=eps, random_state=0)
294 projected = rp.fit_transform(data)
295
296 projected_distances = euclidean_distances(projected, squared=True)
297 projected_distances = projected_distances.ravel()
298
299 # remove 0 distances to avoid division by 0
300 projected_distances = projected_distances[non_identical]
301
302 distances_ratio = projected_distances / original_distances
303
304 # check that the automatically tuned values for the density respect the
305 # contract for eps: pairwise distances are preserved according to the
306 # Johnson-Lindenstrauss lemma
307 assert distances_ratio.max() < 1 + eps
308 assert 1 - eps < distances_ratio.min()
309
310
311@pytest.mark.parametrize("coo_container", COO_CONTAINERS)
312def test_SparseRandomProj_output_representation(coo_container):
313 dense_data = make_sparse_random_data(
314 coo_container,
315 n_samples,
316 n_features,
317 n_nonzeros,
318 random_state=0,
319 sparse_format=None,
320 )
321 sparse_data = make_sparse_random_data(
322 coo_container,
323 n_samples,
324 n_features,
325 n_nonzeros,
326 random_state=0,
327 sparse_format="csr",
328 )
329 for SparseRandomProj in all_SparseRandomProjection:
330 # when using sparse input, the projected data can be forced to be a
331 # dense numpy array
332 rp = SparseRandomProj(n_components=10, dense_output=True, random_state=0)
333 rp.fit(dense_data)
334 assert isinstance(rp.transform(dense_data), np.ndarray)
335 assert isinstance(rp.transform(sparse_data), np.ndarray)
336
337 # the output can be left to a sparse matrix instead
338 rp = SparseRandomProj(n_components=10, dense_output=False, random_state=0)
339 rp = rp.fit(dense_data)
340 # output for dense input will stay dense:
341 assert isinstance(rp.transform(dense_data), np.ndarray)
342
343 # output for sparse output will be sparse:
344 assert sp.issparse(rp.transform(sparse_data))
345
346
347@pytest.mark.parametrize("coo_container", COO_CONTAINERS)
348def test_correct_RandomProjection_dimensions_embedding(
349 coo_container, global_random_seed
350):
351 data = make_sparse_random_data(
352 coo_container,
353 n_samples,
354 n_features,
355 n_nonzeros,
356 random_state=global_random_seed,
357 sparse_format=None,
358 )
359 for RandomProjection in all_RandomProjection:
360 rp = RandomProjection(n_components="auto", random_state=0, eps=0.5).fit(data)
361
362 # the number of components is adjusted from the shape of the training
363 # set
364 assert rp.n_components == "auto"
365 assert rp.n_components_ == 110
366
367 if RandomProjection in all_SparseRandomProjection:
368 assert rp.density == "auto"
369 assert_almost_equal(rp.density_, 0.03, 2)
370
371 assert rp.components_.shape == (110, n_features)
372
373 projected_1 = rp.transform(data)
374 assert projected_1.shape == (n_samples, 110)
375
376 # once the RP is 'fitted' the projection is always the same
377 projected_2 = rp.transform(data)
378 assert_array_equal(projected_1, projected_2)
379
380 # fit transform with same random seed will lead to the same results
381 rp2 = RandomProjection(random_state=0, eps=0.5)
382 projected_3 = rp2.fit_transform(data)
383 assert_array_equal(projected_1, projected_3)
384
385 # Try to transform with an input X of size different from fitted.
386 with pytest.raises(ValueError):
387 rp.transform(data[:, 1:5])
388
389 # it is also possible to fix the number of components and the density
390 # level
391 if RandomProjection in all_SparseRandomProjection:
392 rp = RandomProjection(n_components=100, density=0.001, random_state=0)
393 projected = rp.fit_transform(data)
394 assert projected.shape == (n_samples, 100)
395 assert rp.components_.shape == (100, n_features)
396 assert rp.components_.nnz < 115 # close to 1% density
397 assert 85 < rp.components_.nnz # close to 1% density
398
399
400@pytest.mark.parametrize("coo_container", COO_CONTAINERS)
401def test_warning_n_components_greater_than_n_features(
402 coo_container, global_random_seed
403):
404 n_features = 20
405 n_samples = 5
406 n_nonzeros = int(n_features / 4)
407 data = make_sparse_random_data(
408 coo_container,
409 n_samples,
410 n_features,
411 n_nonzeros,
412 random_state=global_random_seed,
413 sparse_format=None,
414 )
415
416 for RandomProjection in all_RandomProjection:
417 with pytest.warns(DataDimensionalityWarning):
418 RandomProjection(n_components=n_features + 1).fit(data)
419
420
421@pytest.mark.parametrize("coo_container", COO_CONTAINERS)
422def test_works_with_sparse_data(coo_container, global_random_seed):
423 n_features = 20
424 n_samples = 5
425 n_nonzeros = int(n_features / 4)
426 dense_data = make_sparse_random_data(
427 coo_container,
428 n_samples,
429 n_features,
430 n_nonzeros,
431 random_state=global_random_seed,
432 sparse_format=None,
433 )
434 sparse_data = make_sparse_random_data(
435 coo_container,
436 n_samples,
437 n_features,
438 n_nonzeros,
439 random_state=global_random_seed,
440 sparse_format="csr",
441 )
442
443 for RandomProjection in all_RandomProjection:
444 rp_dense = RandomProjection(n_components=3, random_state=1).fit(dense_data)
445 rp_sparse = RandomProjection(n_components=3, random_state=1).fit(sparse_data)
446 assert_array_almost_equal(
447 densify(rp_dense.components_), densify(rp_sparse.components_)
448 )
449
450
451def test_johnson_lindenstrauss_min_dim():
452 """Test Johnson-Lindenstrauss for small eps.
453
454 Regression test for #17111: before #19374, 32-bit systems would fail.
455 """
456 assert johnson_lindenstrauss_min_dim(100, eps=1e-5) == 368416070986
457
458
459@pytest.mark.parametrize("coo_container", COO_CONTAINERS)
460@pytest.mark.parametrize("random_projection_cls", all_RandomProjection)
461def test_random_projection_feature_names_out(
462 coo_container, random_projection_cls, global_random_seed
463):
464 data = make_sparse_random_data(
465 coo_container,
466 n_samples,
467 n_features,
468 n_nonzeros,
469 random_state=global_random_seed,
470 sparse_format=None,
471 )
472 random_projection = random_projection_cls(n_components=2)
473 random_projection.fit(data)
474 names_out = random_projection.get_feature_names_out()
475 class_name_lower = random_projection_cls.__name__.lower()
476 expected_names_out = np.array(
477 [f"{class_name_lower}{i}" for i in range(random_projection.n_components_)],
478 dtype=object,
479 )
480
481 assert_array_equal(names_out, expected_names_out)
482
483
484@pytest.mark.parametrize("coo_container", COO_CONTAINERS)
485@pytest.mark.parametrize("n_samples", (2, 9, 10, 11, 1000))
486@pytest.mark.parametrize("n_features", (2, 9, 10, 11, 1000))
487@pytest.mark.parametrize("random_projection_cls", all_RandomProjection)
488@pytest.mark.parametrize("compute_inverse_components", [True, False])
489def test_inverse_transform(
490 coo_container,
491 n_samples,
492 n_features,
493 random_projection_cls,
494 compute_inverse_components,
495 global_random_seed,
496):
497 n_components = 10
498
499 random_projection = random_projection_cls(
500 n_components=n_components,
501 compute_inverse_components=compute_inverse_components,
502 random_state=global_random_seed,
503 )
504
505 X_dense = make_sparse_random_data(
506 coo_container,
507 n_samples,
508 n_features,
509 n_nonzeros=n_samples * n_features // 100 + 1,
510 random_state=global_random_seed,
511 sparse_format=None,
512 )
513 X_csr = make_sparse_random_data(
514 coo_container,
515 n_samples,
516 n_features,
517 n_nonzeros=n_samples * n_features // 100 + 1,
518 random_state=global_random_seed,
519 sparse_format="csr",
520 )
521
522 for X in [X_dense, X_csr]:
523 with warnings.catch_warnings():
524 warnings.filterwarnings(
525 "ignore",
526 message=(
527 "The number of components is higher than the number of features"
528 ),
529 category=DataDimensionalityWarning,
530 )
531 projected = random_projection.fit_transform(X)
532
533 if compute_inverse_components:
534 assert hasattr(random_projection, "inverse_components_")
535 inv_components = random_projection.inverse_components_
536 assert inv_components.shape == (n_features, n_components)
537
538 projected_back = random_projection.inverse_transform(projected)
539 assert projected_back.shape == X.shape
540
541 projected_again = random_projection.transform(projected_back)
542 if hasattr(projected, "toarray"):
543 projected = projected.toarray()
544 assert_allclose(projected, projected_again, rtol=1e-7, atol=1e-10)
545
546
547@pytest.mark.parametrize("random_projection_cls", all_RandomProjection)
548@pytest.mark.parametrize(
549 "input_dtype, expected_dtype",
550 (
551 (np.float32, np.float32),
552 (np.float64, np.float64),
553 (np.int32, np.float64),
554 (np.int64, np.float64),
555 ),
556)
557def test_random_projection_dtype_match(
558 random_projection_cls, input_dtype, expected_dtype
559):
560 # Verify output matrix dtype
561 rng = np.random.RandomState(42)
562 X = rng.rand(25, 3000)
563 rp = random_projection_cls(random_state=0)
564 transformed = rp.fit_transform(X.astype(input_dtype))
565
566 assert rp.components_.dtype == expected_dtype
567 assert transformed.dtype == expected_dtype
568
569
570@pytest.mark.parametrize("random_projection_cls", all_RandomProjection)
571def test_random_projection_numerical_consistency(random_projection_cls):
572 # Verify numerical consistency among np.float32 and np.float64
573 atol = 1e-5
574 rng = np.random.RandomState(42)
575 X = rng.rand(25, 3000)
576 rp_32 = random_projection_cls(random_state=0)
577 rp_64 = random_projection_cls(random_state=0)
578
579 projection_32 = rp_32.fit_transform(X.astype(np.float32))
580 projection_64 = rp_64.fit_transform(X.astype(np.float64))
581
582 assert_allclose(projection_64, projection_32, atol=atol)
583
584 assert_allclose_dense_sparse(rp_32.components_, rp_64.components_)
585