Aluode/PerceptionLabPortable
0
1"""Linear and quadratic discriminant analysis."""
2
3# Authors: The scikit-learn developers
4# SPDX-License-Identifier: BSD-3-Clause
5
6import warnings
7from numbers import Integral, Real
8
9import numpy as np
10import scipy.linalg
11from scipy import linalg
12
13from .base import (
14 BaseEstimator,
15 ClassifierMixin,
16 ClassNamePrefixFeaturesOutMixin,
17 TransformerMixin,
18 _fit_context,
19)
20from .covariance import empirical_covariance, ledoit_wolf, shrunk_covariance
21from .linear_model._base import LinearClassifierMixin
22from .preprocessing import StandardScaler
23from .utils._array_api import _expit, device, get_namespace, size
24from .utils._param_validation import HasMethods, Interval, StrOptions
25from .utils.extmath import softmax
26from .utils.multiclass import check_classification_targets, unique_labels
27from .utils.validation import check_is_fitted, validate_data
28
29__all__ = ["LinearDiscriminantAnalysis", "QuadraticDiscriminantAnalysis"]
30
31
32def _cov(X, shrinkage=None, covariance_estimator=None):
33 """Estimate covariance matrix (using optional covariance_estimator).
34 Parameters
35 ----------
36 X : array-like of shape (n_samples, n_features)
37 Input data.
38
39 shrinkage : {'empirical', 'auto'} or float, default=None
40 Shrinkage parameter, possible values:
41 - None or 'empirical': no shrinkage (default).
42 - 'auto': automatic shrinkage using the Ledoit-Wolf lemma.
43 - float between 0 and 1: fixed shrinkage parameter.
44
45 Shrinkage parameter is ignored if `covariance_estimator`
46 is not None.
47
48 covariance_estimator : estimator, default=None
49 If not None, `covariance_estimator` is used to estimate
50 the covariance matrices instead of relying on the empirical
51 covariance estimator (with potential shrinkage).
52 The object should have a fit method and a ``covariance_`` attribute
53 like the estimators in :mod:`sklearn.covariance``.
54 if None the shrinkage parameter drives the estimate.
55
56 .. versionadded:: 0.24
57
58 Returns
59 -------
60 s : ndarray of shape (n_features, n_features)
61 Estimated covariance matrix.
62 """
63 if covariance_estimator is None:
64 shrinkage = "empirical" if shrinkage is None else shrinkage
65 if isinstance(shrinkage, str):
66 if shrinkage == "auto":
67 sc = StandardScaler() # standardize features
68 X = sc.fit_transform(X)
69 s = ledoit_wolf(X)[0]
70 # rescale
71 s = sc.scale_[:, np.newaxis] * s * sc.scale_[np.newaxis, :]
72 elif shrinkage == "empirical":
73 s = empirical_covariance(X)
74 elif isinstance(shrinkage, Real):
75 s = shrunk_covariance(empirical_covariance(X), shrinkage)
76 else:
77 if shrinkage is not None and shrinkage != 0:
78 raise ValueError(
79 "covariance_estimator and shrinkage parameters "
80 "are not None. Only one of the two can be set."
81 )
82 covariance_estimator.fit(X)
83 if not hasattr(covariance_estimator, "covariance_"):
84 raise ValueError(
85 "%s does not have a covariance_ attribute"
86 % covariance_estimator.__class__.__name__
87 )
88 s = covariance_estimator.covariance_
89 return s
90
91
92def _class_means(X, y):
93 """Compute class means.
94
95 Parameters
96 ----------
97 X : array-like of shape (n_samples, n_features)
98 Input data.
99
100 y : array-like of shape (n_samples,) or (n_samples, n_targets)
101 Target values.
102
103 Returns
104 -------
105 means : array-like of shape (n_classes, n_features)
106 Class means.
107 """
108 xp, is_array_api_compliant = get_namespace(X)
109 classes, y = xp.unique_inverse(y)
110 means = xp.zeros((classes.shape[0], X.shape[1]), device=device(X), dtype=X.dtype)
111
112 if is_array_api_compliant:
113 for i in range(classes.shape[0]):
114 means[i, :] = xp.mean(X[y == i], axis=0)
115 else:
116 # TODO: Explore the choice of using bincount + add.at as it seems sub optimal
117 # from a performance-wise
118 cnt = np.bincount(y)
119 np.add.at(means, y, X)
120 means /= cnt[:, None]
121 return means
122
123
124def _class_cov(X, y, priors, shrinkage=None, covariance_estimator=None):
125 """Compute weighted within-class covariance matrix.
126
127 The per-class covariance are weighted by the class priors.
128
129 Parameters
130 ----------
131 X : array-like of shape (n_samples, n_features)
132 Input data.
133
134 y : array-like of shape (n_samples,) or (n_samples, n_targets)
135 Target values.
136
137 priors : array-like of shape (n_classes,)
138 Class priors.
139
140 shrinkage : 'auto' or float, default=None
141 Shrinkage parameter, possible values:
142 - None: no shrinkage (default).
143 - 'auto': automatic shrinkage using the Ledoit-Wolf lemma.
144 - float between 0 and 1: fixed shrinkage parameter.
145
146 Shrinkage parameter is ignored if `covariance_estimator` is not None.
147
148 covariance_estimator : estimator, default=None
149 If not None, `covariance_estimator` is used to estimate
150 the covariance matrices instead of relying the empirical
151 covariance estimator (with potential shrinkage).
152 The object should have a fit method and a ``covariance_`` attribute
153 like the estimators in sklearn.covariance.
154 If None, the shrinkage parameter drives the estimate.
155
156 .. versionadded:: 0.24
157
158 Returns
159 -------
160 cov : array-like of shape (n_features, n_features)
161 Weighted within-class covariance matrix
162 """
163 classes = np.unique(y)
164 cov = np.zeros(shape=(X.shape[1], X.shape[1]))
165 for idx, group in enumerate(classes):
166 Xg = X[y == group, :]
167 cov += priors[idx] * np.atleast_2d(_cov(Xg, shrinkage, covariance_estimator))
168 return cov
169
170
171class DiscriminantAnalysisPredictionMixin:
172 """Mixin class for QuadraticDiscriminantAnalysis and NearestCentroid."""
173
174 def decision_function(self, X):
175 """Apply decision function to an array of samples.
176
177 Parameters
178 ----------
179 X : {array-like, sparse matrix} of shape (n_samples, n_features)
180 Array of samples (test vectors).
181
182 Returns
183 -------
184 y_scores : ndarray of shape (n_samples,) or (n_samples, n_classes)
185 Decision function values related to each class, per sample.
186 In the two-class case, the shape is `(n_samples,)`, giving the
187 log likelihood ratio of the positive class.
188 """
189 y_scores = self._decision_function(X)
190 if len(self.classes_) == 2:
191 return y_scores[:, 1] - y_scores[:, 0]
192 return y_scores
193
194 def predict(self, X):
195 """Perform classification on an array of vectors `X`.
196
197 Returns the class label for each sample.
198
199 Parameters
200 ----------
201 X : {array-like, sparse matrix} of shape (n_samples, n_features)
202 Input vectors, where `n_samples` is the number of samples and
203 `n_features` is the number of features.
204
205 Returns
206 -------
207 y_pred : ndarray of shape (n_samples,)
208 Class label for each sample.
209 """
210 scores = self._decision_function(X)
211 return self.classes_.take(scores.argmax(axis=1))
212
213 def predict_proba(self, X):
214 """Estimate class probabilities.
215
216 Parameters
217 ----------
218 X : {array-like, sparse matrix} of shape (n_samples, n_features)
219 Input data.
220
221 Returns
222 -------
223 y_proba : ndarray of shape (n_samples, n_classes)
224 Probability estimate of the sample for each class in the
225 model, where classes are ordered as they are in `self.classes_`.
226 """
227 return np.exp(self.predict_log_proba(X))
228
229 def predict_log_proba(self, X):
230 """Estimate log class probabilities.
231
232 Parameters
233 ----------
234 X : {array-like, sparse matrix} of shape (n_samples, n_features)
235 Input data.
236
237 Returns
238 -------
239 y_log_proba : ndarray of shape (n_samples, n_classes)
240 Estimated log probabilities.
241 """
242 scores = self._decision_function(X)
243 log_likelihood = scores - scores.max(axis=1)[:, np.newaxis]
244 return log_likelihood - np.log(
245 np.exp(log_likelihood).sum(axis=1)[:, np.newaxis]
246 )
247
248
249class LinearDiscriminantAnalysis(
250 ClassNamePrefixFeaturesOutMixin,
251 LinearClassifierMixin,
252 TransformerMixin,
253 BaseEstimator,
254):
255 """Linear Discriminant Analysis.
256
257 A classifier with a linear decision boundary, generated by fitting class
258 conditional densities to the data and using Bayes' rule.
259
260 The model fits a Gaussian density to each class, assuming that all classes
261 share the same covariance matrix.
262
263 The fitted model can also be used to reduce the dimensionality of the input
264 by projecting it to the most discriminative directions, using the
265 `transform` method.
266
267 .. versionadded:: 0.17
268
269 For a comparison between
270 :class:`~sklearn.discriminant_analysis.LinearDiscriminantAnalysis`
271 and :class:`~sklearn.discriminant_analysis.QuadraticDiscriminantAnalysis`, see
272 :ref:`sphx_glr_auto_examples_classification_plot_lda_qda.py`.
273
274 Read more in the :ref:`User Guide <lda_qda>`.
275
276 Parameters
277 ----------
278 solver : {'svd', 'lsqr', 'eigen'}, default='svd'
279 Solver to use, possible values:
280 - 'svd': Singular value decomposition (default).
281 Does not compute the covariance matrix, therefore this solver is
282 recommended for data with a large number of features.
283 - 'lsqr': Least squares solution.
284 Can be combined with shrinkage or custom covariance estimator.
285 - 'eigen': Eigenvalue decomposition.
286 Can be combined with shrinkage or custom covariance estimator.
287
288 .. versionchanged:: 1.2
289 `solver="svd"` now has experimental Array API support. See the
290 :ref:`Array API User Guide <array_api>` for more details.
291
292 shrinkage : 'auto' or float, default=None
293 Shrinkage parameter, possible values:
294 - None: no shrinkage (default).
295 - 'auto': automatic shrinkage using the Ledoit-Wolf lemma.
296 - float between 0 and 1: fixed shrinkage parameter.
297
298 This should be left to None if `covariance_estimator` is used.
299 Note that shrinkage works only with 'lsqr' and 'eigen' solvers.
300
301 For a usage example, see
302 :ref:`sphx_glr_auto_examples_classification_plot_lda.py`.
303
304 priors : array-like of shape (n_classes,), default=None
305 The class prior probabilities. By default, the class proportions are
306 inferred from the training data.
307
308 n_components : int, default=None
309 Number of components (<= min(n_classes - 1, n_features)) for
310 dimensionality reduction. If None, will be set to
311 min(n_classes - 1, n_features). This parameter only affects the
312 `transform` method.
313
314 For a usage example, see
315 :ref:`sphx_glr_auto_examples_decomposition_plot_pca_vs_lda.py`.
316
317 store_covariance : bool, default=False
318 If True, explicitly compute the weighted within-class covariance
319 matrix when solver is 'svd'. The matrix is always computed
320 and stored for the other solvers.
321
322 .. versionadded:: 0.17
323
324 tol : float, default=1.0e-4
325 Absolute threshold for a singular value of X to be considered
326 significant, used to estimate the rank of X. Dimensions whose
327 singular values are non-significant are discarded. Only used if
328 solver is 'svd'.
329
330 .. versionadded:: 0.17
331
332 covariance_estimator : covariance estimator, default=None
333 If not None, `covariance_estimator` is used to estimate
334 the covariance matrices instead of relying on the empirical
335 covariance estimator (with potential shrinkage).
336 The object should have a fit method and a ``covariance_`` attribute
337 like the estimators in :mod:`sklearn.covariance`.
338 if None the shrinkage parameter drives the estimate.
339
340 This should be left to None if `shrinkage` is used.
341 Note that `covariance_estimator` works only with 'lsqr' and 'eigen'
342 solvers.
343
344 .. versionadded:: 0.24
345
346 Attributes
347 ----------
348 coef_ : ndarray of shape (n_features,) or (n_classes, n_features)
349 Weight vector(s).
350
351 intercept_ : ndarray of shape (n_classes,)
352 Intercept term.
353
354 covariance_ : array-like of shape (n_features, n_features)
355 Weighted within-class covariance matrix. It corresponds to
356 `sum_k prior_k * C_k` where `C_k` is the covariance matrix of the
357 samples in class `k`. The `C_k` are estimated using the (potentially
358 shrunk) biased estimator of covariance. If solver is 'svd', only
359 exists when `store_covariance` is True.
360
361 explained_variance_ratio_ : ndarray of shape (n_components,)
362 Percentage of variance explained by each of the selected components.
363 If ``n_components`` is not set then all components are stored and the
364 sum of explained variances is equal to 1.0. Only available when eigen
365 or svd solver is used.
366
367 means_ : array-like of shape (n_classes, n_features)
368 Class-wise means.
369
370 priors_ : array-like of shape (n_classes,)
371 Class priors (sum to 1).
372
373 scalings_ : array-like of shape (rank, n_classes - 1)
374 Scaling of the features in the space spanned by the class centroids.
375 Only available for 'svd' and 'eigen' solvers.
376
377 xbar_ : array-like of shape (n_features,)
378 Overall mean. Only present if solver is 'svd'.
379
380 classes_ : array-like of shape (n_classes,)
381 Unique class labels.
382
383 n_features_in_ : int
384 Number of features seen during :term:`fit`.
385
386 .. versionadded:: 0.24
387
388 feature_names_in_ : ndarray of shape (`n_features_in_`,)
389 Names of features seen during :term:`fit`. Defined only when `X`
390 has feature names that are all strings.
391
392 .. versionadded:: 1.0
393
394 See Also
395 --------
396 QuadraticDiscriminantAnalysis : Quadratic Discriminant Analysis.
397
398 Examples
399 --------
400 >>> import numpy as np
401 >>> from sklearn.discriminant_analysis import LinearDiscriminantAnalysis
402 >>> X = np.array([[-1, -1], [-2, -1], [-3, -2], [1, 1], [2, 1], [3, 2]])
403 >>> y = np.array([1, 1, 1, 2, 2, 2])
404 >>> clf = LinearDiscriminantAnalysis()
405 >>> clf.fit(X, y)
406 LinearDiscriminantAnalysis()
407 >>> print(clf.predict([[-0.8, -1]]))
408 [1]
409 """
410
411 _parameter_constraints: dict = {
412 "solver": [StrOptions({"svd", "lsqr", "eigen"})],
413 "shrinkage": [StrOptions({"auto"}), Interval(Real, 0, 1, closed="both"), None],
414 "n_components": [Interval(Integral, 1, None, closed="left"), None],
415 "priors": ["array-like", None],
416 "store_covariance": ["boolean"],
417 "tol": [Interval(Real, 0, None, closed="left")],
418 "covariance_estimator": [HasMethods("fit"), None],
419 }
420
421 def __init__(
422 self,
423 solver="svd",
424 shrinkage=None,
425 priors=None,
426 n_components=None,
427 store_covariance=False,
428 tol=1e-4,
429 covariance_estimator=None,
430 ):
431 self.solver = solver
432 self.shrinkage = shrinkage
433 self.priors = priors
434 self.n_components = n_components
435 self.store_covariance = store_covariance # used only in svd solver
436 self.tol = tol # used only in svd solver
437 self.covariance_estimator = covariance_estimator
438
439 def _solve_lstsq(self, X, y, shrinkage, covariance_estimator):
440 """Least squares solver.
441
442 The least squares solver computes a straightforward solution of the
443 optimal decision rule based directly on the discriminant functions. It
444 can only be used for classification (with any covariance estimator),
445 because
446 estimation of eigenvectors is not performed. Therefore, dimensionality
447 reduction with the transform is not supported.
448
449 Parameters
450 ----------
451 X : array-like of shape (n_samples, n_features)
452 Training data.
453
454 y : array-like of shape (n_samples,) or (n_samples, n_classes)
455 Target values.
456
457 shrinkage : 'auto', float or None
458 Shrinkage parameter, possible values:
459 - None: no shrinkage.
460 - 'auto': automatic shrinkage using the Ledoit-Wolf lemma.
461 - float between 0 and 1: fixed shrinkage parameter.
462
463 Shrinkage parameter is ignored if `covariance_estimator` i
464 not None
465
466 covariance_estimator : estimator, default=None
467 If not None, `covariance_estimator` is used to estimate
468 the covariance matrices instead of relying the empirical
469 covariance estimator (with potential shrinkage).
470 The object should have a fit method and a ``covariance_`` attribute
471 like the estimators in sklearn.covariance.
472 if None the shrinkage parameter drives the estimate.
473
474 .. versionadded:: 0.24
475
476 Notes
477 -----
478 This solver is based on [1]_, section 2.6.2, pp. 39-41.
479
480 References
481 ----------
482 .. [1] R. O. Duda, P. E. Hart, D. G. Stork. Pattern Classification
483 (Second Edition). John Wiley & Sons, Inc., New York, 2001. ISBN
484 0-471-05669-3.
485 """
486 self.means_ = _class_means(X, y)
487 self.covariance_ = _class_cov(
488 X, y, self.priors_, shrinkage, covariance_estimator
489 )
490 self.coef_ = linalg.lstsq(self.covariance_, self.means_.T)[0].T
491 self.intercept_ = -0.5 * np.diag(np.dot(self.means_, self.coef_.T)) + np.log(
492 self.priors_
493 )
494
495 def _solve_eigen(self, X, y, shrinkage, covariance_estimator):
496 """Eigenvalue solver.
497
498 The eigenvalue solver computes the optimal solution of the Rayleigh
499 coefficient (basically the ratio of between class scatter to within
500 class scatter). This solver supports both classification and
501 dimensionality reduction (with any covariance estimator).
502
503 Parameters
504 ----------
505 X : array-like of shape (n_samples, n_features)
506 Training data.
507
508 y : array-like of shape (n_samples,) or (n_samples, n_targets)
509 Target values.
510
511 shrinkage : 'auto', float or None
512 Shrinkage parameter, possible values:
513 - None: no shrinkage.
514 - 'auto': automatic shrinkage using the Ledoit-Wolf lemma.
515 - float between 0 and 1: fixed shrinkage constant.
516
517 Shrinkage parameter is ignored if `covariance_estimator` i
518 not None
519
520 covariance_estimator : estimator, default=None
521 If not None, `covariance_estimator` is used to estimate
522 the covariance matrices instead of relying the empirical
523 covariance estimator (with potential shrinkage).
524 The object should have a fit method and a ``covariance_`` attribute
525 like the estimators in sklearn.covariance.
526 if None the shrinkage parameter drives the estimate.
527
528 .. versionadded:: 0.24
529
530 Notes
531 -----
532 This solver is based on [1]_, section 3.8.3, pp. 121-124.
533
534 References
535 ----------
536 .. [1] R. O. Duda, P. E. Hart, D. G. Stork. Pattern Classification
537 (Second Edition). John Wiley & Sons, Inc., New York, 2001. ISBN
538 0-471-05669-3.
539 """
540 self.means_ = _class_means(X, y)
541 self.covariance_ = _class_cov(
542 X, y, self.priors_, shrinkage, covariance_estimator
543 )
544
545 Sw = self.covariance_ # within scatter
546 St = _cov(X, shrinkage, covariance_estimator) # total scatter
547 Sb = St - Sw # between scatter
548
549 evals, evecs = linalg.eigh(Sb, Sw)
550 self.explained_variance_ratio_ = np.sort(evals / np.sum(evals))[::-1][
551 : self._max_components
552 ]
553 evecs = evecs[:, np.argsort(evals)[::-1]] # sort eigenvectors
554
555 self.scalings_ = evecs
556 self.coef_ = np.dot(self.means_, evecs).dot(evecs.T)
557 self.intercept_ = -0.5 * np.diag(np.dot(self.means_, self.coef_.T)) + np.log(
558 self.priors_
559 )
560
561 def _solve_svd(self, X, y):
562 """SVD solver.
563
564 Parameters
565 ----------
566 X : array-like of shape (n_samples, n_features)
567 Training data.
568
569 y : array-like of shape (n_samples,) or (n_samples, n_targets)
570 Target values.
571 """
572 xp, is_array_api_compliant = get_namespace(X)
573
574 if is_array_api_compliant:
575 svd = xp.linalg.svd
576 else:
577 svd = scipy.linalg.svd
578
579 n_samples, n_features = X.shape
580 n_classes = self.classes_.shape[0]
581
582 self.means_ = _class_means(X, y)
583 if self.store_covariance:
584 self.covariance_ = _class_cov(X, y, self.priors_)
585
586 Xc = []
587 for idx, group in enumerate(self.classes_):
588 Xg = X[y == group]
589 Xc.append(Xg - self.means_[idx, :])
590
591 self.xbar_ = self.priors_ @ self.means_
592
593 Xc = xp.concat(Xc, axis=0)
594
595 # 1) within (univariate) scaling by with classes std-dev
596 std = xp.std(Xc, axis=0)
597 # avoid division by zero in normalization
598 std[std == 0] = 1.0
599 fac = xp.asarray(1.0 / (n_samples - n_classes), dtype=X.dtype, device=device(X))
600
601 # 2) Within variance scaling
602 X = xp.sqrt(fac) * (Xc / std)
603 # SVD of centered (within)scaled data
604 U, S, Vt = svd(X, full_matrices=False)
605
606 rank = xp.sum(xp.astype(S > self.tol, xp.int32))
607 # Scaling of within covariance is: V' 1/S
608 scalings = (Vt[:rank, :] / std).T / S[:rank]
609 fac = 1.0 if n_classes == 1 else 1.0 / (n_classes - 1)
610
611 # 3) Between variance scaling
612 # Scale weighted centers
613 X = (
614 (xp.sqrt((n_samples * self.priors_) * fac)) * (self.means_ - self.xbar_).T
615 ).T @ scalings
616 # Centers are living in a space with n_classes-1 dim (maximum)
617 # Use SVD to find projection in the space spanned by the
618 # (n_classes) centers
619 _, S, Vt = svd(X, full_matrices=False)
620
621 if self._max_components == 0:
622 self.explained_variance_ratio_ = xp.empty((0,), dtype=S.dtype)
623 else:
624 self.explained_variance_ratio_ = (S**2 / xp.sum(S**2))[
625 : self._max_components
626 ]
627
628 rank = xp.sum(xp.astype(S > self.tol * S[0], xp.int32))
629 self.scalings_ = scalings @ Vt.T[:, :rank]
630 coef = (self.means_ - self.xbar_) @ self.scalings_
631 self.intercept_ = -0.5 * xp.sum(coef**2, axis=1) + xp.log(self.priors_)
632 self.coef_ = coef @ self.scalings_.T
633 self.intercept_ -= self.xbar_ @ self.coef_.T
634
635 @_fit_context(
636 # LinearDiscriminantAnalysis.covariance_estimator is not validated yet
637 prefer_skip_nested_validation=False
638 )
639 def fit(self, X, y):
640 """Fit the Linear Discriminant Analysis model.
641
642 .. versionchanged:: 0.19
643 `store_covariance` and `tol` has been moved to main constructor.
644
645 Parameters
646 ----------
647 X : array-like of shape (n_samples, n_features)
648 Training data.
649
650 y : array-like of shape (n_samples,)
651 Target values.
652
653 Returns
654 -------
655 self : object
656 Fitted estimator.
657 """
658 xp, _ = get_namespace(X)
659
660 X, y = validate_data(
661 self, X, y, ensure_min_samples=2, dtype=[xp.float64, xp.float32]
662 )
663 self.classes_ = unique_labels(y)
664 n_samples, _ = X.shape
665 n_classes = self.classes_.shape[0]
666
667 if n_samples == n_classes:
668 raise ValueError(
669 "The number of samples must be more than the number of classes."
670 )
671
672 if self.priors is None: # estimate priors from sample
673 _, cnts = xp.unique_counts(y) # non-negative ints
674 self.priors_ = xp.astype(cnts, X.dtype) / float(y.shape[0])
675 else:
676 self.priors_ = xp.asarray(self.priors, dtype=X.dtype)
677
678 if xp.any(self.priors_ < 0):
679 raise ValueError("priors must be non-negative")
680
681 if xp.abs(xp.sum(self.priors_) - 1.0) > 1e-5:
682 warnings.warn("The priors do not sum to 1. Renormalizing", UserWarning)
683 self.priors_ = self.priors_ / self.priors_.sum()
684
685 # Maximum number of components no matter what n_components is
686 # specified:
687 max_components = min(n_classes - 1, X.shape[1])
688
689 if self.n_components is None:
690 self._max_components = max_components
691 else:
692 if self.n_components > max_components:
693 raise ValueError(
694 "n_components cannot be larger than min(n_features, n_classes - 1)."
695 )
696 self._max_components = self.n_components
697
698 if self.solver == "svd":
699 if self.shrinkage is not None:
700 raise NotImplementedError("shrinkage not supported with 'svd' solver.")
701 if self.covariance_estimator is not None:
702 raise ValueError(
703 "covariance estimator "
704 "is not supported "
705 "with svd solver. Try another solver"
706 )
707 self._solve_svd(X, y)
708 elif self.solver == "lsqr":
709 self._solve_lstsq(
710 X,
711 y,
712 shrinkage=self.shrinkage,
713 covariance_estimator=self.covariance_estimator,
714 )
715 elif self.solver == "eigen":
716 self._solve_eigen(
717 X,
718 y,
719 shrinkage=self.shrinkage,
720 covariance_estimator=self.covariance_estimator,
721 )
722 if size(self.classes_) == 2: # treat binary case as a special case
723 coef_ = xp.asarray(self.coef_[1, :] - self.coef_[0, :], dtype=X.dtype)
724 self.coef_ = xp.reshape(coef_, (1, -1))
725 intercept_ = xp.asarray(
726 self.intercept_[1] - self.intercept_[0], dtype=X.dtype
727 )
728 self.intercept_ = xp.reshape(intercept_, (1,))
729 self._n_features_out = self._max_components
730 return self
731
732 def transform(self, X):
733 """Project data to maximize class separation.
734
735 Parameters
736 ----------
737 X : array-like of shape (n_samples, n_features)
738 Input data.
739
740 Returns
741 -------
742 X_new : ndarray of shape (n_samples, n_components) or \
743 (n_samples, min(rank, n_components))
744 Transformed data. In the case of the 'svd' solver, the shape
745 is (n_samples, min(rank, n_components)).
746 """
747 if self.solver == "lsqr":
748 raise NotImplementedError(
749 "transform not implemented for 'lsqr' solver (use 'svd' or 'eigen')."
750 )
751 check_is_fitted(self)
752 xp, _ = get_namespace(X)
753 X = validate_data(self, X, reset=False)
754
755 if self.solver == "svd":
756 X_new = (X - self.xbar_) @ self.scalings_
757 elif self.solver == "eigen":
758 X_new = X @ self.scalings_
759
760 return X_new[:, : self._max_components]
761
762 def predict_proba(self, X):
763 """Estimate probability.
764
765 Parameters
766 ----------
767 X : array-like of shape (n_samples, n_features)
768 Input data.
769
770 Returns
771 -------
772 C : ndarray of shape (n_samples, n_classes)
773 Estimated probabilities.
774 """
775 check_is_fitted(self)
776 xp, is_array_api_compliant = get_namespace(X)
777 decision = self.decision_function(X)
778 if size(self.classes_) == 2:
779 proba = _expit(decision, xp)
780 return xp.stack([1 - proba, proba], axis=1)
781 else:
782 return softmax(decision)
783
784 def predict_log_proba(self, X):
785 """Estimate log probability.
786
787 Parameters
788 ----------
789 X : array-like of shape (n_samples, n_features)
790 Input data.
791
792 Returns
793 -------
794 C : ndarray of shape (n_samples, n_classes)
795 Estimated log probabilities.
796 """
797 xp, _ = get_namespace(X)
798 prediction = self.predict_proba(X)
799
800 info = xp.finfo(prediction.dtype)
801 if hasattr(info, "smallest_normal"):
802 smallest_normal = info.smallest_normal
803 else:
804 # smallest_normal was introduced in NumPy 1.22
805 smallest_normal = info.tiny
806
807 prediction[prediction == 0.0] += smallest_normal
808 return xp.log(prediction)
809
810 def decision_function(self, X):
811 """Apply decision function to an array of samples.
812
813 The decision function is equal (up to a constant factor) to the
814 log-posterior of the model, i.e. `log p(y = k | x)`. In a binary
815 classification setting this instead corresponds to the difference
816 `log p(y = 1 | x) - log p(y = 0 | x)`. See :ref:`lda_qda_math`.
817
818 Parameters
819 ----------
820 X : array-like of shape (n_samples, n_features)
821 Array of samples (test vectors).
822
823 Returns
824 -------
825 y_scores : ndarray of shape (n_samples,) or (n_samples, n_classes)
826 Decision function values related to each class, per sample.
827 In the two-class case, the shape is `(n_samples,)`, giving the
828 log likelihood ratio of the positive class.
829 """
830 # Only override for the doc
831 return super().decision_function(X)
832
833 def __sklearn_tags__(self):
834 tags = super().__sklearn_tags__()
835 tags.array_api_support = True
836 return tags
837
838
839class QuadraticDiscriminantAnalysis(
840 DiscriminantAnalysisPredictionMixin, ClassifierMixin, BaseEstimator
841):
842 """Quadratic Discriminant Analysis.
843
844 A classifier with a quadratic decision boundary, generated
845 by fitting class conditional densities to the data
846 and using Bayes' rule.
847
848 The model fits a Gaussian density to each class.
849
850 .. versionadded:: 0.17
851
852 For a comparison between
853 :class:`~sklearn.discriminant_analysis.QuadraticDiscriminantAnalysis`
854 and :class:`~sklearn.discriminant_analysis.LinearDiscriminantAnalysis`, see
855 :ref:`sphx_glr_auto_examples_classification_plot_lda_qda.py`.
856
857 Read more in the :ref:`User Guide <lda_qda>`.
858
859 Parameters
860 ----------
861 priors : array-like of shape (n_classes,), default=None
862 Class priors. By default, the class proportions are inferred from the
863 training data.
864
865 reg_param : float, default=0.0
866 Regularizes the per-class covariance estimates by transforming S2 as
867 ``S2 = (1 - reg_param) * S2 + reg_param * np.eye(n_features)``,
868 where S2 corresponds to the `scaling_` attribute of a given class.
869
870 store_covariance : bool, default=False
871 If True, the class covariance matrices are explicitly computed and
872 stored in the `self.covariance_` attribute.
873
874 .. versionadded:: 0.17
875
876 tol : float, default=1.0e-4
877 Absolute threshold for the covariance matrix to be considered rank
878 deficient after applying some regularization (see `reg_param`) to each
879 `Sk` where `Sk` represents covariance matrix for k-th class. This
880 parameter does not affect the predictions. It controls when a warning
881 is raised if the covariance matrix is not full rank.
882
883 .. versionadded:: 0.17
884
885 Attributes
886 ----------
887 covariance_ : list of len n_classes of ndarray \
888 of shape (n_features, n_features)
889 For each class, gives the covariance matrix estimated using the
890 samples of that class. The estimations are unbiased. Only present if
891 `store_covariance` is True.
892
893 means_ : array-like of shape (n_classes, n_features)
894 Class-wise means.
895
896 priors_ : array-like of shape (n_classes,)
897 Class priors (sum to 1).
898
899 rotations_ : list of len n_classes of ndarray of shape (n_features, n_k)
900 For each class k an array of shape (n_features, n_k), where
901 ``n_k = min(n_features, number of elements in class k)``
902 It is the rotation of the Gaussian distribution, i.e. its
903 principal axis. It corresponds to `V`, the matrix of eigenvectors
904 coming from the SVD of `Xk = U S Vt` where `Xk` is the centered
905 matrix of samples from class k.
906
907 scalings_ : list of len n_classes of ndarray of shape (n_k,)
908 For each class, contains the scaling of
909 the Gaussian distributions along its principal axes, i.e. the
910 variance in the rotated coordinate system. It corresponds to `S^2 /
911 (n_samples - 1)`, where `S` is the diagonal matrix of singular values
912 from the SVD of `Xk`, where `Xk` is the centered matrix of samples
913 from class k.
914
915 classes_ : ndarray of shape (n_classes,)
916 Unique class labels.
917
918 n_features_in_ : int
919 Number of features seen during :term:`fit`.
920
921 .. versionadded:: 0.24
922
923 feature_names_in_ : ndarray of shape (`n_features_in_`,)
924 Names of features seen during :term:`fit`. Defined only when `X`
925 has feature names that are all strings.
926
927 .. versionadded:: 1.0
928
929 See Also
930 --------
931 LinearDiscriminantAnalysis : Linear Discriminant Analysis.
932
933 Examples
934 --------
935 >>> from sklearn.discriminant_analysis import QuadraticDiscriminantAnalysis
936 >>> import numpy as np
937 >>> X = np.array([[-1, -1], [-2, -1], [-3, -2], [1, 1], [2, 1], [3, 2]])
938 >>> y = np.array([1, 1, 1, 2, 2, 2])
939 >>> clf = QuadraticDiscriminantAnalysis()
940 >>> clf.fit(X, y)
941 QuadraticDiscriminantAnalysis()
942 >>> print(clf.predict([[-0.8, -1]]))
943 [1]
944 """
945
946 _parameter_constraints: dict = {
947 "priors": ["array-like", None],
948 "reg_param": [Interval(Real, 0, 1, closed="both")],
949 "store_covariance": ["boolean"],
950 "tol": [Interval(Real, 0, None, closed="left")],
951 }
952
953 def __init__(
954 self, *, priors=None, reg_param=0.0, store_covariance=False, tol=1.0e-4
955 ):
956 self.priors = priors
957 self.reg_param = reg_param
958 self.store_covariance = store_covariance
959 self.tol = tol
960
961 @_fit_context(prefer_skip_nested_validation=True)
962 def fit(self, X, y):
963 """Fit the model according to the given training data and parameters.
964
965 .. versionchanged:: 0.19
966 ``store_covariances`` has been moved to main constructor as
967 ``store_covariance``.
968
969 .. versionchanged:: 0.19
970 ``tol`` has been moved to main constructor.
971
972 Parameters
973 ----------
974 X : array-like of shape (n_samples, n_features)
975 Training vector, where `n_samples` is the number of samples and
976 `n_features` is the number of features.
977
978 y : array-like of shape (n_samples,)
979 Target values (integers).
980
981 Returns
982 -------
983 self : object
984 Fitted estimator.
985 """
986 X, y = validate_data(self, X, y)
987 check_classification_targets(y)
988 self.classes_, y = np.unique(y, return_inverse=True)
989 n_samples, n_features = X.shape
990 n_classes = len(self.classes_)
991 if n_classes < 2:
992 raise ValueError(
993 "The number of classes has to be greater than one; got %d class"
994 % (n_classes)
995 )
996 if self.priors is None:
997 self.priors_ = np.bincount(y) / float(n_samples)
998 else:
999 self.priors_ = np.array(self.priors)
1000
1001 cov = None
1002 store_covariance = self.store_covariance
1003 if store_covariance:
1004 cov = []
1005 means = []
1006 scalings = []
1007 rotations = []
1008 for ind in range(n_classes):
1009 Xg = X[y == ind, :]
1010 meang = Xg.mean(0)
1011 means.append(meang)
1012 if len(Xg) == 1:
1013 raise ValueError(
1014 "y has only 1 sample in class %s, covariance is ill defined."
1015 % str(self.classes_[ind])
1016 )
1017 Xgc = Xg - meang
1018 # Xgc = U * S * V.T
1019 _, S, Vt = np.linalg.svd(Xgc, full_matrices=False)
1020 S2 = (S**2) / (len(Xg) - 1)
1021 S2 = ((1 - self.reg_param) * S2) + self.reg_param
1022 rank = np.sum(S2 > self.tol)
1023 if rank < n_features:
1024 warnings.warn(
1025 f"The covariance matrix of class {ind} is not full rank. "
1026 "Increasing the value of parameter `reg_param` might help"
1027 " reducing the collinearity.",
1028 linalg.LinAlgWarning,
1029 )
1030 if self.store_covariance or store_covariance:
1031 # cov = V * (S^2 / (n-1)) * V.T
1032 cov.append(np.dot(S2 * Vt.T, Vt))
1033 scalings.append(S2)
1034 rotations.append(Vt.T)
1035 if self.store_covariance or store_covariance:
1036 self.covariance_ = cov
1037 self.means_ = np.asarray(means)
1038 self.scalings_ = scalings
1039 self.rotations_ = rotations
1040 return self
1041
1042 def _decision_function(self, X):
1043 # return log posterior, see eq (4.12) p. 110 of the ESL.
1044 check_is_fitted(self)
1045
1046 X = validate_data(self, X, reset=False)
1047 norm2 = []
1048 for i in range(len(self.classes_)):
1049 R = self.rotations_[i]
1050 S = self.scalings_[i]
1051 Xm = X - self.means_[i]
1052 X2 = np.dot(Xm, R * (S ** (-0.5)))
1053 norm2.append(np.sum(X2**2, axis=1))
1054 norm2 = np.array(norm2).T # shape = [len(X), n_classes]
1055 u = np.asarray([np.sum(np.log(s)) for s in self.scalings_])
1056 return -0.5 * (norm2 + u) + np.log(self.priors_)
1057
1058 def decision_function(self, X):
1059 """Apply decision function to an array of samples.
1060
1061 The decision function is equal (up to a constant factor) to the
1062 log-posterior of the model, i.e. `log p(y = k | x)`. In a binary
1063 classification setting this instead corresponds to the difference
1064 `log p(y = 1 | x) - log p(y = 0 | x)`. See :ref:`lda_qda_math`.
1065
1066 Parameters
1067 ----------
1068 X : array-like of shape (n_samples, n_features)
1069 Array of samples (test vectors).
1070
1071 Returns
1072 -------
1073 C : ndarray of shape (n_samples,) or (n_samples, n_classes)
1074 Decision function values related to each class, per sample.
1075 In the two-class case, the shape is `(n_samples,)`, giving the
1076 log likelihood ratio of the positive class.
1077 """
1078 return super().decision_function(X)
1079
1080 def predict(self, X):
1081 """Perform classification on an array of test vectors X.
1082
1083 The predicted class C for each sample in X is returned.
1084
1085 Parameters
1086 ----------
1087 X : array-like of shape (n_samples, n_features)
1088 Vector to be scored, where `n_samples` is the number of samples and
1089 `n_features` is the number of features.
1090
1091 Returns
1092 -------
1093 C : ndarray of shape (n_samples,)
1094 Estimated probabilities.
1095 """
1096 return super().predict(X)
1097
1098 def predict_proba(self, X):
1099 """Return posterior probabilities of classification.
1100
1101 Parameters
1102 ----------
1103 X : array-like of shape (n_samples, n_features)
1104 Array of samples/test vectors.
1105
1106 Returns
1107 -------
1108 C : ndarray of shape (n_samples, n_classes)
1109 Posterior probabilities of classification per class.
1110 """
1111 # compute the likelihood of the underlying gaussian models
1112 # up to a multiplicative constant.
1113 return super().predict_proba(X)
1114
1115 def predict_log_proba(self, X):
1116 """Return log of posterior probabilities of classification.
1117
1118 Parameters
1119 ----------
1120 X : array-like of shape (n_samples, n_features)
1121 Array of samples/test vectors.
1122
1123 Returns
1124 -------
1125 C : ndarray of shape (n_samples, n_classes)
1126 Posterior log-probabilities of classification per class.
1127 """
1128 # XXX : can do better to avoid precision overflows
1129 return super().predict_log_proba(X)
1130