Matlab Indexing and Repmat Replacement

Quick Reference

Repmat with One-Dimensional Arrays

x = [1:5]'

repmat(x, [1 3])

x(:, ones(3, 1))

repmat(x, [3 1])

idx = [1 : size(x,1)]';
x(idx(:, ones(3, 1)), :)

x([1:size(x,1)]' * ones(1,3), :)

Repmat with Multi-Dimensional Arrays

A = zeros(3, 4);
A(:) = 1:12;

repmat(A, [2 3])

rowIdx = [1 : size(A,1)]';
colIdx = [1 : size(A,2)]';

A(rowIdx(:, ones(2,1)), colIdx(:, ones(3,1)))

A([1:size(A,1)]' * ones(1,2),
  [1:size(A,2)]' * ones(1,3))

repmat(A, [1 1 3])

A(:, :, ones(3, 1))

Indexing Arrays

Indexing into Matlab arrays is very often glossed over, but it is actually quite a powerful (and fast) technique.  There are a surprising number of things that can be done by just indexing into a Matlab array.  These indexing tricks come in handy to avoid "for" loops and in vectorizing code.

Consider the following column vector:

x = [1:5]'

1
2
3
4
5

The following lines of code give the same result and they simply list all the elements of the array:

x(:)
x(1:end)
x([1 2 3 4 5])

1
2
3
4
5

Instead of indexing all the elements in the array in order, you can index into the array in any arbitrary way, here are a couple examples:

% reverse the array
x(end:-1:1)

5
4
3
2
1

% every second element
x(1:2:end)

1
3
5

% specific elements in the array
x([1 1 4 2 5 3 5])

1
1
4
2
5
3
5

% random elements from the array
x(floor(rand(4, 1) * 5) + 1)

3
1
5
3

Indexing can also be done on "singleton" dimensions (dimensions with size 1).  For example, x is a column vector and its dimensions are [5 1], which means that its second dimension is a singleton.

Consequently, the following two statements are equivalent.  In the second statement, all the rows of x are displayed in addition to the first (and only) column of x.

x(:)
x(:, 1)

1
2
3
4
5

Interesting results are obtained when indexing the singleton dimension multiple times.  The following two lines of code are equivalent and they index all the rows of x and show the first column 3 times.

x(:, [1 1 1])
x(:, ones(3,1))

1   1   1
2   2   2
3   3   3
4   4   4
5   5   5

Note that this is the same as

repmat(x, [1 3])

 

Consider the following array:

A = zeros(3, 4);
A(:) = 1:12;

A

1  4  7  10
2  5  8  11
3  6  9  12

Indexing into the array in different ways is quite powerful:

% reverse the rows and columns
A(end:-1:1, end:-1:1)

12  9  6  3
11  8  5  2
10  7  4  1

% get only the 1st and 3rd rows
% and the 2nd and 4th columns
A([1 3], [2 4])

4  10
6  12

% get the diagonal elements
A(1 : size(A,1)+1 : end)

1  5  9

% tile the array over a 7x6 area,
% circularly wrapping the array
A([1 2 3 1 2 3 1], [1 2 3 4 1 2])

% alternatively,
A(mod([0:6], size(A,1))+1, mod([0:5], size(A,2))+1)

1  4  7  10  1  4
2  5  8  11  2  5
3  6  9  12  3  6
1  4  7  10  1  4
2  5  8  11  2  5
3  6  9  12  3  6
1  4  7  10  1  4

% remove the 2nd and 4th columns
A(:, [2 4]) = []

1  7
2  8
3  9

Repmat without using Repmat

Repmat replicates and tiles an array.  It is a powerful and commonly used function as it can be used to "vectorize" code.  Vectorizing code generally means to take code with looping constructs and re-code them to utilize Matlab's inherent ability to work with arrays.  The reasoning for this is that Matlab has traditionally been very slow in executing code with loops, but is quite fast with matrix operations.

Consider the following problem.  You have a matrix A and you want to normalize each row independently so that each row sums to 1.  If you were coding like you would in C/C++/C#/Java, you'd write something like this, where you sum up each element in a row in A, then divide each element in that row with that value.

A = rand(3, 4);

% for each row
for i = 1 : size(A, 1)
    z = 0;

    % sum all elements in the "i"th row
    for j = 1 : size(A, 2)
        z = z + A(i, j);
    end

    % divide all elements in the "i"th row by the sum, z
    for j = 1 : size(A, 2)
        A(i, j) = A(i, j) / z;
    end
end

This is a dumb way write code in Matlab.  It's much too verbose and slow.  The vectorized version (no loops) would look something like this:

z = sum(A, 2);
A = A ./ repmat(z, [1 size(A,2)]);

This requires fewer lines of code, is more concise, is arguably easier to read, and most importantly, runs faster.

Repmat is slow (faster than loops, but not fast enough!)

So, repmat can do some cool stuff.  The problem is that repmat is not a built-in function (compare what you get when you type "edit repmat" and "edit reshape" into the Matlab command window).  Repmat is actually a Matlab script.  If you call repmat many times, the overhead can be quite large.  Not just the fact that it's a script, but because it has to check the arguments, and has "if" statements and "for" loops so that it can be general enough to handle arrays with arbitrary number of dimensions and sizes.

Consider the following statement that was used to normalize the row vectors:

A = A ./ repmat(z, [1 size(A,2)]);

In the previous section, we saw the following, where a column vector, x, could be replicated across multiple columns by indexing the first (and only) column multiple times.

repmat(x, [1 3])

x(:, ones(3,1))

This is exactly the operation needed here and the normalization can be done without using repmat by indexing into the array intelligently:

A = A ./ z(:, ones(size(A,2), 1));

The following table shows how you can do different kinds of normalizations:

% normalize the rows

A = rand(3, 4);

z = sum(A, 2);
A = A ./ z(:, ones(size(A,2), 1));

% normalize the columns

A = rand(3, 4);

z = sum(A, 1);
A = A ./ z(ones(size(A,1), 1), :)

% normalize a multi-dimensional array over one of the dimensions

C = rand(4, 2, 3);

z = sum(C, 3);
C = C ./ z(:, :, ones(3,1));

% normalize over multiple dimensions

D = rand(4, 2, 3, 5);

z = sum(sum(D, 3), 4);
D = D ./ z(:, :, ones(3,1), ones(5,1));

 

Performing repmat over a singleton dimension is quite easy as you simply need to index "1" in that dimension multiple times. Slightly different code is required to do the following repmat operation, which duplicates the rows of x, a column vector, 3 times.

repmat(x, [3 1])

1
2
3
4
5
1
2
3
4
5
1
2
3
4
5

The key thing to note here is that the new array is obtained simply by repeating the rows 3 times.  We can achieve this effect by indexing the desired rows directly from x, namely the 1st through 5th rows followed by the 1st through 5th rows again and then one more time.  The follow lines of code both give the same result as the previous repmat code:

x([1 2 3 4 5 1 2 3 4 5 1 2 3 4 5])
x([1:end 1:end 1:end])

A fancier way to do this involves more intelligently creating the index into x.

idx = [1 : size(x,1)]'

1
2
3
4
5

repmatIdx = idx(:, ones(3, 1))

1   1   1
2   2   2
3   3   3
4   4   4
5   5   5

repmatIdx(:)

1
2
3
4
5
1
2
3
4
5
1
2
3
4
5

% These give the same result
% as repmat(x, [3 1])

x(repmatIdx(:))
x(repmatIdx, 1)
x(repmatIdx, :)

1
2
3
4
5
1
2
3
4
5
1
2
3
4
5

Note that when you give a multi-dimensional array as the index into an array (as in the last line of the table above), Matlab will usually convert it to a column vector (it will automatically convert repmatIdx into repmatIdx(:)), except when only one index is provided.  When only one index is provided, it will index according to the structure of the multi-dimensional array.  As a result, x(repmatIdx)does not give us the desired result and we also need to index the first (and only) column, x(repmatIdx, 1), or equivalently,  x(repmatIdx, :).

More concisely:

repmat(x, [3 1])

idx = [1 : size(x,1)]';
x(idx(:, ones(3, 1)), :)

x([1:size(x,1)]' * ones(1,3), :)

Where the alternative formulation (second row on the right) uses matrix multiplication to obtain the desired index map, i.e.

idx(:, ones(3, 1)) == [1:size(x,1)]' * ones(1,3)

The extension to multi-dimensional arrays is quite straight forward.  Consider the following matrix:

A = zeros(3, 4);
A(:) = 1:12;

A

1   4   7   10
2   5   8   11
3   6   9   12

Repmat can be done by repeated indexing of the rows and columns.  The following lines of code give the same result:

repmat(A, [2 3])

A([1 2 3 1 2 3],
  [1 2 3 4 1 2 3 4 1 2 3 4])

A([1:end 1:end],
  [1:end 1:end 1:end])

1  4  7  10  1  4  7  10  1  4  7  10
2  5  8  11  2  5  8  11  2  5  8  11
3  6  9  12  3  6  9  12  3  6  9  12
1  4  7  10  1  4  7  10  1  4  7  10
2  5  8  11  2  5  8  11  2  5  8  10
3  6  9  12  3  6  9  12  3  6  9  12

Again, the index into the array can be done more intelligently:

rowIdx = [1:size(A,1)]';
colIdx = [1:size(A,2)]';

repmatRowIdx = rowIdx(:, ones(2,1));
repmatColIdx = colIdx(:, ones(3,1));

Now, repmatRowIdx(:) and repmatColIdx(:) give us the same indices that we used in the table above and the following gives us what we want:

A(repmatRowIdx(:), repmatColIdx(:))

However, since we are indexing into multiple dimensions, Matlab will automatically change repmatRowIdx into repmatRowIdx(:) and repmatColIdx into repmatColIdx(:) when we do the following and this also gives us what we want:

A(repmatRowIdx, repmatColIdx)

More concisely:

repmat(A, [2 3])

rowIdx = [1 : size(A,1)]';
colIdx = [1 : size(A,2)]';

A(rowIdx(:, ones(2,1)), colIdx(:, ones(3,1)))

A([1:size(A,1)]' * ones(1,2),
  [1:size(A,2)]' * ones(1,3))

Where, as before, you can use matrix multiplication to obtain the array indices in order to do everything in one line of code.

References

  1. Matlab array manipulation tips and tricks
    Excellent reference on how to do incredible things with arrays.  A must read.  It's my (Matlab) bible.

  2. MathWorks' Matlab matrix indexing
    Shows how to do some fancy indexing into arrays.

  3. MathWorks' Matlab code vectorization guide
    The basics in code vectorization.

  4. MathWorks' Matlab online documentation
    So you can read about the hundreds of functionalities of Matlab.