代码之家 › 专栏 › 技术社区 › shamalaia

回归图计算的矢量化

covariance vectorization regression performance matlab

shamalaia · 技术社区 · 7 年前

我计算时间序列的回归图 A(t) 在田野上 B(x,y,t) 以如下方式:

A=1:10; %time
B=rand(100,100,10); %x,y,time

rc=nan(size(B,1),size(B,2));
for ii=size(B,1)
  for jj=1:size(B,2)
     tmp = cov(A,squeeze(B(ii,jj,:))); %covariance matrix
     rc(ii,jj) = tmp(1,2); %covariance A and B
  end
end
rc = rc/var(A); %regression coefficient

有没有办法使代码矢量化/加速或者是一些我不知道的内置函数来达到同样的效果?

1 回复 | 直到 7 年前

Dev-iL 7 年前

为了使这个算法矢量化,你必须“弄脏你的手”,自己计算协方差。如果你看看里面 cov 您将看到它有许多输入检查行和很少的实际计算行来总结关键步骤:

y = varargin{1};
x = x(:);
y = y(:);
x = [x y];
[m,~] = size(x);
denom = m - 1;
xc = x - sum(x,1)./m;  % Remove mean
c = (xc' * xc) ./ denom;

为了稍微简化上述内容:

x = [x(:) y(:)];
m = size(x,1);
xc = x - sum(x,1)./m;
c = (xc' * xc) ./ (m - 1);

现在这是一个相当简单的矢量化…

function q51466884
A = 1:10; %time
B = rand(200,200,10); %x,y,time
%% Test Equivalence:
assert( norm(sol1-sol2) < 1E-10);
%% Benchmark:
disp([timeit(@sol1), timeit(@sol2)]);

%%
function rc = sol1()
rc=nan(size(B,1),size(B,2));
for ii=1:size(B,1)
  for jj=1:size(B,2)
     tmp = cov(A,squeeze(B(ii,jj,:))); %covariance matrix
     rc(ii,jj) = tmp(1,2); %covariance A and B
  end
end
rc = rc/var(A); %regression coefficient
end

function rC = sol2()  
m = numel(A);
rB = reshape(B,[],10).'; % reshape
% Center:
cA = A(:) - sum(A)./m;
cB = rB - sum(rB,1)./m;
% Multiply:
rC = reshape( (cA.' * cB) ./ (m-1), size(B(:,:,1)) ) ./ var(A);
end

end

我得到这些时间: [0.5381 0.0025] 这意味着我们在运行时节省了两个数量级:)

请注意,优化算法的很大一部分是假设数据中没有任何“奇怪之处”,比如 NaN 价值观等等,看看里面 cov.m 看看我们跳过的所有支票。

推荐文章

Costinean Sebastian · 如何优化使用for循环处理大型数组的MATLAB脚本的执行时间?"

1 年前

Gyro · 为什么我的%xmm3寄存器在vbroadcastsd中使用第一个参数,而不是第四个参数?

1 年前

aeiou · 对日志表进行矢量化处理,以确定最新的可用性

1 年前

IchKenneDeinenNamen · 当将向量存储为类属性时,为什么NumPy的矢量化评估速度较慢?

2 年前

roshoka · 如何在不使用for循环的情况下将函数应用于林空间的每个元素

2 年前

Xavier · 变长矢量的矢量化和

2 年前

bbbb · 矢量化或在大型数据帧上使用多重处理[重复]

2 年前

Sterling Butters · NumPy 1D近邻中循环的避免

2 年前

AnthonyML · 修改2d numpy数组中分隔符后的行中的所有元素

2 年前

Tim · np.vectorize和relativelta返回“relativeltaonlydiffs datetime/date”

2 年前