function [var_share] = vardcmp_rec(Q,sxi,BD,RD)
global T;global F;global ind65;
nt=48;qnt=48*4;
IND=ind65;
%load indfile;IND=ind2;

Q1 = zeros(qnt,qnt);
Q1(1,1)=Q(1,1);
Q2 = zeros(qnt,qnt);
Q2(2,2)=Q(2,2);
Q3 = zeros(qnt,qnt);
Q3(3,3)=Q(3,3);
Q4 = zeros(qnt,qnt);
Q4(4,4)=Q(4,4);

tau = zeros(nt,nt);
 for i=1:nt
 for j=2:nt
     k=i-j+1;
     if k>0
     tau(i,j-1)=k;
     else 
     tau(i,j-1)=0;
     end
 end
 end

for nn=1:nt
    tmp = 4*(nn-1);
    ssxi(tmp+1,nn:nn+T-1) = sxi(1,:);
    ssxi(tmp+2,nn:nn+T-1) = sxi(2,:);
    ssxi(tmp+3,nn:nn+T-1) = sxi(3,:);
    ssxi(tmp+4,nn:nn+T-1) = sxi(4,:);
end

 dssxi = ssxi(:,nt:T);
 dsxi = sxi(:,nt:T);
 TT = T-nt+1;
 nn = 49;
 Xi=dssxi;

vVY_total=zeros(nn,TT);
vVY_1=zeros(nn,TT);
vVY_2=zeros(nn,TT);
vVY_3=zeros(nn,TT);
vVY_2_q1=zeros(nn,TT);
vVY_2_q2=zeros(nn,TT);
vVY_2_q3=zeros(nn,TT);
vVY_2_q4=zeros(nn,TT);

for n=1:48 
W1 = zeros(nt,TT);
W2 = zeros(nt,TT);
X1 = zeros(nt,TT);
X2 = zeros(nt,TT);
MW1 = ones(nt,nt,TT);
MW2 = ones(nt,nt,TT);
MC1 = zeros(nt,nt,TT);
MC2 = zeros(nt,nt,TT);
MP1 = zeros(nt,nt,TT);
MP2 = zeros(nt,nt,TT);
D_MW1 = zeros(nt,nt,TT);
D_MW2 = zeros(nt,nt,TT);
D_MP1 = zeros(nt,nt,TT);
D_MP2 = zeros(nt,nt,TT);
M_JH1 = zeros(nt,qnt,TT);

UL=zeros(nt,TT);
UH=zeros(nt,TT);
TU=zeros(nt,TT);
yhat = zeros(5,TT);
fJH1 = zeros(5,qnt,TT);
fJH = zeros(qnt,5,TT);

Xq=zeros(qnt,qnt,n-1);
Xq1=zeros(qnt,qnt,n-1);
Xq2=zeros(qnt,qnt,n-1);
Xq3=zeros(qnt,qnt,n-1);
Xq4=zeros(qnt,qnt,n-1);

SXq=zeros(qnt,qnt);
SXq1=zeros(qnt,qnt);
SXq2=zeros(qnt,qnt);
SXq3=zeros(qnt,qnt);
SXq4=zeros(qnt,qnt);

fVY_total=zeros(5,5,TT);
fVY_1=zeros(5,5,TT);
fVY_2=zeros(5,5,TT);
fVY_2_q1=zeros(5,5,TT);
fVY_2_q2=zeros(5,5,TT);
fVY_2_q3=zeros(5,5,TT);
fVY_2_q4=zeros(5,5,TT);
fVY_3=zeros(5,5,TT);

fxi0(:,n+1) = (F^n)*Xi(:,1);

for t=n+1:TT 

if IND(t+47,1)==0;
 D=BD;
else IND(t+47,1)==1;
 D=RD;
end
   
 for i=1:nt
 W1(i,t) = fxi0(1+4*(i-1),t);
 W2(i,t) = fxi0(2+4*(i-1),t);
 X1(i,t) = fxi0(3+4*(i-1),t);
 X2(i,t) = fxi0(4+4*(i-1),t);
 end

 for i=2:nt 
 for j=1:nt            
     if tau(i,j)>0
 MW1(i,j,t) = exp(-exp(X1(j,t)+D(tau(i,j),1)));
 MW2(i,j,t) = exp(-exp(X2(j,t)+D(tau(i,j),1)));
     else
 MW1(i,j,t) = 1;
 MW2(i,j,t) = 1;
 end
 end
 end

 D_MW1(1,1,t)= 1;
 D_MW2(1,1,t)= 1;
 
 for j=1:nt
 D_MW1(j,j,t)= prod(MW1(j,:,t),2);
 D_MW2(j,j,t)= prod(MW2(j,:,t),2);
 end

 MP1(1,:,t)=zeros(1,nt);
 MP2(1,:,t)=zeros(1,nt);

 for i=2:nt
 MP1(i,1:i-1,t)=W1(i,t)*D_MW1(i,i,t);
 MP2(i,1:i-1,t)=W2(i,t)*D_MW2(i,i,t);
 end

 MC1(1,:,t)=zeros(1,nt);
 MC2(1,:,t)=zeros(1,nt);

 for i=2:nt
 for j=1:nt 
     if tau(i,j)>0
 MC1(i,j,t) = -exp(X1(j,t)+D(tau(i,j),1));
 MC2(i,j,t) = -exp(X2(j,t)+D(tau(i,j),1));
     else
 MC1(i,j,t) = 0;
 MC2(i,j,t) = 0;
 end
 end
 end
 
 for i=1:nt
 for j=1:nt 
 D_MP1(i,j,t)= MP1(i,j,t).*MC1(i,j,t);
 D_MP2(i,j,t)= MP2(i,j,t).*MC2(i,j,t);
 end
 end

for i=1:nt
M_JH1(:,1+4*(i-1),t)=D_MW1(:,i,t);
M_JH1(:,2+4*(i-1),t)=D_MW2(:,i,t);
M_JH1(:,3+4*(i-1),t)=D_MP1(:,i,t);
M_JH1(:,4+4*(i-1),t)=D_MP2(:,i,t);
end

fJH1(1,:,t)= M_JH1(1,:,t);
fJH1(2,:,t)= M_JH1(2,:,t)+M_JH1(3,:,t);
fJH1(3,:,t)= M_JH1(4,:,t)+M_JH1(5,:,t)+M_JH1(6,:,t);
fJH1(4,:,t)= sum(M_JH1(7:12,:,t),1);
fJH1(5,:,t)= sum(M_JH1(13:48,:,t),1);

FJH1(t,:) = ones(1,5)*fJH1(:,:,t);
FJH(:,:,t) = reshape(FJH1(t,:)',4,nt);

if n<nt
FJH(:,n+1:nt,t) = zeros(4,nt-n);
CFJH(:,1:n,t)=cumsum(FJH(:,1:n,t),2);
sq_CFJH(:,1:n,t)=CFJH(:,1:n,t).*CFJH(:,1:n,t);
csum_CFJH(:,t,n)=sum(sq_CFJH(:,1:n,t),2);
else
CFJH(:,1:nt,t)=cumsum(FJH(:,1:nt,t),2);
sq_CFJH(:,1:nt,t)=CFJH(:,1:nt,t).*CFJH(:,1:nt,t);
csum_CFJH(:,t,n)=sum(sq_CFJH(:,1:nt,t),2);
end

vq=[Q(1,1),Q(2,2),Q(3,3),Q(4,4)];
vVY_2(t,n)= vq*csum_CFJH(:,t,n);
vVY_2_q1(t,n)= vq(1,1)*csum_CFJH(1,t,n);
vVY_2_q2(t,n)= vq(1,2)*csum_CFJH(2,t,n);
vVY_2_q3(t,n)= vq(1,3)*csum_CFJH(3,t,n);
vVY_2_q4(t,n)= vq(1,4)*csum_CFJH(4,t,n);

fxi0(:,t+1) = (F^n)*Xi(:,t+1-n);
end
sum_TU(n,:) = sum(TU,1);
sum_vVY_q(n,1) = sum(vVY_2(:,n),1)/(TT-n);
sum_vVY_q(n,2) = sum(vVY_2_q1(:,n),1)/(TT-n);
sum_vVY_q(n,3) = sum(vVY_2_q2(:,n),1)/(TT-n);
sum_vVY_q(n,4) = sum(vVY_2_q3(:,n),1)/(TT-n);
sum_vVY_q(n,5) = sum(vVY_2_q4(:,n),1)/(TT-n);

end
for s=1:4
var_share(:,s) = sum_vVY_q(:,s+1)./sum_vVY_q(:,1);
end
end

