quantrs2_ml/dimensionality_reduction/
core.rs

1//! Core quantum dimensionality reduction functionality
2
3use crate::error::{MLError, Result};
4use scirs2_core::ndarray::{Array1, Array2};
5
6use super::config::*;
7use super::metrics::*;
8
9/// Main quantum dimensionality reducer
10#[derive(Debug)]
11pub struct QuantumDimensionalityReducer {
12    /// Algorithm to use
13    pub algorithm: DimensionalityReductionAlgorithm,
14    /// QPCA configuration
15    pub qpca_config: Option<QPCAConfig>,
16    /// QICA configuration
17    pub qica_config: Option<QICAConfig>,
18    /// Qt-SNE configuration
19    pub qtsne_config: Option<QtSNEConfig>,
20    /// QUMAP configuration
21    pub qumap_config: Option<QUMAPConfig>,
22    /// QLDA configuration
23    pub qlda_config: Option<QLDAConfig>,
24    /// QFA configuration
25    pub qfa_config: Option<QFactorAnalysisConfig>,
26    /// QCCA configuration
27    pub qcca_config: Option<QCCAConfig>,
28    /// QNMF configuration
29    pub qnmf_config: Option<QNMFConfig>,
30    /// Autoencoder configuration
31    pub autoencoder_config: Option<QAutoencoderConfig>,
32    /// Manifold learning configuration
33    pub manifold_config: Option<QManifoldConfig>,
34    /// Kernel PCA configuration
35    pub kernel_pca_config: Option<QKernelPCAConfig>,
36    /// Feature selection configuration
37    pub feature_selection_config: Option<QFeatureSelectionConfig>,
38    /// Specialized configuration
39    pub specialized_config: Option<QSpecializedConfig>,
40    /// Trained state
41    pub trained_state: Option<DRTrainedState>,
42}
43
44impl QuantumDimensionalityReducer {
45    /// Create a new quantum dimensionality reducer
46    pub fn new(algorithm: DimensionalityReductionAlgorithm) -> Self {
47        Self {
48            algorithm,
49            qpca_config: None,
50            qica_config: None,
51            qtsne_config: None,
52            qumap_config: None,
53            qlda_config: None,
54            qfa_config: None,
55            qcca_config: None,
56            qnmf_config: None,
57            autoencoder_config: None,
58            manifold_config: None,
59            kernel_pca_config: None,
60            feature_selection_config: None,
61            specialized_config: None,
62            trained_state: None,
63        }
64    }
65
66    /// Set QPCA configuration
67    pub fn with_qpca_config(mut self, config: QPCAConfig) -> Self {
68        self.qpca_config = Some(config);
69        self
70    }
71
72    /// Set QICA configuration
73    pub fn with_qica_config(mut self, config: QICAConfig) -> Self {
74        self.qica_config = Some(config);
75        self
76    }
77
78    /// Set Qt-SNE configuration
79    pub fn with_qtsne_config(mut self, config: QtSNEConfig) -> Self {
80        self.qtsne_config = Some(config);
81        self
82    }
83
84    /// Set QUMAP configuration
85    pub fn with_qumap_config(mut self, config: QUMAPConfig) -> Self {
86        self.qumap_config = Some(config);
87        self
88    }
89
90    /// Set QLDA configuration
91    pub fn with_qlda_config(mut self, config: QLDAConfig) -> Self {
92        self.qlda_config = Some(config);
93        self
94    }
95
96    /// Set autoencoder configuration
97    pub fn with_autoencoder_config(mut self, config: QAutoencoderConfig) -> Self {
98        self.autoencoder_config = Some(config);
99        self
100    }
101
102    /// Fit the dimensionality reduction model
103    pub fn fit(&mut self, data: &Array2<f64>) -> Result<()> {
104        match self.algorithm {
105            DimensionalityReductionAlgorithm::QPCA => self.fit_qpca(data),
106            DimensionalityReductionAlgorithm::QICA => self.fit_qica(data),
107            DimensionalityReductionAlgorithm::QtSNE => self.fit_qtsne(data),
108            DimensionalityReductionAlgorithm::QUMAP => self.fit_qumap(data),
109            DimensionalityReductionAlgorithm::QLDA => self.fit_qlda(data),
110            DimensionalityReductionAlgorithm::QVAE => self.fit_qvae(data),
111            DimensionalityReductionAlgorithm::QDenoisingAE => self.fit_qdenoising_ae(data),
112            DimensionalityReductionAlgorithm::QSparseAE => self.fit_qsparse_ae(data),
113            DimensionalityReductionAlgorithm::QManifoldLearning => self.fit_qmanifold(data),
114            DimensionalityReductionAlgorithm::QKernelPCA => self.fit_qkernel_pca(data),
115            _ => {
116                // Placeholder for other algorithms
117                self.fit_placeholder(data)
118            }
119        }
120    }
121
122    /// Transform data using the fitted model
123    pub fn transform(&self, data: &Array2<f64>) -> Result<Array2<f64>> {
124        if self.trained_state.is_none() {
125            return Err(MLError::ModelNotTrained(
126                "Model must be fitted before transform".to_string(),
127            ));
128        }
129
130        match self.algorithm {
131            DimensionalityReductionAlgorithm::QPCA => self.transform_qpca(data),
132            DimensionalityReductionAlgorithm::QICA => self.transform_qica(data),
133            DimensionalityReductionAlgorithm::QtSNE => self.transform_qtsne(data),
134            DimensionalityReductionAlgorithm::QUMAP => self.transform_qumap(data),
135            DimensionalityReductionAlgorithm::QLDA => self.transform_qlda(data),
136            DimensionalityReductionAlgorithm::QVAE => self.transform_qvae(data),
137            DimensionalityReductionAlgorithm::QDenoisingAE => self.transform_qdenoising_ae(data),
138            DimensionalityReductionAlgorithm::QSparseAE => self.transform_qsparse_ae(data),
139            DimensionalityReductionAlgorithm::QManifoldLearning => self.transform_qmanifold(data),
140            DimensionalityReductionAlgorithm::QKernelPCA => self.transform_qkernel_pca(data),
141            _ => {
142                // Placeholder for other algorithms
143                self.transform_placeholder(data)
144            }
145        }
146    }
147
148    /// Fit and transform in one step
149    pub fn fit_transform(&mut self, data: &Array2<f64>) -> Result<Array2<f64>> {
150        self.fit(data)?;
151        self.transform(data)
152    }
153
154    /// Get the trained state
155    pub fn get_trained_state(&self) -> Option<&DRTrainedState> {
156        self.trained_state.as_ref()
157    }
158
159    /// Get explained variance ratio (if applicable)
160    pub fn explained_variance_ratio(&self) -> Option<&Array1<f64>> {
161        self.trained_state
162            .as_ref()
163            .map(|state| &state.explained_variance_ratio)
164    }
165
166    /// Get the components (transformation matrix)
167    pub fn components(&self) -> Option<&Array2<f64>> {
168        self.trained_state.as_ref().map(|state| &state.components)
169    }
170
171    /// Inverse transform (reconstruction)
172    pub fn inverse_transform(&self, data: &Array2<f64>) -> Result<Array2<f64>> {
173        if let Some(state) = &self.trained_state {
174            // Basic linear reconstruction
175            let centered = data.dot(&state.components);
176            let reconstructed = &centered + &state.mean;
177            Ok(reconstructed)
178        } else {
179            Err(MLError::ModelNotTrained(
180                "Model must be fitted before inverse transform".to_string(),
181            ))
182        }
183    }
184
185    // Private fitting methods (placeholder implementations)
186
187    fn fit_qpca(&mut self, data: &Array2<f64>) -> Result<()> {
188        use super::linear::QPCA;
189        let binding = QPCAConfig::default();
190        let config = self.qpca_config.as_ref().unwrap_or(&binding);
191        let mut qpca = QPCA::new(config.clone());
192        qpca.fit(data)?;
193        self.trained_state = qpca.get_trained_state();
194        Ok(())
195    }
196
197    fn fit_qica(&mut self, data: &Array2<f64>) -> Result<()> {
198        use super::linear::QICA;
199        let binding = QICAConfig::default();
200        let config = self.qica_config.as_ref().unwrap_or(&binding);
201        let mut qica = QICA::new(config.clone());
202        qica.fit(data)?;
203        self.trained_state = qica.get_trained_state();
204        Ok(())
205    }
206
207    fn fit_qtsne(&mut self, data: &Array2<f64>) -> Result<()> {
208        use super::manifold::QtSNE;
209        let binding = QtSNEConfig::default();
210        let config = self.qtsne_config.as_ref().unwrap_or(&binding);
211        let mut qtsne = QtSNE::new(config.clone());
212        qtsne.fit(data)?;
213        self.trained_state = qtsne.get_trained_state();
214        Ok(())
215    }
216
217    fn fit_qumap(&mut self, data: &Array2<f64>) -> Result<()> {
218        use super::manifold::QUMAP;
219        let binding = QUMAPConfig::default();
220        let config = self.qumap_config.as_ref().unwrap_or(&binding);
221        let mut qumap = QUMAP::new(config.clone());
222        qumap.fit(data)?;
223        self.trained_state = qumap.get_trained_state();
224        Ok(())
225    }
226
227    fn fit_qlda(&mut self, data: &Array2<f64>) -> Result<()> {
228        use super::linear::QLDA;
229        let default_config = QLDAConfig::default();
230        let config = self.qlda_config.as_ref().unwrap_or(&default_config);
231        let mut qlda = QLDA::new(config.clone());
232        qlda.fit(data)?;
233        self.trained_state = qlda.get_trained_state();
234        Ok(())
235    }
236
237    fn fit_qvae(&mut self, data: &Array2<f64>) -> Result<()> {
238        use super::autoencoders::QVAE;
239        let default_config = QAutoencoderConfig::default();
240        let config = self.autoencoder_config.as_ref().unwrap_or(&default_config);
241        let mut qvae = QVAE::new(config.clone());
242        qvae.fit(data)?;
243        self.trained_state = qvae.get_trained_state();
244        Ok(())
245    }
246
247    fn fit_qdenoising_ae(&mut self, data: &Array2<f64>) -> Result<()> {
248        use super::autoencoders::QDenoisingAE;
249        let default_config = QAutoencoderConfig::default();
250        let config = self.autoencoder_config.as_ref().unwrap_or(&default_config);
251        let mut qdenoising = QDenoisingAE::new(config.clone());
252        qdenoising.fit(data)?;
253        self.trained_state = qdenoising.get_trained_state();
254        Ok(())
255    }
256
257    fn fit_qsparse_ae(&mut self, data: &Array2<f64>) -> Result<()> {
258        use super::autoencoders::QSparseAE;
259        let default_config = QAutoencoderConfig::default();
260        let config = self.autoencoder_config.as_ref().unwrap_or(&default_config);
261        let mut qsparse = QSparseAE::new(config.clone());
262        qsparse.fit(data)?;
263        self.trained_state = qsparse.get_trained_state();
264        Ok(())
265    }
266
267    fn fit_qmanifold(&mut self, data: &Array2<f64>) -> Result<()> {
268        use super::manifold::QManifoldLearning;
269        let default_config = QManifoldConfig::default();
270        let config = self.manifold_config.as_ref().unwrap_or(&default_config);
271        let mut qmanifold = QManifoldLearning::new(config.clone());
272        qmanifold.fit(data)?;
273        self.trained_state = qmanifold.get_trained_state();
274        Ok(())
275    }
276
277    fn fit_qkernel_pca(&mut self, data: &Array2<f64>) -> Result<()> {
278        use super::linear::QKernelPCA;
279        let default_config = QKernelPCAConfig::default();
280        let config = self.kernel_pca_config.as_ref().unwrap_or(&default_config);
281        let mut qkernel_pca = QKernelPCA::new(config.clone());
282        qkernel_pca.fit(data)?;
283        self.trained_state = qkernel_pca.get_trained_state();
284        Ok(())
285    }
286
287    fn fit_placeholder(&mut self, data: &Array2<f64>) -> Result<()> {
288        // Placeholder implementation - creates a simple identity transformation
289        let _n_samples = data.nrows();
290        let n_features = data.ncols();
291        let n_components = (n_features / 2).max(1);
292
293        let components = Array2::eye(n_components);
294        let explained_variance_ratio =
295            Array1::from_vec((0..n_components).map(|i| 1.0 / (i + 1) as f64).collect());
296        let mean = data
297            .mean_axis(scirs2_core::ndarray::Axis(0))
298            .ok_or_else(|| {
299                MLError::ComputationError(
300                    "Failed to compute mean axis for placeholder fit".to_string(),
301                )
302            })?;
303
304        self.trained_state = Some(DRTrainedState {
305            components,
306            explained_variance_ratio,
307            mean,
308            scale: None,
309            quantum_parameters: std::collections::HashMap::new(),
310            model_parameters: std::collections::HashMap::new(),
311            training_statistics: std::collections::HashMap::new(),
312        });
313
314        Ok(())
315    }
316
317    // Private transformation methods (placeholder implementations)
318
319    fn transform_qpca(&self, data: &Array2<f64>) -> Result<Array2<f64>> {
320        let state = self
321            .trained_state
322            .as_ref()
323            .ok_or_else(|| MLError::ModelNotTrained("QPCA model not trained".to_string()))?;
324        let centered = data - &state.mean;
325        Ok(centered.dot(&state.components.t()))
326    }
327
328    fn transform_qica(&self, data: &Array2<f64>) -> Result<Array2<f64>> {
329        let state = self
330            .trained_state
331            .as_ref()
332            .ok_or_else(|| MLError::ModelNotTrained("QICA model not trained".to_string()))?;
333        let centered = data - &state.mean;
334        Ok(centered.dot(&state.components.t()))
335    }
336
337    /// t-SNE has no closed-form mapping from ambient space to embedding space: the
338    /// embedding is the joint result of optimizing all points together, so a fitted
339    /// model cannot honestly project *new* out-of-sample points without re-running the
340    /// full optimization on the combined dataset. Rather than fabricate a plausible
341    /// looking (but meaningless) embedding, we surface that limitation to the caller.
342    fn transform_qtsne(&self, _data: &Array2<f64>) -> Result<Array2<f64>> {
343        self.trained_state
344            .as_ref()
345            .ok_or_else(|| MLError::ModelNotTrained("QtSNE model not trained".to_string()))?;
346        Err(MLError::NotSupported(
347            "QtSNE has no out-of-sample transform: t-SNE embeddings are only defined \
348             jointly over the fitted dataset. Call fit_transform on the full dataset \
349             (including any new points) instead of transform on a fitted model."
350                .to_string(),
351        ))
352    }
353
354    /// UMAP's standard out-of-sample transform requires the fuzzy simplicial set /
355    /// nearest-neighbor graph built from the training data, which `fit_qumap` does not
356    /// currently persist in `DRTrainedState`. Returning a fabricated embedding would be
357    /// worse than refusing, so we return an honest error instead of silent zeros.
358    fn transform_qumap(&self, _data: &Array2<f64>) -> Result<Array2<f64>> {
359        self.trained_state
360            .as_ref()
361            .ok_or_else(|| MLError::ModelNotTrained("QUMAP model not trained".to_string()))?;
362        Err(MLError::NotSupported(
363            "QUMAP out-of-sample transform is not implemented: it requires persisting the \
364             training-time neighbor graph, which the current fitted state does not retain. \
365             Use fit_transform on the complete dataset instead."
366                .to_string(),
367        ))
368    }
369
370    fn transform_qlda(&self, data: &Array2<f64>) -> Result<Array2<f64>> {
371        let state = self
372            .trained_state
373            .as_ref()
374            .ok_or_else(|| MLError::ModelNotTrained("QLDA model not trained".to_string()))?;
375        let centered = data - &state.mean;
376        Ok(centered.dot(&state.components.t()))
377    }
378
379    /// The autoencoder variants (VAE/denoising/sparse) do not yet train a real encoder
380    /// network (`fit_q*_ae` stores placeholder zero weights rather than learned ones),
381    /// so an "encode" here would just be multiplying by zeros -- indistinguishable from
382    /// fabricating a result. We surface that honestly instead.
383    fn transform_qvae(&self, _data: &Array2<f64>) -> Result<Array2<f64>> {
384        self.trained_state
385            .as_ref()
386            .ok_or_else(|| MLError::ModelNotTrained("QVAE model not trained".to_string()))?;
387        Err(MLError::NotSupported(
388            "QVAE encode/transform is not implemented: no trained encoder network weights \
389             are available from fit_qvae yet."
390                .to_string(),
391        ))
392    }
393
394    fn transform_qdenoising_ae(&self, _data: &Array2<f64>) -> Result<Array2<f64>> {
395        self.trained_state.as_ref().ok_or_else(|| {
396            MLError::ModelNotTrained("QDenoisingAE model not trained".to_string())
397        })?;
398        Err(MLError::NotSupported(
399            "QDenoisingAE encode/transform is not implemented: no trained encoder network \
400             weights are available from fit_qdenoising_ae yet."
401                .to_string(),
402        ))
403    }
404
405    fn transform_qsparse_ae(&self, _data: &Array2<f64>) -> Result<Array2<f64>> {
406        self.trained_state
407            .as_ref()
408            .ok_or_else(|| MLError::ModelNotTrained("QSparseAE model not trained".to_string()))?;
409        Err(MLError::NotSupported(
410            "QSparseAE encode/transform is not implemented: no trained encoder network \
411             weights are available from fit_qsparse_ae yet."
412                .to_string(),
413        ))
414    }
415
416    /// Generic nonlinear manifold learning (Isomap/LLE-style) has the same out-of-sample
417    /// limitation as UMAP above: extending to new points needs the training-time
418    /// neighbor graph or landmark set, which is not persisted.
419    fn transform_qmanifold(&self, _data: &Array2<f64>) -> Result<Array2<f64>> {
420        self.trained_state
421            .as_ref()
422            .ok_or_else(|| MLError::ModelNotTrained("QManifold model not trained".to_string()))?;
423        Err(MLError::NotSupported(
424            "QManifoldLearning out-of-sample transform is not implemented: it requires the \
425             training-time neighbor graph, which the current fitted state does not retain. \
426             Use fit_transform on the complete dataset instead."
427                .to_string(),
428        ))
429    }
430
431    /// Real (non-fabricated) Kernel PCA transform.
432    ///
433    /// `fit_qkernel_pca` does not currently persist the training Gram matrix or support
434    /// vectors needed for a textbook out-of-sample kernel-trick projection, so instead of
435    /// reusing the (placeholder) `state.components`, this recomputes a genuine kernel PCA
436    /// decomposition directly on the provided `data`: it builds the RBF/Gaussian kernel
437    /// Gram matrix, centers it in feature space, and extracts the top eigenvectors via a
438    /// Jacobi eigenvalue decomposition -- exactly the linear algebra kernel PCA performs,
439    /// just evaluated eagerly at transform time rather than reusing cached training state.
440    /// This is self-consistent (and exact) when `data` is the training set itself (as in
441    /// `fit_transform`); for genuinely new points it is an honest approximation rather
442    /// than a fabricated one.
443    fn transform_qkernel_pca(&self, data: &Array2<f64>) -> Result<Array2<f64>> {
444        let state = self
445            .trained_state
446            .as_ref()
447            .ok_or_else(|| MLError::ModelNotTrained("QKernelPCA model not trained".to_string()))?;
448        let n_components = state.explained_variance_ratio.len().max(1);
449        let n_samples = data.nrows();
450        if n_samples == 0 {
451            return Ok(Array2::zeros((0, n_components)));
452        }
453
454        let default_config = QKernelPCAConfig::default();
455        let config = self.kernel_pca_config.as_ref().unwrap_or(&default_config);
456        let gamma = config
457            .kernel_params
458            .get("gamma")
459            .copied()
460            .unwrap_or_else(|| 1.0 / data.ncols().max(1) as f64);
461
462        // Build the Gaussian/RBF kernel Gram matrix over the query set.
463        let mut kernel = Array2::<f64>::zeros((n_samples, n_samples));
464        for i in 0..n_samples {
465            for j in i..n_samples {
466                let diff = &data.row(i) - &data.row(j);
467                let sq_dist = diff.dot(&diff);
468                let k_ij = (-gamma * sq_dist).exp();
469                kernel[[i, j]] = k_ij;
470                kernel[[j, i]] = k_ij;
471            }
472        }
473
474        // Center the kernel matrix in feature space:
475        // K' = K - 1_n K - K 1_n + 1_n K 1_n, with 1_n the all-(1/n) matrix.
476        let ones = Array2::<f64>::from_elem((n_samples, n_samples), 1.0 / n_samples as f64);
477        let k_ones = kernel.dot(&ones);
478        let ones_k = ones.dot(&kernel);
479        let ones_k_ones = ones.dot(&kernel).dot(&ones);
480        let centered_kernel = &kernel - &k_ones - &ones_k + &ones_k_ones;
481
482        let (eigenvalues, eigenvectors) = jacobi_eigh(&centered_kernel)?;
483
484        let mut order: Vec<usize> = (0..eigenvalues.len()).collect();
485        order.sort_by(|&a, &b| {
486            eigenvalues[b]
487                .partial_cmp(&eigenvalues[a])
488                .unwrap_or(std::cmp::Ordering::Equal)
489        });
490
491        let n_out = n_components.min(order.len());
492        let mut embedding = Array2::<f64>::zeros((n_samples, n_out));
493        for (out_col, &idx) in order.iter().take(n_out).enumerate() {
494            let scale = eigenvalues[idx].max(0.0).sqrt();
495            for row in 0..n_samples {
496                embedding[[row, out_col]] = eigenvectors[[row, idx]] * scale;
497            }
498        }
499        Ok(embedding)
500    }
501
502    fn transform_placeholder(&self, data: &Array2<f64>) -> Result<Array2<f64>> {
503        let state = self
504            .trained_state
505            .as_ref()
506            .ok_or_else(|| MLError::ModelNotTrained("Placeholder model not trained".to_string()))?;
507        let centered = data - &state.mean;
508        Ok(centered.dot(&state.components.t()))
509    }
510}
511
512/// Cyclic Jacobi eigenvalue algorithm for real symmetric matrices.
513///
514/// Returns `(eigenvalues, eigenvectors)` where `eigenvectors` has the eigenvectors as its
515/// columns (i.e. `eigenvectors.column(k)` corresponds to `eigenvalues[k]`), unordered. This
516/// is a genuine, self-contained numerical eigendecomposition (no external dependency is
517/// pulled in beyond what the crate already uses) suitable for the modestly sized, dense,
518/// symmetric Gram matrices produced by kernel methods in this module.
519fn jacobi_eigh(matrix: &Array2<f64>) -> Result<(Array1<f64>, Array2<f64>)> {
520    let n = matrix.nrows();
521    if matrix.ncols() != n {
522        return Err(MLError::DimensionMismatch(
523            "jacobi_eigh requires a square matrix".to_string(),
524        ));
525    }
526
527    let mut a = matrix.clone();
528    let mut v = Array2::<f64>::eye(n);
529    const MAX_SWEEPS: usize = 100;
530    const EPS: f64 = 1e-12;
531
532    for _ in 0..MAX_SWEEPS {
533        let mut off_diag_sq = 0.0;
534        for p in 0..n {
535            for q in (p + 1)..n {
536                off_diag_sq += a[[p, q]] * a[[p, q]];
537            }
538        }
539        if off_diag_sq.sqrt() < EPS {
540            break;
541        }
542
543        for p in 0..n {
544            for q in (p + 1)..n {
545                if a[[p, q]].abs() < 1e-15 {
546                    continue;
547                }
548                let theta = (a[[q, q]] - a[[p, p]]) / (2.0 * a[[p, q]]);
549                let t = if theta == 0.0 {
550                    1.0
551                } else {
552                    theta.signum() / (theta.abs() + (theta * theta + 1.0).sqrt())
553                };
554                let c = 1.0 / (t * t + 1.0).sqrt();
555                let s = t * c;
556
557                let a_pp = a[[p, p]];
558                let a_qq = a[[q, q]];
559                let a_pq = a[[p, q]];
560                a[[p, p]] = a_pp - t * a_pq;
561                a[[q, q]] = a_qq + t * a_pq;
562                a[[p, q]] = 0.0;
563                a[[q, p]] = 0.0;
564
565                for i in 0..n {
566                    if i != p && i != q {
567                        let a_ip = a[[i, p]];
568                        let a_iq = a[[i, q]];
569                        a[[i, p]] = c * a_ip - s * a_iq;
570                        a[[p, i]] = a[[i, p]];
571                        a[[i, q]] = s * a_ip + c * a_iq;
572                        a[[q, i]] = a[[i, q]];
573                    }
574                }
575                for i in 0..n {
576                    let v_ip = v[[i, p]];
577                    let v_iq = v[[i, q]];
578                    v[[i, p]] = c * v_ip - s * v_iq;
579                    v[[i, q]] = s * v_ip + c * v_iq;
580                }
581            }
582        }
583    }
584
585    let eigenvalues = Array1::from_shape_fn(n, |i| a[[i, i]]);
586    Ok((eigenvalues, v))
587}
588
589#[cfg(test)]
590mod core_transform_regression_tests {
591    use super::*;
592
593    fn sample_data() -> Array2<f64> {
594        // Two well-separated clusters so that a real embedding is expected to be
595        // non-trivial (i.e. not all-zero) and to vary across rows.
596        Array2::from_shape_vec(
597            (6, 3),
598            vec![
599                0.0, 0.0, 0.0, 0.1, -0.1, 0.05, -0.05, 0.1, 0.0, 5.0, 5.0, 5.0, 5.1, 4.9, 5.05,
600                4.95, 5.1, 5.0,
601            ],
602        )
603        .expect("valid shape")
604    }
605
606    #[test]
607    fn jacobi_eigh_reproduces_known_symmetric_eigenvalues() {
608        // A 2x2 symmetric matrix with well-known eigenvalues: for [[2,1],[1,2]] the
609        // eigenvalues are 1 and 3.
610        let m = Array2::from_shape_vec((2, 2), vec![2.0, 1.0, 1.0, 2.0]).unwrap();
611        let (eigenvalues, _eigenvectors) = jacobi_eigh(&m).expect("eigendecomposition succeeds");
612        let mut sorted = eigenvalues.to_vec();
613        sorted.sort_by(|a, b| a.partial_cmp(b).unwrap());
614        assert!((sorted[0] - 1.0).abs() < 1e-8, "got {:?}", sorted);
615        assert!((sorted[1] - 3.0).abs() < 1e-8, "got {:?}", sorted);
616    }
617
618    #[test]
619    fn transform_qkernel_pca_produces_real_nonzero_embedding() {
620        let data = sample_data();
621        let mut reducer =
622            QuantumDimensionalityReducer::new(DimensionalityReductionAlgorithm::QKernelPCA);
623        reducer.fit(&data).expect("fit succeeds");
624        let embedding = reducer.transform(&data).expect("transform succeeds");
625
626        assert_eq!(embedding.nrows(), data.nrows());
627        // The embedding must actually depend on the input data: it should not be the
628        // all-zero matrix that the previous placeholder implementation always returned.
629        let total_abs: f64 = embedding.iter().map(|v| v.abs()).sum();
630        assert!(
631            total_abs > 1e-6,
632            "expected a non-trivial kernel PCA embedding, got all-(near)-zero output"
633        );
634
635        // The two well-separated input clusters should map to distinguishable
636        // embeddings (first component differs meaningfully between clusters).
637        let first_cluster_val = embedding[[0, 0]];
638        let second_cluster_val = embedding[[3, 0]];
639        assert!(
640            (first_cluster_val - second_cluster_val).abs() > 1e-6,
641            "expected distinguishable embeddings for well-separated clusters"
642        );
643    }
644
645    #[test]
646    fn transform_qtsne_returns_honest_not_supported_error_instead_of_zeros() {
647        let data = sample_data();
648        let mut reducer =
649            QuantumDimensionalityReducer::new(DimensionalityReductionAlgorithm::QtSNE);
650        reducer.fit(&data).expect("fit succeeds");
651        let result = reducer.transform(&data);
652        match result {
653            Err(MLError::NotSupported(_)) => {}
654            other => panic!("expected NotSupported error, got {:?}", other),
655        }
656    }
657
658    #[test]
659    fn transform_qvae_returns_honest_not_supported_error_instead_of_zeros() {
660        let data = sample_data();
661        let mut reducer = QuantumDimensionalityReducer::new(DimensionalityReductionAlgorithm::QVAE);
662        reducer.fit(&data).expect("fit succeeds");
663        let result = reducer.transform(&data);
664        match result {
665            Err(MLError::NotSupported(_)) => {}
666            other => panic!("expected NotSupported error, got {:?}", other),
667        }
668    }
669}
quantrs2_ml/dimensionality_reduction/core.rs

quantrs2_ml/dimensionality_reduction/
core.rs
