%!!!open in "current folder" the folder where the data is stored!!!
%otherwise, the code will not work!!!
%each section of code is independent and could be run separately
%% 1986-2018 + diagnostics 
clear all
%Import return data
opts = spreadsheetImportOptions("NumVariables", 103);
opts.Sheet = "ri (stock data)";
opts.DataRange = "A2:CY397";
opts.VariableNames = ["Date", "MktRF", "RF", "ADAMSEXPRESSCO", "ABBOTTLABORATORIES", "FIRSTFINANCIALBANCORPOHIO", "FIRSTMIDWESTBANCORPDE", "ROGERSCORP", "ROLLINSINC", "SEABOARDCORP", "STSOURCECORP", "SHERWINWILLIAMSCO", "FLEXSTEELINDUSTRIESINC", "FOSTERLBCO", "FRANKLINELECTRICINC", "FRANKLINRESOURCESINC", "FULLERHBCO", "GENTEXCORP", "DONNELLEYRRSONSCO", "AMREPCORP", "STANDARDMOTORPRODUCTSINC", "APACHECORP", "MATTELINC", "BECTONDICKINSONCO", "GRACOINC", "BARNESGROUPINC", "WEYERHAEUSERCO", "INTERNATIONALFLAVORSFRAGINC", "AVONPRODUCTSINC", "HAWAIIANELECTRICINDUSTRIESINC", "HAVERTYFURNITURECOSINC", "HEALTHCARESERVICESGROUPINC", "PARKERHANNIFINCORP", "WEISMARKETSINC", "WOLVERINEWORLDWIDEINC", "BARNWELLINDUSTRIESINC", "SMUCKERJMCO", "MOOGINC", "MEREDITHCORP", "HUNTJBTRANSPORTSERVICESINC", "HUNTINGTONBANCSHARESINC", "HURCOCOMPANY", "SYNALLOYCORP", "MCDONALDSCORP", "STANDEXINTERNATIONALCORP", "VFCORP", "IMMUNOMEDICSINC", "BEMISCOINC", "KENNAMETALINC", "PFINDUSTRIESINC", "TELEFLEXINC", "INTEGRATEDDEVICETECHNOLOGY", "AUTOMATICDATAPROCESSINGINC", "DIODESINC", "INVACARECORP", "MARSHMCLENNANCOSINC", "ISRAMCOINC", "SAFEGUARDSCIENTIFICSINC", "CLOROXCO", "GENUINEPARTSCO", "CABOTCORP", "RITEAIDCORP", "KAMANCORP", "KELLYSERVICESINC", "NEWYORKTIMESCO", "KEYTRONICSCORP", "CNAFINANCIALCORP", "FEDERALSIGNALCORP", "KIMBALLINTERNATIONALINC", "KULICKESOFFAINDSINC", "LAMRESHCORP", "LANCASTERCOLONYCORP", "HUMANAINC", "UNIONPACIFICCORP", "LAWSONPRODUCTSINC", "LINCOLNNATIONALCORP", "SIFCOINDUSTRIESINC", "ALEXANDERSINC", "BLOCKHRINC", "DANAHERCORP", "SERVOTRONICSINC", "LEEENTERPRISESINC", "MTSSYSTEMSCORP", "PENNSYLVANIAREALESTATEINVTTR", "LILLYELICO", "SEMTECHCORP", "WINNEBAGOINDUSTRIESINC", "SONYCORP", "MANITOWOCCOINC", "TERADYNEINC", "MARCUSCORP", "MARINEPETROLEUMTRUST", "SERVICECORPINTERNATIONAL", "SYSCOCORP", "MCCORMICKCOINC", "JACOBSENGINEERINGGROUPINC", "SUPERIORINDUSTRIESINTLINC", "EQUIFAXINC", "COCACOLABOTTLINGCOCONS", "DYCOMINDUSTRIESINC", "GATXCORP", "GENERALDYNAMICSCORP", "GENERALELECTRICCO"];
opts.VariableTypes = ["string", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double"];
opts = setvaropts(opts, 1, "WhitespaceRule", "preserve");
opts = setvaropts(opts, 1, "EmptyFieldRule", "auto");
data = readtable("data.xlsx", opts, "UseExcel", false);
clear opts

%test for normality
stockdatasorted100=data
return_asset=stockdatasorted100(:,4:103);
dataS1 =stockdatasorted100(:,2);
n=width(return_asset);
jb=zeros(1,n);
pval=zeros(1,n);
skew=zeros(1,n);
kurt=zeros(1,n);
t_stat=zeros(1,n);
T=396
for i=(1:n)
a=return_asset(:,i);
normality=[dataS1 a];
%change name in column 2
normality.Properties.VariableNames = {'ExcessMarket' 'StockReturn'};
[h, pvalue, jbstat, critval] = jbtest(normality.StockReturn)
jb(i)=jbstat
pval(i)=pvalue
skew(i)=skewness(normality.StockReturn);
kurt(i)=kurtosis(normality.StockReturn);
t_stat(i)=T *( skew(i)^2 / 6 + (kurt(i) - 3)^2 / 24 );
end
reject_rule= pval<0.05;
percent_of_abnormal_data=sum(reject_rule)/n

%import excess return data
opts = spreadsheetImportOptions("NumVariables", 102);
opts.Sheet = "ri-rf";
opts.DataRange = "A2:CX397";
opts.VariableNames = ["Date", "MktRF", "ADAMSEXPRESSCO", "ABBOTTLABORATORIES", "FIRSTFINANCIALBANCORPOHIO", "FIRSTMIDWESTBANCORPDE", "ROGERSCORP", "ROLLINSINC", "SEABOARDCORP", "STSOURCECORP", "SHERWINWILLIAMSCO", "FLEXSTEELINDUSTRIESINC", "FOSTERLBCO", "FRANKLINELECTRICINC", "FRANKLINRESOURCESINC", "FULLERHBCO", "GENTEXCORP", "DONNELLEYRRSONSCO", "AMREPCORP", "STANDARDMOTORPRODUCTSINC", "APACHECORP", "MATTELINC", "BECTONDICKINSONCO", "GRACOINC", "BARNESGROUPINC", "WEYERHAEUSERCO", "INTERNATIONALFLAVORSFRAGINC", "AVONPRODUCTSINC", "HAWAIIANELECTRICINDUSTRIESINC", "HAVERTYFURNITURECOSINC", "HEALTHCARESERVICESGROUPINC", "PARKERHANNIFINCORP", "WEISMARKETSINC", "WOLVERINEWORLDWIDEINC", "BARNWELLINDUSTRIESINC", "SMUCKERJMCO", "MOOGINC", "MEREDITHCORP", "HUNTJBTRANSPORTSERVICESINC", "HUNTINGTONBANCSHARESINC", "HURCOCOMPANY", "SYNALLOYCORP", "MCDONALDSCORP", "STANDEXINTERNATIONALCORP", "VFCORP", "IMMUNOMEDICSINC", "BEMISCOINC", "KENNAMETALINC", "PFINDUSTRIESINC", "TELEFLEXINC", "INTEGRATEDDEVICETECHNOLOGY", "AUTOMATICDATAPROCESSINGINC", "DIODESINC", "INVACARECORP", "MARSHMCLENNANCOSINC", "ISRAMCOINC", "SAFEGUARDSCIENTIFICSINC", "CLOROXCO", "GENUINEPARTSCO", "CABOTCORP", "RITEAIDCORP", "KAMANCORP", "KELLYSERVICESINC", "NEWYORKTIMESCO", "KEYTRONICSCORP", "CNAFINANCIALCORP", "FEDERALSIGNALCORP", "KIMBALLINTERNATIONALINC", "KULICKESOFFAINDSINC", "LAMRESHCORP", "LANCASTERCOLONYCORP", "HUMANAINC", "UNIONPACIFICCORP", "LAWSONPRODUCTSINC", "LINCOLNNATIONALCORP", "SIFCOINDUSTRIESINC", "ALEXANDERSINC", "BLOCKHRINC", "DANAHERCORP", "SERVOTRONICSINC", "LEEENTERPRISESINC", "MTSSYSTEMSCORP", "PENNSYLVANIAREALESTATEINVTTR", "LILLYELICO", "SEMTECHCORP", "WINNEBAGOINDUSTRIESINC", "SONYCORP", "MANITOWOCCOINC", "TERADYNEINC", "MARCUSCORP", "MARINEPETROLEUMTRUST", "SERVICECORPINTERNATIONAL", "SYSCOCORP", "MCCORMICKCOINC", "JACOBSENGINEERINGGROUPINC", "SUPERIORINDUSTRIESINTLINC", "EQUIFAXINC", "COCACOLABOTTLINGCOCONS", "DYCOMINDUSTRIESINC", "GATXCORP", "GENERALDYNAMICSCORP", "GENERALELECTRICCO"];
opts.VariableTypes = ["string", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double"];
opts = setvaropts(opts, 1, "WhitespaceRule", "preserve");
opts = setvaropts(opts, 1, "EmptyFieldRule", "auto");
dataS1 = readtable("data.xlsx", opts, "UseExcel", false);
clear opts

%regression
%step 1
stockdatasorted100S1=dataS1
excess_asset=stockdatasorted100S1(:,3:102);
market_excess=stockdatasorted100S1(:,2);
n=width(excess_asset);
beta=zeros(1,n);
p=zeros(1,n);
alpha=zeros(1,n);
alpha_p=zeros(1,n);
for i=(1:n);
a=excess_asset(:,i);
regression=[market_excess a];
regression.Properties.VariableNames = {'ExcessMarket' 'ExcessStock'};
modelFit = fitlm(regression,'ExcessStock ~ ExcessMarket')
beta(i) = modelFit.Coefficients.Estimate(2);
p(i)= modelFit.Coefficients.pValue(2);
alpha(i)=modelFit.Coefficients.Estimate(1);
alpha_p(i)=modelFit.Coefficients.pValue(1);
end
output =[beta;p]
reject_rule1= p<0.05;
percent_of_significant_coeff=sum(reject_rule1)/n

%import average excess return of companies
opts = spreadsheetImportOptions("NumVariables", 1);
opts.Sheet = "average_excess_return";
opts.DataRange = "B2:B101";
opts.VariableNames = "average_excess_return";
opts.VariableTypes = "double";
dataS2 = readtable("data.xlsx", opts, "UseExcel", false);
clear opts

%step 2
beta=beta';
dataS2.beta = beta;
regression2=dataS2;
modelFit1 = fitlm(regression2,'average_excess_return ~ beta')

%Expended CAPM
%step 1 
market_excess.market_excess_sq=market_excess.MktRF.^2
market_excess.market_excess_cubus=market_excess.MktRF.^3
beta1=zeros(1,n);
beta2=zeros(1,n);
beta3=zeros(1,n);

for i=(1:n);
a=excess_asset(:,i);
regression3=[market_excess a];
regression3.Properties.VariableNames = {'ExcessMarket' 'ExcessMarket2' 'ExcessMarket3' 'ExcessStock'};
modelFit3 = fitlm(regression3,'ExcessStock ~ ExcessMarket+ExcessMarket2+ExcessMarket3')
beta1(i) = modelFit3.Coefficients.Estimate(2);
beta2(i)=modelFit3.Coefficients.Estimate(3);
beta3(i)=modelFit3.Coefficients.Estimate(4);
end

%step 2
beta1=beta1';
beta2=beta2';
beta3=beta3';
dataS2.beta1 = beta1;
dataS2.beta2 = beta2;
dataS2.beta3 = beta3;
regression4=dataS2;
model = fitlm(regression4,'average_excess_return ~ beta1+beta2+beta3')

% Tests
% 1.expected mean of residuals=0
uhat = model.Residuals.Raw;
T=model.NumObservations;
k=model.NumEstimatedCoefficients;
resid_mean = mean(uhat)
resid_std = std(uhat)
%H0:resid_mean=0 vs H1: resid_mean<>0
tstat_resid=(resid_mean-0)/resid_std
%tcrit with T-2 d.o.f. for two-sided test
alpha=.05
tcrit=tinv(1-1/2*alpha, T-2)
reject = tstat_resid>tcrit
%our case: cannot reject
% assumption not violated
p = 2*(1-tcdf(abs(tstat_resid),T-k))
%p-value is 1 - cannot reject H0, assumption not violated

% 2. White test
% Export residuals from the first regression
u11 = model.Residuals.Raw;
% Calculate squared residuals
r = u11.^2
% Create data for auxiliary regression
data2=table(r, regression4.beta1, regression4.beta2, regression4.beta3, regression4.beta1.^2, regression4.beta2.^2, regression4.beta3.^2, 'VariableNames',{'u2','beta1','beta2','beta3','beta1sq','beta2sq','beta3sq'})
% Run auxiliary regression with level and squared
unr_regr=fitlm(data2, 'u2~beta1+beta2+beta3+beta1sq+beta2sq+beta3sq')
% Export residual sum of squares from unrestricted regression
URSS = unr_regr.SSE;
% Export the number of observations
T = unr_regr.NumObservations;
% Export the number of estimated coefficients
k = unr_regr.NumEstimatedCoefficients
% Test H0: homoscedasticity (H0: a2=0 and a3=0 ... and a7=0 vs. H1: a1
% not equal to 0 or ... a7 not equal to 0)
c=zeros(7,1);
H = [ 0 0 0 0 0 0 0 ;
0 1 0 0 0 0 0 ;
0 0 1 0 0 0 0 ;
0 0 0 1 0 0 0 ;
0 0 0 0 1 0 0 ;
0 0 0 0 0 1 0 ;
0 0 0 0 0 0 1 ];
[p1, F1] = coefTest(unr_regr, H, c)
alpha = [0.05];
data2.ones = ones(height(data2), 1);
% Run the restricted regression
res_regr = fitlm(data2, 'u2 ~ -1 + ones');
% Calculate the number of regressors in the auxiliary regression excluding the constant term
m = k - 1;
% Export residual sum of squares from restricted regression
RRSS = res_regr.SSE;
% Calculate F-stat
FStat = (RRSS - URSS)/URSS * (T-k)/m
% Calculate F-critical
FCrit = finv(1-alpha, m, T - k)
% Calculate p-value
p = 1 - fcdf(FStat, m, T - k)
% The condition of rejection H0
reject = FStat>FCrit
% 1: condition is true, i.e., reject H0
% Calculate T*r^2~chi^2(m)
t_stat = unr_regr.Rsquared.Ordinary * T
% Calcualte critical value
crit_val = chi2inv(1-alpha, m)
% Calculate p-value
p = 1 - cdf('Chisquare', t_stat, m)
% The condition of rejection H0
reject1 = t_stat > crit_val

% Correction model
% Create the matrix of independent variables
x = [regression4.beta1, regression4.beta2, regression4.beta3];
% Create dependent variable
y = regression4.average_excess_return;
% Export the number of observations
T = model.NumObservations;
% Estimate regression model with White's heteroscedasticity consistent standard errors
% assumption: estimated residuals are heteroscedastic, but NOT autocorrelated
white = table();
const = 1; % 1: constant is included, 0: constant is not included
nwlags = 0; % Newey-West lags
% you have to install the MFE toolbox!!!
[beta, tstat, ~, vcvm, ~,Rbar, ~ ] = olsnw(y, x, const, nwlags)
white.beta = beta;
white.SE = sqrt(diag(vcvm));
white.tstat = tstat;
white.p = 2*(1-tcdf(abs(tstat),T-2));
white.Properties.RowNames = {'(Intercept)', 'beta1', 'beta2', 'beta3' };
white

% 3 Autocorrelation between error terms
%H0: zero autocorrelation between error terms Cov(uCC;i; uCC;j) = 0 for i \neq j
%H1: not zero autocorrelation between error terms
u= model.Residuals.Raw
%Breusch-Godfrey Serial Correlation LM Test with r=4 lags
regression4.u = u;
regression4.u1 = lagmatrix(u, 1);
regression4.u2 = lagmatrix(u, 2);
regression4.u3 = lagmatrix(u, 3);
regression4.u4 = lagmatrix(u, 4);
% Setting all nan's to zero
regression4.u1(isnan(regression4.u1))=0
regression4.u2(isnan(regression4.u2))=0
regression4.u3(isnan(regression4.u3))=0
regression4.u4(isnan(regression4.u4))=0
% running auxiliary regression
aux_regr = fitlm(regression4, ['u ~ beta1 + beta2 + beta3+ u1 + u2 + u3 + u4 '])
alpha = [0.10, 0.05, 0.01];
r = 4;
T= model.NumObservations;
% Calculate T*r^2~chi^2(m)
% test_statistic
test_stat = aux_regr.Rsquared.Ordinary * (T - r)
crit_val = chi2inv(1-alpha, r)
p = 1 - cdf('Chisquare', test_stat, r)
reject = test_stat > crit_val
% 0: condition is not true, i.e., do not reject H0
% Calculate F-stat
c = zeros(8,1); % c and H specify the null hypothesis in the form H*b = c, where b ...
H = zeros(8); % ... is the vector of unknown parameters estimated by beta
H(5,5) = 1;
H(6,6) = 1;
H(7,7) = 1;
H(8,8) = 1;
[p1, F1] = coefTest(aux_regr, H, c)
%there is no autocorrelation of residuals

% 4. there is no relation between the residuals and the independent variables,
corr_beta1= corr( model.Residuals.Raw , model.Variables.beta1);
corr_beta2=corr( model.Residuals.Raw , model.Variables.beta2);
corr_beta3=corr( model.Residuals.Raw , model.Variables.beta3);
[corr_beta1 corr_beta2 corr_beta3 ]

% 5 Testing the normality assumption
[h, p, jbstat, critval] = jbtest(u)
alpha = [0.10, 0.05, 0.01];
skew = mean(u.^3) / mean(u.^2).^(3/2)
kurt = mean(u.^4) / mean(u.^2).^2
test_stat = T * ( skew^2 / 6 + (kurt - 3)^2 / 24 )
crit_val = chi2inv(1-alpha, 2)
p = 1 - cdf('Chisquare', test_stat, 2)
reject = test_stat > crit_val
% 1: condition is true, i.e., reject H0
% 0: condition is not true, i.e., do not reject H0

% graphs
%Time series plot of estimated residuals
figure(1);
plot( regression4.u, 'b');
hold on;
xlabel('Period'); % adds a label to x-axis
ylabel('Estimated residuals'); % adds a label to y-axis
title ('Time series plot of the estimated residuals')
%Histogram of estimated residuals
figure(2);
histogram(regression4.u);
hold on;
ylabel('Estimated residuals'); % adds a label to y-axis
title ('Histogram of the estimated residuals')

% Multicolinearity problem
X=[ model.Variables.beta1 model.Variables.beta2 model.Variables.beta3]
%correlation matrix:
R=corrcoef(X)

%% 1986-2018 with interaction term + diagnostics 
clear all
%import excess return data
opts = spreadsheetImportOptions("NumVariables", 102);
opts.Sheet = "ri-rf";
opts.DataRange = "A2:CX397";
opts.VariableNames = ["Date", "MktRF", "ADAMSEXPRESSCO", "ABBOTTLABORATORIES", "FIRSTFINANCIALBANCORPOHIO", "FIRSTMIDWESTBANCORPDE", "ROGERSCORP", "ROLLINSINC", "SEABOARDCORP", "STSOURCECORP", "SHERWINWILLIAMSCO", "FLEXSTEELINDUSTRIESINC", "FOSTERLBCO", "FRANKLINELECTRICINC", "FRANKLINRESOURCESINC", "FULLERHBCO", "GENTEXCORP", "DONNELLEYRRSONSCO", "AMREPCORP", "STANDARDMOTORPRODUCTSINC", "APACHECORP", "MATTELINC", "BECTONDICKINSONCO", "GRACOINC", "BARNESGROUPINC", "WEYERHAEUSERCO", "INTERNATIONALFLAVORSFRAGINC", "AVONPRODUCTSINC", "HAWAIIANELECTRICINDUSTRIESINC", "HAVERTYFURNITURECOSINC", "HEALTHCARESERVICESGROUPINC", "PARKERHANNIFINCORP", "WEISMARKETSINC", "WOLVERINEWORLDWIDEINC", "BARNWELLINDUSTRIESINC", "SMUCKERJMCO", "MOOGINC", "MEREDITHCORP", "HUNTJBTRANSPORTSERVICESINC", "HUNTINGTONBANCSHARESINC", "HURCOCOMPANY", "SYNALLOYCORP", "MCDONALDSCORP", "STANDEXINTERNATIONALCORP", "VFCORP", "IMMUNOMEDICSINC", "BEMISCOINC", "KENNAMETALINC", "PFINDUSTRIESINC", "TELEFLEXINC", "INTEGRATEDDEVICETECHNOLOGY", "AUTOMATICDATAPROCESSINGINC", "DIODESINC", "INVACARECORP", "MARSHMCLENNANCOSINC", "ISRAMCOINC", "SAFEGUARDSCIENTIFICSINC", "CLOROXCO", "GENUINEPARTSCO", "CABOTCORP", "RITEAIDCORP", "KAMANCORP", "KELLYSERVICESINC", "NEWYORKTIMESCO", "KEYTRONICSCORP", "CNAFINANCIALCORP", "FEDERALSIGNALCORP", "KIMBALLINTERNATIONALINC", "KULICKESOFFAINDSINC", "LAMRESHCORP", "LANCASTERCOLONYCORP", "HUMANAINC", "UNIONPACIFICCORP", "LAWSONPRODUCTSINC", "LINCOLNNATIONALCORP", "SIFCOINDUSTRIESINC", "ALEXANDERSINC", "BLOCKHRINC", "DANAHERCORP", "SERVOTRONICSINC", "LEEENTERPRISESINC", "MTSSYSTEMSCORP", "PENNSYLVANIAREALESTATEINVTTR", "LILLYELICO", "SEMTECHCORP", "WINNEBAGOINDUSTRIESINC", "SONYCORP", "MANITOWOCCOINC", "TERADYNEINC", "MARCUSCORP", "MARINEPETROLEUMTRUST", "SERVICECORPINTERNATIONAL", "SYSCOCORP", "MCCORMICKCOINC", "JACOBSENGINEERINGGROUPINC", "SUPERIORINDUSTRIESINTLINC", "EQUIFAXINC", "COCACOLABOTTLINGCOCONS", "DYCOMINDUSTRIESINC", "GATXCORP", "GENERALDYNAMICSCORP", "GENERALELECTRICCO"];
opts.VariableTypes = ["string", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double"];
opts = setvaropts(opts, 1, "WhitespaceRule", "preserve");
opts = setvaropts(opts, 1, "EmptyFieldRule", "auto");
dataS1 = readtable("data.xlsx", opts, "UseExcel", false);
clear opts

%regression
%step 1
stockdatasorted100S1=dataS1
excess_asset=stockdatasorted100S1(:,3:102);
market_excess=stockdatasorted100S1(:,2);
n=width(excess_asset);

%import average excess return of companies
opts = spreadsheetImportOptions("NumVariables", 1);
opts.Sheet = "average_excess_return";
opts.DataRange = "B2:B101";
opts.VariableNames = "average_excess_return";
opts.VariableTypes = "double";
dataS2 = readtable("data.xlsx", opts, "UseExcel", false);
clear opts

%Expended CAPM
%step 1 
market_excess.market_excess_sq=market_excess.MktRF.^2
market_excess.market_excess_cubus=market_excess.MktRF.^3
market_excess.market_excess_interaction=market_excess.market_excess_sq.*market_excess.market_excess_cubus
beta1=zeros(1,n);
beta2=zeros(1,n);
beta3=zeros(1,n);
beta4=zeros(1,n);
for i=(1:n);
a=excess_asset(:,i);
regression3=[market_excess a];
regression3.Properties.VariableNames = {'ExcessMarket' 'ExcessMarket2' 'ExcessMarket3' 'InteractionTerm' 'ExcessStock'};

modelFit3 = fitlm(regression3,'ExcessStock ~ ExcessMarket+ExcessMarket2+ExcessMarket3+InteractionTerm')
beta1(i) = modelFit3.Coefficients.Estimate(2);
beta2(i)=modelFit3.Coefficients.Estimate(3);
beta3(i)=modelFit3.Coefficients.Estimate(4);
beta4(i)=modelFit3.Coefficients.Estimate(5);
end

%step 2
beta1=beta1';
beta2=beta2';
beta3=beta3';
beta4=beta4';
dataS2.beta1 = beta1;
dataS2.beta2 = beta2;
dataS2.beta3 = beta3;
dataS2.beta4 = beta4;
regression4=dataS2;
model = fitlm(regression4,'average_excess_return ~ beta1+beta2+beta3+beta4')

% Tests
% 1.expected mean of residuals=0
uhat = model.Residuals.Raw;
T=model.NumObservations;
k=model.NumEstimatedCoefficients;
resid_mean = mean(uhat)
resid_std = std(uhat)
%H0:resid_mean=0 vs H1: resid_mean<>0
tstat_resid=(resid_mean-0)/resid_std
%tcrit with T-2 d.o.f. for two-sided test
alpha=.05
tcrit=tinv(1-1/2*alpha, T-2)
reject = tstat_resid>tcrit
%our case: cannot reject
% assumption not violated
p = 2*(1-tcdf(abs(tstat_resid),T-k))
%p-value is 1 - cannot reject H0, assumption not violated

%2  White test
% Export residuals from the first regression
u11 = model.Residuals.Raw;
% Calculate squared residuals
res.u2 = u11.^2
% Create data for auxiliary regression
data2=table(res.u2, regression4.beta1, regression4.beta2, regression4.beta3, regression4.beta4, regression4.beta1.^2, regression4.beta2.^2, regression4.beta3.^2, regression4.beta4.^2, 'VariableNames',{'u2','beta1','beta2','beta3','beta4','beta1sq','beta2sq','beta3sq', 'beta4sq'})
% Run auxiliary regression with level and squared
unr_regr=fitlm(data2, 'u2~beta1+beta2+beta3+beta4+beta1sq+beta2sq+beta3sq+beta4sq')
% Export residual sum of squares from unrestricted regression
URSS = unr_regr.SSE;
% Export the number of observations
T = unr_regr.NumObservations;
% Export the number of estimated coefficients
k = unr_regr.NumEstimatedCoefficients
% Test H0: homoscedasticity (H0: a2=0 and a3=0 ... and a7=0 vs. H1: a1
% not equal to 0 or ... a7 not equal to 0)
c=zeros(9,1);
H = [ 0 0 0 0 0 0 0 0 0;
0 1 0 0 0 0 0 0 0;
0 0 1 0 0 0 0 0 0;
0 0 0 1 0 0 0 0 0;
0 0 0 0 1 0 0 0 0;
0 0 0 0 0 1 0 0 0;
0 0 0 0 0 0 1 0 0;
0 0 0 0 0 0 0 1 0;
0 0 0 0 0 0 0 0 1];
[p1, F1] = coefTest(unr_regr, H, c)
alpha = [0.05];
data2.ones = ones(height(data2), 1);
% Run the restricted regression
res_regr = fitlm(data2, 'u2 ~ -1 + ones');
% Calculate the number of regressors in the auxiliary regression excluding the constant term
m = k - 1;
% Export residual sum of squares from restricted regression
RRSS = res_regr.SSE;
% Calculate F-stat
FStat = (RRSS - URSS)/URSS * (T-k)/m
% Calculate F-critical
FCrit = finv(1-alpha, m, T - k)
% Calculate p-value
p = 1 - fcdf(FStat, m, T - k)
% The condition of rejection H0
reject = FStat>FCrit
% 1: condition is true, i.e., reject H0
% Calculate T*r^2~chi^2(m)
t_stat = unr_regr.Rsquared.Ordinary * T
% Calcualte critical value
crit_val = chi2inv(1-alpha, m)
% Calculate p-value
p = 1 - cdf('Chisquare', t_stat, m)
% The condition of rejection H0
reject1 = t_stat > crit_val

% Correction of heteroscedasticity
% Create the matrix of independent variables
x = [regression4.beta1, regression4.beta2, regression4.beta3 regression4.beta4];
% Create dependent variable
y = regression4.average_excess_return;
% Export the number of observations
T = model.NumObservations;
% Estimate regression model with White's heteroscedasticity consistent standard errors
% assumption: estimated residuals are heteroscedastic, but NOT autocorrelated
white = table();
const = 1; % 1: constant is included, 0: constant is not included
nwlags = 0; % Newey-West lags
% you have to install the MFE toolbox!!!
[beta, tstat, ~, vcvm, ~,Rbar, ~ ] = olsnw(y, x, const, nwlags)
white.beta = beta;
white.SE = sqrt(diag(vcvm));
white.tstat = tstat;
white.p = 2*(1-tcdf(abs(tstat),T-2));
white.Properties.RowNames = {'(Intercept)', 'beta1', 'beta2', 'beta3','beta4' };
white

% 3 Autocorrelation between error terms 
%H0: zero autocorrelation between error terms Cov(uCC;i; uCC;j) = 0 for i \neq j
%H1: not zero autocorrelation between error terms
u= model.Residuals.Raw
%Breusch-Godfrey Serial Correlation LM Test with r=4 lags
regression4.u = u;
regression4.u1 = lagmatrix(u, 1);
regression4.u2 = lagmatrix(u, 2);
regression4.u3 = lagmatrix(u, 3);
regression4.u4 = lagmatrix(u, 4);
% Setting all nan's to zero
regression4.u1(isnan(regression4.u1))=0
regression4.u2(isnan(regression4.u2))=0
regression4.u3(isnan(regression4.u3))=0
regression4.u4(isnan(regression4.u4))=0
% running auxiliary regression
aux_regr = fitlm(regression4, ['u ~ beta1 + beta2 + beta3+beta4+ u1 + u2 + u3 + u4 '])
alpha = [0.10, 0.05, 0.01];
r = 4;
T= model.NumObservations;
% Calculate T*r^2~chi^2(m)
% test_statistic
test_stat = aux_regr.Rsquared.Ordinary * (T - r)
crit_val = chi2inv(1-alpha, r)
p = 1 - cdf('Chisquare', test_stat, r)
reject = test_stat > crit_val
% 0: condition is not true, i.e., do not reject H0
% Calculate F-stat
c = zeros(9,1); % c and H specify the null hypothesis in the form H*b = c, where b ...
H = zeros(9); % ... is the vector of unknown parameters estimated by beta
H(6,6) = 1;
H(7,7) = 1;
H(8,8) = 1;
H(9,9) = 1;
[p1, F1] = coefTest(aux_regr, H, c)
%there is no autocorrelation of residuals

% 4. there is no relation between the residuals and the independent variables,
corr_beta1= corr( model.Residuals.Raw , model.Variables.beta1);
corr_beta2=corr( model.Residuals.Raw , model.Variables.beta2);
corr_beta3=corr( model.Residuals.Raw , model.Variables.beta3);
corr_beta4=corr( model.Residuals.Raw , model.Variables.beta4);
result = [corr_beta1 corr_beta2 corr_beta3 corr_beta4]

% 5 Testing the normality assumption
[h, p, jbstat, critval] = jbtest(u)
alpha = [0.10, 0.05, 0.01];
skew = mean(u.^3) / mean(u.^2).^(3/2)
kurt = mean(u.^4) / mean(u.^2).^2
test_stat = T * ( skew^2 / 6 + (kurt - 3)^2 / 24 )
crit_val = chi2inv(1-alpha, 2)
p = 1 - cdf('Chisquare', test_stat, 2)
reject = test_stat > crit_val
% 1: condition is true, i.e., reject H0
% 0: condition is not true, i.e., do not reject H0

% graphs
%Time series plot of estimated residuals
figure(1);
plot( regression4.u, 'b');
hold on;
xlabel('Period'); % adds a label to x-axis
ylabel('Estimated residuals'); % adds a label to y-axis
title ('Time series plot of the estimated residuals')
%Histogram of estimated residuals
figure(2);
histogram(regression4.u);
hold on;
ylabel('Estimated residuals'); % adds a label to y-axis
title ('Histogram of the estimated residuals')

% Multicolinearity
X=[ model.Variables.beta1 model.Variables.beta2 model.Variables.beta3 model.Variables.beta4]
%correlation matrix:
R=corrcoef(X)

%% 1993-1997
clear all
%Import return data
opts = spreadsheetImportOptions("NumVariables", 103);
opts.Sheet = "ri (stock data)";
opts.DataRange = "A2:CY397";
opts.VariableNames = ["Date", "MktRF", "RF", "ADAMSEXPRESSCO", "ABBOTTLABORATORIES", "FIRSTFINANCIALBANCORPOHIO", "FIRSTMIDWESTBANCORPDE", "ROGERSCORP", "ROLLINSINC", "SEABOARDCORP", "STSOURCECORP", "SHERWINWILLIAMSCO", "FLEXSTEELINDUSTRIESINC", "FOSTERLBCO", "FRANKLINELECTRICINC", "FRANKLINRESOURCESINC", "FULLERHBCO", "GENTEXCORP", "DONNELLEYRRSONSCO", "AMREPCORP", "STANDARDMOTORPRODUCTSINC", "APACHECORP", "MATTELINC", "BECTONDICKINSONCO", "GRACOINC", "BARNESGROUPINC", "WEYERHAEUSERCO", "INTERNATIONALFLAVORSFRAGINC", "AVONPRODUCTSINC", "HAWAIIANELECTRICINDUSTRIESINC", "HAVERTYFURNITURECOSINC", "HEALTHCARESERVICESGROUPINC", "PARKERHANNIFINCORP", "WEISMARKETSINC", "WOLVERINEWORLDWIDEINC", "BARNWELLINDUSTRIESINC", "SMUCKERJMCO", "MOOGINC", "MEREDITHCORP", "HUNTJBTRANSPORTSERVICESINC", "HUNTINGTONBANCSHARESINC", "HURCOCOMPANY", "SYNALLOYCORP", "MCDONALDSCORP", "STANDEXINTERNATIONALCORP", "VFCORP", "IMMUNOMEDICSINC", "BEMISCOINC", "KENNAMETALINC", "PFINDUSTRIESINC", "TELEFLEXINC", "INTEGRATEDDEVICETECHNOLOGY", "AUTOMATICDATAPROCESSINGINC", "DIODESINC", "INVACARECORP", "MARSHMCLENNANCOSINC", "ISRAMCOINC", "SAFEGUARDSCIENTIFICSINC", "CLOROXCO", "GENUINEPARTSCO", "CABOTCORP", "RITEAIDCORP", "KAMANCORP", "KELLYSERVICESINC", "NEWYORKTIMESCO", "KEYTRONICSCORP", "CNAFINANCIALCORP", "FEDERALSIGNALCORP", "KIMBALLINTERNATIONALINC", "KULICKESOFFAINDSINC", "LAMRESHCORP", "LANCASTERCOLONYCORP", "HUMANAINC", "UNIONPACIFICCORP", "LAWSONPRODUCTSINC", "LINCOLNNATIONALCORP", "SIFCOINDUSTRIESINC", "ALEXANDERSINC", "BLOCKHRINC", "DANAHERCORP", "SERVOTRONICSINC", "LEEENTERPRISESINC", "MTSSYSTEMSCORP", "PENNSYLVANIAREALESTATEINVTTR", "LILLYELICO", "SEMTECHCORP", "WINNEBAGOINDUSTRIESINC", "SONYCORP", "MANITOWOCCOINC", "TERADYNEINC", "MARCUSCORP", "MARINEPETROLEUMTRUST", "SERVICECORPINTERNATIONAL", "SYSCOCORP", "MCCORMICKCOINC", "JACOBSENGINEERINGGROUPINC", "SUPERIORINDUSTRIESINTLINC", "EQUIFAXINC", "COCACOLABOTTLINGCOCONS", "DYCOMINDUSTRIESINC", "GATXCORP", "GENERALDYNAMICSCORP", "GENERALELECTRICCO"];
opts.VariableTypes = ["string", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double"];
opts = setvaropts(opts, 1, "WhitespaceRule", "preserve");
opts = setvaropts(opts, 1, "EmptyFieldRule", "auto");
data = readtable("data.xlsx", opts, "UseExcel", false);
clear opts

%test for normality
stockdatasorted100=data(85:144,:);
return_asset=stockdatasorted100(:,4:103);
dataS1 =stockdatasorted100(:,2);
n=width(return_asset);
jb=zeros(1,n);
pval=zeros(1,n);
skew=zeros(1,n);
kurt=zeros(1,n);
t_stat=zeros(1,n);
T=396
for i=(1:n)
a=return_asset(:,i);
normality=[dataS1 a];
%change name in column 2
normality.Properties.VariableNames = {'ExcessMarket' 'StockReturn'};
[h, pvalue, jbstat, critval] = jbtest(normality.StockReturn)
jb(i)=jbstat
pval(i)=pvalue
skew(i)=skewness(normality.StockReturn);
kurt(i)=kurtosis(normality.StockReturn);
t_stat(i)=T *( skew(i)^2 / 6 + (kurt(i) - 3)^2 / 24 );
end
reject_rule= pval<0.05;
percent_of_abnormal_data=sum(reject_rule)/n

%import excess return data
opts = spreadsheetImportOptions("NumVariables", 102);
opts.Sheet = "ri-rf";
opts.DataRange = "A2:CX397";
opts.VariableNames = ["Date", "MktRF", "ADAMSEXPRESSCO", "ABBOTTLABORATORIES", "FIRSTFINANCIALBANCORPOHIO", "FIRSTMIDWESTBANCORPDE", "ROGERSCORP", "ROLLINSINC", "SEABOARDCORP", "STSOURCECORP", "SHERWINWILLIAMSCO", "FLEXSTEELINDUSTRIESINC", "FOSTERLBCO", "FRANKLINELECTRICINC", "FRANKLINRESOURCESINC", "FULLERHBCO", "GENTEXCORP", "DONNELLEYRRSONSCO", "AMREPCORP", "STANDARDMOTORPRODUCTSINC", "APACHECORP", "MATTELINC", "BECTONDICKINSONCO", "GRACOINC", "BARNESGROUPINC", "WEYERHAEUSERCO", "INTERNATIONALFLAVORSFRAGINC", "AVONPRODUCTSINC", "HAWAIIANELECTRICINDUSTRIESINC", "HAVERTYFURNITURECOSINC", "HEALTHCARESERVICESGROUPINC", "PARKERHANNIFINCORP", "WEISMARKETSINC", "WOLVERINEWORLDWIDEINC", "BARNWELLINDUSTRIESINC", "SMUCKERJMCO", "MOOGINC", "MEREDITHCORP", "HUNTJBTRANSPORTSERVICESINC", "HUNTINGTONBANCSHARESINC", "HURCOCOMPANY", "SYNALLOYCORP", "MCDONALDSCORP", "STANDEXINTERNATIONALCORP", "VFCORP", "IMMUNOMEDICSINC", "BEMISCOINC", "KENNAMETALINC", "PFINDUSTRIESINC", "TELEFLEXINC", "INTEGRATEDDEVICETECHNOLOGY", "AUTOMATICDATAPROCESSINGINC", "DIODESINC", "INVACARECORP", "MARSHMCLENNANCOSINC", "ISRAMCOINC", "SAFEGUARDSCIENTIFICSINC", "CLOROXCO", "GENUINEPARTSCO", "CABOTCORP", "RITEAIDCORP", "KAMANCORP", "KELLYSERVICESINC", "NEWYORKTIMESCO", "KEYTRONICSCORP", "CNAFINANCIALCORP", "FEDERALSIGNALCORP", "KIMBALLINTERNATIONALINC", "KULICKESOFFAINDSINC", "LAMRESHCORP", "LANCASTERCOLONYCORP", "HUMANAINC", "UNIONPACIFICCORP", "LAWSONPRODUCTSINC", "LINCOLNNATIONALCORP", "SIFCOINDUSTRIESINC", "ALEXANDERSINC", "BLOCKHRINC", "DANAHERCORP", "SERVOTRONICSINC", "LEEENTERPRISESINC", "MTSSYSTEMSCORP", "PENNSYLVANIAREALESTATEINVTTR", "LILLYELICO", "SEMTECHCORP", "WINNEBAGOINDUSTRIESINC", "SONYCORP", "MANITOWOCCOINC", "TERADYNEINC", "MARCUSCORP", "MARINEPETROLEUMTRUST", "SERVICECORPINTERNATIONAL", "SYSCOCORP", "MCCORMICKCOINC", "JACOBSENGINEERINGGROUPINC", "SUPERIORINDUSTRIESINTLINC", "EQUIFAXINC", "COCACOLABOTTLINGCOCONS", "DYCOMINDUSTRIESINC", "GATXCORP", "GENERALDYNAMICSCORP", "GENERALELECTRICCO"];
opts.VariableTypes = ["string", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double"];
opts = setvaropts(opts, 1, "WhitespaceRule", "preserve");
opts = setvaropts(opts, 1, "EmptyFieldRule", "auto");
dataS1 = readtable("data.xlsx", opts, "UseExcel", false);
clear opts

%regression
%step 1
stockdatasorted100S1=dataS1(85:144,:);
excess_asset=stockdatasorted100S1(:,3:102);
market_excess=stockdatasorted100S1(:,2);
n=width(excess_asset);
beta=zeros(1,n);
p=zeros(1,n);
alpha=zeros(1,n);
alpha_p=zeros(1,n);
for i=(1:n);
a=excess_asset(:,i);
regression=[market_excess a];
regression.Properties.VariableNames = {'ExcessMarket' 'ExcessStock'};
modelFit = fitlm(regression,'ExcessStock ~ ExcessMarket')
beta(i) = modelFit.Coefficients.Estimate(2);
p(i)= modelFit.Coefficients.pValue(2);
alpha(i)=modelFit.Coefficients.Estimate(1);
alpha_p(i)=modelFit.Coefficients.pValue(1);
end
output =[beta;p]
reject_rule1= p<0.05;
percent_of_significant_coeff=sum(reject_rule1)/n

%import average excess return of companies
opts = spreadsheetImportOptions("NumVariables", 1);
opts.Sheet = "average_excess_return";
opts.DataRange = "E2:E101";
opts.VariableNames = "average_excess_return";
opts.VariableTypes = "double";
dataS2 = readtable("data.xlsx", opts, "UseExcel", false);
clear opts

%step 2
beta=beta';
dataS2.beta = beta;
regression2=dataS2;
modelFit1 = fitlm(regression2,'average_excess_return ~ beta')

%Expended CAPM
%step 1 
market_excess.market_excess_sq=market_excess.MktRF.^2
market_excess.market_excess_cubus=market_excess.MktRF.^3
market_excess.market_excess_interaction=market_excess.market_excess_sq.*market_excess.market_excess_cubus
beta1=zeros(1,n);
beta2=zeros(1,n);
beta3=zeros(1,n);
beta4=zeros(1,n);
for i=(1:n);
a=excess_asset(:,i);
regression3=[market_excess a];
regression3.Properties.VariableNames = {'ExcessMarket' 'ExcessMarket2' 'ExcessMarket3' 'InteractionTerm' 'ExcessStock'};

modelFit3 = fitlm(regression3,'ExcessStock ~ ExcessMarket+ExcessMarket2+ExcessMarket3+InteractionTerm')
beta1(i) = modelFit3.Coefficients.Estimate(2);
beta2(i)=modelFit3.Coefficients.Estimate(3);
beta3(i)=modelFit3.Coefficients.Estimate(4);
beta4(i)=modelFit3.Coefficients.Estimate(5);
end

%step 2
beta1=beta1';
beta2=beta2';
beta3=beta3';
beta4=beta4';
dataS2.beta1 = beta1;
dataS2.beta2 = beta2;
dataS2.beta3 = beta3;
dataS2.beta4 = beta4;
regression4=dataS2;
model = fitlm(regression4,'average_excess_return ~ beta1+beta2+beta3+beta4')

%% 1998-2002
clear all
%Import return data
opts = spreadsheetImportOptions("NumVariables", 103);
opts.Sheet = "ri (stock data)";
opts.DataRange = "A2:CY397";
opts.VariableNames = ["Date", "MktRF", "RF", "ADAMSEXPRESSCO", "ABBOTTLABORATORIES", "FIRSTFINANCIALBANCORPOHIO", "FIRSTMIDWESTBANCORPDE", "ROGERSCORP", "ROLLINSINC", "SEABOARDCORP", "STSOURCECORP", "SHERWINWILLIAMSCO", "FLEXSTEELINDUSTRIESINC", "FOSTERLBCO", "FRANKLINELECTRICINC", "FRANKLINRESOURCESINC", "FULLERHBCO", "GENTEXCORP", "DONNELLEYRRSONSCO", "AMREPCORP", "STANDARDMOTORPRODUCTSINC", "APACHECORP", "MATTELINC", "BECTONDICKINSONCO", "GRACOINC", "BARNESGROUPINC", "WEYERHAEUSERCO", "INTERNATIONALFLAVORSFRAGINC", "AVONPRODUCTSINC", "HAWAIIANELECTRICINDUSTRIESINC", "HAVERTYFURNITURECOSINC", "HEALTHCARESERVICESGROUPINC", "PARKERHANNIFINCORP", "WEISMARKETSINC", "WOLVERINEWORLDWIDEINC", "BARNWELLINDUSTRIESINC", "SMUCKERJMCO", "MOOGINC", "MEREDITHCORP", "HUNTJBTRANSPORTSERVICESINC", "HUNTINGTONBANCSHARESINC", "HURCOCOMPANY", "SYNALLOYCORP", "MCDONALDSCORP", "STANDEXINTERNATIONALCORP", "VFCORP", "IMMUNOMEDICSINC", "BEMISCOINC", "KENNAMETALINC", "PFINDUSTRIESINC", "TELEFLEXINC", "INTEGRATEDDEVICETECHNOLOGY", "AUTOMATICDATAPROCESSINGINC", "DIODESINC", "INVACARECORP", "MARSHMCLENNANCOSINC", "ISRAMCOINC", "SAFEGUARDSCIENTIFICSINC", "CLOROXCO", "GENUINEPARTSCO", "CABOTCORP", "RITEAIDCORP", "KAMANCORP", "KELLYSERVICESINC", "NEWYORKTIMESCO", "KEYTRONICSCORP", "CNAFINANCIALCORP", "FEDERALSIGNALCORP", "KIMBALLINTERNATIONALINC", "KULICKESOFFAINDSINC", "LAMRESHCORP", "LANCASTERCOLONYCORP", "HUMANAINC", "UNIONPACIFICCORP", "LAWSONPRODUCTSINC", "LINCOLNNATIONALCORP", "SIFCOINDUSTRIESINC", "ALEXANDERSINC", "BLOCKHRINC", "DANAHERCORP", "SERVOTRONICSINC", "LEEENTERPRISESINC", "MTSSYSTEMSCORP", "PENNSYLVANIAREALESTATEINVTTR", "LILLYELICO", "SEMTECHCORP", "WINNEBAGOINDUSTRIESINC", "SONYCORP", "MANITOWOCCOINC", "TERADYNEINC", "MARCUSCORP", "MARINEPETROLEUMTRUST", "SERVICECORPINTERNATIONAL", "SYSCOCORP", "MCCORMICKCOINC", "JACOBSENGINEERINGGROUPINC", "SUPERIORINDUSTRIESINTLINC", "EQUIFAXINC", "COCACOLABOTTLINGCOCONS", "DYCOMINDUSTRIESINC", "GATXCORP", "GENERALDYNAMICSCORP", "GENERALELECTRICCO"];
opts.VariableTypes = ["string", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double"];
opts = setvaropts(opts, 1, "WhitespaceRule", "preserve");
opts = setvaropts(opts, 1, "EmptyFieldRule", "auto");
data = readtable("data.xlsx", opts, "UseExcel", false);
clear opts

%test for normality
stockdatasorted100=data(145:204,:);
return_asset=stockdatasorted100(:,4:103);
dataS1 =stockdatasorted100(:,2);
n=width(return_asset);
jb=zeros(1,n);
pval=zeros(1,n);
skew=zeros(1,n);
kurt=zeros(1,n);
t_stat=zeros(1,n);
T=396
for i=(1:n)
a=return_asset(:,i);
normality=[dataS1 a];
%change name in column 2
normality.Properties.VariableNames = {'ExcessMarket' 'StockReturn'};
[h, pvalue, jbstat, critval] = jbtest(normality.StockReturn)
jb(i)=jbstat
pval(i)=pvalue
skew(i)=skewness(normality.StockReturn);
kurt(i)=kurtosis(normality.StockReturn);
t_stat(i)=T *( skew(i)^2 / 6 + (kurt(i) - 3)^2 / 24 );
end
reject_rule= pval<0.05;
percent_of_abnormal_data=sum(reject_rule)/n

%import excess return data
opts = spreadsheetImportOptions("NumVariables", 102);
opts.Sheet = "ri-rf";
opts.DataRange = "A2:CX397";
opts.VariableNames = ["Date", "MktRF", "ADAMSEXPRESSCO", "ABBOTTLABORATORIES", "FIRSTFINANCIALBANCORPOHIO", "FIRSTMIDWESTBANCORPDE", "ROGERSCORP", "ROLLINSINC", "SEABOARDCORP", "STSOURCECORP", "SHERWINWILLIAMSCO", "FLEXSTEELINDUSTRIESINC", "FOSTERLBCO", "FRANKLINELECTRICINC", "FRANKLINRESOURCESINC", "FULLERHBCO", "GENTEXCORP", "DONNELLEYRRSONSCO", "AMREPCORP", "STANDARDMOTORPRODUCTSINC", "APACHECORP", "MATTELINC", "BECTONDICKINSONCO", "GRACOINC", "BARNESGROUPINC", "WEYERHAEUSERCO", "INTERNATIONALFLAVORSFRAGINC", "AVONPRODUCTSINC", "HAWAIIANELECTRICINDUSTRIESINC", "HAVERTYFURNITURECOSINC", "HEALTHCARESERVICESGROUPINC", "PARKERHANNIFINCORP", "WEISMARKETSINC", "WOLVERINEWORLDWIDEINC", "BARNWELLINDUSTRIESINC", "SMUCKERJMCO", "MOOGINC", "MEREDITHCORP", "HUNTJBTRANSPORTSERVICESINC", "HUNTINGTONBANCSHARESINC", "HURCOCOMPANY", "SYNALLOYCORP", "MCDONALDSCORP", "STANDEXINTERNATIONALCORP", "VFCORP", "IMMUNOMEDICSINC", "BEMISCOINC", "KENNAMETALINC", "PFINDUSTRIESINC", "TELEFLEXINC", "INTEGRATEDDEVICETECHNOLOGY", "AUTOMATICDATAPROCESSINGINC", "DIODESINC", "INVACARECORP", "MARSHMCLENNANCOSINC", "ISRAMCOINC", "SAFEGUARDSCIENTIFICSINC", "CLOROXCO", "GENUINEPARTSCO", "CABOTCORP", "RITEAIDCORP", "KAMANCORP", "KELLYSERVICESINC", "NEWYORKTIMESCO", "KEYTRONICSCORP", "CNAFINANCIALCORP", "FEDERALSIGNALCORP", "KIMBALLINTERNATIONALINC", "KULICKESOFFAINDSINC", "LAMRESHCORP", "LANCASTERCOLONYCORP", "HUMANAINC", "UNIONPACIFICCORP", "LAWSONPRODUCTSINC", "LINCOLNNATIONALCORP", "SIFCOINDUSTRIESINC", "ALEXANDERSINC", "BLOCKHRINC", "DANAHERCORP", "SERVOTRONICSINC", "LEEENTERPRISESINC", "MTSSYSTEMSCORP", "PENNSYLVANIAREALESTATEINVTTR", "LILLYELICO", "SEMTECHCORP", "WINNEBAGOINDUSTRIESINC", "SONYCORP", "MANITOWOCCOINC", "TERADYNEINC", "MARCUSCORP", "MARINEPETROLEUMTRUST", "SERVICECORPINTERNATIONAL", "SYSCOCORP", "MCCORMICKCOINC", "JACOBSENGINEERINGGROUPINC", "SUPERIORINDUSTRIESINTLINC", "EQUIFAXINC", "COCACOLABOTTLINGCOCONS", "DYCOMINDUSTRIESINC", "GATXCORP", "GENERALDYNAMICSCORP", "GENERALELECTRICCO"];
opts.VariableTypes = ["string", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double"];
opts = setvaropts(opts, 1, "WhitespaceRule", "preserve");
opts = setvaropts(opts, 1, "EmptyFieldRule", "auto");
dataS1 = readtable("data.xlsx", opts, "UseExcel", false);
clear opts

%regression
%step 1
stockdatasorted100S1=dataS1(145:204,:);
excess_asset=stockdatasorted100S1(:,3:102);
market_excess=stockdatasorted100S1(:,2);
n=width(excess_asset);
beta=zeros(1,n);
p=zeros(1,n);
alpha=zeros(1,n);
alpha_p=zeros(1,n);
for i=(1:n);
a=excess_asset(:,i);
regression=[market_excess a];
regression.Properties.VariableNames = {'ExcessMarket' 'ExcessStock'};
modelFit = fitlm(regression,'ExcessStock ~ ExcessMarket')
beta(i) = modelFit.Coefficients.Estimate(2);
p(i)= modelFit.Coefficients.pValue(2);
alpha(i)=modelFit.Coefficients.Estimate(1);
alpha_p(i)=modelFit.Coefficients.pValue(1);
end
output =[beta;p]
reject_rule1= p<0.05;
percent_of_significant_coeff=sum(reject_rule1)/n

%import average excess return of companies
opts = spreadsheetImportOptions("NumVariables", 1);
opts.Sheet = "average_excess_return";
opts.DataRange = "F2:F101";
opts.VariableNames = "average_excess_return";
opts.VariableTypes = "double";
dataS2 = readtable("data.xlsx", opts, "UseExcel", false);
clear opts

%step 2
beta=beta';
dataS2.beta = beta;
regression2=dataS2;
modelFit1 = fitlm(regression2,'average_excess_return ~ beta')

%Expended CAPM
%step 1 
market_excess.market_excess_sq=market_excess.MktRF.^2
market_excess.market_excess_cubus=market_excess.MktRF.^3
market_excess.market_excess_interaction=market_excess.market_excess_sq.*market_excess.market_excess_cubus
beta1=zeros(1,n);
beta2=zeros(1,n);
beta3=zeros(1,n);
beta4=zeros(1,n);
for i=(1:n);
a=excess_asset(:,i);
regression3=[market_excess a];
regression3.Properties.VariableNames = {'ExcessMarket' 'ExcessMarket2' 'ExcessMarket3' 'InteractionTerm' 'ExcessStock'};

modelFit3 = fitlm(regression3,'ExcessStock ~ ExcessMarket+ExcessMarket2+ExcessMarket3+InteractionTerm')
beta1(i) = modelFit3.Coefficients.Estimate(2);
beta2(i)=modelFit3.Coefficients.Estimate(3);
beta3(i)=modelFit3.Coefficients.Estimate(4);
beta4(i)=modelFit3.Coefficients.Estimate(5);
end

%step 2
beta1=beta1';
beta2=beta2';
beta3=beta3';
beta4=beta4';
dataS2.beta1 = beta1;
dataS2.beta2 = beta2;
dataS2.beta3 = beta3;
dataS2.beta4 = beta4;
regression4=dataS2;
model = fitlm(regression4,'average_excess_return ~ beta1+beta2+beta3+beta4')

%% 2003-2007
clear all
%Import return data
opts = spreadsheetImportOptions("NumVariables", 103);
opts.Sheet = "ri (stock data)";
opts.DataRange = "A2:CY397";
opts.VariableNames = ["Date", "MktRF", "RF", "ADAMSEXPRESSCO", "ABBOTTLABORATORIES", "FIRSTFINANCIALBANCORPOHIO", "FIRSTMIDWESTBANCORPDE", "ROGERSCORP", "ROLLINSINC", "SEABOARDCORP", "STSOURCECORP", "SHERWINWILLIAMSCO", "FLEXSTEELINDUSTRIESINC", "FOSTERLBCO", "FRANKLINELECTRICINC", "FRANKLINRESOURCESINC", "FULLERHBCO", "GENTEXCORP", "DONNELLEYRRSONSCO", "AMREPCORP", "STANDARDMOTORPRODUCTSINC", "APACHECORP", "MATTELINC", "BECTONDICKINSONCO", "GRACOINC", "BARNESGROUPINC", "WEYERHAEUSERCO", "INTERNATIONALFLAVORSFRAGINC", "AVONPRODUCTSINC", "HAWAIIANELECTRICINDUSTRIESINC", "HAVERTYFURNITURECOSINC", "HEALTHCARESERVICESGROUPINC", "PARKERHANNIFINCORP", "WEISMARKETSINC", "WOLVERINEWORLDWIDEINC", "BARNWELLINDUSTRIESINC", "SMUCKERJMCO", "MOOGINC", "MEREDITHCORP", "HUNTJBTRANSPORTSERVICESINC", "HUNTINGTONBANCSHARESINC", "HURCOCOMPANY", "SYNALLOYCORP", "MCDONALDSCORP", "STANDEXINTERNATIONALCORP", "VFCORP", "IMMUNOMEDICSINC", "BEMISCOINC", "KENNAMETALINC", "PFINDUSTRIESINC", "TELEFLEXINC", "INTEGRATEDDEVICETECHNOLOGY", "AUTOMATICDATAPROCESSINGINC", "DIODESINC", "INVACARECORP", "MARSHMCLENNANCOSINC", "ISRAMCOINC", "SAFEGUARDSCIENTIFICSINC", "CLOROXCO", "GENUINEPARTSCO", "CABOTCORP", "RITEAIDCORP", "KAMANCORP", "KELLYSERVICESINC", "NEWYORKTIMESCO", "KEYTRONICSCORP", "CNAFINANCIALCORP", "FEDERALSIGNALCORP", "KIMBALLINTERNATIONALINC", "KULICKESOFFAINDSINC", "LAMRESHCORP", "LANCASTERCOLONYCORP", "HUMANAINC", "UNIONPACIFICCORP", "LAWSONPRODUCTSINC", "LINCOLNNATIONALCORP", "SIFCOINDUSTRIESINC", "ALEXANDERSINC", "BLOCKHRINC", "DANAHERCORP", "SERVOTRONICSINC", "LEEENTERPRISESINC", "MTSSYSTEMSCORP", "PENNSYLVANIAREALESTATEINVTTR", "LILLYELICO", "SEMTECHCORP", "WINNEBAGOINDUSTRIESINC", "SONYCORP", "MANITOWOCCOINC", "TERADYNEINC", "MARCUSCORP", "MARINEPETROLEUMTRUST", "SERVICECORPINTERNATIONAL", "SYSCOCORP", "MCCORMICKCOINC", "JACOBSENGINEERINGGROUPINC", "SUPERIORINDUSTRIESINTLINC", "EQUIFAXINC", "COCACOLABOTTLINGCOCONS", "DYCOMINDUSTRIESINC", "GATXCORP", "GENERALDYNAMICSCORP", "GENERALELECTRICCO"];
opts.VariableTypes = ["string", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double"];
opts = setvaropts(opts, 1, "WhitespaceRule", "preserve");
opts = setvaropts(opts, 1, "EmptyFieldRule", "auto");
data = readtable("data.xlsx", opts, "UseExcel", false);
clear opts

%test for normality
stockdatasorted100=data(205:264,:);
return_asset=stockdatasorted100(:,4:103);
dataS1 =stockdatasorted100(:,2);
n=width(return_asset);
jb=zeros(1,n);
pval=zeros(1,n);
skew=zeros(1,n);
kurt=zeros(1,n);
t_stat=zeros(1,n);
T=396
for i=(1:n)
a=return_asset(:,i);
normality=[dataS1 a];
%change name in column 2
normality.Properties.VariableNames = {'ExcessMarket' 'StockReturn'};
[h, pvalue, jbstat, critval] = jbtest(normality.StockReturn)
jb(i)=jbstat
pval(i)=pvalue
skew(i)=skewness(normality.StockReturn);
kurt(i)=kurtosis(normality.StockReturn);
t_stat(i)=T *( skew(i)^2 / 6 + (kurt(i) - 3)^2 / 24 );
end
reject_rule= pval<0.05;
percent_of_abnormal_data=sum(reject_rule)/n

%import excess return data
opts = spreadsheetImportOptions("NumVariables", 102);
opts.Sheet = "ri-rf";
opts.DataRange = "A2:CX397";
opts.VariableNames = ["Date", "MktRF", "ADAMSEXPRESSCO", "ABBOTTLABORATORIES", "FIRSTFINANCIALBANCORPOHIO", "FIRSTMIDWESTBANCORPDE", "ROGERSCORP", "ROLLINSINC", "SEABOARDCORP", "STSOURCECORP", "SHERWINWILLIAMSCO", "FLEXSTEELINDUSTRIESINC", "FOSTERLBCO", "FRANKLINELECTRICINC", "FRANKLINRESOURCESINC", "FULLERHBCO", "GENTEXCORP", "DONNELLEYRRSONSCO", "AMREPCORP", "STANDARDMOTORPRODUCTSINC", "APACHECORP", "MATTELINC", "BECTONDICKINSONCO", "GRACOINC", "BARNESGROUPINC", "WEYERHAEUSERCO", "INTERNATIONALFLAVORSFRAGINC", "AVONPRODUCTSINC", "HAWAIIANELECTRICINDUSTRIESINC", "HAVERTYFURNITURECOSINC", "HEALTHCARESERVICESGROUPINC", "PARKERHANNIFINCORP", "WEISMARKETSINC", "WOLVERINEWORLDWIDEINC", "BARNWELLINDUSTRIESINC", "SMUCKERJMCO", "MOOGINC", "MEREDITHCORP", "HUNTJBTRANSPORTSERVICESINC", "HUNTINGTONBANCSHARESINC", "HURCOCOMPANY", "SYNALLOYCORP", "MCDONALDSCORP", "STANDEXINTERNATIONALCORP", "VFCORP", "IMMUNOMEDICSINC", "BEMISCOINC", "KENNAMETALINC", "PFINDUSTRIESINC", "TELEFLEXINC", "INTEGRATEDDEVICETECHNOLOGY", "AUTOMATICDATAPROCESSINGINC", "DIODESINC", "INVACARECORP", "MARSHMCLENNANCOSINC", "ISRAMCOINC", "SAFEGUARDSCIENTIFICSINC", "CLOROXCO", "GENUINEPARTSCO", "CABOTCORP", "RITEAIDCORP", "KAMANCORP", "KELLYSERVICESINC", "NEWYORKTIMESCO", "KEYTRONICSCORP", "CNAFINANCIALCORP", "FEDERALSIGNALCORP", "KIMBALLINTERNATIONALINC", "KULICKESOFFAINDSINC", "LAMRESHCORP", "LANCASTERCOLONYCORP", "HUMANAINC", "UNIONPACIFICCORP", "LAWSONPRODUCTSINC", "LINCOLNNATIONALCORP", "SIFCOINDUSTRIESINC", "ALEXANDERSINC", "BLOCKHRINC", "DANAHERCORP", "SERVOTRONICSINC", "LEEENTERPRISESINC", "MTSSYSTEMSCORP", "PENNSYLVANIAREALESTATEINVTTR", "LILLYELICO", "SEMTECHCORP", "WINNEBAGOINDUSTRIESINC", "SONYCORP", "MANITOWOCCOINC", "TERADYNEINC", "MARCUSCORP", "MARINEPETROLEUMTRUST", "SERVICECORPINTERNATIONAL", "SYSCOCORP", "MCCORMICKCOINC", "JACOBSENGINEERINGGROUPINC", "SUPERIORINDUSTRIESINTLINC", "EQUIFAXINC", "COCACOLABOTTLINGCOCONS", "DYCOMINDUSTRIESINC", "GATXCORP", "GENERALDYNAMICSCORP", "GENERALELECTRICCO"];
opts.VariableTypes = ["string", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double"];
opts = setvaropts(opts, 1, "WhitespaceRule", "preserve");
opts = setvaropts(opts, 1, "EmptyFieldRule", "auto");
dataS1 = readtable("data.xlsx", opts, "UseExcel", false);
clear opts

%regression
%step 1
stockdatasorted100S1=dataS1(205:264,:);
excess_asset=stockdatasorted100S1(:,3:102);
market_excess=stockdatasorted100S1(:,2);
n=width(excess_asset);
beta=zeros(1,n);
p=zeros(1,n);
alpha=zeros(1,n);
alpha_p=zeros(1,n);
for i=(1:n);
a=excess_asset(:,i);
regression=[market_excess a];
regression.Properties.VariableNames = {'ExcessMarket' 'ExcessStock'};
modelFit = fitlm(regression,'ExcessStock ~ ExcessMarket')
beta(i) = modelFit.Coefficients.Estimate(2);
p(i)= modelFit.Coefficients.pValue(2);
alpha(i)=modelFit.Coefficients.Estimate(1);
alpha_p(i)=modelFit.Coefficients.pValue(1);
end
output =[beta;p]
reject_rule1= p<0.05;
percent_of_significant_coeff=sum(reject_rule1)/n

%import average excess return of companies
opts = spreadsheetImportOptions("NumVariables", 1);
opts.Sheet = "average_excess_return";
opts.DataRange = "G2:G101";
opts.VariableNames = "average_excess_return";
opts.VariableTypes = "double";
dataS2 = readtable("data.xlsx", opts, "UseExcel", false);
clear opts

%step 2
beta=beta';
dataS2.beta = beta;
regression2=dataS2;
modelFit1 = fitlm(regression2,'average_excess_return ~ beta')

%Expended CAPM
%step 1 
market_excess.market_excess_sq=market_excess.MktRF.^2
market_excess.market_excess_cubus=market_excess.MktRF.^3
market_excess.market_excess_interaction=market_excess.market_excess_sq.*market_excess.market_excess_cubus
beta1=zeros(1,n);
beta2=zeros(1,n);
beta3=zeros(1,n);
beta4=zeros(1,n);
for i=(1:n);
a=excess_asset(:,i);
regression3=[market_excess a];
regression3.Properties.VariableNames = {'ExcessMarket' 'ExcessMarket2' 'ExcessMarket3' 'InteractionTerm' 'ExcessStock'};

modelFit3 = fitlm(regression3,'ExcessStock ~ ExcessMarket+ExcessMarket2+ExcessMarket3+InteractionTerm')
beta1(i) = modelFit3.Coefficients.Estimate(2);
beta2(i)=modelFit3.Coefficients.Estimate(3);
beta3(i)=modelFit3.Coefficients.Estimate(4);
beta4(i)=modelFit3.Coefficients.Estimate(5);
end

%step 2
beta1=beta1';
beta2=beta2';
beta3=beta3';
beta4=beta4';
dataS2.beta1 = beta1;
dataS2.beta2 = beta2;
dataS2.beta3 = beta3;
dataS2.beta4 = beta4;
regression4=dataS2;
model = fitlm(regression4,'average_excess_return ~ beta1+beta2+beta3+beta4')

%% 2008-2012
clear all
%Import return data
opts = spreadsheetImportOptions("NumVariables", 103);
opts.Sheet = "ri (stock data)";
opts.DataRange = "A2:CY397";
opts.VariableNames = ["Date", "MktRF", "RF", "ADAMSEXPRESSCO", "ABBOTTLABORATORIES", "FIRSTFINANCIALBANCORPOHIO", "FIRSTMIDWESTBANCORPDE", "ROGERSCORP", "ROLLINSINC", "SEABOARDCORP", "STSOURCECORP", "SHERWINWILLIAMSCO", "FLEXSTEELINDUSTRIESINC", "FOSTERLBCO", "FRANKLINELECTRICINC", "FRANKLINRESOURCESINC", "FULLERHBCO", "GENTEXCORP", "DONNELLEYRRSONSCO", "AMREPCORP", "STANDARDMOTORPRODUCTSINC", "APACHECORP", "MATTELINC", "BECTONDICKINSONCO", "GRACOINC", "BARNESGROUPINC", "WEYERHAEUSERCO", "INTERNATIONALFLAVORSFRAGINC", "AVONPRODUCTSINC", "HAWAIIANELECTRICINDUSTRIESINC", "HAVERTYFURNITURECOSINC", "HEALTHCARESERVICESGROUPINC", "PARKERHANNIFINCORP", "WEISMARKETSINC", "WOLVERINEWORLDWIDEINC", "BARNWELLINDUSTRIESINC", "SMUCKERJMCO", "MOOGINC", "MEREDITHCORP", "HUNTJBTRANSPORTSERVICESINC", "HUNTINGTONBANCSHARESINC", "HURCOCOMPANY", "SYNALLOYCORP", "MCDONALDSCORP", "STANDEXINTERNATIONALCORP", "VFCORP", "IMMUNOMEDICSINC", "BEMISCOINC", "KENNAMETALINC", "PFINDUSTRIESINC", "TELEFLEXINC", "INTEGRATEDDEVICETECHNOLOGY", "AUTOMATICDATAPROCESSINGINC", "DIODESINC", "INVACARECORP", "MARSHMCLENNANCOSINC", "ISRAMCOINC", "SAFEGUARDSCIENTIFICSINC", "CLOROXCO", "GENUINEPARTSCO", "CABOTCORP", "RITEAIDCORP", "KAMANCORP", "KELLYSERVICESINC", "NEWYORKTIMESCO", "KEYTRONICSCORP", "CNAFINANCIALCORP", "FEDERALSIGNALCORP", "KIMBALLINTERNATIONALINC", "KULICKESOFFAINDSINC", "LAMRESHCORP", "LANCASTERCOLONYCORP", "HUMANAINC", "UNIONPACIFICCORP", "LAWSONPRODUCTSINC", "LINCOLNNATIONALCORP", "SIFCOINDUSTRIESINC", "ALEXANDERSINC", "BLOCKHRINC", "DANAHERCORP", "SERVOTRONICSINC", "LEEENTERPRISESINC", "MTSSYSTEMSCORP", "PENNSYLVANIAREALESTATEINVTTR", "LILLYELICO", "SEMTECHCORP", "WINNEBAGOINDUSTRIESINC", "SONYCORP", "MANITOWOCCOINC", "TERADYNEINC", "MARCUSCORP", "MARINEPETROLEUMTRUST", "SERVICECORPINTERNATIONAL", "SYSCOCORP", "MCCORMICKCOINC", "JACOBSENGINEERINGGROUPINC", "SUPERIORINDUSTRIESINTLINC", "EQUIFAXINC", "COCACOLABOTTLINGCOCONS", "DYCOMINDUSTRIESINC", "GATXCORP", "GENERALDYNAMICSCORP", "GENERALELECTRICCO"];
opts.VariableTypes = ["string", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double"];
opts = setvaropts(opts, 1, "WhitespaceRule", "preserve");
opts = setvaropts(opts, 1, "EmptyFieldRule", "auto");
data = readtable("data.xlsx", opts, "UseExcel", false);
clear opts

%test for normality
stockdatasorted100=data(265:324,:);
return_asset=stockdatasorted100(:,4:103);
dataS1 =stockdatasorted100(:,2);
n=width(return_asset);
jb=zeros(1,n);
pval=zeros(1,n);
skew=zeros(1,n);
kurt=zeros(1,n);
t_stat=zeros(1,n);
T=396
for i=(1:n)
a=return_asset(:,i);
normality=[dataS1 a];
%change name in column 2
normality.Properties.VariableNames = {'ExcessMarket' 'StockReturn'};
[h, pvalue, jbstat, critval] = jbtest(normality.StockReturn)
jb(i)=jbstat
pval(i)=pvalue
skew(i)=skewness(normality.StockReturn);
kurt(i)=kurtosis(normality.StockReturn);
t_stat(i)=T *( skew(i)^2 / 6 + (kurt(i) - 3)^2 / 24 );
end
reject_rule= pval<0.05;
percent_of_abnormal_data=sum(reject_rule)/n

%import excess return data
opts = spreadsheetImportOptions("NumVariables", 102);
opts.Sheet = "ri-rf";
opts.DataRange = "A2:CX397";
opts.VariableNames = ["Date", "MktRF", "ADAMSEXPRESSCO", "ABBOTTLABORATORIES", "FIRSTFINANCIALBANCORPOHIO", "FIRSTMIDWESTBANCORPDE", "ROGERSCORP", "ROLLINSINC", "SEABOARDCORP", "STSOURCECORP", "SHERWINWILLIAMSCO", "FLEXSTEELINDUSTRIESINC", "FOSTERLBCO", "FRANKLINELECTRICINC", "FRANKLINRESOURCESINC", "FULLERHBCO", "GENTEXCORP", "DONNELLEYRRSONSCO", "AMREPCORP", "STANDARDMOTORPRODUCTSINC", "APACHECORP", "MATTELINC", "BECTONDICKINSONCO", "GRACOINC", "BARNESGROUPINC", "WEYERHAEUSERCO", "INTERNATIONALFLAVORSFRAGINC", "AVONPRODUCTSINC", "HAWAIIANELECTRICINDUSTRIESINC", "HAVERTYFURNITURECOSINC", "HEALTHCARESERVICESGROUPINC", "PARKERHANNIFINCORP", "WEISMARKETSINC", "WOLVERINEWORLDWIDEINC", "BARNWELLINDUSTRIESINC", "SMUCKERJMCO", "MOOGINC", "MEREDITHCORP", "HUNTJBTRANSPORTSERVICESINC", "HUNTINGTONBANCSHARESINC", "HURCOCOMPANY", "SYNALLOYCORP", "MCDONALDSCORP", "STANDEXINTERNATIONALCORP", "VFCORP", "IMMUNOMEDICSINC", "BEMISCOINC", "KENNAMETALINC", "PFINDUSTRIESINC", "TELEFLEXINC", "INTEGRATEDDEVICETECHNOLOGY", "AUTOMATICDATAPROCESSINGINC", "DIODESINC", "INVACARECORP", "MARSHMCLENNANCOSINC", "ISRAMCOINC", "SAFEGUARDSCIENTIFICSINC", "CLOROXCO", "GENUINEPARTSCO", "CABOTCORP", "RITEAIDCORP", "KAMANCORP", "KELLYSERVICESINC", "NEWYORKTIMESCO", "KEYTRONICSCORP", "CNAFINANCIALCORP", "FEDERALSIGNALCORP", "KIMBALLINTERNATIONALINC", "KULICKESOFFAINDSINC", "LAMRESHCORP", "LANCASTERCOLONYCORP", "HUMANAINC", "UNIONPACIFICCORP", "LAWSONPRODUCTSINC", "LINCOLNNATIONALCORP", "SIFCOINDUSTRIESINC", "ALEXANDERSINC", "BLOCKHRINC", "DANAHERCORP", "SERVOTRONICSINC", "LEEENTERPRISESINC", "MTSSYSTEMSCORP", "PENNSYLVANIAREALESTATEINVTTR", "LILLYELICO", "SEMTECHCORP", "WINNEBAGOINDUSTRIESINC", "SONYCORP", "MANITOWOCCOINC", "TERADYNEINC", "MARCUSCORP", "MARINEPETROLEUMTRUST", "SERVICECORPINTERNATIONAL", "SYSCOCORP", "MCCORMICKCOINC", "JACOBSENGINEERINGGROUPINC", "SUPERIORINDUSTRIESINTLINC", "EQUIFAXINC", "COCACOLABOTTLINGCOCONS", "DYCOMINDUSTRIESINC", "GATXCORP", "GENERALDYNAMICSCORP", "GENERALELECTRICCO"];
opts.VariableTypes = ["string", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double"];
opts = setvaropts(opts, 1, "WhitespaceRule", "preserve");
opts = setvaropts(opts, 1, "EmptyFieldRule", "auto");
dataS1 = readtable("data.xlsx", opts, "UseExcel", false);
clear opts

%regression
%step 1
stockdatasorted100S1=dataS1(265:324,:);
excess_asset=stockdatasorted100S1(:,3:102);
market_excess=stockdatasorted100S1(:,2);
n=width(excess_asset);
beta=zeros(1,n);
p=zeros(1,n);
alpha=zeros(1,n);
alpha_p=zeros(1,n);
for i=(1:n);
a=excess_asset(:,i);
regression=[market_excess a];
regression.Properties.VariableNames = {'ExcessMarket' 'ExcessStock'};
modelFit = fitlm(regression,'ExcessStock ~ ExcessMarket')
beta(i) = modelFit.Coefficients.Estimate(2);
p(i)= modelFit.Coefficients.pValue(2);
alpha(i)=modelFit.Coefficients.Estimate(1);
alpha_p(i)=modelFit.Coefficients.pValue(1);
end
output =[beta;p]
reject_rule1= p<0.05;
percent_of_significant_coeff=sum(reject_rule1)/n

%import average excess return of companies
opts = spreadsheetImportOptions("NumVariables", 1);
opts.Sheet = "average_excess_return";
opts.DataRange = "h2:h101";
opts.VariableNames = "average_excess_return";
opts.VariableTypes = "double";
dataS2 = readtable("data.xlsx", opts, "UseExcel", false);
clear opts

%step 2
beta=beta';
dataS2.beta = beta;
regression2=dataS2;
modelFit1 = fitlm(regression2,'average_excess_return ~ beta')

%% 2013-2017
clear all
%Import return data
opts = spreadsheetImportOptions("NumVariables", 103);
opts.Sheet = "ri (stock data)";
opts.DataRange = "A2:CY397";
opts.VariableNames = ["Date", "MktRF", "RF", "ADAMSEXPRESSCO", "ABBOTTLABORATORIES", "FIRSTFINANCIALBANCORPOHIO", "FIRSTMIDWESTBANCORPDE", "ROGERSCORP", "ROLLINSINC", "SEABOARDCORP", "STSOURCECORP", "SHERWINWILLIAMSCO", "FLEXSTEELINDUSTRIESINC", "FOSTERLBCO", "FRANKLINELECTRICINC", "FRANKLINRESOURCESINC", "FULLERHBCO", "GENTEXCORP", "DONNELLEYRRSONSCO", "AMREPCORP", "STANDARDMOTORPRODUCTSINC", "APACHECORP", "MATTELINC", "BECTONDICKINSONCO", "GRACOINC", "BARNESGROUPINC", "WEYERHAEUSERCO", "INTERNATIONALFLAVORSFRAGINC", "AVONPRODUCTSINC", "HAWAIIANELECTRICINDUSTRIESINC", "HAVERTYFURNITURECOSINC", "HEALTHCARESERVICESGROUPINC", "PARKERHANNIFINCORP", "WEISMARKETSINC", "WOLVERINEWORLDWIDEINC", "BARNWELLINDUSTRIESINC", "SMUCKERJMCO", "MOOGINC", "MEREDITHCORP", "HUNTJBTRANSPORTSERVICESINC", "HUNTINGTONBANCSHARESINC", "HURCOCOMPANY", "SYNALLOYCORP", "MCDONALDSCORP", "STANDEXINTERNATIONALCORP", "VFCORP", "IMMUNOMEDICSINC", "BEMISCOINC", "KENNAMETALINC", "PFINDUSTRIESINC", "TELEFLEXINC", "INTEGRATEDDEVICETECHNOLOGY", "AUTOMATICDATAPROCESSINGINC", "DIODESINC", "INVACARECORP", "MARSHMCLENNANCOSINC", "ISRAMCOINC", "SAFEGUARDSCIENTIFICSINC", "CLOROXCO", "GENUINEPARTSCO", "CABOTCORP", "RITEAIDCORP", "KAMANCORP", "KELLYSERVICESINC", "NEWYORKTIMESCO", "KEYTRONICSCORP", "CNAFINANCIALCORP", "FEDERALSIGNALCORP", "KIMBALLINTERNATIONALINC", "KULICKESOFFAINDSINC", "LAMRESHCORP", "LANCASTERCOLONYCORP", "HUMANAINC", "UNIONPACIFICCORP", "LAWSONPRODUCTSINC", "LINCOLNNATIONALCORP", "SIFCOINDUSTRIESINC", "ALEXANDERSINC", "BLOCKHRINC", "DANAHERCORP", "SERVOTRONICSINC", "LEEENTERPRISESINC", "MTSSYSTEMSCORP", "PENNSYLVANIAREALESTATEINVTTR", "LILLYELICO", "SEMTECHCORP", "WINNEBAGOINDUSTRIESINC", "SONYCORP", "MANITOWOCCOINC", "TERADYNEINC", "MARCUSCORP", "MARINEPETROLEUMTRUST", "SERVICECORPINTERNATIONAL", "SYSCOCORP", "MCCORMICKCOINC", "JACOBSENGINEERINGGROUPINC", "SUPERIORINDUSTRIESINTLINC", "EQUIFAXINC", "COCACOLABOTTLINGCOCONS", "DYCOMINDUSTRIESINC", "GATXCORP", "GENERALDYNAMICSCORP", "GENERALELECTRICCO"];
opts.VariableTypes = ["string", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double"];
opts = setvaropts(opts, 1, "WhitespaceRule", "preserve");
opts = setvaropts(opts, 1, "EmptyFieldRule", "auto");
data = readtable("data.xlsx", opts, "UseExcel", false);
clear opts

%test for normality
stockdatasorted100=data(325:384,:);
return_asset=stockdatasorted100(:,4:103);
dataS1 =stockdatasorted100(:,2);
n=width(return_asset);
jb=zeros(1,n);
pval=zeros(1,n);
skew=zeros(1,n);
kurt=zeros(1,n);
t_stat=zeros(1,n);
T=396
for i=(1:n)
a=return_asset(:,i);
normality=[dataS1 a];
%change name in column 2
normality.Properties.VariableNames = {'ExcessMarket' 'StockReturn'};
[h, pvalue, jbstat, critval] = jbtest(normality.StockReturn)
jb(i)=jbstat
pval(i)=pvalue
skew(i)=skewness(normality.StockReturn);
kurt(i)=kurtosis(normality.StockReturn);
t_stat(i)=T *( skew(i)^2 / 6 + (kurt(i) - 3)^2 / 24 );
end
reject_rule= pval<0.05;
percent_of_abnormal_data=sum(reject_rule)/n

%import excess return data
opts = spreadsheetImportOptions("NumVariables", 102);
opts.Sheet = "ri-rf";
opts.DataRange = "A2:CX397";
opts.VariableNames = ["Date", "MktRF", "ADAMSEXPRESSCO", "ABBOTTLABORATORIES", "FIRSTFINANCIALBANCORPOHIO", "FIRSTMIDWESTBANCORPDE", "ROGERSCORP", "ROLLINSINC", "SEABOARDCORP", "STSOURCECORP", "SHERWINWILLIAMSCO", "FLEXSTEELINDUSTRIESINC", "FOSTERLBCO", "FRANKLINELECTRICINC", "FRANKLINRESOURCESINC", "FULLERHBCO", "GENTEXCORP", "DONNELLEYRRSONSCO", "AMREPCORP", "STANDARDMOTORPRODUCTSINC", "APACHECORP", "MATTELINC", "BECTONDICKINSONCO", "GRACOINC", "BARNESGROUPINC", "WEYERHAEUSERCO", "INTERNATIONALFLAVORSFRAGINC", "AVONPRODUCTSINC", "HAWAIIANELECTRICINDUSTRIESINC", "HAVERTYFURNITURECOSINC", "HEALTHCARESERVICESGROUPINC", "PARKERHANNIFINCORP", "WEISMARKETSINC", "WOLVERINEWORLDWIDEINC", "BARNWELLINDUSTRIESINC", "SMUCKERJMCO", "MOOGINC", "MEREDITHCORP", "HUNTJBTRANSPORTSERVICESINC", "HUNTINGTONBANCSHARESINC", "HURCOCOMPANY", "SYNALLOYCORP", "MCDONALDSCORP", "STANDEXINTERNATIONALCORP", "VFCORP", "IMMUNOMEDICSINC", "BEMISCOINC", "KENNAMETALINC", "PFINDUSTRIESINC", "TELEFLEXINC", "INTEGRATEDDEVICETECHNOLOGY", "AUTOMATICDATAPROCESSINGINC", "DIODESINC", "INVACARECORP", "MARSHMCLENNANCOSINC", "ISRAMCOINC", "SAFEGUARDSCIENTIFICSINC", "CLOROXCO", "GENUINEPARTSCO", "CABOTCORP", "RITEAIDCORP", "KAMANCORP", "KELLYSERVICESINC", "NEWYORKTIMESCO", "KEYTRONICSCORP", "CNAFINANCIALCORP", "FEDERALSIGNALCORP", "KIMBALLINTERNATIONALINC", "KULICKESOFFAINDSINC", "LAMRESHCORP", "LANCASTERCOLONYCORP", "HUMANAINC", "UNIONPACIFICCORP", "LAWSONPRODUCTSINC", "LINCOLNNATIONALCORP", "SIFCOINDUSTRIESINC", "ALEXANDERSINC", "BLOCKHRINC", "DANAHERCORP", "SERVOTRONICSINC", "LEEENTERPRISESINC", "MTSSYSTEMSCORP", "PENNSYLVANIAREALESTATEINVTTR", "LILLYELICO", "SEMTECHCORP", "WINNEBAGOINDUSTRIESINC", "SONYCORP", "MANITOWOCCOINC", "TERADYNEINC", "MARCUSCORP", "MARINEPETROLEUMTRUST", "SERVICECORPINTERNATIONAL", "SYSCOCORP", "MCCORMICKCOINC", "JACOBSENGINEERINGGROUPINC", "SUPERIORINDUSTRIESINTLINC", "EQUIFAXINC", "COCACOLABOTTLINGCOCONS", "DYCOMINDUSTRIESINC", "GATXCORP", "GENERALDYNAMICSCORP", "GENERALELECTRICCO"];
opts.VariableTypes = ["string", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double"];
opts = setvaropts(opts, 1, "WhitespaceRule", "preserve");
opts = setvaropts(opts, 1, "EmptyFieldRule", "auto");
dataS1 = readtable("data.xlsx", opts, "UseExcel", false);
clear opts

%regression
%step 1
stockdatasorted100S1=dataS1(325:384,:);
excess_asset=stockdatasorted100S1(:,3:102);
market_excess=stockdatasorted100S1(:,2);
n=width(excess_asset);
beta=zeros(1,n);
p=zeros(1,n);
alpha=zeros(1,n);
alpha_p=zeros(1,n);
for i=(1:n);
a=excess_asset(:,i);
regression=[market_excess a];
regression.Properties.VariableNames = {'ExcessMarket' 'ExcessStock'};
modelFit = fitlm(regression,'ExcessStock ~ ExcessMarket')
beta(i) = modelFit.Coefficients.Estimate(2);
p(i)= modelFit.Coefficients.pValue(2);
alpha(i)=modelFit.Coefficients.Estimate(1);
alpha_p(i)=modelFit.Coefficients.pValue(1);
end
output =[beta;p]
reject_rule1= p<0.05;
percent_of_significant_coeff=sum(reject_rule1)/n

%import average excess return of companies
opts = spreadsheetImportOptions("NumVariables", 1);
opts.Sheet = "average_excess_return";
opts.DataRange = "i2:i101";
opts.VariableNames = "average_excess_return";
opts.VariableTypes = "double";
dataS2 = readtable("data.xlsx", opts, "UseExcel", false);
clear opts

%step 2
beta=beta';
dataS2.beta = beta;
regression2=dataS2;
modelFit1 = fitlm(regression2,'average_excess_return ~ beta')

%% 2008-2018
clear all
%Import return data
opts = spreadsheetImportOptions("NumVariables", 103);
opts.Sheet = "ri (stock data)";
opts.DataRange = "A2:CY397";
opts.VariableNames = ["Date", "MktRF", "RF", "ADAMSEXPRESSCO", "ABBOTTLABORATORIES", "FIRSTFINANCIALBANCORPOHIO", "FIRSTMIDWESTBANCORPDE", "ROGERSCORP", "ROLLINSINC", "SEABOARDCORP", "STSOURCECORP", "SHERWINWILLIAMSCO", "FLEXSTEELINDUSTRIESINC", "FOSTERLBCO", "FRANKLINELECTRICINC", "FRANKLINRESOURCESINC", "FULLERHBCO", "GENTEXCORP", "DONNELLEYRRSONSCO", "AMREPCORP", "STANDARDMOTORPRODUCTSINC", "APACHECORP", "MATTELINC", "BECTONDICKINSONCO", "GRACOINC", "BARNESGROUPINC", "WEYERHAEUSERCO", "INTERNATIONALFLAVORSFRAGINC", "AVONPRODUCTSINC", "HAWAIIANELECTRICINDUSTRIESINC", "HAVERTYFURNITURECOSINC", "HEALTHCARESERVICESGROUPINC", "PARKERHANNIFINCORP", "WEISMARKETSINC", "WOLVERINEWORLDWIDEINC", "BARNWELLINDUSTRIESINC", "SMUCKERJMCO", "MOOGINC", "MEREDITHCORP", "HUNTJBTRANSPORTSERVICESINC", "HUNTINGTONBANCSHARESINC", "HURCOCOMPANY", "SYNALLOYCORP", "MCDONALDSCORP", "STANDEXINTERNATIONALCORP", "VFCORP", "IMMUNOMEDICSINC", "BEMISCOINC", "KENNAMETALINC", "PFINDUSTRIESINC", "TELEFLEXINC", "INTEGRATEDDEVICETECHNOLOGY", "AUTOMATICDATAPROCESSINGINC", "DIODESINC", "INVACARECORP", "MARSHMCLENNANCOSINC", "ISRAMCOINC", "SAFEGUARDSCIENTIFICSINC", "CLOROXCO", "GENUINEPARTSCO", "CABOTCORP", "RITEAIDCORP", "KAMANCORP", "KELLYSERVICESINC", "NEWYORKTIMESCO", "KEYTRONICSCORP", "CNAFINANCIALCORP", "FEDERALSIGNALCORP", "KIMBALLINTERNATIONALINC", "KULICKESOFFAINDSINC", "LAMRESHCORP", "LANCASTERCOLONYCORP", "HUMANAINC", "UNIONPACIFICCORP", "LAWSONPRODUCTSINC", "LINCOLNNATIONALCORP", "SIFCOINDUSTRIESINC", "ALEXANDERSINC", "BLOCKHRINC", "DANAHERCORP", "SERVOTRONICSINC", "LEEENTERPRISESINC", "MTSSYSTEMSCORP", "PENNSYLVANIAREALESTATEINVTTR", "LILLYELICO", "SEMTECHCORP", "WINNEBAGOINDUSTRIESINC", "SONYCORP", "MANITOWOCCOINC", "TERADYNEINC", "MARCUSCORP", "MARINEPETROLEUMTRUST", "SERVICECORPINTERNATIONAL", "SYSCOCORP", "MCCORMICKCOINC", "JACOBSENGINEERINGGROUPINC", "SUPERIORINDUSTRIESINTLINC", "EQUIFAXINC", "COCACOLABOTTLINGCOCONS", "DYCOMINDUSTRIESINC", "GATXCORP", "GENERALDYNAMICSCORP", "GENERALELECTRICCO"];
opts.VariableTypes = ["string", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double"];
opts = setvaropts(opts, 1, "WhitespaceRule", "preserve");
opts = setvaropts(opts, 1, "EmptyFieldRule", "auto");
data = readtable("data.xlsx", opts, "UseExcel", false);
clear opts

%test for normality
stockdatasorted100=data(265:396,:);
return_asset=stockdatasorted100(:,4:103);
dataS1 =stockdatasorted100(:,2);
n=width(return_asset);
jb=zeros(1,n);
pval=zeros(1,n);
skew=zeros(1,n);
kurt=zeros(1,n);
t_stat=zeros(1,n);
T=396
for i=(1:n)
a=return_asset(:,i);
normality=[dataS1 a];
%change name in column 2
normality.Properties.VariableNames = {'ExcessMarket' 'StockReturn'};
[h, pvalue, jbstat, critval] = jbtest(normality.StockReturn)
jb(i)=jbstat
pval(i)=pvalue
skew(i)=skewness(normality.StockReturn);
kurt(i)=kurtosis(normality.StockReturn);
t_stat(i)=T *( skew(i)^2 / 6 + (kurt(i) - 3)^2 / 24 );
end
reject_rule= pval<0.05;
percent_of_abnormal_data=sum(reject_rule)/n

%import excess return data
opts = spreadsheetImportOptions("NumVariables", 102);
opts.Sheet = "ri-rf";
opts.DataRange = "A2:CX397";
opts.VariableNames = ["Date", "MktRF", "ADAMSEXPRESSCO", "ABBOTTLABORATORIES", "FIRSTFINANCIALBANCORPOHIO", "FIRSTMIDWESTBANCORPDE", "ROGERSCORP", "ROLLINSINC", "SEABOARDCORP", "STSOURCECORP", "SHERWINWILLIAMSCO", "FLEXSTEELINDUSTRIESINC", "FOSTERLBCO", "FRANKLINELECTRICINC", "FRANKLINRESOURCESINC", "FULLERHBCO", "GENTEXCORP", "DONNELLEYRRSONSCO", "AMREPCORP", "STANDARDMOTORPRODUCTSINC", "APACHECORP", "MATTELINC", "BECTONDICKINSONCO", "GRACOINC", "BARNESGROUPINC", "WEYERHAEUSERCO", "INTERNATIONALFLAVORSFRAGINC", "AVONPRODUCTSINC", "HAWAIIANELECTRICINDUSTRIESINC", "HAVERTYFURNITURECOSINC", "HEALTHCARESERVICESGROUPINC", "PARKERHANNIFINCORP", "WEISMARKETSINC", "WOLVERINEWORLDWIDEINC", "BARNWELLINDUSTRIESINC", "SMUCKERJMCO", "MOOGINC", "MEREDITHCORP", "HUNTJBTRANSPORTSERVICESINC", "HUNTINGTONBANCSHARESINC", "HURCOCOMPANY", "SYNALLOYCORP", "MCDONALDSCORP", "STANDEXINTERNATIONALCORP", "VFCORP", "IMMUNOMEDICSINC", "BEMISCOINC", "KENNAMETALINC", "PFINDUSTRIESINC", "TELEFLEXINC", "INTEGRATEDDEVICETECHNOLOGY", "AUTOMATICDATAPROCESSINGINC", "DIODESINC", "INVACARECORP", "MARSHMCLENNANCOSINC", "ISRAMCOINC", "SAFEGUARDSCIENTIFICSINC", "CLOROXCO", "GENUINEPARTSCO", "CABOTCORP", "RITEAIDCORP", "KAMANCORP", "KELLYSERVICESINC", "NEWYORKTIMESCO", "KEYTRONICSCORP", "CNAFINANCIALCORP", "FEDERALSIGNALCORP", "KIMBALLINTERNATIONALINC", "KULICKESOFFAINDSINC", "LAMRESHCORP", "LANCASTERCOLONYCORP", "HUMANAINC", "UNIONPACIFICCORP", "LAWSONPRODUCTSINC", "LINCOLNNATIONALCORP", "SIFCOINDUSTRIESINC", "ALEXANDERSINC", "BLOCKHRINC", "DANAHERCORP", "SERVOTRONICSINC", "LEEENTERPRISESINC", "MTSSYSTEMSCORP", "PENNSYLVANIAREALESTATEINVTTR", "LILLYELICO", "SEMTECHCORP", "WINNEBAGOINDUSTRIESINC", "SONYCORP", "MANITOWOCCOINC", "TERADYNEINC", "MARCUSCORP", "MARINEPETROLEUMTRUST", "SERVICECORPINTERNATIONAL", "SYSCOCORP", "MCCORMICKCOINC", "JACOBSENGINEERINGGROUPINC", "SUPERIORINDUSTRIESINTLINC", "EQUIFAXINC", "COCACOLABOTTLINGCOCONS", "DYCOMINDUSTRIESINC", "GATXCORP", "GENERALDYNAMICSCORP", "GENERALELECTRICCO"];
opts.VariableTypes = ["string", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double"];
opts = setvaropts(opts, 1, "WhitespaceRule", "preserve");
opts = setvaropts(opts, 1, "EmptyFieldRule", "auto");
dataS1 = readtable("data.xlsx", opts, "UseExcel", false);
clear opts

%regression
%step 1
stockdatasorted100S1=dataS1(265:396,:);
excess_asset=stockdatasorted100S1(:,3:102);
market_excess=stockdatasorted100S1(:,2);
n=width(excess_asset);
beta=zeros(1,n);
p=zeros(1,n);
alpha=zeros(1,n);
alpha_p=zeros(1,n);
for i=(1:n);
a=excess_asset(:,i);
regression=[market_excess a];
regression.Properties.VariableNames = {'ExcessMarket' 'ExcessStock'};
modelFit = fitlm(regression,'ExcessStock ~ ExcessMarket')
beta(i) = modelFit.Coefficients.Estimate(2);
p(i)= modelFit.Coefficients.pValue(2);
alpha(i)=modelFit.Coefficients.Estimate(1);
alpha_p(i)=modelFit.Coefficients.pValue(1);
end
output =[beta;p]
reject_rule1= p<0.05;
percent_of_significant_coeff=sum(reject_rule1)/n

%import average excess return of companies
opts = spreadsheetImportOptions("NumVariables", 1);
opts.Sheet = "average_excess_return";
opts.DataRange = "c2:c101";
opts.VariableNames = "average_excess_return";
opts.VariableTypes = "double";
dataS2 = readtable("data.xlsx", opts, "UseExcel", false);
clear opts

%step 2
beta=beta';
dataS2.beta = beta;
regression2=dataS2;
modelFit1 = fitlm(regression2,'average_excess_return ~ beta')

%% 1963-2018 (portfolio data) with interaction term 
clear all
%Import return data
opts = spreadsheetImportOptions("NumVariables", 34);
opts.Sheet = "rp (portfolio)";
opts.DataRange = "A2:AH668";
opts.VariableNames = ["VarName1", "RF", "SMALLLoOPLoINV", "ME1OP1INV2", "ME1OP1INV3", "SMALLLoOPHiINV", "ME1OP2INV1", "ME1OP2INV2", "ME1OP2INV3", "ME1OP2INV4", "ME1OP3INV1", "ME1OP3INV2", "ME1OP3INV3", "ME1OP3INV4", "SMALLHiOPLoINV", "ME1OP4INV2", "ME1OP4INV3", "SMALLHiOPHiINV", "BIGLoOPLoINV", "ME2OP1INV2", "ME2OP1INV3", "BIGLoOPHiINV", "ME2OP2INV1", "ME2OP2INV2", "ME2OP2INV3", "ME2OP2INV4", "ME2OP3INV1", "ME2OP3INV2", "ME2OP3INV3", "ME2OP3INV4", "BIGHiOPLoINV", "ME2OP4INV2", "ME2OP4INV3", "BIGHiOPHiINV"];
opts.VariableTypes = ["double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double"];
data32portfolios = readtable("data.xlsx", opts, "UseExcel", false);
clear opts

%test for normality
return_asset=data32portfolios(:,4:34);
market_excess =data32portfolios(:,2);
n=width(return_asset);
jb=zeros(1,n);
pval=zeros(1,n);
skew=zeros(1,n);
kurt=zeros(1,n);
t_stat=zeros(1,n);
T=height(return_asset);
for i=(1:n)
a=return_asset(:,i);
normality=[market_excess a];
normality.Properties.VariableNames = {'ExcessMarket' 'StockReturn'};
[h, pvalue, jbstat, critval] = jbtest(normality.StockReturn)
jb(i)=jbstat
pval(i)=pvalue
skew(i)=skewness(normality.StockReturn);
kurt(i)=kurtosis(normality.StockReturn);
t_stat(i)=T *( skew(i)^2 / 6 + (kurt(i) - 3)^2 / 24 );
end
reject_rule= pval<0.05;
percent_of_abnormal_data=sum(reject_rule)/n

%regression
%import data
opts = spreadsheetImportOptions("NumVariables", 34);
opts.Sheet = "rp-rf";
opts.DataRange = "A2:AH668";
opts.VariableNames = ["VarName1", "MktRF", "SMALLLoOPLoINV", "ME1OP1INV2", "ME1OP1INV3", "SMALLLoOPHiINV", "ME1OP2INV1", "ME1OP2INV2", "ME1OP2INV3", "ME1OP2INV4", "ME1OP3INV1", "ME1OP3INV2", "ME1OP3INV3", "ME1OP3INV4", "SMALLHiOPLoINV", "ME1OP4INV2", "ME1OP4INV3", "SMALLHiOPHiINV", "BIGLoOPLoINV", "ME2OP1INV2", "ME2OP1INV3", "BIGLoOPHiINV", "ME2OP2INV1", "ME2OP2INV2", "ME2OP2INV3", "ME2OP2INV4", "ME2OP3INV1", "ME2OP3INV2", "ME2OP3INV3", "ME2OP3INV4", "BIGHiOPLoINV", "ME2OP4INV2", "ME2OP4INV3", "BIGHiOPHiINV"];
opts.VariableTypes = ["double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double"];
data32portfoliosS1 = readtable("data.xlsx", opts, "UseExcel", false);
clear opts

%step 1
excess_asset=data32portfoliosS1(:,3:34);
market_excess=data32portfoliosS1(:,2);
n=width(excess_asset);
beta=zeros(1,n);
p=zeros(1,n);
alpha=zeros(1,n);
alpha_p=zeros(1,n);
for i=(1:n);
a=excess_asset(:,i);
%regression.a=a;
regression=[market_excess a];
%change name in column 2
regression.Properties.VariableNames = {'ExcessMarket' 'ExcessStock'};
modelFit = fitlm(regression,'ExcessStock ~ ExcessMarket')
beta(i) = modelFit.Coefficients.Estimate(2);
p(i)= modelFit.Coefficients.pValue(2);
alpha(i)=modelFit.Coefficients.Estimate(1);
alpha_p(i)=modelFit.Coefficients.pValue(1);
end
output =[beta;p]
reject_rule1= p<0.05;
percent_of_significant_coeff=sum(reject_rule1)/n

%import data
opts = spreadsheetImportOptions("NumVariables", 1);
opts.Sheet = "average_excess_return_p";
opts.DataRange = "B2:B33";
opts.VariableNames = "average_excess";
opts.VariableTypes = "double";
data32portfoliosS3 = readtable("data.xlsx", opts, "UseExcel", false);
clear opts

%step 2
beta=beta';
data32portfoliosS3.beta = beta;
regression2=data32portfoliosS3;
modelFit1 = fitlm(regression2,'average_excess~ beta')
%gamma0 is significantly different from 0.
%When evidence indicates that the excess return is not exclusively explained by beta, an alternative model should be build, which incorporates other factors.
%b is not significant, CAPM does not work