Non-continguous matrix slices

[kernelmachine]

@bluss we discussed this on the IRC, but thought I'd record the desired functionality here.

For a project I'm working on, I've been in need of non-contiguous row/column slices of ndarray matrices.

So far, I've been performing the operation via a copy -- stacking row/column views. Note that because stack requires each to-be-stacked array to be 2-dimensional, I use into_shape on each row view in the closure.

/// Rectangular matrix
pub type Mat<A> = OwnedArray<A, (Ix, Ix)>;

pub fn noncontig_2d_row_slice(mat : &Mat<f64>, indices : &Vec<usize>)  -> Mat<f64> {
    let mat = indices.iter()
                        .map(|&x| mat.row(x).into_shape((1,mat.shape()[1]))
                        .ok().expect("indexing Error")).collect::<Vec<_>>();
    stack(Axis(0), mat.as_slice()).ok().expect("stacking error")
}

For example:

fn main(){
    let x = arr2(&[[0.0, 1.0], [1.0,0.0],[1.0,0.0],[1.0,0.0],[1.0,0.0],[0.0, 1.0],[0.0, 1.0]]);
    let s = noncontig_2d_slice(&x,&vec![1,3,5]);
    let target = arr2(&[[1.0,0.0],[1.0,0.0],[0.0, 1.0]]);
    assert!(s.all_close(&target,1e-8))
}

Is this functionality worth integrating into the codebase?

[bluss]

Yes it is worth integrating. It's something that was wished for in #178 too. We can't call it exactly slice since it's not a view type anymore (this is why I said it needs a new type of array implementation, I thought you wanted a non-contiguous view type).

There's another missing feature here, being able to stack multiple 1-d arrays into a 2-d array. That's a function we want either way.

[kernelmachine]

Cool, I can submit a PR for stacking multiple 1-d arrays into a 2-d array, and extracting submatrices via this stacking.

[bluss]

My aim with such operations is always to make them as general as possible, so we want to support the general n-dimensional case.

Btw, I can't help but want to tell you to write float literals as 1. and 0., you save so much typing by omitting those zeros.. 😄

[kernelmachine]

sounds good! do you believe that this should be a separate function from the current stack, or can we make stack generic enough to support both 2d/1d view stacking?

and yeah i'll follow that style :)

[bluss]

I don't think stack can do it-- the number of axes is encoded in the types. Right now it always returns an array of the same number of axes as the input array views.

I see how at least the select along one axis can be implemented in a simple way

[sgururangan-r7]

Yeah, I guess I could write a new stacking function that takes an array with dim D::Smaller, and have it output an array with dim D, where D: RemoveAxis. However, unless I'm missing something, it doesn't seem like the D::Smaller -> D transformation is supported. Do we need another impl on the RemoveAxis trait to do this?

I'm thinking the generic version of my implementation above can go something like the following:

pub fn noncontig_subview(&mut self, axis : Axis,  indices : &Vec<usize>)  
-> Result<OwnedArray<A, D>, ShapeError>
            where A: Copy,
                  D: RemoveAxis
        {
        let v = indices.iter().map(|&x| self.subview(axis,x)).collect::<Vec<_>>();
        new_stack(axis, v.as_slice())
    }

pub fn new_stack<'a, A, D>(axis: Axis, arrays: &[ArrayView<'a, A, D::Smaller>])
    -> Result<OwnedArray<A, D>, ShapeError>
    where A: Copy,
               D: RemoveAxis,
{
.
.
.
}

This isn't a fully generic stacking function b/c you still can't stack arrays that are of fewer dimensions than D::Smaller to get an array of dim D. But I'm not sure if there's an obvious use case for that, so maybe it doesn't matter.

[bluss]

@sgururangan-r7 You're completely right. For this particular thing we can actually use the old stack and isubview, which does the same thing as subview, except it doesn't remove the collapsed axis.

Since this kind of method only has indexing-related errors, I guess it's ok that it uses panic for errors, like other methods using indexes.

This seems to work for that:

    pub fn select_along(&self, axis: Axis, indices: &[Ix]) -> OwnedArray<A, D>
        where A: Copy,
              D: RemoveAxis,
    {
        let v = self.view();
        let mut subs = vec![v; indices.len()];
        for (&i, sub) in zipsl(indices, &mut subs[..]) {
            sub.isubview(axis, i);
        }
        if subs.is_empty() {
            let mut dim = self.dim();
            dim.set_axis(axis, 0);
            unsafe {
                OwnedArray::from_vec_dim_unchecked(dim, vec![])
            }
        } else {
            stack(axis, &subs).unwrap()
        }
    }

[kernelmachine]

Cool, isubview is the trick! Submitted PR around this.