

load d_ibmvwewsp6203.txt -ascii;
r = d_ibmvwewsp6203;
plot(r(:,5));
rsp = r(1:3000,5);
plot(rsp)

%ESTIMACION GMM
p = sacf(rsp,12);
theta = (-1+sqrt(1-4*(p(1)^2)))/(2*p(1));
v = var(rsp);
display('Estimación de la constante')
mu = mean(rsp)
display('Estimación del parámetro del MA(1)')
sigma = v/(1+theta^2)



%   SIMULAR RENTABILIDAD
%   Hacemos un supuesto sobre la innovacion: suponemos normalidad
M = randn(10,500)*sqrt(sigma);
dt_eps = sqrt(sigma);

mean(M');         % Promedios para cada periodo. 
sqrt(var(M,0,2)); % Para cada período.
                  %  Debe aproximarse a 0.66 que es la raiz de la estimación de sigma

%generas la primera posicion de las 500 sendas
s(1,:) = mu*ones(1,500) + M(1,:);
for i= 2:10
s(i,:) = mu*ones(1,500) + M(i,:) - theta.*M((i-1),:);
end

plot(s(:,1:5))

SP = [];
SP(1,:) = 1784.*(1+s(1,:));
for i =2:10
    SP(i,:) = SP(i-1,:).*(1+s(i,:));
end
size(SP);
mean(SP,2);
SP = [1784.*ones(1,500);SP];
plot(SP(:,1:5));
plot(SP) %representacion de la 500 sendas

%valorar call
strike = 1850;
result = SP(11,:)-strike.*ones(1,500);
result = max(result,0);
r = 0.03;
payoffdescon = result./(1.03^(10/360));
valorcall = mean(payoffdescon)
