2017-01-27 133 views
1

矢量化我想「矢量化」這個循環在Matlab的計算效率矩陣二次方程式在MATLAB

for t=1:T 
    j=1; 
    for m=1:M 
     for n=1:N 
     y(t,j) = v{m,n} + data(t,:)*b{m,n} + data(t,:)*f{m,n}*data(t,:)'; 
     j=j+1; 
     end 
    end 
end 

哪裏v是標量的(M X N)細胞。 b是(Kx1)個載體的(M xN)單元。 f是(K x K)矩陣的(M x N)單元格。 data是(T x K)數組。

爲了讓我的意思我曾經向量化同一迴路中無需二次項的代碼是一個示例:

B = [reshape(cell2mat(v)',1,N*M);cell2mat(reshape(b'),1,M*N)]; 
X = [ones(T,1),data]; 
y = X*B; 

謝謝!

回答

2

對於那些有興趣在這裏被我發現

f = f'; 
tMat = blkdiag(f{:})+(blkdiag(f{:}))'; 
y2BB = [reshape(cell2mat(v)',1,N*M);... 
     cell2mat(reshape(b',1,M*N));... 
     reshape(diag(blkdiag(f{:})),K,N*M);... 
     reshape(tMat((tril(tMat,-1)~=0)),sum(1:K-1),M*N)]; 
y2YBar = [ones(T,1),data,data.^2]; 

jj=1; 
kk=1; 
ll=1; 
for k=1:sum(1:K-1) 
    y2YBar = [y2YBar,data(:,jj).*data(:,kk+jj)]; 
    if kk<(K-ll) 
     kk=kk+1; 
    else 
     kk=1; 
     jj=jj+1; 
     ll=ll+1; 
    end 
end 
y = y2YBar*y2BB; 
+0

最後一個循環是草率的,但我的大腦傷害。如果您想到更有效的方式來添加這些元素,請告訴我。 – hipHopMetropolisHastings

0

解決方案下面是針對性能的最量化的形式 -

% Extract as multi-dim arrays 
vA = reshape([v{:}],M,N); 
bA = reshape([b{:}],K,M,N); 
fA = reshape([f{:}],K,K,M,N); 

% Perform : data(t,:)*f{m,n} for all iterations 
data_f_mult = reshape(data*reshape(fA,K,[]),T,K,M,N); 

% Now there are three parts : 
% v{m,n} 
% data(t,:)*b{m,n} 
% data(t,:)*f{m,n}*data(t,:)'; 

% Compute those parts one by one 
parte1 = vA(:).'; 
parte2 = data*reshape(bA,[],M*N); 

parte3 = zeros(T,M*N); 
for t = 1:T 
    parte3(t,:) = data(t,:)*reshape(data_f_mult(t,:,:),K,[]); 
end 

% Finally sum those up and to present in desired format permute dims 
sums = bsxfun(@plus, parte1, parte2 + parte3); 
out = reshape(permute(reshape(sums,T,M,N),[1,3,2]),[],M*N); 
+0

不錯的工作!我會測試哪一種效率最高,並選擇答案。考慮到我的應用程序T = 152,N = K = 3和M = 8,我的循環次數要少得多。 – hipHopMetropolisHastings

+0

@hipHopMetropolisHastings有沒有更新? – Divakar

+0

是的,你的代碼運行在0.025和我的0.008左右。再次感謝您的回答。 – hipHopMetropolisHastings