Add Barnes-Hut implementation of t-SNE (#101)

* Initial commit

* Add generic type definition to wrapper

* Add example

* Improve error handling and documentation

* Add example two website and run rustfmt

* Fix example in documentation

* Make t-SNE tests deterministic random

* Address review

Author: bytesnake

Cargo.toml

@@ -76,6 +76,7 @@ members = [
     "algorithms/linfa-bayes",
     "algorithms/linfa-elasticnet",
     "algorithms/linfa-pls",
+    "algorithms/linfa-tsne",
     "datasets",
 ]
 

README.md

@@ -35,13 +35,14 @@ Where does `linfa` stand right now? [Are we learning yet?](http://www.arewelearn
 | [linear](algorithms/linfa-linear/) | Linear regression | Tested  | Partial fit | Contains Ordinary Least Squares (OLS), Generalized Linear Models (GLM) | 
 | [elasticnet](algorithms/linfa-elasticnet/) | Elastic Net | Tested | Supervised learning | Linear regression with elastic net constraints |
 | [logistic](algorithms/linfa-logistic/) | Logistic regression | Tested  | Partial fit | Builds two-class logistic regression models
-| [reduction](algorithms/linfa-reduction/) | Dimensionality reduction | Tested  | Pre-processing | Diffusion mapping and Principal Component Analysis (PCA) |
+| [reduction](algorithms/linfa-reduction/) | Dimensionality reduction | Tested | Pre-processing | Diffusion mapping and Principal Component Analysis (PCA) |
 | [trees](algorithms/linfa-trees/) | Decision trees | Experimental  | Supervised learning | Linear decision trees
 | [svm](algorithms/linfa-svm/) | Support Vector Machines | Tested  | Supervised learning | Classification or regression analysis of labeled datasets | 
 | [hierarchical](algorithms/linfa-hierarchical/) | Agglomerative hierarchical clustering | Tested | Unsupervised learning | Cluster and build hierarchy of clusters |
 | [bayes](algorithms/linfa-bayes/) | Naive Bayes | Tested | Supervised learning | Contains Gaussian Naive Bayes |
 | [ica](algorithms/linfa-ica/) | Independent component analysis | Tested | Unsupervised learning | Contains FastICA implementation |
 | [pls](algorithms/linfa-pls/) | Partial Least Squares | Tested | Supervised learning | Contains PLS estimators for dimensionality reduction and regression |
+| [tsne](algorithms/linfa-tsne/) | Dimensionality reduction| Tested | Unsupervised learning | Contains exact solution and Barnes-Hut approximation t-SNE |
 
 We believe that only a significant community effort can nurture, build, and sustain a machine learning ecosystem in Rust - there is no other way forward.
 

algorithms/linfa-logistic/src/lib.rs

@@ -384,12 +384,12 @@ fn convert_params<F: Float>(n_features: usize, w: &Array1<F>) -> (Array1<F>, F)
     } else if n_features + 1 == w.len() {
         (w.slice(s![..w.len() - 1]).to_owned(), w[w.len() - 1])
     } else {
-        panic!(format!(
+        panic!(
             "Unexpected length of parameter vector `w`, exected {} or {}, found {}",
             n_features,
             n_features + 1,
             w.len()
-        ));
+        );
     }
 }
 

algorithms/linfa-reduction/src/pca.rs

@@ -28,7 +28,7 @@ use serde_crate::{Deserialize, Serialize};
 
 use linfa::{
     dataset::Records,
-    traits::{Fit, PredictRef},
+    traits::{Fit, PredictRef, Transformer},
     DatasetBase, Float,
 };
 
@@ -173,6 +173,22 @@ impl<F: Float, D: Data<Elem = F>> PredictRef<ArrayBase<D, Ix2>, Array2<F>> for P
     }
 }
 
+impl<F: Float, D: Data<Elem = F>, T>
+    Transformer<DatasetBase<ArrayBase<D, Ix2>, T>, DatasetBase<Array2<F>, T>> for Pca<F>
+{
+    fn transform(&self, ds: DatasetBase<ArrayBase<D, Ix2>, T>) -> DatasetBase<Array2<F>, T> {
+        let DatasetBase {
+            records,
+            targets,
+            weights,
+            ..
+        } = ds;
+
+        let new_records = self.predict_ref(&records);
+
+        DatasetBase::new(new_records, targets).with_weights(weights)
+    }
+}
 #[cfg(test)]
 mod tests {
     use super::*;

algorithms/linfa-tsne/Cargo.toml

@@ -0,0 +1,29 @@
+[package]
+name = "linfa-tsne"
+version = "0.3.1"
+authors = ["Lorenz Schmidt <[email protected]>"]
+edition = "2018"
+
+description = "Barnes-Hut t-distributed stochastic neighbor embedding"
+license = "MIT/Apache-2.0"
+
+repository = "https://github.com/rust-ml/linfa"
+readme = "README.md"
+
+keywords = ["tsne", "data visualization", "clustering", "machine-learning", "linfa"]
+categories = ["algorithms", "mathematics", "science"]
+
+[dependencies]
+thiserror = "1"
+ndarray = { version = "0.13", default-features = false }
+ndarray-rand = "0.11"
+bhtsne = "0.4.0"
+
+linfa = { version = "0.3.1", path = "../.." }
+
+[dev-dependencies]
+rand = "0.7"
+approx = "0.3"
+
+linfa-datasets = { version = "0.3.1", path = "../../datasets", features = ["iris"] }
+linfa-reduction = { version = "0.3.1", path = "../linfa-reduction" }

algorithms/linfa-tsne/README.md

@@ -0,0 +1,32 @@
+# t-SNE
+
+`linfa-tsne` provides a pure Rust implementation of exact and Barnes-Hut t-SNE.
+
+## The Big Picture
+
+`linfa-tsne` is a crate in the [`linfa`](https://crates.io/crates/linfa) ecosystem, an effort to create a toolkit for classical Machine Learning implemented in pure Rust, akin to Python's `scikit-learn`.
+
+## Current state
+
+`linfa-tsne` currently provides an implementation of the following methods: 
+
+- exact solution t-SNE
+- Barnes-Hut t-SNE
+
+It wraps the [bhtsne](https://github.com/frjnn/bhtsne) crate, all kudos to them.
+
+## Examples
+
+There is an usage example in the `examples/` directory. The example uses a BLAS backend, to run it and use the `intel-mkl` library do:
+
+```bash
+$ cargo run --example tsne --features linfa/intel-mkl-system
+```
+
+You have to install the `gnuplot` library for plotting. Also take a look at the [README](https://github.com/rust-ml/linfa#blaslapack-backend) to see what BLAS/LAPACK backends are possible.
+
+## License
+Dual-licensed to be compatible with the Rust project.
+
+Licensed under the Apache License, Version 2.0 <http://www.apache.org/licenses/LICENSE-2.0> or the MIT license <http://opensource.org/licenses/MIT>, at your option. This file may not be copied, modified, or distributed except according to those terms.
+

algorithms/linfa-tsne/examples/iris_plot.plt

@@ -0,0 +1,7 @@
+set style increment user
+set style line 1 lc rgb 'red'
+set style line 2 lc rgb 'blue'
+set style line 3 lc rgb 'green'
+
+set style data points
+plot 'iris.dat' using 1:2:3 linecolor variable pt 7 ps 2 t ''

algorithms/linfa-tsne/examples/tsne.rs

@@ -0,0 +1,31 @@
+use linfa::traits::{Fit, Transformer};
+use linfa_reduction::Pca;
+use linfa_tsne::{Result, TSne};
+use std::{io::Write, process::Command};
+
+fn main() -> Result<()> {
+    let ds = linfa_datasets::iris();
+    let ds = Pca::params(3).whiten(true).fit(&ds).transform(ds);
+
+    let ds = TSne::embedding_size(2)
+        .perplexity(10.0)
+        .approx_threshold(0.1)
+        .transform(ds)?;
+
+    let mut f = std::fs::File::create("examples/iris.dat").unwrap();
+
+    for (x, y) in ds.sample_iter() {
+        f.write(format!("{} {} {}\n", x[0], x[1], y[0]).as_bytes())
+            .unwrap();
+    }
+
+    Command::new("gnuplot")
+        .arg("-p")
+        .arg("examples/iris_plot.plt")
+        .spawn()
+        .expect(
+            "Failed to launch gnuplot. Pleasure ensure that gnuplot is installed and on the $PATH.",
+        );
+
+    Ok(())
+}

algorithms/linfa-tsne/src/error.rs

@@ -0,0 +1,21 @@
+use thiserror::Error;
+
+pub type Result<T> = std::result::Result<T, TSneError>;
+
+#[derive(Error, Debug)]
+pub enum TSneError {
+    #[error("negative perplexity")]
+    NegativePerplexity,
+    #[error("perplexity too large for number of samples")]
+    PerplexityTooLarge,
+    #[error("negative approximation threshold")]
+    NegativeApproximationThreshold,
+    #[error("embedding size larger than original dimensionality")]
+    EmbeddingSizeTooLarge,
+    #[error("number of preliminary iterations larger than total iterations")]
+    PreliminaryIterationsTooLarge,
+    #[error("invalid shaped array {0}")]
+    InvalidShape(#[from] ndarray::ShapeError),
+    #[error(transparent)]
+    BaseCrate(#[from] linfa::Error),
+}