y=csvread('Task_2.csv');

N=length(y);
u=zeros(1,N);
e=randn(1,N);

initial_theta_1=0.2;
initial_theta_2=0.5;
initial_theta_3=0.8;

l=100;

for i = 4:N
    for j=4:N-l
        y_noise(i)= y(i-1)*initial_theta_1+y(i-2)*initial_theta_2+y(i-3)*initial_theta_3+e(j);
    end
end

for i = 1:N
    r = rand();
    if r <= .33
        u(i) = -1;
    elseif r <= .66
        u(i) = 0;
    else
        u(1) = 1;
    end
end

Y = y_noise(1,4:end)';
X=[y_noise(1,3:end-1)', y_noise(1,2:end-2)', y_noise(1,1:end-3)'];
theta_ls = (X'*X) \ X'*Y;

z = zeros(1,N);

for i = 2:N
    z(i) = -theta_ls(1)*z(i-1) + theta_ls(2)*u(i-1);
end
KSZI1 = [z(1:end-3)', u(1:end-3)'];

theta_iv_1 = (KSZI1'*X) \ KSZI1'*Y;

z2 = zeros(1,N);

for i = 2:N
    z2(i) = -theta_iv_1(1)*z2(i-1) + theta_iv_1(2)*u(i-1);
end

KSZI2 = [z2(1:end-3)', u(1:end-3)'];

theta_iv_2 = (KSZI2'*X) \ KSZI2'*Y;