4

I currently have an array A which has dimensions N x t x t. I want to create a 2D matrix, N x t, of the form:

B = [ A[:,1,1] A[:,2,2],...,A[:,t,t]]

Obviously, 2 ways I could do this are writing it out in full (impractical since t is large) and loops (potentially slow). Is there a way to do this without loops. I figured it would work if I did:

B = A[:,[1:end],[1:end]]

but that just gave me back the original matrix.

Amro
  • 123,847
  • 25
  • 243
  • 454
Greg
  • 161
  • 6

3 Answers3

8

Basically, you need to start thinking about how to reorganize your matrix.

From 

A = randn([5 3 3]);

Let's look at 

A(:,:)

Basically you want columns 1, 5, 9. Thinking about it, knowing t = 3, from the present column you want to increment by t + 1. The formula is basically:

((1:3)-1)*(3+1)+1 %(or (0:2)*(3+1) + 1)

Which plugged in A yields your solution

A(:,((1:3)-1)*(3+1)+1)

In a general format, you can do:

A(:,((1:t)-1)*(t+1)+1)

EDIT:

Amro basically just put me to shame. The idea is still the same, it just becomes more readable thanks to end

Therefore use:

A(:,1:t+1:end)
Rasman
  • 5,349
  • 1
  • 25
  • 38
  • +1: yup, definitely the best way to go. Wouldn't win any prizes for readability, but there's comments for that :) – Rody Oldenhuis Jun 26 '13 at 13:53
  • @RodyOldenhuis Unfortunately readability and math seldom go hand in hand. – Rasman Jun 26 '13 at 14:06
  • @Rasman: all the credit goes to you :) I'm trying to find a mention of this in [the documentation](http://www.mathworks.com/help/matlab/math/multidimensional-arrays.html#f1-86846) but I cant.. – Amro Jun 26 '13 at 14:26
  • @Amro: Correct me if I'm misreading, but in fact the documentation states for this case: "Two subscripts always access the **first page** of a multidimensional array", implying that the behavior we see here for `(:, [...])` is a bug...? – Rody Oldenhuis Jun 26 '13 at 14:52
  • @RodyOldenhuis: weird right? fortunately Octave behaves the same way here, so it must be a feature not a bug – Amro Jun 26 '13 at 15:47
  • @Amro: well, a "useful bug" any way :) – Rody Oldenhuis Jun 26 '13 at 15:57
3

Here is my take on that subject:

mask = repmat(logical(permute(eye(size(A,2)), [3 1 2])), size(A,1), 1);
newA = reshape(A(mask), size(A,1), []);

Just generate mask, apply it and reshape result to the right shape.

Edit:

Or even simpler:

newA = A(:, logical(eye(size(A,2))))

or along the same lines,

newA = A(:, diag(true(size(A,2),1)));

which is a tad faster.

Rody Oldenhuis
  • 37,726
  • 7
  • 50
  • 96
j_kubik
  • 6,062
  • 1
  • 23
  • 42
  • +1: not bad (although I fear @Rasman cannot be beaten in this case :) – Rody Oldenhuis Jun 26 '13 at 14:48
  • @RodyOldenhuis I also think so, but with a bit of tinkering my version can get much neater too ;) See my edit. – j_kubik Jun 26 '13 at 15:44
  • +1 that second one is also very interesting, and I'm surprised it even works! Logical indexing usually requires a mask of the same size as the matrix: http://www.mathworks.com/help/matlab/math/matrix-indexing.html#bq7egb6-1 – Amro Jun 26 '13 at 15:49
  • @Amro: I have to agree. I've stumbled upon this notation before, used it a few times thought nothing of it (and then forgot all about it). Now that we're talking about it, it actually makes less and less sense... – Rody Oldenhuis Jun 26 '13 at 15:55
1

MATLAB is pretty good with vectors and matrices, but when it comes to "general arrays" you'll often have to switch to "general methods". When you're used to the ease with which you manipulate matrices/vectors, this will seem really clumsy, backward and not very convenient at all (there's actually very defendable reasons for that, but that's a discussion for another time :).

The version below loops through each page via arrayfun (which is only faster than a normal loop for large matrices) and calls diag on each page:

% "flip" the array, so that 3rd dimension becomes the 1st (rows), and 
% the 1st dimension becomes the 3rd, so that each "page" is a regular 
% matrix. You want the diagonals of all these matrices.
b = permute(a, [3 2 1]);

% Now "loop" through the pages, and collect the matrix diagonal for each page
c = arrayfun(@(ii)diag(b(:,:,ii)), 1:size(b,3), 'UniformOutput', false);

% The above will output a cell-array, which you can cast back to 
% a numeric matrix through this operation: 
A = [c{:}];
Rody Oldenhuis
  • 37,726
  • 7
  • 50
  • 96
  • Quite the contrast from "A(:,1:t+1:end)" ;) What do you mean general arrays though? A cell array? Because you can more or less index into cell arrays just like matrices I believe. – Shaun314 Jun 26 '13 at 14:46
  • @Shaun314: no, I meant an array that is not 1- or 2 dimensional (vector or matrix), but `N`-dimensional. – Rody Oldenhuis Jun 26 '13 at 14:50
  • Ah gotcha, makes a little bit more sense now :) – Shaun314 Jun 26 '13 at 14:52