%  Shiller_data.m 
%  Properties of returns and consumption growth using Shiller's long data set.   
%  Backus, Chernov, and Martin, "Equity index options"
%  Data read from data sheet constructed from Shiller, headings at top
%  NB:  in his data, prices and rates are start of period, returns and growth rates
%  are forward looking (see monthly version for comparison) 
%  Written:  February 2009 and after
format compact
clear all
close all

disp(' ')
disp('Properties of Shiller data') 
disp('---------------------------------------------------------------')

%  1. Data input 

data = xlsread('Shiller_data.xls','DKB_data_new');

%return

%  2. Statistics 

[nobs,nvar] = size(data);
year = data(:,1);
logr1 = log(data(:,8));
logre = data(:,14);
logxr = logre - logr1;
xr = exp(logre) - exp(logr1);
cons = data(:,9);
consgr = 0*cons;
consgr(1:nobs-1) = cons(2:nobs)./cons(1:nobs-1);
logg = log(consgr);

date1 = find(year==1889);
date2 = find(year==1986);
date3 = find(year==2009);
longdates = [date1:date3]';
shortdates = [date2:date3]';

newdata = [logg logr1 logre logxr xr year/10000];
data_long = newdata(longdates,:);
data_short = newdata(shortdates,:);
year_long = year(longdates,:);

%disp('logg logr1 logre logxr xr year')
%newdata

%return 

disp('Variables (logg logr1 logre logxr xr year)')
disp('Stats (mean, std, skew, kurt, auto)')
disp(' ')
disp('Long sample') 
stats_long = [mean(data_long); std(data_long); skewness(data_long); kurtosis(data_long)-3; acf(data_long,1)]' 

disp('Short sample') 
stats_short = [mean(data_short); std(data_short); skewness(data_short); kurtosis(data_short)-3; acf(data_short,1)]' 

%return 

%  3. Relation between excess return and consumption growth 

FontSize = 14;
FontName = 'Helvetica';  % or 'Times' 
LineWidth = 2;

figure(1) 
plot(data_long(:,1),data_long(:,4),'LineWidth',1,'Color','k','LineStyle','o')
hold on
plot(data_long(:,1),data_long(:,4),'LineWidth',LineWidth,'Color','k','LineStyle','*')
set(gca,'LineWidth',LineWidth,'FontSize',FontSize,'FontName',FontName)
axis([-0.125 0.125 -0.7 0.7]) 
line([0; 0],[-0.7, 0.7])
line([-0.125;0.125],[0; 0])
%text(-0.092,-0.58,'1931','FontSize',FontSize,'FontName',FontName)
xlabel('Log Consumption Growth','FontSize',FontSize,'FontName',FontName)
ylabel('Log Excess Return on Equity','FontSize',FontSize,'FontName',FontName)

%print -depsc scatter_gxr.eps
%return 

y = data_long(:,4);
X = [ones([size(data_long(:,1))]) data_long(:,1)];

disp('Regression parameters (xs return on c gr)') 
lambda = inv(X'*X)*X'*y
%Px = X*(X'*X)*X';
res = y - X*lambda;
disp('Moments of residuals') 
[mean(res) std(res) skewness(res) kurtosis(res)-3] 
disp('Correlation') 
corr = corrcoef(data_long(:,1),data_long(:,4));
corr(1,2)

%disp('Data (year, logg, logxr)')
%[year_long/10000 data_long(:,1) data_long(:,4)]


return