2014-08-31 88 views
1

我對二叉樹有一個向後遞歸。在每個節點上,一個未知的C以這樣一種方式進入,即在起始節點處我們得到一個公式,A(1,1),取決於C。代碼如下:在Matlab中加速符號遞歸

A=sym(zeros(1,Steps)); 
B=zeros(1,Steps); 
syms C; % The unknown that enters A at every node 
tic 
for t=Steps-1:-1:1 

% Values needed in A and B 
Lambda=1-exp(-(1./S(t,1:t).^b).*h); 
Q=((1./D(t))./(1-Lambda)-d)/(u-d); 
R=normcdf(a0+a1*Lambda); 

% the backward recursion for A and B 
A(1:t)=D(t)*C+D(t)*... 
    (Q.*(1-Lambda).*A(1:t) ... 
     + (1-Q).*(1-Lambda).*A(2:t+1)); 

B(1:t)=Lambda.*(1-R)+D(t)*... 
    (Q.*(1-Lambda).*B(1:t)... 
     + (1-Q.*(1-Lambda).*B(2:t+1))); 
end 

C = solve(A(1,1)==sym(B(1,1)),C); 

此代碼大約需要4秒,如果Steps = 104。但是,如果我們刪除C並將矩陣A設置爲常規雙矩陣,則只需要大約0.02秒。因此使用syms可將計算時間增加200倍。這對我來說似乎太多了。任何建議加快這一點?

我用Matlab 2013b上的一臺MacBook Air 13英寸的2013年春季。此外,如果你有興趣在上面的部分(不知道是否是相關的)之前的代碼:

a0 = 0.9; 
a1 = -3.2557; 
b = 1.2594; 
S0=18.57; 
sigma=0.6579; 
h=1/104; 
T=1; 
Steps=T/h; 
f=transpose(normrnd(0.04, 0.001 [1 pl])); 
D=exp(-h*f); % discount values 
pl=T/h; % pathlength - amount of steps in maturity 

u=exp(sigma*sqrt(h)); 
d=1/u; 

u_row = repmat(cumprod([1 u*ones(1,pl-1)]),pl,1); 
d_row = cumprod(tril(d*ones(pl),-1)+triu(ones(pl)),1); 
path = tril(u_row.*d_row); 
S=S0*path; 
+0

您能否提供包含所有必需變量的代碼?我不認爲有很多需要改進的地方,但我會仔細看看是否可以運行代碼。 – Daniel 2014-08-31 17:15:03

+0

'f'是未定義的,它是什麼?如果你意識到什麼是符號數學相對於浮點運算的話,200的因子也不算太壞。但是,是的,有辦法加快速度。是否有理由使用符號數學呢?有一件事你可以做到這一點,可能會加速版本是從'for'循環中刪除'Lambda','Q'和'R'的所有或部分 - 創建向量和索引。 – horchler 2014-08-31 20:02:47

+0

@Daniel:我的確忘記了包含變量f。但現在已經完成。 – 2014-09-01 07:12:58

回答

0

除非我錯過了一些東西,沒有必要使用符號數學或使用未知變量。你可以在你的遞歸關係中有效地假設C = 1並在最後解決實際值。下面是一些其他的改進全碼:

rng(1); % Always seed your random number generator 

a0 = 0.9; 
a1 = -3.2557; 
b = 1.2594; 
S0 = 18.57; 
sigma = 0.6579; 
h = 1/104; 
T = 1; 
Steps = T/h; 
pl = T/h; 
f = 0.04+0.001*randn(pl,1); 
D = exp(-h*f); 
u = exp(sigma*sqrt(h)); 
d = 1/u; 

u_row = repmat(cumprod([1 u*ones(1,pl-1)]),pl,1); 
d_row = cumprod(tril(d*ones(pl),-1)+triu(ones(pl)),1); 
pth = tril(u_row.*d_row); 
S = S0*pth; 

A = zeros(1,Steps); 
B = zeros(1,Steps); 
tic 
for t = Steps-1:-1:1 
    Lambda = 1-exp(-h./S(t,1:t).^b); 
    Q = ((1./D(t))./(1-Lambda)-d)/(u-d); 
    R = 0.5*erfc((-a0-a1*Lambda)/sqrt(2)); % Faster than normcdf 

    % Backward recursion for A and B 
    A = D(t)+D(t)*(Q.*(1-Lambda).*A(1:end-1) + ... 
        (1-Q).*(1-Lambda).*A(2:end)); 
    B = Lambda.*(1-R)+D(t)*(Q.*(1-Lambda).*B(1:end-1) + ... 
          (1-Q.*(1-Lambda).*B(2:end))); 
end 
C = B/A 
toc 

這需要大約0.005秒,在我的MacBook Pro運行。當然,你可以做出其他改進。在多個地方有很多變量組合(例如,1-LambdaD(t)*(1-Lambda)),可以計算一次。 Matlab可能會嘗試優化這一點。您可以嘗試將Lambda,QR移出循環 - 或者至少在外部計算它們的一部分並將結果保存在數組中。

+0

謝謝,我會盡快檢查出來。 – 2014-09-02 13:11:45

+0

你完全正確。 'C'可以簡單地從遞歸中取出。謝謝! – 2014-09-03 06:45:25