vectorize lsqnonneg?

2011 11 23
2 个回答
回答已采纳
4 次查看（30 天）

0 个投票

Hi,
I'm currently working with sequences of images and I'm interested in the temporal dynamics of individual pixels. Typically I have a series of approx. 50 images with 250 000 pixels of interest and I want to fit them to 5-10 time-activity curves of reference. I need to do a non-negative fitting because negative coefficients are physically nonsensical in my problem.
At the moment, I'm basically doing this
nPixels = 250 000;
nRefTACs = 5;
nImages = 50;
% imTS := image time series data [nImages x nPixels]
% refTACs := reference time-activity curves [nImages x nRefTACs]
% nnCoeff := fitting coefficients [nRefTACs x nPixels]
nnCoeff = zeros(nRefTACs, nPixels);
for k = 1:nPixels
nnCoeff(:,k) = lsqnonneg(refTACs, imTS(:,k));
end
I'm wondering if there would be a way to speed this up by vectorization or by using another optimization function. I've been browsing the optimization toolbox documentation for some time, but I'm not very familiar with optimization in general.
Any help would be appreciated.

请先登录，再回答此问题。

 采纳的回答

Matthias Klemm
Matthias Klemm 2011-12-9

1 个投票

lsqnoneg is written in Matlab code. You can open it with the editor and vectorize it yourself to gain some speed. You can't avoid the loops in there I guess. Still you should get some improvement out of it. The only show stopper is mldivide (\) as it is not vectorized. In your case it should work as you already can use mldivide to compute a (not ideal) solution for your problem (because your problem is not overdetermined I think).

2 个评论
显示 无 隐藏 无

David Provencher
David Provencher 2011-12-13
For some reason, I did not think about doing this, thanks for the suggestion. I managed to cut down computation time by four for my specific problem, so I guess that's something.
About mldivide, it is vectorized, but since the variables in the positive space (P) change from one data vector to the next, there needs to be a loop. In my case, I ended up looping over the unique columns of P (the nth column is the equivalent of P in the original code for the nth data vector) and solving for all variables at once which have identical combination of models (column of P) on a given iteration. However, I'm guessing the efficiency of this can vary on the number of models and data vectors.
Silas Leavesley
Silas Leavesley 2012-12-1
Hi David,
Would you mind posting your modified lsqnonneg, or at least the modified portion of it? I am currently using lsqnonneg for spectral image analysis and would like to see if this would speed things up. Thanks.

请先登录，再进行评论。

更多回答（1 个）

Sean de Wolski
Sean de Wolski 2011-11-23

0 个投票

preallocate nnCoeff so that it doesn't change size on each iteration.
nnCoeff = zeros(nRedTACs,nPixels);
for k ...
nnCoeff(:,k) = lsqnonneg(refTACs, imTS);
end

3 个评论
显示 1更早的评论 隐藏 1更早的评论

David Provencher
David Provencher 2011-11-23
It is already preallocated, although it might not be clear from the snippet of (pseudo-) code I wrote. I'm rather looking for something along the lines of vectorizing the code or using another optimization function to get rid of the loop. For instance, I could use
nnCoeff = refTACs \ imTS;
but this would not satisfy the non-negativity constraint. However, since it runs in approx. 0.5s compared to 50s for my loop, I was hoping I could speed things up a bit.
Steve Grikschat
Steve Grikschat 2011-12-13
Can you expand this into one large problem? That may or may not make things slower, but it could be worth a try. Additionally, if you have Parallel Computing Toolbox, you could solve these problems in parallel easily.
The reason that (\) is faster than lsqnonneg is that lsqnonneg has to solve a series of linear systems using (\) in order to satisfy the non-negativity constraint.
FYI: I assume that you are actually changing the RHS (imTS) in the loop as well. Your code snippet doesn't show that.
David Provencher
David Provencher 2011-12-14
I'm not sure I understand what you mean by expanding it into one large problem. I would think it means to solve everything at once, as in
Coeff = refTACs \ imTS;
but that is what I was trying to do, except with lsqnonneg instead of mldivide (\).
What I ended up doing was to vectorize all I could inside lsqnonneg and solve as many variables as I could with each call of (\). In short, each call of (\) in lsqnonneg fits the data to a subset of the models, so when I vectorized the code, I solve for all data vectors that call for the same subset of models at once. So basically, now it looks like this :
nnCoeff = lsqnonnegvect(refTACs, imTS);
instead of
for k = 1:nPixels
nnCoeff(:,k) = lsqnonneg(refTACs, imTS(:,k));
end
Also, you are right about imTS, only one pixel at a time is used. Thanks for pointing it out, I just updated it

请先登录，再进行评论。

类别

在 帮助中心 和 File Exchange 中查找有关 Choose a Solver 的更多信息

产品

标签

Community Treasure Hunt

Find the treasures in MATLAB Central and discover how the community can help you!

Start Hunting!

Translated by