Octave:这些FOR循环如何被矢量化？

Question

我正在写一个Octave脚本来计算欧洲期权的价格.

第一部分使用蒙特卡罗模拟n个时间段内的标的资产价格.这重复了几次.

Octave使得设置初始矩阵变得非常容易.但我还没有找到以矢量化方式完成任务的方法,避免使用FOR循环:

%% Octave simplifies creation of 'e', 'dlns', and 'Prices'
e = norminv(rand(nIter,n));
dlns = cat(2, ones(nIter,1), exp((adj_r+0.5*sigma^2)*dt+sigma*e.*sqrt(dt)));
Prices = zeros(nIter, n+1);

for i = 1:nIter              % IS THERE A WAY TO VECTORIZE THESE FOR LOOPS?
  for j = 1:n+1
    if j == 1
      Prices(i,j)=S0;
    else
      Prices(i,j)=Prices(i,j-1)*dlns(i,j);
    end
  endfor
endfor

请注意,n中的价格等于n-1倍因子的价格,因此以下情况不起作用......

Prices(i,:) = S0 * dlns(i,:)

...因为它需要S0并乘以所有因子,产生与预期的随机游走不同的结果.

Answer 1

由于迭代之间依赖于获取每个新列相对于前一列的结果,似乎您需要至少一个循环,但是以矢量化方式在列中执行所有操作,这可能会为您加速.两个嵌套循环的矢量化替换看起来像这样 -

Prices(:,1)=S0;
for j = 2:n+1
    Prices(:,j) = Prices(:,j-1).*dlns(:,j);
endfor

我刚刚意识到依赖关系可以cumprod得到cumulative product解决,这将使我们基本上在这里完成,从而导致无循环解决方案!这是实施 -

Prices = [repmat(S0,nIter,1) cumprod(dlns(:,2:end),2)*S0]

在MATLAB上进行基准测试

基准代码 -

%// Parameters as told by OP and then create the inputs
nIter= 100000;
n = 100;
adj_r = 0.03;
sigma = 0.2;
dt = 1/n;
S0 = 60;
e = norminv(rand(nIter,n));
dlns = cat(2, ones(nIter,1), exp((adj_r+0.5*sigma^2)*dt+sigma*e.*sqrt(dt)));

disp('-------------------------------------- With Original Approach')
tic
Prices = zeros(nIter, n+1);
for i = 1:nIter
    for j = 1:n+1
        if j == 1
            Prices(i,j)=S0;
        else
            Prices(i,j)=Prices(i,j-1)*dlns(i,j);
        end
    end
end
toc, clear Prices

disp('-------------------------------------- With Proposed Approach - I')
tic
Prices2(nIter, n+1)=0; %// faster pre-allocation scheme
Prices2(:,1)=S0;
for j = 2:n+1
    Prices2(:,j)=Prices2(:,j-1).*dlns(:,j);
end
toc, clear Prices2

disp('-------------------------------------- With Proposed Approach - II')
tic
Prices3 = [repmat(S0,nIter,1) cumprod(dlns(:,2:end),2)*S0];
toc, clear Prices3

运行时结果 -

-------------------------------------- With Original Approach
Elapsed time is 0.259054 seconds.
-------------------------------------- With Proposed Approach - I
Elapsed time is 0.020566 seconds.
-------------------------------------- With Proposed Approach - II
Elapsed time is 0.067292 seconds.

现在,运行时确实表明第一个提出的方法可能更适合这里!