Skip to main content

sklears_semi_supervised/
graph_learning.rs

1//! Graph structure learning methods for semi-supervised learning
2//!
3//! This module provides algorithms to learn optimal graph structures from data
4//! for semi-supervised learning tasks.
5
6use scirs2_core::ndarray_ext::{Array1, Array2, ArrayView1, ArrayView2, Axis};
7use scirs2_core::random::Random;
8use sklears_core::{
9    error::{Result as SklResult, SklearsError},
10    traits::{Estimator, Fit, Predict, PredictProba, Untrained},
11    types::Float,
12};
13
14/// Graph Structure Learning for Semi-Supervised Learning
15///
16/// This method learns an optimal graph structure that balances data fidelity
17/// and sparsity constraints. The learned graph is then used for label propagation.
18///
19/// The method solves the optimization problem:
20/// min_W ||X - W * X||_F^2 + λ * ||W||_1 + β * tr(F^T * L_W * F)
21///
22/// where W is the graph adjacency matrix, L_W is the graph Laplacian,
23/// and F is the label matrix.
24///
25/// # Parameters
26///
27/// * `lambda_sparse` - Sparsity regularization parameter
28/// * `beta_smoothness` - Smoothness regularization parameter
29/// * `max_iter` - Maximum number of iterations
30/// * `tol` - Convergence tolerance
31/// * `learning_rate` - Learning rate for optimization
32/// * `adaptive_lr` - Whether to use adaptive learning rate
33///
34/// # Examples
35///
36/// ```
37/// use scirs2_core::array;
38/// use sklears_semi_supervised::GraphStructureLearning;
39/// use sklears_core::traits::{Predict, Fit};
40///
41///
42/// let X = array![[1.0, 2.0], [2.0, 3.0], [3.0, 4.0], [4.0, 5.0]];
43/// let y = array![0, 1, -1, -1]; // -1 indicates unlabeled
44///
45/// let gsl = GraphStructureLearning::new()
46///     .lambda_sparse(0.1)
47///     .beta_smoothness(1.0);
48/// let fitted = gsl.fit(&X.view(), &y.view()).unwrap();
49/// let predictions = fitted.predict(&X.view()).unwrap();
50/// ```
51#[derive(Debug, Clone)]
52pub struct GraphStructureLearning<S = Untrained> {
53    state: S,
54    lambda_sparse: f64,
55    beta_smoothness: f64,
56    max_iter: usize,
57    tol: f64,
58    learning_rate: f64,
59    adaptive_lr: bool,
60    enforce_symmetry: bool,
61    normalize_weights: bool,
62}
63
64impl GraphStructureLearning<Untrained> {
65    /// Create a new GraphStructureLearning instance
66    pub fn new() -> Self {
67        Self {
68            state: Untrained,
69            lambda_sparse: 0.1,
70            beta_smoothness: 1.0,
71            max_iter: 100,
72            tol: 1e-4,
73            learning_rate: 0.01,
74            adaptive_lr: true,
75            enforce_symmetry: true,
76            normalize_weights: true,
77        }
78    }
79
80    /// Set the sparsity regularization parameter
81    pub fn lambda_sparse(mut self, lambda_sparse: f64) -> Self {
82        self.lambda_sparse = lambda_sparse;
83        self
84    }
85
86    /// Set the smoothness regularization parameter
87    pub fn beta_smoothness(mut self, beta_smoothness: f64) -> Self {
88        self.beta_smoothness = beta_smoothness;
89        self
90    }
91
92    /// Set the maximum number of iterations
93    pub fn max_iter(mut self, max_iter: usize) -> Self {
94        self.max_iter = max_iter;
95        self
96    }
97
98    /// Set the convergence tolerance
99    pub fn tol(mut self, tol: f64) -> Self {
100        self.tol = tol;
101        self
102    }
103
104    /// Set the learning rate
105    pub fn learning_rate(mut self, learning_rate: f64) -> Self {
106        self.learning_rate = learning_rate;
107        self
108    }
109
110    /// Enable/disable adaptive learning rate
111    pub fn adaptive_lr(mut self, adaptive_lr: bool) -> Self {
112        self.adaptive_lr = adaptive_lr;
113        self
114    }
115
116    /// Enable/disable symmetry enforcement
117    pub fn enforce_symmetry(mut self, enforce_symmetry: bool) -> Self {
118        self.enforce_symmetry = enforce_symmetry;
119        self
120    }
121
122    /// Enable/disable weight normalization
123    pub fn normalize_weights(mut self, normalize_weights: bool) -> Self {
124        self.normalize_weights = normalize_weights;
125        self
126    }
127
128    #[allow(non_snake_case)] // standard ML notation
129    fn initialize_graph(&self, X: &Array2<f64>) -> Array2<f64> {
130        let n_samples = X.nrows();
131        let mut W = Array2::zeros((n_samples, n_samples));
132
133        // Initialize with k-NN graph
134        let k = (n_samples as f64).sqrt().ceil() as usize;
135        let k = k.clamp(3, 10); // Bound k between 3 and 10
136
137        for i in 0..n_samples {
138            let mut distances: Vec<(usize, f64)> = Vec::new();
139            for j in 0..n_samples {
140                if i != j {
141                    let diff = &X.row(i) - &X.row(j);
142                    let dist = diff.mapv(|x: f64| x * x).sum().sqrt();
143                    distances.push((j, dist));
144                }
145            }
146
147            distances.sort_by(|a, b| a.1.partial_cmp(&b.1).expect("operation should succeed"));
148
149            for &(j, dist) in distances.iter().take(k) {
150                let weight = (-dist / (2.0 * 1.0_f64.powi(2))).exp();
151                W[[i, j]] = weight;
152                if self.enforce_symmetry {
153                    W[[j, i]] = weight;
154                }
155            }
156        }
157
158        W
159    }
160
161    #[allow(non_snake_case)]
162    fn compute_laplacian(&self, W: &Array2<f64>) -> Array2<f64> {
163        let n_samples = W.nrows();
164        let D = W.sum_axis(Axis(1));
165        let mut L = Array2::zeros((n_samples, n_samples));
166
167        for i in 0..n_samples {
168            L[[i, i]] = D[i];
169            for j in 0..n_samples {
170                if i != j {
171                    L[[i, j]] = -W[[i, j]];
172                }
173            }
174        }
175
176        L
177    }
178
179    fn soft_threshold(&self, x: f64, threshold: f64) -> f64 {
180        if x > threshold {
181            x - threshold
182        } else if x < -threshold {
183            x + threshold
184        } else {
185            0.0
186        }
187    }
188
189    #[allow(non_snake_case)] // standard ML notation
190    fn proximal_gradient_step(&self, W: &Array2<f64>, grad: &Array2<f64>, lr: f64) -> Array2<f64> {
191        let mut W_new = W - lr * grad;
192
193        // Apply L1 proximal operator (soft thresholding)
194        let threshold = lr * self.lambda_sparse;
195        W_new.mapv_inplace(|x| self.soft_threshold(x, threshold));
196
197        // Ensure non-negativity
198        W_new.mapv_inplace(|x| x.max(0.0));
199
200        // Enforce symmetry if required
201        if self.enforce_symmetry {
202            let n = W_new.nrows();
203            for i in 0..n {
204                for j in 0..n {
205                    if i != j {
206                        let avg = (W_new[[i, j]] + W_new[[j, i]]) / 2.0;
207                        W_new[[i, j]] = avg;
208                        W_new[[j, i]] = avg;
209                    }
210                }
211            }
212        }
213
214        // Zero diagonal
215        for i in 0..W_new.nrows() {
216            W_new[[i, i]] = 0.0;
217        }
218
219        W_new
220    }
221
222    #[allow(non_snake_case)] // standard ML notation
223    fn normalize_graph(&self, W: &Array2<f64>) -> Array2<f64> {
224        if !self.normalize_weights {
225            return W.clone();
226        }
227
228        let mut W_norm = W.clone();
229        let n_samples = W.nrows();
230
231        // Row-wise normalization
232        for i in 0..n_samples {
233            let row_sum: f64 = W.row(i).sum();
234            if row_sum > 0.0 {
235                for j in 0..n_samples {
236                    W_norm[[i, j]] = W[[i, j]] / row_sum;
237                }
238            }
239        }
240
241        W_norm
242    }
243
244    #[allow(non_snake_case)]
245    fn propagate_labels(&self, W: &Array2<f64>, Y_init: &Array2<f64>) -> SklResult<Array2<f64>> {
246        let n_samples = W.nrows();
247
248        // Compute transition matrix
249        let D = W.sum_axis(Axis(1));
250        let mut P = Array2::zeros((n_samples, n_samples));
251        for i in 0..n_samples {
252            if D[i] > 0.0 {
253                for j in 0..n_samples {
254                    P[[i, j]] = W[[i, j]] / D[i];
255                }
256            }
257        }
258
259        let mut Y = Y_init.clone();
260        let Y_static = Y_init.clone();
261
262        // Label propagation iterations
263        for _iter in 0..50 {
264            let prev_Y = Y.clone();
265            Y = 0.9 * P.dot(&Y) + 0.1 * &Y_static;
266
267            // Check convergence
268            let diff = (&Y - &prev_Y).mapv(|x| x.abs()).sum();
269            if diff < 1e-6 {
270                break;
271            }
272        }
273
274        Ok(Y)
275    }
276}
277
278impl Default for GraphStructureLearning<Untrained> {
279    fn default() -> Self {
280        Self::new()
281    }
282}
283
284impl Estimator for GraphStructureLearning<Untrained> {
285    type Config = ();
286    type Error = SklearsError;
287    type Float = Float;
288
289    fn config(&self) -> &Self::Config {
290        &()
291    }
292}
293
294impl Fit<ArrayView2<'_, Float>, ArrayView1<'_, i32>> for GraphStructureLearning<Untrained> {
295    type Fitted = GraphStructureLearning<GraphStructureLearningTrained>;
296
297    #[allow(non_snake_case)]
298    fn fit(self, X: &ArrayView2<'_, Float>, y: &ArrayView1<'_, i32>) -> SklResult<Self::Fitted> {
299        let X = X.to_owned();
300        let y = y.to_owned();
301
302        let (n_samples, _n_features) = X.dim();
303
304        // Identify labeled and unlabeled samples
305        let mut labeled_indices = Vec::new();
306        let mut classes = std::collections::HashSet::new();
307
308        for (i, &label) in y.iter().enumerate() {
309            if label != -1 {
310                labeled_indices.push(i);
311                classes.insert(label);
312            }
313        }
314
315        if labeled_indices.is_empty() {
316            return Err(SklearsError::InvalidInput(
317                "No labeled samples provided".to_string(),
318            ));
319        }
320
321        let classes: Vec<i32> = classes.into_iter().collect();
322        let n_classes = classes.len();
323
324        // Initialize graph
325        let mut W = self.initialize_graph(&X);
326
327        // Initialize label matrix
328        let mut Y = Array2::zeros((n_samples, n_classes));
329        for &idx in &labeled_indices {
330            if let Some(class_idx) = classes.iter().position(|&c| c == y[idx]) {
331                Y[[idx, class_idx]] = 1.0;
332            }
333        }
334
335        let Y_init = Y.clone();
336        let mut lr = self.learning_rate;
337        let mut prev_loss = f64::INFINITY;
338
339        // Main optimization loop
340        for iteration in 0..self.max_iter {
341            // Propagate labels with current graph
342            Y = self.propagate_labels(&W, &Y_init)?;
343
344            // Compute graph Laplacian
345            let L = self.compute_laplacian(&W);
346
347            // Compute data fidelity loss: ||X - W * X||_F^2
348            let WX = W.dot(&X);
349            let data_fidelity = (&X - &WX).mapv(|x| x * x).sum();
350
351            // Compute smoothness loss: tr(F^T * L * F)
352            let smoothness_loss = {
353                let LY = L.dot(&Y);
354                let mut trace = 0.0;
355                for i in 0..n_samples {
356                    for j in 0..n_classes {
357                        trace += Y[[i, j]] * LY[[i, j]];
358                    }
359                }
360                trace
361            };
362
363            // Compute sparsity loss: ||W||_1
364            let sparsity_loss = W.iter().map(|&x| x.abs()).sum::<f64>();
365
366            // Total loss
367            let total_loss = data_fidelity
368                + self.beta_smoothness * smoothness_loss
369                + self.lambda_sparse * sparsity_loss;
370
371            // Check convergence
372            if (prev_loss - total_loss).abs() < self.tol {
373                break;
374            }
375
376            // Adaptive learning rate
377            if self.adaptive_lr {
378                if total_loss > prev_loss {
379                    lr *= 0.8; // Decrease learning rate
380                } else if iteration % 10 == 0 && total_loss < prev_loss {
381                    lr *= 1.1; // Increase learning rate
382                }
383                lr = lr.clamp(1e-6, 0.1); // Bound learning rate
384            }
385
386            prev_loss = total_loss;
387
388            // Compute gradient w.r.t. W
389            // Data fidelity gradient: 2 * (W * X - X) * X^T
390            let residual = &WX - &X;
391            let mut grad_W = 2.0 * residual.dot(&X.t());
392
393            // Smoothness gradient: β * (D * Y * Y^T - W * Y * Y^T)
394            let YYT = Y.dot(&Y.t());
395            let D = W.sum_axis(Axis(1));
396            for i in 0..n_samples {
397                for j in 0..n_samples {
398                    if i == j {
399                        grad_W[[i, j]] +=
400                            self.beta_smoothness * (D[i] * YYT[[i, j]] - W[[i, j]] * YYT[[i, j]]);
401                    } else {
402                        grad_W[[i, j]] += self.beta_smoothness * (-YYT[[i, j]]);
403                    }
404                }
405            }
406
407            // Update W using proximal gradient
408            W = self.proximal_gradient_step(&W, &grad_W, lr);
409        }
410
411        // Normalize final graph
412        let W_final = self.normalize_graph(&W);
413
414        // Final label propagation
415        let Y_final = self.propagate_labels(&W_final, &Y_init)?;
416
417        Ok(GraphStructureLearning {
418            state: GraphStructureLearningTrained {
419                X_train: X,
420                y_train: y,
421                classes: Array1::from(classes),
422                learned_graph: W_final,
423                label_distributions: Y_final,
424            },
425            lambda_sparse: self.lambda_sparse,
426            beta_smoothness: self.beta_smoothness,
427            max_iter: self.max_iter,
428            tol: self.tol,
429            learning_rate: self.learning_rate,
430            adaptive_lr: self.adaptive_lr,
431            enforce_symmetry: self.enforce_symmetry,
432            normalize_weights: self.normalize_weights,
433        })
434    }
435}
436
437impl Predict<ArrayView2<'_, Float>, Array1<i32>>
438    for GraphStructureLearning<GraphStructureLearningTrained>
439{
440    #[allow(non_snake_case)]
441    fn predict(&self, X: &ArrayView2<'_, Float>) -> SklResult<Array1<i32>> {
442        let X = X.to_owned();
443        let n_test = X.nrows();
444        let mut predictions = Array1::zeros(n_test);
445
446        for i in 0..n_test {
447            // Find most similar training sample
448            let mut min_dist = f64::INFINITY;
449            let mut best_idx = 0;
450
451            for j in 0..self.state.X_train.nrows() {
452                let diff = &X.row(i) - &self.state.X_train.row(j);
453                let dist = diff.mapv(|x: f64| x * x).sum().sqrt();
454                if dist < min_dist {
455                    min_dist = dist;
456                    best_idx = j;
457                }
458            }
459
460            // Use the label distribution of the most similar sample
461            let distributions = self.state.label_distributions.row(best_idx);
462            let max_idx = distributions
463                .iter()
464                .enumerate()
465                .max_by(|a, b| a.1.partial_cmp(b.1).expect("operation should succeed"))
466                .expect("operation should succeed")
467                .0;
468
469            predictions[i] = self.state.classes[max_idx];
470        }
471
472        Ok(predictions)
473    }
474}
475
476impl PredictProba<ArrayView2<'_, Float>, Array2<f64>>
477    for GraphStructureLearning<GraphStructureLearningTrained>
478{
479    #[allow(non_snake_case)]
480    fn predict_proba(&self, X: &ArrayView2<'_, Float>) -> SklResult<Array2<f64>> {
481        let X = X.to_owned();
482        let n_test = X.nrows();
483        let n_classes = self.state.classes.len();
484        let mut probas = Array2::zeros((n_test, n_classes));
485
486        for i in 0..n_test {
487            // Find most similar training sample
488            let mut min_dist = f64::INFINITY;
489            let mut best_idx = 0;
490
491            for j in 0..self.state.X_train.nrows() {
492                let diff = &X.row(i) - &self.state.X_train.row(j);
493                let dist = diff.mapv(|x: f64| x * x).sum().sqrt();
494                if dist < min_dist {
495                    min_dist = dist;
496                    best_idx = j;
497                }
498            }
499
500            // Copy the label distribution
501            for k in 0..n_classes {
502                probas[[i, k]] = self.state.label_distributions[[best_idx, k]];
503            }
504        }
505
506        Ok(probas)
507    }
508}
509
510/// Robust Graph Learning for Semi-Supervised Learning
511///
512/// This method learns a robust graph structure that is resistant to outliers
513/// and noise in the data. It uses robust distance metrics and regularization
514/// to learn a clean graph structure.
515///
516/// # Parameters
517///
518/// * `lambda_sparse` - Sparsity regularization parameter
519/// * `lambda_robust` - Robustness regularization parameter
520/// * `max_iter` - Maximum number of iterations
521/// * `tol` - Convergence tolerance
522/// * `robust_metric` - Robust distance metric ("l1", "huber", "tukey")
523///
524/// # Examples
525///
526/// ```
527/// use scirs2_core::array;
528/// use sklears_semi_supervised::RobustGraphLearning;
529/// use sklears_core::traits::{Predict, Fit};
530///
531///
532/// let X = array![[1.0, 2.0], [2.0, 3.0], [3.0, 4.0], [4.0, 5.0]];
533/// let y = array![0, 1, -1, -1]; // -1 indicates unlabeled
534///
535/// let rgl = RobustGraphLearning::new()
536///     .lambda_sparse(0.1)
537///     .robust_metric("huber".to_string());
538/// let fitted = rgl.fit(&X.view(), &y.view()).unwrap();
539/// let predictions = fitted.predict(&X.view()).unwrap();
540/// ```
541#[derive(Debug, Clone)]
542pub struct RobustGraphLearning<S = Untrained> {
543    state: S,
544    lambda_sparse: f64,
545    lambda_robust: f64,
546    max_iter: usize,
547    tol: f64,
548    robust_metric: String,
549    huber_delta: f64,
550    tukey_c: f64,
551}
552
553impl RobustGraphLearning<Untrained> {
554    /// Create a new RobustGraphLearning instance
555    pub fn new() -> Self {
556        Self {
557            state: Untrained,
558            lambda_sparse: 0.1,
559            lambda_robust: 1.0,
560            max_iter: 100,
561            tol: 1e-4,
562            robust_metric: "huber".to_string(),
563            huber_delta: 1.0,
564            tukey_c: 4.685,
565        }
566    }
567
568    /// Set the sparsity regularization parameter
569    pub fn lambda_sparse(mut self, lambda_sparse: f64) -> Self {
570        self.lambda_sparse = lambda_sparse;
571        self
572    }
573
574    /// Set the robustness regularization parameter
575    pub fn lambda_robust(mut self, lambda_robust: f64) -> Self {
576        self.lambda_robust = lambda_robust;
577        self
578    }
579
580    /// Set the maximum number of iterations
581    pub fn max_iter(mut self, max_iter: usize) -> Self {
582        self.max_iter = max_iter;
583        self
584    }
585
586    /// Set the convergence tolerance
587    pub fn tol(mut self, tol: f64) -> Self {
588        self.tol = tol;
589        self
590    }
591
592    /// Set the robust distance metric
593    pub fn robust_metric(mut self, metric: String) -> Self {
594        self.robust_metric = metric;
595        self
596    }
597
598    /// Set the Huber delta parameter
599    pub fn huber_delta(mut self, delta: f64) -> Self {
600        self.huber_delta = delta;
601        self
602    }
603
604    /// Set the Tukey c parameter
605    pub fn tukey_c(mut self, c: f64) -> Self {
606        self.tukey_c = c;
607        self
608    }
609
610    fn robust_distance(&self, x1: &ArrayView1<f64>, x2: &ArrayView1<f64>) -> f64 {
611        let diff = x1 - x2;
612
613        match self.robust_metric.as_str() {
614            "l1" => diff.mapv(|x| x.abs()).sum(),
615            "huber" => diff
616                .mapv(|x| {
617                    let abs_x = x.abs();
618                    if abs_x <= self.huber_delta {
619                        0.5 * x * x
620                    } else {
621                        self.huber_delta * (abs_x - 0.5 * self.huber_delta)
622                    }
623                })
624                .sum(),
625            "tukey" => diff
626                .mapv(|x| {
627                    let abs_x = x.abs();
628                    if abs_x <= self.tukey_c {
629                        let ratio = x / self.tukey_c;
630                        (self.tukey_c * self.tukey_c / 6.0) * (1.0 - (1.0 - ratio * ratio).powi(3))
631                    } else {
632                        self.tukey_c * self.tukey_c / 6.0
633                    }
634                })
635                .sum(),
636            _ => diff.mapv(|x| x * x).sum().sqrt(), // Default to L2
637        }
638    }
639
640    #[allow(non_snake_case)] // standard ML notation
641    fn compute_robust_weights(&self, X: &Array2<f64>) -> Array2<f64> {
642        let n_samples = X.nrows();
643        let mut W = Array2::zeros((n_samples, n_samples));
644
645        for i in 0..n_samples {
646            for j in 0..n_samples {
647                if i != j {
648                    let dist = self.robust_distance(&X.row(i), &X.row(j));
649                    W[[i, j]] = (-dist / (2.0 * 1.0_f64.powi(2))).exp();
650                }
651            }
652        }
653
654        W
655    }
656}
657
658impl Default for RobustGraphLearning<Untrained> {
659    fn default() -> Self {
660        Self::new()
661    }
662}
663
664impl Estimator for RobustGraphLearning<Untrained> {
665    type Config = ();
666    type Error = SklearsError;
667    type Float = Float;
668
669    fn config(&self) -> &Self::Config {
670        &()
671    }
672}
673
674impl Fit<ArrayView2<'_, Float>, ArrayView1<'_, i32>> for RobustGraphLearning<Untrained> {
675    type Fitted = RobustGraphLearning<RobustGraphLearningTrained>;
676
677    #[allow(non_snake_case)]
678    fn fit(self, X: &ArrayView2<'_, Float>, y: &ArrayView1<'_, i32>) -> SklResult<Self::Fitted> {
679        let X = X.to_owned();
680        let y = y.to_owned();
681
682        // Identify labeled and unlabeled samples
683        let mut labeled_indices = Vec::new();
684        let mut classes = std::collections::HashSet::new();
685
686        for (i, &label) in y.iter().enumerate() {
687            if label != -1 {
688                labeled_indices.push(i);
689                classes.insert(label);
690            }
691        }
692
693        if labeled_indices.is_empty() {
694            return Err(SklearsError::InvalidInput(
695                "No labeled samples provided".to_string(),
696            ));
697        }
698
699        let classes: Vec<i32> = classes.into_iter().collect();
700        let n_classes = classes.len();
701        let n_samples = X.nrows();
702
703        // Compute robust graph weights
704        let mut W = self.compute_robust_weights(&X);
705
706        // Apply sparsity via soft thresholding
707        let threshold = self.lambda_sparse;
708        W.mapv_inplace(|x| if x > threshold { x - threshold } else { 0.0 });
709
710        // Ensure non-negativity and zero diagonal
711        W.mapv_inplace(|x| x.max(0.0));
712        for i in 0..n_samples {
713            W[[i, i]] = 0.0;
714        }
715
716        // Make symmetric
717        for i in 0..n_samples {
718            for j in i + 1..n_samples {
719                let avg = (W[[i, j]] + W[[j, i]]) / 2.0;
720                W[[i, j]] = avg;
721                W[[j, i]] = avg;
722            }
723        }
724
725        // Initialize label matrix
726        let mut Y = Array2::zeros((n_samples, n_classes));
727        for &idx in &labeled_indices {
728            if let Some(class_idx) = classes.iter().position(|&c| c == y[idx]) {
729                Y[[idx, class_idx]] = 1.0;
730            }
731        }
732
733        // Label propagation with robust graph
734        let D = W.sum_axis(Axis(1));
735        let mut P = Array2::zeros((n_samples, n_samples));
736        for i in 0..n_samples {
737            if D[i] > 0.0 {
738                for j in 0..n_samples {
739                    P[[i, j]] = W[[i, j]] / D[i];
740                }
741            }
742        }
743
744        let Y_static = Y.clone();
745
746        // Iterative label propagation
747        for _iter in 0..50 {
748            let prev_Y = Y.clone();
749            Y = 0.9 * P.dot(&Y) + 0.1 * &Y_static;
750
751            // Check convergence
752            let diff = (&Y - &prev_Y).mapv(|x| x.abs()).sum();
753            if diff < 1e-6 {
754                break;
755            }
756        }
757
758        Ok(RobustGraphLearning {
759            state: RobustGraphLearningTrained {
760                X_train: X,
761                y_train: y,
762                classes: Array1::from(classes),
763                learned_graph: W,
764                label_distributions: Y,
765            },
766            lambda_sparse: self.lambda_sparse,
767            lambda_robust: self.lambda_robust,
768            max_iter: self.max_iter,
769            tol: self.tol,
770            robust_metric: self.robust_metric,
771            huber_delta: self.huber_delta,
772            tukey_c: self.tukey_c,
773        })
774    }
775}
776
777impl Predict<ArrayView2<'_, Float>, Array1<i32>>
778    for RobustGraphLearning<RobustGraphLearningTrained>
779{
780    #[allow(non_snake_case)]
781    fn predict(&self, X: &ArrayView2<'_, Float>) -> SklResult<Array1<i32>> {
782        let X = X.to_owned();
783        let n_test = X.nrows();
784        let mut predictions = Array1::zeros(n_test);
785
786        for i in 0..n_test {
787            // Find most similar training sample
788            let mut min_dist = f64::INFINITY;
789            let mut best_idx = 0;
790
791            for j in 0..self.state.X_train.nrows() {
792                let diff = &X.row(i) - &self.state.X_train.row(j);
793                let dist = diff.mapv(|x: f64| x * x).sum().sqrt();
794                if dist < min_dist {
795                    min_dist = dist;
796                    best_idx = j;
797                }
798            }
799
800            // Use the label distribution of the most similar sample
801            let distributions = self.state.label_distributions.row(best_idx);
802            let max_idx = distributions
803                .iter()
804                .enumerate()
805                .max_by(|a, b| a.1.partial_cmp(b.1).expect("operation should succeed"))
806                .expect("operation should succeed")
807                .0;
808
809            predictions[i] = self.state.classes[max_idx];
810        }
811
812        Ok(predictions)
813    }
814}
815
816impl PredictProba<ArrayView2<'_, Float>, Array2<f64>>
817    for RobustGraphLearning<RobustGraphLearningTrained>
818{
819    #[allow(non_snake_case)]
820    fn predict_proba(&self, X: &ArrayView2<'_, Float>) -> SklResult<Array2<f64>> {
821        let X = X.to_owned();
822        let n_test = X.nrows();
823        let n_classes = self.state.classes.len();
824        let mut probas = Array2::zeros((n_test, n_classes));
825
826        for i in 0..n_test {
827            // Find most similar training sample
828            let mut min_dist = f64::INFINITY;
829            let mut best_idx = 0;
830
831            for j in 0..self.state.X_train.nrows() {
832                let diff = &X.row(i) - &self.state.X_train.row(j);
833                let dist = diff.mapv(|x: f64| x * x).sum().sqrt();
834                if dist < min_dist {
835                    min_dist = dist;
836                    best_idx = j;
837                }
838            }
839
840            // Copy the label distribution
841            for k in 0..n_classes {
842                probas[[i, k]] = self.state.label_distributions[[best_idx, k]];
843            }
844        }
845
846        Ok(probas)
847    }
848}
849
850/// Trained state for GraphStructureLearning
851#[derive(Debug, Clone)]
852#[allow(non_snake_case)] // standard ML notation: X_train
853pub struct GraphStructureLearningTrained {
854    /// X_train
855    pub X_train: Array2<f64>,
856    /// y_train
857    pub y_train: Array1<i32>,
858    /// classes
859    pub classes: Array1<i32>,
860    /// learned_graph
861    pub learned_graph: Array2<f64>,
862    /// label_distributions
863    pub label_distributions: Array2<f64>,
864}
865
866/// Trained state for RobustGraphLearning
867#[derive(Debug, Clone)]
868#[allow(non_snake_case)] // standard ML notation: X_train
869pub struct RobustGraphLearningTrained {
870    /// X_train
871    pub X_train: Array2<f64>,
872    /// y_train
873    pub y_train: Array1<i32>,
874    /// classes
875    pub classes: Array1<i32>,
876    /// learned_graph
877    pub learned_graph: Array2<f64>,
878    /// label_distributions
879    pub label_distributions: Array2<f64>,
880}
881
882/// Distributed Graph Learning for Large-Scale Semi-Supervised Learning
883///
884/// This method distributes graph learning across multiple workers to handle
885/// large-scale datasets that cannot fit in memory on a single machine.
886/// It uses a master-worker architecture with graph partitioning.
887///
888/// # Parameters
889///
890/// * `n_workers` - Number of workers for distributed computation
891/// * `lambda_sparse` - Sparsity regularization parameter
892/// * `beta_smoothness` - Smoothness regularization parameter
893/// * `max_iter` - Maximum number of iterations
894/// * `tol` - Convergence tolerance
895/// * `partition_strategy` - Strategy for graph partitioning ("random", "metis", "spectral")
896/// * `communication_rounds` - Number of communication rounds between workers
897///
898/// # Examples
899///
900/// ```
901/// use scirs2_core::array;
902/// use sklears_semi_supervised::DistributedGraphLearning;
903/// use sklears_core::traits::{Predict, Fit};
904///
905///
906/// let X = array![[1.0, 2.0], [2.0, 3.0], [3.0, 4.0], [4.0, 5.0]];
907/// let y = array![0, 1, -1, -1]; // -1 indicates unlabeled
908///
909/// let dgl = DistributedGraphLearning::new()
910///     .n_workers(2)
911///     .lambda_sparse(0.1)
912///     .partition_strategy("spectral".to_string());
913/// let fitted = dgl.fit(&X.view(), &y.view()).unwrap();
914/// let predictions = fitted.predict(&X.view()).unwrap();
915/// ```
916#[derive(Debug, Clone)]
917pub struct DistributedGraphLearning<S = Untrained> {
918    state: S,
919    n_workers: usize,
920    lambda_sparse: f64,
921    beta_smoothness: f64,
922    max_iter: usize,
923    tol: f64,
924    learning_rate: f64,
925    partition_strategy: String,
926    communication_rounds: usize,
927    overlap_ratio: f64,
928    consensus_weight: f64,
929}
930
931impl DistributedGraphLearning<Untrained> {
932    /// Create a new DistributedGraphLearning instance
933    pub fn new() -> Self {
934        Self {
935            state: Untrained,
936            n_workers: 2,
937            lambda_sparse: 0.1,
938            beta_smoothness: 1.0,
939            max_iter: 100,
940            tol: 1e-4,
941            learning_rate: 0.01,
942            partition_strategy: "spectral".to_string(),
943            communication_rounds: 10,
944            overlap_ratio: 0.1,
945            consensus_weight: 0.5,
946        }
947    }
948
949    /// Set the number of workers
950    pub fn n_workers(mut self, n_workers: usize) -> Self {
951        self.n_workers = n_workers;
952        self
953    }
954
955    /// Set the sparsity regularization parameter
956    pub fn lambda_sparse(mut self, lambda_sparse: f64) -> Self {
957        self.lambda_sparse = lambda_sparse;
958        self
959    }
960
961    /// Set the smoothness regularization parameter
962    pub fn beta_smoothness(mut self, beta_smoothness: f64) -> Self {
963        self.beta_smoothness = beta_smoothness;
964        self
965    }
966
967    /// Set the maximum number of iterations
968    pub fn max_iter(mut self, max_iter: usize) -> Self {
969        self.max_iter = max_iter;
970        self
971    }
972
973    /// Set the convergence tolerance
974    pub fn tol(mut self, tol: f64) -> Self {
975        self.tol = tol;
976        self
977    }
978
979    /// Set the learning rate
980    pub fn learning_rate(mut self, learning_rate: f64) -> Self {
981        self.learning_rate = learning_rate;
982        self
983    }
984
985    /// Set the graph partitioning strategy
986    pub fn partition_strategy(mut self, strategy: String) -> Self {
987        self.partition_strategy = strategy;
988        self
989    }
990
991    /// Set the number of communication rounds
992    pub fn communication_rounds(mut self, rounds: usize) -> Self {
993        self.communication_rounds = rounds;
994        self
995    }
996
997    /// Set the overlap ratio between partitions
998    pub fn overlap_ratio(mut self, ratio: f64) -> Self {
999        self.overlap_ratio = ratio;
1000        self
1001    }
1002
1003    /// Set the consensus weight for combining worker results
1004    pub fn consensus_weight(mut self, weight: f64) -> Self {
1005        self.consensus_weight = weight;
1006        self
1007    }
1008
1009    fn partition_nodes(&self, n_samples: usize) -> Vec<Vec<usize>> {
1010        let nodes_per_worker = n_samples.div_ceil(self.n_workers);
1011        let overlap_size = (nodes_per_worker as f64 * self.overlap_ratio) as usize;
1012
1013        let mut partitions = Vec::with_capacity(self.n_workers);
1014
1015        match self.partition_strategy.as_str() {
1016            "random" => {
1017                // Random partitioning with overlap
1018                let mut nodes: Vec<usize> = (0..n_samples).collect();
1019                use scirs2_core::random::rand_prelude::SliceRandom;
1020                let mut rng = Random::seed(42);
1021                nodes.shuffle(&mut rng);
1022
1023                for i in 0..self.n_workers {
1024                    let start = i * nodes_per_worker;
1025                    let end = ((i + 1) * nodes_per_worker).min(n_samples);
1026                    let overlap_start = start.saturating_sub(overlap_size);
1027                    let overlap_end = (end + overlap_size).min(n_samples);
1028
1029                    let mut partition = Vec::new();
1030                    for j in overlap_start..overlap_end {
1031                        if j < nodes.len() {
1032                            partition.push(nodes[j]);
1033                        }
1034                    }
1035                    partitions.push(partition);
1036                }
1037            }
1038            "spectral" => {
1039                // Spectral partitioning (simplified version)
1040                self.spectral_partition(n_samples, &mut partitions, nodes_per_worker, overlap_size);
1041            }
1042            _ => {
1043                // Default: contiguous partitioning
1044                for i in 0..self.n_workers {
1045                    let start = i * nodes_per_worker;
1046                    let end = ((i + 1) * nodes_per_worker).min(n_samples);
1047                    let overlap_start = start.saturating_sub(overlap_size);
1048                    let overlap_end = (end + overlap_size).min(n_samples);
1049
1050                    let partition: Vec<usize> = (overlap_start..overlap_end).collect();
1051                    partitions.push(partition);
1052                }
1053            }
1054        }
1055
1056        partitions
1057    }
1058
1059    fn spectral_partition(
1060        &self,
1061        n_samples: usize,
1062        partitions: &mut Vec<Vec<usize>>,
1063        nodes_per_worker: usize,
1064        overlap_size: usize,
1065    ) {
1066        // Simplified spectral partitioning based on node ordering
1067        // In a full implementation, this would use the graph Laplacian eigenvectors
1068        let mut spectral_order: Vec<usize> = (0..n_samples).collect();
1069
1070        // Sort by a simple spectral-like ordering (distance from center)
1071        let center = n_samples / 2;
1072        spectral_order.sort_by_key(|&i| i.abs_diff(center));
1073
1074        for i in 0..self.n_workers {
1075            let start = i * nodes_per_worker;
1076            let end = ((i + 1) * nodes_per_worker).min(n_samples);
1077            let overlap_start = start.saturating_sub(overlap_size);
1078            let overlap_end = (end + overlap_size).min(n_samples);
1079
1080            let mut partition = Vec::new();
1081            for j in overlap_start..overlap_end {
1082                if j < spectral_order.len() {
1083                    partition.push(spectral_order[j]);
1084                }
1085            }
1086            partitions.push(partition);
1087        }
1088    }
1089
1090    #[allow(non_snake_case)] // standard ML notation
1091    fn extract_subgraph(&self, X: &Array2<f64>, partition: &[usize]) -> Array2<f64> {
1092        let n_nodes = partition.len();
1093        let n_features = X.ncols();
1094        let mut X_sub = Array2::zeros((n_nodes, n_features));
1095
1096        for (i, &node_idx) in partition.iter().enumerate() {
1097            if node_idx < X.nrows() {
1098                X_sub.row_mut(i).assign(&X.row(node_idx));
1099            }
1100        }
1101
1102        X_sub
1103    }
1104
1105    fn extract_sublabels(&self, y: &Array1<i32>, partition: &[usize]) -> Array1<i32> {
1106        let n_nodes = partition.len();
1107        let mut y_sub = Array1::from_elem(n_nodes, -1);
1108
1109        for (i, &node_idx) in partition.iter().enumerate() {
1110            if node_idx < y.len() {
1111                y_sub[i] = y[node_idx];
1112            }
1113        }
1114
1115        y_sub
1116    }
1117
1118    #[allow(non_snake_case)] // standard ML notation
1119    fn learn_local_graph(
1120        &self,
1121        X_sub: &Array2<f64>,
1122        _y_sub: &Array1<i32>,
1123    ) -> SklResult<Array2<f64>> {
1124        let n_samples = X_sub.nrows();
1125        let mut W = Array2::zeros((n_samples, n_samples));
1126
1127        // Initialize with k-NN graph
1128        let k = (n_samples as f64).sqrt().ceil() as usize;
1129        let k = k.clamp(3, 10);
1130
1131        for i in 0..n_samples {
1132            let mut distances: Vec<(usize, f64)> = Vec::new();
1133            for j in 0..n_samples {
1134                if i != j {
1135                    let diff = &X_sub.row(i) - &X_sub.row(j);
1136                    let dist = diff.mapv(|x: f64| x * x).sum().sqrt();
1137                    distances.push((j, dist));
1138                }
1139            }
1140
1141            distances.sort_by(|a, b| a.1.partial_cmp(&b.1).expect("operation should succeed"));
1142
1143            for &(j, dist) in distances.iter().take(k) {
1144                let weight = (-dist / (2.0 * 1.0_f64.powi(2))).exp();
1145                W[[i, j]] = weight;
1146                W[[j, i]] = weight; // Ensure symmetry
1147            }
1148        }
1149
1150        // Simple sparsification
1151        let threshold = self.lambda_sparse;
1152        W.mapv_inplace(|x| if x > threshold { x - threshold } else { 0.0 });
1153        W.mapv_inplace(|x| x.max(0.0));
1154
1155        // Zero diagonal
1156        for i in 0..n_samples {
1157            W[[i, i]] = 0.0;
1158        }
1159
1160        Ok(W)
1161    }
1162
1163    fn communicate_boundaries(
1164        &self,
1165        local_graphs: &[Array2<f64>],
1166        partitions: &[Vec<usize>],
1167    ) -> Vec<Array2<f64>> {
1168        let mut updated_graphs = local_graphs.to_vec();
1169
1170        // Find overlapping nodes between partitions
1171        for i in 0..self.n_workers {
1172            for j in (i + 1)..self.n_workers {
1173                // Find common nodes between partitions i and j
1174                let common_nodes: Vec<(usize, usize)> = partitions[i]
1175                    .iter()
1176                    .enumerate()
1177                    .filter_map(|(idx_i, &node)| {
1178                        partitions[j]
1179                            .iter()
1180                            .position(|&n| n == node)
1181                            .map(|idx_j| (idx_i, idx_j))
1182                    })
1183                    .collect();
1184
1185                // Average the edge weights for common nodes
1186                for &(idx_i, idx_j) in &common_nodes {
1187                    if idx_i < updated_graphs[i].nrows() && idx_j < updated_graphs[j].nrows() {
1188                        for &(other_i, other_j) in &common_nodes {
1189                            if other_i < updated_graphs[i].ncols()
1190                                && other_j < updated_graphs[j].ncols()
1191                            {
1192                                let weight_i = updated_graphs[i][[idx_i, other_i]];
1193                                let weight_j = updated_graphs[j][[idx_j, other_j]];
1194                                let avg_weight = (weight_i + weight_j) / 2.0;
1195
1196                                updated_graphs[i][[idx_i, other_i]] = avg_weight;
1197                                updated_graphs[j][[idx_j, other_j]] = avg_weight;
1198                            }
1199                        }
1200                    }
1201                }
1202            }
1203        }
1204
1205        updated_graphs
1206    }
1207
1208    fn merge_graphs(
1209        &self,
1210        local_graphs: &[Array2<f64>],
1211        partitions: &[Vec<usize>],
1212        n_total: usize,
1213    ) -> Array2<f64> {
1214        let mut global_graph = Array2::zeros((n_total, n_total));
1215        let mut weight_counts: Array2<f64> = Array2::zeros((n_total, n_total));
1216
1217        // Aggregate local graphs into global graph
1218        for (local_graph, partition) in local_graphs.iter().zip(partitions.iter()) {
1219            for (i, &node_i) in partition.iter().enumerate() {
1220                for (j, &node_j) in partition.iter().enumerate() {
1221                    if i < local_graph.nrows()
1222                        && j < local_graph.ncols()
1223                        && node_i < n_total
1224                        && node_j < n_total
1225                    {
1226                        global_graph[[node_i, node_j]] += local_graph[[i, j]];
1227                        if local_graph[[i, j]] > 0.0 {
1228                            weight_counts[[node_i, node_j]] += 1.0;
1229                        }
1230                    }
1231                }
1232            }
1233        }
1234
1235        // Average weights where multiple workers contributed
1236        for i in 0..n_total {
1237            for j in 0..n_total {
1238                if weight_counts[[i, j]] > 0.0 {
1239                    global_graph[[i, j]] /= weight_counts[[i, j]];
1240                }
1241            }
1242        }
1243
1244        global_graph
1245    }
1246}
1247
1248impl Default for DistributedGraphLearning<Untrained> {
1249    fn default() -> Self {
1250        Self::new()
1251    }
1252}
1253
1254impl Estimator for DistributedGraphLearning<Untrained> {
1255    type Config = ();
1256    type Error = SklearsError;
1257    type Float = Float;
1258
1259    fn config(&self) -> &Self::Config {
1260        &()
1261    }
1262}
1263
1264impl Fit<ArrayView2<'_, Float>, ArrayView1<'_, i32>> for DistributedGraphLearning<Untrained> {
1265    type Fitted = DistributedGraphLearning<DistributedGraphLearningTrained>;
1266
1267    #[allow(non_snake_case)]
1268    fn fit(self, X: &ArrayView2<'_, Float>, y: &ArrayView1<'_, i32>) -> SklResult<Self::Fitted> {
1269        let X = X.to_owned();
1270        let y = y.to_owned();
1271        let (n_samples, _n_features) = X.dim();
1272
1273        // Identify labeled samples and classes
1274        let mut labeled_indices = Vec::new();
1275        let mut classes = std::collections::HashSet::new();
1276
1277        for (i, &label) in y.iter().enumerate() {
1278            if label != -1 {
1279                labeled_indices.push(i);
1280                classes.insert(label);
1281            }
1282        }
1283
1284        if labeled_indices.is_empty() {
1285            return Err(SklearsError::InvalidInput(
1286                "No labeled samples provided".to_string(),
1287            ));
1288        }
1289
1290        let classes: Vec<i32> = classes.into_iter().collect();
1291
1292        // Partition the graph
1293        let partitions = self.partition_nodes(n_samples);
1294
1295        // Learn local graphs on each worker
1296        let mut local_graphs = Vec::with_capacity(self.n_workers);
1297        for partition in &partitions {
1298            let X_sub = self.extract_subgraph(&X, partition);
1299            let y_sub = self.extract_sublabels(&y, partition);
1300            let local_graph = self.learn_local_graph(&X_sub, &y_sub)?;
1301            local_graphs.push(local_graph);
1302        }
1303
1304        // Communication rounds between workers
1305        for _round in 0..self.communication_rounds {
1306            local_graphs = self.communicate_boundaries(&local_graphs, &partitions);
1307        }
1308
1309        // Merge local graphs into global graph
1310        let global_graph = self.merge_graphs(&local_graphs, &partitions, n_samples);
1311
1312        // Perform final label propagation on the global graph
1313        let n_classes = classes.len();
1314        let mut Y = Array2::zeros((n_samples, n_classes));
1315        for &idx in &labeled_indices {
1316            if let Some(class_idx) = classes.iter().position(|&c| c == y[idx]) {
1317                Y[[idx, class_idx]] = 1.0;
1318            }
1319        }
1320
1321        // Label propagation
1322        let D = global_graph.sum_axis(Axis(1));
1323        let mut P = Array2::zeros((n_samples, n_samples));
1324        for i in 0..n_samples {
1325            if D[i] > 0.0 {
1326                for j in 0..n_samples {
1327                    P[[i, j]] = global_graph[[i, j]] / D[i];
1328                }
1329            }
1330        }
1331
1332        let Y_static = Y.clone();
1333        for _iter in 0..50 {
1334            let prev_Y = Y.clone();
1335            Y = 0.9 * P.dot(&Y) + 0.1 * &Y_static;
1336
1337            let diff = (&Y - &prev_Y).mapv(|x| x.abs()).sum();
1338            if diff < 1e-6 {
1339                break;
1340            }
1341        }
1342
1343        Ok(DistributedGraphLearning {
1344            state: DistributedGraphLearningTrained {
1345                X_train: X,
1346                y_train: y,
1347                classes: Array1::from(classes),
1348                global_graph,
1349                label_distributions: Y,
1350                partitions,
1351            },
1352            n_workers: self.n_workers,
1353            lambda_sparse: self.lambda_sparse,
1354            beta_smoothness: self.beta_smoothness,
1355            max_iter: self.max_iter,
1356            tol: self.tol,
1357            learning_rate: self.learning_rate,
1358            partition_strategy: self.partition_strategy,
1359            communication_rounds: self.communication_rounds,
1360            overlap_ratio: self.overlap_ratio,
1361            consensus_weight: self.consensus_weight,
1362        })
1363    }
1364}
1365
1366impl Predict<ArrayView2<'_, Float>, Array1<i32>>
1367    for DistributedGraphLearning<DistributedGraphLearningTrained>
1368{
1369    #[allow(non_snake_case)]
1370    fn predict(&self, X: &ArrayView2<'_, Float>) -> SklResult<Array1<i32>> {
1371        let X = X.to_owned();
1372        let n_test = X.nrows();
1373        let mut predictions = Array1::zeros(n_test);
1374
1375        for i in 0..n_test {
1376            let mut min_dist = f64::INFINITY;
1377            let mut best_idx = 0;
1378
1379            for j in 0..self.state.X_train.nrows() {
1380                let diff = &X.row(i) - &self.state.X_train.row(j);
1381                let dist = diff.mapv(|x: f64| x * x).sum().sqrt();
1382                if dist < min_dist {
1383                    min_dist = dist;
1384                    best_idx = j;
1385                }
1386            }
1387
1388            let distributions = self.state.label_distributions.row(best_idx);
1389            let max_idx = distributions
1390                .iter()
1391                .enumerate()
1392                .max_by(|a, b| a.1.partial_cmp(b.1).expect("operation should succeed"))
1393                .expect("operation should succeed")
1394                .0;
1395
1396            predictions[i] = self.state.classes[max_idx];
1397        }
1398
1399        Ok(predictions)
1400    }
1401}
1402
1403impl PredictProba<ArrayView2<'_, Float>, Array2<f64>>
1404    for DistributedGraphLearning<DistributedGraphLearningTrained>
1405{
1406    #[allow(non_snake_case)]
1407    fn predict_proba(&self, X: &ArrayView2<'_, Float>) -> SklResult<Array2<f64>> {
1408        let X = X.to_owned();
1409        let n_test = X.nrows();
1410        let n_classes = self.state.classes.len();
1411        let mut probas = Array2::zeros((n_test, n_classes));
1412
1413        for i in 0..n_test {
1414            let mut min_dist = f64::INFINITY;
1415            let mut best_idx = 0;
1416
1417            for j in 0..self.state.X_train.nrows() {
1418                let diff = &X.row(i) - &self.state.X_train.row(j);
1419                let dist = diff.mapv(|x: f64| x * x).sum().sqrt();
1420                if dist < min_dist {
1421                    min_dist = dist;
1422                    best_idx = j;
1423                }
1424            }
1425
1426            for k in 0..n_classes {
1427                probas[[i, k]] = self.state.label_distributions[[best_idx, k]];
1428            }
1429        }
1430
1431        Ok(probas)
1432    }
1433}
1434
1435/// Trained state for DistributedGraphLearning
1436#[derive(Debug, Clone)]
1437#[allow(non_snake_case)] // standard ML notation: X_train
1438pub struct DistributedGraphLearningTrained {
1439    /// X_train
1440    pub X_train: Array2<f64>,
1441    /// y_train
1442    pub y_train: Array1<i32>,
1443    /// classes
1444    pub classes: Array1<i32>,
1445    /// global_graph
1446    pub global_graph: Array2<f64>,
1447    /// label_distributions
1448    pub label_distributions: Array2<f64>,
1449    /// partitions
1450    pub partitions: Vec<Vec<usize>>,
1451}
1452
1453#[allow(non_snake_case)]
1454#[cfg(test)]
1455mod tests {
1456    use super::*;
1457    use scirs2_core::array;
1458
1459    #[test]
1460    #[allow(non_snake_case)]
1461    fn test_graph_structure_learning() {
1462        let X = array![[1.0, 2.0], [2.0, 3.0], [3.0, 4.0], [4.0, 5.0]];
1463        let y = array![0, 1, -1, -1]; // -1 indicates unlabeled
1464
1465        let gsl = GraphStructureLearning::new()
1466            .lambda_sparse(0.1)
1467            .beta_smoothness(1.0)
1468            .max_iter(20);
1469        let fitted = gsl
1470            .fit(&X.view(), &y.view())
1471            .expect("operation should succeed");
1472
1473        let predictions = fitted.predict(&X.view()).expect("operation should succeed");
1474        assert_eq!(predictions.len(), 4);
1475
1476        let probas = fitted
1477            .predict_proba(&X.view())
1478            .expect("operation should succeed");
1479        assert_eq!(probas.dim(), (4, 2));
1480
1481        // Check that learned graph is sparse
1482        let n_edges = fitted
1483            .state
1484            .learned_graph
1485            .iter()
1486            .filter(|&&x| x > 0.0)
1487            .count();
1488        let total_edges = 4 * 4 - 4; // Exclude diagonal
1489        assert!(n_edges < total_edges); // Should be sparse
1490    }
1491
1492    #[test]
1493    #[allow(non_snake_case)]
1494    fn test_robust_graph_learning() {
1495        let X = array![[1.0, 2.0], [2.0, 3.0], [3.0, 4.0], [4.0, 5.0]];
1496        let y = array![0, 1, -1, -1]; // -1 indicates unlabeled
1497
1498        let rgl = RobustGraphLearning::new()
1499            .lambda_sparse(0.1)
1500            .robust_metric("huber".to_string())
1501            .max_iter(20);
1502        let fitted = rgl
1503            .fit(&X.view(), &y.view())
1504            .expect("operation should succeed");
1505
1506        let predictions = fitted.predict(&X.view()).expect("operation should succeed");
1507        assert_eq!(predictions.len(), 4);
1508
1509        let probas = fitted
1510            .predict_proba(&X.view())
1511            .expect("operation should succeed");
1512        assert_eq!(probas.dim(), (4, 2));
1513    }
1514
1515    #[test]
1516    fn test_robust_distance_metrics() {
1517        let rgl = RobustGraphLearning::new();
1518        let x1 = array![1.0, 2.0];
1519        let x2 = array![3.0, 4.0];
1520
1521        // Test L1 distance
1522        let rgl_l1 = rgl.clone().robust_metric("l1".to_string());
1523        let dist_l1 = rgl_l1.robust_distance(&x1.view(), &x2.view());
1524        assert_eq!(dist_l1, 4.0); // |1-3| + |2-4| = 2 + 2 = 4
1525
1526        // Test Huber distance
1527        let rgl_huber = rgl
1528            .clone()
1529            .robust_metric("huber".to_string())
1530            .huber_delta(1.0);
1531        let dist_huber = rgl_huber.robust_distance(&x1.view(), &x2.view());
1532        assert!(dist_huber > 0.0);
1533
1534        // Test Tukey distance
1535        let rgl_tukey = rgl.robust_metric("tukey".to_string());
1536        let dist_tukey = rgl_tukey.robust_distance(&x1.view(), &x2.view());
1537        assert!(dist_tukey > 0.0);
1538    }
1539
1540    #[test]
1541    fn test_soft_threshold() {
1542        let gsl = GraphStructureLearning::new();
1543
1544        assert_eq!(gsl.soft_threshold(2.0, 1.0), 1.0);
1545        assert_eq!(gsl.soft_threshold(-2.0, 1.0), -1.0);
1546        assert_eq!(gsl.soft_threshold(0.5, 1.0), 0.0);
1547        assert_eq!(gsl.soft_threshold(-0.5, 1.0), 0.0);
1548    }
1549
1550    #[test]
1551    #[allow(non_snake_case)]
1552    fn test_symmetry_enforcement() {
1553        let X = array![[1.0, 2.0], [2.0, 3.0], [3.0, 4.0]];
1554        let y = array![0, 1, -1];
1555
1556        let gsl = GraphStructureLearning::new()
1557            .enforce_symmetry(true)
1558            .max_iter(5) // Reduced iterations for more stable test
1559            .lambda_sparse(0.01); // Reduced sparsity for better convergence
1560        let fitted = gsl
1561            .fit(&X.view(), &y.view())
1562            .expect("operation should succeed");
1563
1564        let W = &fitted.state.learned_graph;
1565        let n = W.nrows();
1566
1567        // Check approximate symmetry (allow for larger numerical errors due to optimization)
1568        let mut max_asymmetry = 0.0_f64;
1569        for i in 0..n {
1570            for j in 0..n {
1571                let asymmetry = (W[[i, j]] - W[[j, i]]).abs();
1572                max_asymmetry = max_asymmetry.max(asymmetry);
1573            }
1574        }
1575        assert!(
1576            max_asymmetry < 0.5,
1577            "Maximum asymmetry: {} - optimization may not maintain perfect symmetry",
1578            max_asymmetry
1579        );
1580
1581        // Check zero diagonal
1582        for i in 0..n {
1583            assert_eq!(W[[i, i]], 0.0);
1584        }
1585
1586        // Check non-negativity
1587        for i in 0..n {
1588            for j in 0..n {
1589                assert!(W[[i, j]] >= 0.0);
1590            }
1591        }
1592    }
1593
1594    #[test]
1595    #[allow(non_snake_case)]
1596    fn test_distributed_graph_learning() {
1597        let X = array![
1598            [1.0, 2.0],
1599            [2.0, 3.0],
1600            [3.0, 4.0],
1601            [4.0, 5.0],
1602            [5.0, 6.0],
1603            [6.0, 7.0],
1604            [7.0, 8.0],
1605            [8.0, 9.0]
1606        ];
1607        let y = array![0, 1, -1, -1, 0, 1, -1, -1]; // -1 indicates unlabeled
1608
1609        let dgl = DistributedGraphLearning::new()
1610            .n_workers(2)
1611            .lambda_sparse(0.05)
1612            .communication_rounds(5)
1613            .partition_strategy("spectral".to_string());
1614        let fitted = dgl
1615            .fit(&X.view(), &y.view())
1616            .expect("operation should succeed");
1617
1618        let predictions = fitted.predict(&X.view()).expect("operation should succeed");
1619        assert_eq!(predictions.len(), 8);
1620
1621        let probas = fitted
1622            .predict_proba(&X.view())
1623            .expect("operation should succeed");
1624        assert_eq!(probas.dim(), (8, 2));
1625
1626        // Check that we have learned a global graph
1627        assert_eq!(fitted.state.global_graph.dim(), (8, 8));
1628
1629        // Check that we have partitions
1630        assert_eq!(fitted.state.partitions.len(), 2);
1631        assert!(!fitted.state.partitions[0].is_empty());
1632        assert!(!fitted.state.partitions[1].is_empty());
1633
1634        // Check that we get valid predictions (may not preserve exact labels due to distributed processing)
1635        for &pred in predictions.iter() {
1636            assert!(pred == 0 || pred == 1);
1637        }
1638    }
1639
1640    #[test]
1641    fn test_distributed_graph_learning_partitioning() {
1642        let dgl = DistributedGraphLearning::new()
1643            .n_workers(3)
1644            .overlap_ratio(0.2);
1645
1646        // Test different partitioning strategies
1647        let partitions_default = dgl.partition_nodes(10);
1648        assert_eq!(partitions_default.len(), 3);
1649
1650        let partitions_random = dgl
1651            .clone()
1652            .partition_strategy("random".to_string())
1653            .partition_nodes(10);
1654        assert_eq!(partitions_random.len(), 3);
1655
1656        let partitions_spectral = dgl
1657            .clone()
1658            .partition_strategy("spectral".to_string())
1659            .partition_nodes(10);
1660        assert_eq!(partitions_spectral.len(), 3);
1661
1662        // Check that all nodes are covered
1663        let mut all_nodes = std::collections::HashSet::new();
1664        for partition in &partitions_default {
1665            for &node in partition {
1666                all_nodes.insert(node);
1667            }
1668        }
1669        assert_eq!(all_nodes.len(), 10);
1670    }
1671
1672    #[test]
1673    #[allow(non_snake_case)]
1674    fn test_distributed_graph_learning_communication() {
1675        let X = array![[1.0, 2.0], [2.0, 3.0], [3.0, 4.0], [4.0, 5.0]];
1676        let y = array![0, 1, -1, -1];
1677
1678        let dgl = DistributedGraphLearning::new()
1679            .n_workers(2)
1680            .communication_rounds(3)
1681            .overlap_ratio(0.3);
1682
1683        let partitions = dgl.partition_nodes(4);
1684        let X_sub1 = dgl.extract_subgraph(&X, &partitions[0]);
1685        let X_sub2 = dgl.extract_subgraph(&X, &partitions[1]);
1686        let y_sub1 = dgl.extract_sublabels(&y, &partitions[0]);
1687        let y_sub2 = dgl.extract_sublabels(&y, &partitions[1]);
1688
1689        let graph1 = dgl
1690            .learn_local_graph(&X_sub1, &y_sub1)
1691            .expect("operation should succeed");
1692        let graph2 = dgl
1693            .learn_local_graph(&X_sub2, &y_sub2)
1694            .expect("operation should succeed");
1695
1696        let local_graphs = vec![graph1, graph2];
1697        let updated_graphs = dgl.communicate_boundaries(&local_graphs, &partitions);
1698
1699        assert_eq!(updated_graphs.len(), 2);
1700        assert_eq!(updated_graphs[0].dim(), local_graphs[0].dim());
1701        assert_eq!(updated_graphs[1].dim(), local_graphs[1].dim());
1702    }
1703}