|
function []=C()sigma=0.2465
x0=[0,0,0]
lambda=4
S0=129.03
deltat=1/252
T=132/252
n=132
r=0.0041
w=1000
K=70for s = 1:w
St(s,:) = ones(1,n+1);
St(s,1) = S0;
for k = 1:n
X = random('Poisson',deltat*lambda);%compute the number of jumps in the period of deltat
Y=1.1+rand(1,X)*(1.3-1.1);%generate the intensity of each jump
M=prod(Y);%compute the cumulative intensity of and jumps in the period of deltat
deltaw = random('Normal',0,sqrt(deltat),1,1);
St(s,k+1) = St(s,k)*exp((r-sigma^2/2)*deltat+sigma*deltaw)*M; % solution by Equation
end
endfor j=n:-1:1
p=find(St(:,j)<K);
i=1:length(p);
a(i,j)=St(p(i),j);
c(i,j)=max(K-St(p(i),j),0)*exp(-r*deltat);
xdate=a(i,j);
ydate=c(i,j);
a=lsqcurvefit(@Regfun,x0,xdate,ydate);
C1(i,j)=a(1)+a(2)*a(i,j)+a(3)*a(i,j).^2;
if C1(i,j)>K-a(i,j);
St(i,j)=K+1;
else
St(i,j+1)=K+1;
end
endfor i=1:w
j=1:n
C=sum(max(K-St(i,j),0)*exp(-r*deltat*j))/w
end哪里错了 |