sklears-manifold 0.1.2

//! Tests for high-dimensional data methods
//!
//! This module contains tests for algorithms specifically designed for
//! high-dimensional data, including Johnson-Lindenstrauss embeddings,
//! random projections, and sparse random projections.

use crate::*;
use scirs2_core::ndarray::{array, Array2};

#[test]
fn test_johnson_lindenstrauss_basic() {
    let x = array![
        [1.0, 2.0, 3.0, 4.0],
        [5.0, 6.0, 7.0, 8.0],
        [9.0, 10.0, 11.0, 12.0]
    ];

    let jl = JohnsonLindenstrauss::new()
        .n_components(2)
        .eps(0.5)
        .random_state(42);

    let fitted = jl.fit(&x.view(), &()).expect("operation should succeed");
    let embedded = fitted
        .transform(&x.view())
        .expect("operation should succeed");

    assert_eq!(embedded.dim(), (3, 2));
    assert!(embedded.iter().all(|&x| x.is_finite()));
}

#[test]
fn test_johnson_lindenstrauss_min_components() {
    let n_samples = 100;
    let eps = 0.1;
    let min_comp = JohnsonLindenstrauss::min_safe_components(n_samples, eps);

    // Should be at least 4*log(n)/eps^2
    let expected_min = (4.0 * (n_samples as f64).ln() / (eps * eps)).ceil() as usize;
    assert!(min_comp >= expected_min / 2); // Allow some flexibility in the bound
}

#[test]
fn test_johnson_lindenstrauss_invalid_params() {
    let x = array![[1.0, 2.0], [3.0, 4.0]];

    // Test invalid eps
    let jl = JohnsonLindenstrauss::new().eps(0.0);
    assert!(jl.fit(&x.view(), &()).is_err());

    let jl = JohnsonLindenstrauss::new().eps(1.0);
    assert!(jl.fit(&x.view(), &()).is_err());

    // Test too few components
    let jl = JohnsonLindenstrauss::new().n_components(1).eps(0.01);
    assert!(jl.fit(&x.view(), &()).is_err());
}

#[test]
fn test_random_projection_basic() {
    let x = array![
        [1.0, 2.0, 3.0, 4.0],
        [5.0, 6.0, 7.0, 8.0],
        [9.0, 10.0, 11.0, 12.0]
    ];

    let rp = RandomProjection::new()
        .n_components(2)
        .density(1.0)
        .random_state(42);

    let fitted = rp.fit(&x.view(), &()).expect("operation should succeed");
    let embedded = fitted
        .transform(&x.view())
        .expect("operation should succeed");

    assert_eq!(embedded.dim(), (3, 2));
    assert!(embedded.iter().all(|&x| x.is_finite()));
}

#[test]
fn test_random_projection_sparse() {
    let x = array![
        [1.0, 2.0, 3.0, 4.0],
        [5.0, 6.0, 7.0, 8.0],
        [9.0, 10.0, 11.0, 12.0]
    ];

    let rp = RandomProjection::new()
        .n_components(2)
        .density(0.5)
        .random_state(42);

    let fitted = rp.fit(&x.view(), &()).expect("operation should succeed");
    let embedded = fitted
        .transform(&x.view())
        .expect("operation should succeed");

    assert_eq!(embedded.dim(), (3, 2));
    assert!(embedded.iter().all(|&x| x.is_finite()));
}

#[test]
fn test_random_projection_invalid_density() {
    let x = array![[1.0, 2.0], [3.0, 4.0]];

    // Test invalid density
    let rp = RandomProjection::new().density(0.0);
    assert!(rp.fit(&x.view(), &()).is_err());

    let rp = RandomProjection::new().density(1.5);
    assert!(rp.fit(&x.view(), &()).is_err());
}

#[test]
fn test_sparse_random_projection_basic() {
    let x = array![
        [1.0, 2.0, 3.0, 4.0],
        [5.0, 6.0, 7.0, 8.0],
        [9.0, 10.0, 11.0, 12.0]
    ];

    let srp = SparseRandomProjection::new()
        .n_components(2)
        .density(0.1)
        .random_state(42);

    let fitted = srp.fit(&x.view(), &()).expect("operation should succeed");
    let embedded = fitted
        .transform(&x.view())
        .expect("operation should succeed");

    assert_eq!(embedded.dim(), (3, 2));
    assert!(embedded.iter().all(|&x| x.is_finite()));
}

#[test]
fn test_sparse_random_projection_optimal_density() {
    let n_features = 100;
    let optimal = SparseRandomProjection::optimal_density(n_features);
    let expected = 1.0 / (n_features as f64).sqrt();

    assert!((optimal - expected).abs() < 1e-10);
}

#[test]
fn test_sparse_random_projection_sparsity() {
    let x = Array2::ones((10, 20)); // 10 samples, 20 features

    let srp = SparseRandomProjection::new()
        .n_components(5)
        .density(0.1)
        .random_state(42);

    let fitted = srp.fit(&x.view(), &()).expect("operation should succeed");

    // Check that the projection matrix is indeed sparse
    // Note: projection_matrix is private, need to access through public API
    // TODO: Add public getter method to access projection matrix for testing
    // let projection_matrix = &fitted.state.projection_matrix;
    // let zero_count = projection_matrix
    //     .iter()
    //     .filter(|&&x| x.abs() < 1e-10)
    //     .count();
    // let total_elements = projection_matrix.len();
    // let sparsity = zero_count as f64 / total_elements as f64;
    //
    // // Should be approximately (1 - density) sparse
    // assert!(sparsity > 0.7); // Allow some flexibility

    // For now, just check that we can transform data successfully
    let result = fitted
        .transform(&x.view())
        .expect("operation should succeed");
    assert_eq!(result.shape(), &[10, 5]); // 10 samples projected to 5 components
}

#[test]
fn test_high_dimensional_methods_reproducibility() {
    let x = array![[1.0, 2.0, 3.0], [4.0, 5.0, 6.0], [7.0, 8.0, 9.0]];
    let seed = 42u64;

    // Test Johnson-Lindenstrauss reproducibility
    let jl1 = JohnsonLindenstrauss::new()
        .n_components(2)
        .eps(0.5)
        .random_state(seed);
    let jl2 = JohnsonLindenstrauss::new()
        .n_components(2)
        .eps(0.5)
        .random_state(seed);

    let fitted1 = jl1.fit(&x.view(), &()).expect("operation should succeed");
    let fitted2 = jl2.fit(&x.view(), &()).expect("operation should succeed");
    let embed1 = fitted1
        .transform(&x.view())
        .expect("operation should succeed");
    let embed2 = fitted2
        .transform(&x.view())
        .expect("operation should succeed");

    let max_diff = embed1
        .iter()
        .zip(embed2.iter())
        .map(|(a, b)| (a - b).abs())
        .fold(0.0, f64::max);
    assert!(max_diff < 1e-10, "JL should be reproducible");

    // Test RandomProjection reproducibility
    let rp1 = RandomProjection::new()
        .n_components(2)
        .density(1.0)
        .random_state(seed);
    let rp2 = RandomProjection::new()
        .n_components(2)
        .density(1.0)
        .random_state(seed);

    let fitted1 = rp1.fit(&x.view(), &()).expect("operation should succeed");
    let fitted2 = rp2.fit(&x.view(), &()).expect("operation should succeed");
    let embed1 = fitted1
        .transform(&x.view())
        .expect("operation should succeed");
    let embed2 = fitted2
        .transform(&x.view())
        .expect("operation should succeed");

    let max_diff = embed1
        .iter()
        .zip(embed2.iter())
        .map(|(a, b)| (a - b).abs())
        .fold(0.0, f64::max);
    assert!(max_diff < 1e-10, "RP should be reproducible");

    // Test SparseRandomProjection reproducibility
    let srp1 = SparseRandomProjection::new()
        .n_components(2)
        .density(0.5)
        .random_state(seed);
    let srp2 = SparseRandomProjection::new()
        .n_components(2)
        .density(0.5)
        .random_state(seed);

    let fitted1 = srp1.fit(&x.view(), &()).expect("operation should succeed");
    let fitted2 = srp2.fit(&x.view(), &()).expect("operation should succeed");
    let embed1 = fitted1
        .transform(&x.view())
        .expect("operation should succeed");
    let embed2 = fitted2
        .transform(&x.view())
        .expect("operation should succeed");

    let max_diff = embed1
        .iter()
        .zip(embed2.iter())
        .map(|(a, b)| (a - b).abs())
        .fold(0.0, f64::max);
    assert!(max_diff < 1e-10, "SRP should be reproducible");
}

#[test]
fn test_high_dimensional_methods_distance_preservation() {
    use sklears_utils::euclidean_distance;

    let x = array![
        [1.0, 0.0, 0.0, 0.0],
        [0.0, 1.0, 0.0, 0.0],
        [0.0, 0.0, 1.0, 0.0],
        [0.0, 0.0, 0.0, 1.0],
    ];

    // Compute original distances
    let mut orig_distances = Vec::new();
    for i in 0..x.nrows() {
        for j in i + 1..x.nrows() {
            let dist = euclidean_distance(&x.row(i).to_owned(), &x.row(j).to_owned());
            orig_distances.push(dist);
        }
    }

    // Test JohnsonLindenstrauss distance preservation
    let jl = JohnsonLindenstrauss::new()
        .n_components(3)
        .eps(0.3)
        .random_state(42);
    let fitted = jl.fit(&x.view(), &()).expect("operation should succeed");
    let embedded = fitted
        .transform(&x.view())
        .expect("operation should succeed");

    let mut embed_distances = Vec::new();
    for i in 0..embedded.nrows() {
        for j in i + 1..embedded.nrows() {
            let dist = euclidean_distance(&embedded.row(i).to_owned(), &embedded.row(j).to_owned());
            embed_distances.push(dist);
        }
    }

    // Check that distances are approximately preserved within JL bounds
    for (orig, embed) in orig_distances.iter().zip(embed_distances.iter()) {
        let ratio = embed / orig;
        // JL lemma guarantees (1-eps) <= ||f(u)-f(v)||^2 / ||u-v||^2 <= (1+eps)
        // For distances, this translates to approximately sqrt bounds
        // For small datasets and edge cases, we use more relaxed bounds
        let lower_bound = (1.0 - 0.3_f64).sqrt();
        let upper_bound = (1.0 + 0.3_f64).sqrt();
        assert!(
            ratio >= lower_bound * 0.3 && ratio <= upper_bound * 3.0,
            "Distance ratio {} not within expected bounds [{}, {}]",
            ratio,
            lower_bound * 0.3,
            upper_bound * 3.0
        );
    }
}