Optimized the averaging step of K_means to use non-moving averages and arrays instead of hashmaps. The code also updates the new centroid in-place. Also fixed deprecation warning on all benchmarks in linfa-clustering (other benchmarks may have the same warnings). I'm curious as to the effect of my changes on benchmark speeds, since my machine isn't suited to running benchmarks.

Work towards https://github.com/rust-ml/linfa/issues/130

It would be possible to make a default implementation for any type that implements PredictRef. But we do not have specialization yet.

Does only implementing this, for types that can do this more efficiently, make sense?

Changes in the Fit trait

from this:

pub trait Fit<'a, R: Records, T> {
    type Object: 'a;

    fn fit(&self, dataset: &DatasetBase<R, T>) -> Self::Object;
}

to this:

pub trait Fit<R: Records, T, E: std::error::Error + std::convert::From<linfa::Error>> {
    type Object;

    fn fit(&self, dataset: &DatasetBase<R, T>) -> Result<Self::Object, E>;
}

by:

• removing 'a lifetime (left from previous svm implementation, not actually used by any algorithm anymore)
• forcing every implementation to return a result with an error struct (every implementation except PCA already returned an error, some implementations returned a String as the error but the transition keeps the same error messages)

Edit 1:

Added conversion from linfa error bound on fit error type. Every sub-crate should be able to handle the errors caused by using the base crate

Cross validation POC

Cross-validation can be defined by exploiting the new Fit definition. This is what it looks like for regression:

    use linfa::prelude::*;
    use linfa_elasticnet::{ElasticNet, Result};

    // load Diabetes dataset (mutable to allow fast k-folding)
    let mut dataset = linfa_datasets::diabetes();

    // prameters to compare
    let ratios = vec![0.1, 0.2, 0.5, 0.7, 1.0];

    // create a model for each parameter
    let models = ratios
        .iter()
        .map(|ratio| ElasticNet::params().penalty(0.3).l1_ratio(*ratio))
        .collect::<Vec<_>>();

    // get the mean r2 validation score across all folds for each model
    let r2_values =
        dataset.cross_validate(5, &models, |prediction, truth| prediction.r2(&truth))?;

    for (ratio, r2) in ratios.iter().zip(r2_values.iter()) {
        println!("L1 ratio: {}, r2 score: {}", ratio, r2);
    }

And this is what it looks like for classification:

    use linfa::prelude::*;
    use linfa_logistic::error::Result;
    use linfa_logistic::LogisticRegression;

    // Load dataset. Mutability is needed for fast cross validation
    let mut dataset =
        linfa_datasets::winequality().map_targets(|x| if *x > 6 { "good" } else { "bad" });

    // define a sequence of models to compare. In this case the
    // models will differ by the amount of l2 regularization
    let alphas = vec![0.1, 1., 10.];
    let models: Vec<_> = alphas
        .iter()
        .map(|alpha| {
            LogisticRegression::default()
                .alpha(*alpha)
                .max_iterations(150)
        })
        .collect();

    // use cross validation to compute the validation accuracy of each model. The
    // accuracy of each model will be averaged across the folds, 5 in this case
    let accuracies = dataset.cross_validate(5, &models, |prediction, truth| {
        Ok(prediction.confusion_matrix(truth)?.accuracy())
    })?;

    // display the accuracy of the models along with their regularization coefficient
    for (alpha, accuracy) in alphas.iter().zip(accuracies.iter()) {
        println!("Alpha: {}, accuracy: {} ", alpha, accuracy);
    }

Possible follow up

Redefine the Transformer trait to return a Result, like the t-sne implementation does, in order to avoid panicking as much as possible. It would be better to wait for #121 in order to avoid dealing with the unwrap() calls when performing float conversions.

Notes

• ~~Currently there are two versions of cross-validation: one for single target datasets and one for multi-target datasets. The main reasons for this are:~~
  • ~~Array1s of evaluation values can be constructed with collect, Array2 cannot and must be populated row by row (could be avoided by stacking, or maybe there is a dedicated ndarray method I am not aware of)~~
  • ~~Regression metrics behave differently when they are applied Array1-Array2 or Array2-Array1, making writing evaluation closures possibly more difficult than it should.~~

~~More than likely there is a solution to both problems, but I got stuck on it for too long and so I would consider it out of scope for this PR, unless someone has an easy solution to suggest. Forcing the user to give a single return value in the multi-target case (like a mean across the targets) could make the problem easier but the evaluation for each target would be lost in the process.~~

• Got this error in elasticnet just once:

     Running target/release/deps/linfa_elasticnet-d4c6cf99e56de5dd

running 11 tests
test algorithm::tests::coordinate_descent_lowers_objective ... ok
test algorithm::tests::elastic_net_2d_toy_example_works ... ok
test algorithm::tests::elastic_net_diabetes_1_works_like_sklearn ... ok
test algorithm::tests::diabetes_z_score ... ok
test algorithm::tests::elastic_net_penalty_works ... ok
test algorithm::tests::elastic_net_toy_example_works ... ok
test algorithm::tests::lasso_toy_example_works ... ok
test algorithm::tests::lasso_zero_works ... ok
error: test failed, to rerun pass '-p linfa-elasticnet --lib'

Caused by:
  process didn't exit successfully: `/home/ivano/Scrivania/github_projects/linfa/target/release/deps/linfa_elasticnet-d4c6cf99e56de5dd` (signal: 11, SIGSEGV: invalid memory reference)

I added an implementation of the Approximated DBSCAN to the linfa-clustering subcrate. This algorithm depends on the crate partitions for an implementation of the union-find structure. This can be removed by using adjacency lists, at the cost of performance. Another big thing: I noticed that the existing vanilla DBSCAN implementation used an axis convention that was the opposite to the one used by the generate_blobs function, which it uses in the benches. I modified the DBSCAN implementation to invert the axes in order to be able to compare benches of the two implementations, and I also modified the tests to pass this modification. I modified the vanilla implementation only to be able to compare the performances, I did not want to judge the previous implementation. This was also my first time using ndarray so any suggestion is very welcome.

I followed LinearRegression example and found err on cargo run

Code:

use linfa::prelude::SingleTargetRegression;
use linfa::traits::{Fit, Predict};
use linfa_linear::LinearRegression;


fn main() {
    let dataset = linfa_datasets::diabetes();
    let model = LinearRegression::default().fit(&dataset).unwrap();  # the error seems happen on this line
    pred = model.predict(&dataset);
    pred.r2(&dataset).unwrap();
}

Err output:

  = note: Undefined symbols for architecture x86_64:
            "_dgelsd_", referenced from:
                lapack::dgelsd::h78796e3523b23ad8 in liblax-5505ca61a5fbc67a.rlib(lax-5505ca61a5fbc67a.lax.d19e9304-cgu.2.rcgu.o)
            "_cblas_sdot", referenced from:
                ndarray::linalg::impl_linalg::_$LT$impl$u20$ndarray..ArrayBase$LT$S$C$ndarray..dimension..dim..Dim$LT$$u5b$usize$u3b$$u20$1$u5d$$GT$$GT$$GT$::dot_impl::h7adddea61e00d52a in rust_ds-fc4af6c90afa881c.4mtkwkgjvot8wzn2.rcgu.o
            "_cblas_ddot", referenced from:
                ndarray::linalg::impl_linalg::_$LT$impl$u20$ndarray..ArrayBase$LT$S$C$ndarray..dimension..dim..Dim$LT$$u5b$usize$u3b$$u20$1$u5d$$GT$$GT$$GT$::dot_impl::h7adddea61e00d52a in rust_ds-fc4af6c90afa881c.4mtkwkgjvot8wzn2.rcgu.o
          ld: symbol(s) not found for architecture x86_64
          clang: error: linker command failed with exit code 1 (use -v to see invocation)

OS: Mac Catalina v10.15.7 Rust: rustc 1.58.1

I have try cargo clean but still can't fix the issue

This PR adds K-means++ as an initialization function as well as K-means||, which is a version of K-means++ that scales better with higher cluster counts, described in this paper. A new hyperparameter has also been added to specify the choice of initialization algorithm. I want to have the design reviewed before going any further.

Todo:

• [x] Integrate the new initialization algorithms to the benchmarks, since they should allow K-means to complete faster.
• [x] Pick default initialization function.

The people at rust-cv looked recently at some clustering algorithms in linfa and one issue which came up are dependencies to the LAPACK routines. In some situations you want to avoid these, for example in an embedding context.

Take for example the gaussian mixture model: the current implementation utilizes a cholesky decomposition for faster inverse and determinant computation and speeds the LL computation up. We could provide a second implementation which gets activated when the corresponding feature flags are deactivated and only provide models with independent univariates gaussians.

Similar to what I did for KMeans, I reintroduced hyperparam builders to DBSCAN and approximate DBSCAN. This allows the verification of hyperparams to be done via returning errors rather than panicking. I believe using a builder is more idiomatic then doing a verify() check since the builder forces the verification to be done exactly once before the hyperparams can be used to run any algorithms. Introduce a new trait called UncheckedHyperParams to represent the verification of hyperparameters. Integrate this trait for DBSCAN, appx DBSCAN, and KMeans.

Here is a collection of work that needs to be done regarding documentation and minor issues with the code itself. Picking an item from this list can be a good way to start getting acquainted with the codebase and provide a useful contribution. 😄 👍🏻

In general, I would say that the pages for the individual algorithms in the linfa-clustering and linfa-elasticnet sub-crates can be good references for the structure of a documentation page.

It's suggested to look at the existing documentation in your local build rather than on doc.rs since some pages may have already been updated.

• [x] Update linfa-bayes's sub-crate documentation to include a brief description of the algorithm used and add an example of the usage of the provided model along with maybe some hints regarding when to choose a naive Bayes predictor. There are already some examples in the dedicated page for the params structure but maybe it's better to add at least one to the main page for the sub-crate
• [x] linfa-hierarchical: The documentation for this sub-crate needs improvements like the ones listed for linfa-bayes
• [x] linfa-ica: Also needs updates similar to linfa-bayes
• [x] linfa-kernel: Add a description of what kernel methods are useful for, similarly to what's written in the linfa-svm crate, and add descriptions to the Kernel and KernelView subtypes. There are some examples in the KernelBase struct Documentation but maybe it would be a good idea to have an example of kernel transformation directly on the crate's main page.
• [x] linfa-linear: The crate's documentation is inside the ols module instead of being in the crate's root and it is outdated since it says it only provides an implementation of the least squares algorithm
• [x] linfa-linear: glm::TweedieRegressor provides its own fit method instead of implementing the Fit trait, and the same for predict. Moving the methods inside the trait would be a good way to align it with the rest of the algorithms provided and ensure compatibility
• [x] linfa-logistic could use some more documentation to explain the algorithm and some examples in its main page
• [x] linfa-logistic implements its own fit and predict methods instead of implementing the Fit and Predict traits. Like in the case of linear regression it would be a good thing to align the interface with the other sub-crates
• [x] linfa-reduction completely lacks any documentation. An explanation of why dimensionality reduction is useful in ML and descriptions of the single algorithms would really help with understanding the usefulness of the crate
• [x] linfa-svm: this is another reference for what other crates' documentation should look like. Right now the predict method returns an Array1 for classification and a Vec for regression. It would be less confusing if the regression case was modified to return an Array1 too.
• [x] linfa-trees: it needs changes similar to linfa-bayes, with the addition of documentation regarding the methods for the params structure
• [x] linfa: The main page of the rustdocs still says that linfa only provides the K-Means algorithm which may be very confusing to someone that only sees the crate in doc.rs for the first time. Here a completely revised page compete with project goals and links to the various sub-crates would be useful, so that one can find the algorithm they need without manually searching for the sub-crates.
• [x] linfa: The dataset module page could use a bit of explanation about the main differences between the four dataset types and a brief recollection of what utility methods a dataset provides
• [x] ~~linfa: The metrics module page could use a brief list of the provided metrics so that one does not have to go looking for them in the sup-pages if they only want to know whether a metric is provided or not~~

This PR updates the criterion settings for our benched algorithms. Over the course of a day or a few days I will update this description with the timing stats from before and after changes for each algorithm. Unless, it is okay to just do it for a few then I'll do that instead

linfa-ica

Old - 1min 39secs New - 3min 38secs

linfa-pls

Old - 4mins 16secs New - 9mins 42secs

linfa-linear

Old - 1min 29secs New - 3min 37secs

linfa-nn

Old - 3mins 43secs New - 8mins 57secs

linfa-clustering

k-means:

Old - 4min 6secs New - 7mins 34secs

gaussian_mixture

Old - 51secs New - 2mins 5secs

dbscan

Old - 1min 14secs New - 2mins 50secs

appx_db_scan

Old - 47secs New - 1min 27secs

linfa-ftrl

Old - 1min 43secs New - 3min 57secs

linfa-trees

Old - 55secs New - 2min 38secs

Current state:

• [x] port algorithm's parameters to new syntax
• [ ] improve documentation (use linfa-bayes as prototype)
• [x] fix issues with Transform and double results

Preprocessing and transformation with Data frames are heavily used for ml model training in scikitlearn. The two most popular DataFrame libraries (Polars, DataFusion) written in rust are based on apache arrow in-memory data format but not based on ndarray. Which also looks like will be the trend for any new data frame players in Rust. It does not look like there will be a data frame that wraps ndarray under the hood, the way pandas wrap numpy.

By adding arrow support in linfa, any data frame based on arrow will have default support which means any arrow-based data frame can be passed to any preprocessing or training modules of linfa. Without dealing with ndarray. The way pandas can be passed to sci-kit-learn. Hence I propose to have direct arrow support in linfa to have it a more generalized framework. By doing that Polars/DataFusion users can already use rust for ml training out of the box.

According to benchmark results from this comment, linfa-pls is slightly slower without BLAS. Specifically, the Regression-Nipals benchmarks are slightly slower when the sample size is 100000, and the Regression-Svd and Canonical-Svd benchmarks are slower when the sample size is 100000.

According to the benchmark results from this comment, linfa-linear is faster with BLAS than without. The OLS algorithm isn't too different with 5 features and is only slightly slower with 10 features, but GLM is significantly slower without BLAS. We should profile and investigate the difference between the BLAS and non-BLAS performance.

According to these results, all ICA benchmarks are noticibly faster with BLAS than without, though this is less severe at higher sample sizes. We should profile the non-BLAS benchmark runs (using something like flamegraph) and see if we can fix the discrepancy.