import numpy as np
import pandas as pd
import scipy as sp
import scipy.stats as ss
import statsmodels
import statsmodels.api as sm
import statsmodels.formula.api as smf
from statsmodels.iolib.summary2 import summary_col
from wooldridge import *
df = dataWoo("phillips")
df = df[(df['year']<1997)]
df['uhat1'] = smf.ols('df.inf ~ unem + 1', data=df).fit().resid
print(smf.ols('uhat1 ~ uhat1.shift(1)', data=df).fit().summary())
df['uhat2'] = smf.ols('df.cinf ~ unem', data=df).fit().resid
print(smf.ols('uhat2 ~ uhat2.shift(1)', data=df).fit().summary())
df = dataWoo("prminwge")
df['uhat'] = smf.ols('lprepop ~ lmincov + lprgnp + lusgnp + t', data=df).fit().resid
AR1c = smf.ols('uhat ~ lmincov + lprgnp + lusgnp + t + uhat.shift(1)' , data=df).fit()
AR1s = smf.ols('uhat ~ uhat.shift(1)', data=df).fit()
print(summary_col([AR1c, AR1s],stars=True,float_format='%0.3f',
model_names=['AR1c\n(b/se)','AR1s\n(b/se)'],
info_dict={'N':lambda x: "{0:d}".format(int(x.nobs)),
'R2':lambda x: "{:.3f}".format(x.rsquared),
'Adj.R2':lambda x: "{:.3f}".format(x.rsquared_adj)}))
df = dataWoo("barium")
df['u'] = smf.ols('lchnimp ~ lchempi + lgas + lrtwex + befile6 + affile6 + afdec6', data=df).fit().resid
AR3 = smf.ols('u ~ lchempi + lgas + lrtwex + befile6 + affile6 + afdec6 + u.shift(1) + u.shift(2) + u.shift(3) + 1', data = df).fit()
print(AR3.summary())
hypotheses = '(u.shift(1) = u.shift(2) = u.shift(3) = 0)'
f_test = AR3.f_test(hypotheses)
print(f_test)
df = dataWoo("barium")
OLS = smf.ols('lchnimp ~ lchempi + lgas + lrtwex + befile6 + affile6 + afdec6', data=df).fit()
def ols_ar1(model,rho,drop1=True):
x = model.model.exog
y = model.model.endog
ystar = y[1:]-rho*y[:-1]
xstar = x[1:,]-rho*x[:-1,]
if drop1 == False:
ystar = np.append(np.sqrt(1-rho**2)*y[0],ystar)
xstar = np.append([np.sqrt(1-rho**2)*x[0,]],xstar,axis=0)
model_ar1 = sm.OLS(ystar,xstar).fit()
return(model_ar1)
def OLSAR1(model,drop1=True):
x = model.model.exog
y = model.model.endog
e = y-x@model.params
e1 = e[:-1]; e0 = e[1:]
rho0 = np.dot(e1,e[1:])/np.dot(e1,e1)
rdiff = 1.0
while(rdiff>1.0e-5):
model1 = ols_ar1(model,rho0,drop1)
e = y - (x @ model1.params)
e1 = e[:-1]; e0 = e[1:]
rho1 = np.dot(e1,e[1:])/np.dot(e1,e1)
rdiff = np.sqrt((rho1-rho0)**2)
rho0 = rho1
print('Rho = ', rho0)
model1 = ols_ar1(model,rho0,drop1)
return(model1)
ar1_pw = OLSAR1(OLS ,drop1=False)
print(ar1_pw.summary())
df = dataWoo("phillips")
df = df[(df['year']<1997)]
ols = smf.ols('df.inf ~ unem', data=df).fit()
print(ols.summary())
ar1_pw = OLSAR1(ols ,drop1=False)
print(ar1_pw.summary())
df = dataWoo("intdef")
def1 = df['def']
ureg = smf.ols('i3 ~ df.inf + def1', data=df).fit()
df['u'] = ureg.resid
AR1u = smf.ols('u ~ u.shift(1)', data=df).fit()
ereg = smf.ols('ci3 ~ df.cinf + cdef', data=df).fit()
df['e'] = ereg.resid
AR1e = smf.ols('e ~ e.shift(1)', data=df).fit()
print(summary_col([AR1u, AR1e, ureg, ereg],stars=True,float_format='%0.3f',
model_names=['AR1u\n(b/se)','AR1e\n(b/se)', 'ureg\n(b/se)','ereg\n(b/se)'],
info_dict={'N':lambda x: "{0:d}".format(int(x.nobs)),
'R2':lambda x: "{:.3f}".format(x.rsquared),
'Adj.R2':lambda x: "{:.3f}".format(x.rsquared_adj)}))
df = dataWoo("prminwge")
OLS2 =smf.ols('lprepop ~ lmincov + lprgnp + lusgnp + t', data=df).fit()
Newey = OLS2.get_robustcov_results(cov_type='HAC',maxlags=1)
print(summary_col([OLS, Newey],stars=True,float_format='%0.3f',
model_names=['OLS\n(b/se)','Newey\n(b/se)'],
info_dict={'N':lambda x: "{0:d}".format(int(x.nobs)),
'R2':lambda x: "{:.3f}".format(x.rsquared),
'Adj.R2':lambda x: "{:.3f}".format(x.rsquared_adj)}))
print(OLSAR1(OLS2 ,drop1=False).summary()) #PW
df = dataWoo("nyse")
reg1 = smf.ols('df[("return")] ~ return_1', data=df).fit()
print(reg1.summary())
df['u2'] = np.square(reg1.resid)
print(smf.ols('u2 ~ return_1', data=df).fit().summary())
df = dataWoo("nyse")
df['u2']= np.square(smf.ols('df[("return")] ~ return_1', data=df).fit().resid)
print(smf.ols('u2 ~ u2.shift(1)', data=df).fit().summary())
df['u']= smf.ols('df[("return")] ~ return_1', data=df).fit().resid
print(smf.ols('u ~ u.shift(1)', data=df).fit().summary())