 Hello,

I’ve implemented matrix centering and scaling (as in the R function `scale`) as follow:

``````    let column_means = mat.column_mean();

let centered_mat = mat.map_with_location(|row, col, n| {
n - column_means.get(row).unwrap()
});

let variance = mat.map(|n| n * n).sum() / mat.len() as f64;
let stddev = variance.sqrt();
let scaled_mat = centered_mat.map(|n| n / stddev);
``````

Now I would like to implement both `center(matrix) -> matrix` and a `scale(matrix) -> matrix` generic over any matrix.

So far I was not able to find which type my functions should be generic over, any help on how to define those would be really appreciated!

Thanks!

Hi!

Here is what this function can look like:

``````use nalgebra::{Matrix, MatrixMN, RealField, Dim, DefaultAllocator};
use nalgebra::allocator::Allocator;
use nalgebra::storage::Storage;

fn center_scale<N, R, C, S>(mat: &Matrix<N, R, C, S>) -> MatrixMN<N, R, C>
where N: RealField,
R: Dim,
C: Dim,
S: Storage<N, R, C>,
DefaultAllocator: Allocator<N, R, C> + Allocator<N, R> {
let column_means = mat.column_mean();

let centered_mat = mat.map_with_location(|row, col, n| {
n - column_means[row]
});

let variance = mat.map(|n| n * n).sum() / nalgebra::convert(mat.len() as f64);
let stddev = variance.sqrt();
centered_mat.map(|n| n / stddev)
}
``````
• The type parameter `N` means that the matrix is generic wrt. the scalar type. Here it can be any float type (`f32` or `f64`).
• The type parameter `R` means that the matrix is generic wrt. the number of rows.
• The type parameter `C` means that the matrix is generic wrt. the number of columns.
• The type parameter `S` means that the matrix is generic wrt. its internal storage type. This allows this function to be used on both matrixces that own their data as well as matrix slices.
• The return type is `MatrixMN<N, R, C>` does not include a parameter for its storage type. This is because its storage type must be deduced automatically by the `DefaultAllocator`.
• The `DefaultAllocator: Allocator<N, R, C>` is here to deduce the storage type of the result matrix.
• The `DefaultAllocator: Allocator<N, R>` is required by `mat.column_mean()`.
• I replaced `column_means.get(row).unwrap()` by direct indexing `column_means[row]` because you unwrap it anyways.
• I replaced `mat.len() as f64` by `nalgebra::convert(mat.len() as f64)` which will convert `usize -> f64 -> N`.

For even more genericity, you can do this:

``````use nalgebra::{Matrix, MatrixMN, SimdComplexField, Dim, DefaultAllocator};
use nalgebra::allocator::Allocator;
use nalgebra::storage::Storage;

fn center_scale<N, R, C, S>(mat: &Matrix<N, R, C, S>) -> MatrixMN<N, R, C>
where N: SimdComplexField,
R: Dim,
C: Dim,
S: Storage<N, R, C>,
DefaultAllocator: Allocator<N, R, C> + Allocator<N, R> {
let column_means = mat.column_mean();

let centered_mat = mat.map_with_location(|row, _, n| {
n - column_means[row]
});

let zero: N::SimdRealField = nalgebra::zero();
let variance: N::SimdRealField = mat.fold(zero, |acc, n| acc + n.simd_modulus_squared()) / nalgebra::convert(mat.len() as f64);
let stddev = variance.simd_sqrt();
centered_mat.unscale(stddev)
}
``````

This will allow your function to work on matrices containing complex numbers too, as well as SIMD types (useful if you want to center four matrices simultaneously for example). Note that I replaced your `mat.map(...).sum()` by a `mat.fold`. You can do this on the first piece of code I mentioned too and it will avoid the allocation made by the `map`.

1 Like

Thank you for the quick and outstanding quality response!

This is very instructive and exactly what I was looking for, maybe this kind of tutorial could be documented so it can be easier to find.

What I failed to grasp was that I could define the return type without specifying the allocator explicitly.