clear all
close all

x = load('us.txt');
t = x(:,1);
gdp = x(:,2);
c = x(:,3);
cpi = x(:,4);
ir = x(:,5);
w = x(:,6);
inv = x(:,7);

lgdp = log(gdp);
lc = log(c);
lcpi = log(cpi);
ir = ir/100;
lw = log(w);
linv = log(inv);

[lgdp_trend,lgdp_gap] = hp(lgdp,1600);
[lc_trend,lc_gap] = hp(lc,1600);
[linv_trend,linv_gap] = hp(linv,1600);

figure
plot(t,[lgdp_gap, lc_gap])
legend('mezera HDP','mezera C')

figure
plot(t,[lgdp_gap, linv_gap])
legend('mezera HDP','mezera INV')

c1 = corrfn(lgdp_gap,lc_gap,5);
c2 = corrfn(lgdp_gap,linv_gap,5);

k = (5:-1:-5)';

[k c1];
[k c2];

figure
bar(k,c1)
xlabel('k')
ylabel('korelacni koef.')
title('HDP_t vs. C_{t+k}')

figure
bar(k,c2)
xlabel('k')
ylabel('korelacni koef.')
title('HDP_t vs. INV_{t+k}')

cy_ratio = c ./ gdp;
iy_ratio = inv./ gdp;
prum_cy = mean(cy_ratio);
prum_iy = mean(iy_ratio);

std_y = std(lgdp_gap);
std_c = std(lc_gap);
std_inv = std(linv_gap);

disp('Podil na HDP: C          I')
disp([prum_cy   prum_iy])
disp('Smerodatna odchylka: Y     C      I')
disp([std_y   std_c  std_inv]*100)

%%  rw, pie4
lrw = lw - lcpi;
%  rw = w ./ cpi; lrw = log(rw);   alternativni vypocet

pie4 = lcpi(5:end) - lcpi(1:end-4);
pie4 = [nan(4,1);pie4];

[lrw_trend,lrw_gap] = hp(lrw,1600);
[pie4_trend,pie4_gap] = hp(pie4,1600);

figure
plot(t,[lgdp_gap, lrw_gap])
legend('mezera HDP','mezera RW')

figure
plot(t,[lgdp_gap, pie4_gap])
legend('mezera HDP','mezera \pi')

c3 = corrfn(lgdp_gap,lrw_gap,5);
c4 = corrfn(lgdp_gap,pie4_gap,5);

figure
bar(k,c3)
xlabel('k')
ylabel('korelacni koef.')
title('HDP_t x RW_{t+k}')

figure
bar(k,c4)
xlabel('k')
ylabel('korelacni koef.')
title('HDP_t x \pi_{t+k}')

%% cy_ratio, iy_ratio

[cy_ratio_trend,cy_ratio_gap] = hp(cy_ratio,1600);
[iy_ratio_trend,iy_ratio_gap] = hp(iy_ratio,1600);

c5 = corrfn(lgdp_gap,cy_ratio_gap,0);
c6 = corrfn(lgdp_gap,iy_ratio_gap,0);

figure
plot(t,[lgdp_gap, cy_ratio_gap])
legend('mezera HDP','mezera C/Y')

figure
plot(t,[lgdp_gap, iy_ratio_gap])
legend('mezera HDP','mezera I/Y')

%% rir
rir = ir(1:end-4) - pie4(5:end);
rir = [rir;nan(4,1)];

[rir_trend, rir_gap] = hp(rir,1600);
c7 = corrfn(lgdp_gap,rir_gap,5);

figure
plot(t,[lgdp_gap, rir_gap])
legend('mezera HDP', 'mezera RIR')

figure
bar(k,c7)
xlabel('k')
ylabel('korelacni koeficient')
title('HDP_t x RIR_{t+k}')

%% test vyznamnosti korelacniho koeficientu

n = length(t);

test.c = [c1 k corrtest(c1,n)];
test.i = [c2 k corrtest(c2,n)];
test.rw =[c3 k corrtest(c3,n)];
test.pie = [c4 k corrtest(c4,n-4)];
test.rir =[c7 k corrtest(c7,n-4)];



























