iterate over single array dimension

I am having trouble iterating over a single array dimension. For instance, let's say I have a dozen sets of data, and for each dataset I want to solve the quadratic equation. The polynomial coefficient list for all datasets is stored in a single 2D array with dimensions [12,3].

data = randn(12,3);

For a single set of 3 coefficients, I can use the roots function to get the two roots to the equation. How can I iterate over all 12 datasets? Do I need to write a 'for' loop?

In the past, I used the bsxfun something like the following:

quadraticRoots = bsxfun(@(x,y)roots(x), data, data);

In the days before implicit expansion, this would give me an output array with dimensions [12,2] which is what I want. There are two roots for each of the 12 datasets. With the updates to bsxfun, I get a message that the input must be a vector. Instead of single-dimension expansion, bsxfun implicitly expands everything! Arrayfun is no help because it also implicitly expands each element in the 2D array.

The above case is relatively simple, but I do this type of operation very frequently. I have a list of vectors stored as a 2D array, and want to apply a function to each row or column in the list. Is there a general way to iterate over a single array dimension without a 'for' loop?

[Guillaume]

As far as I know the behaviour of bsxfun has not changed at all since implicit expansion was introduced. Your usage of bsxfun doesn't make much sense, you've got an anonymous function with two inputs, x and y, but it doesn't use y. And since you're passing twice the same input, there's obviously no expansion to be done. So your equation reduces to (regardless of version):

quadraticRoots = roots(data);

which obviously doesn't work with matrices.

[G Smith]

I wish I still had an old version of MATLAB installed, but in the old days, they listed bsxfun as having singleton expansion instead of implicit expansion mentioned today. You can see that on the internet archive:https://web.archive.org/web/20160103023616/https://www.mathworks.com/help/matlab/ref/bsxfun.html

What I did with the function above was make use of that singleton expansion to iterate over the first array index only. It was compact to write, very fast to calculate, and used little memory. Repeating the data a second time was just a trick to fit into bsxfun requirement for two input variables. As you said, the second copy would not be used.

Of course this is all just relying on my memory since I don't still have an old version. It is possible the function would not work as written, but it was a useful trick so I believe the concept is correct. I never before used the trick with the roots function, but I frequently used it with other anonymous functions that do similar kind of processing.

Answer 1

[Jon]

You can do it with arrayfun for example

coef = randn(12,3)
qroots = cell2mat(arrayfun(@(n) roots(coef(n,:)),1:size(coef,1),'UniformOutput',false))

This returns a 2 by 12 array of roots, one column for each row of coef. You could of course transpose it if you wanted.

This seems quite arcane though. I don't know if there is any performance advantage of doing it this way rather than just putting it in a loop. Certainly a little loop would be more readable. Maybe someone know how to do this in a cleaner way.

[G Smith]

Nice idea. Iterating over an index is something I had not thought of. It is basically rolling the loop into a single statement. As you say, the loop is probably more readable though.

I profiled your solution vs a 'for' loop with a pre-allocated array, and the loop was actually a bit faster. Since I do this type of operation frequently I was hoping for a compact, vectorized solution. I hate to start cluttering my code with lots of loops, but perhaps it is the best solution.

[Guillaume]

The advantage of arrayfun or cellfun over explicit loops is that you don't have to worry about preallocation, getting the correct bounds, etc. It's also a form of functional programming preferred by some. In my opinion, it's also clearer than a loop, it clearly says: do this for each element. However, it will rarely outperform the speed of a loop due to the extra function call and here the cost of the anynomous function. If your only criteria is speed (as opposed to clarity and maintainability which are often more important) then don't use arrayfun.

I'd use cellfun personally so you don't even have to compute any size:

qroots = cell2mat(cellfun(@roots, num2cell(coef, 1), 'UniformOutput', false))

This also loses the anonymous function (less costly) but adds a conversion from matrix to cell array (more costly).

[G Smith]

Nice. Swapping back-and-forth between cell array is another good solution. It is a bit cleaner than using the index for the arrayfun.

I had to modify your code a little to make it work.

qroots = cell2mat(cellfun(@roots, num2cell(coef, 2), 'UniformOutput', false).').';

As you pointed out, the timing is still a bit slower than a simple loop (likely because of the overhead in converting to/from a cell array), but at least it is easier to read.

[Jon]

I'm not sure about performance, but for readability you could always make yourself a little function like this (would need a little more of a header to fit standard documentation)

function B = rowiter(fun,A)
% iterate input fun over rows of A
% fun is a function handle to a function that  accepts a single vector
% input and returns a vector output

% get dimensions
numRows = size(A,1);

% determine length of vector returned by fun
% (assume it is the same for all rows and just check the first)
numCols = numel(fun(A(1,:)));

% preallocate output array
B = zeros(numRows,numCols);

% iterate through rows of A
for k = 1:numRows
    out = fun(A(k,:));
    B(k,:) = out(:); % use (:) to make sure it is a column
end

and then for the specific example you had you could just apply it for example

coef = randn(12,3);
qroots = rowiter(@roots,coef)