事件研究法-python-萬礦

moser_freshman發表於2019-07-05
事件Python
import pandas as pd
from WindPy import *
import numpy as np
import scipy as sp
import matplotlib.pyplot as plt
from scipy.optimize import leastsq

w.start()
data_df_QTJY=w.wsd("300359.SZ", "close", "2018-12-01", "2019-02-28", "Fill=Previous;PriceAdj=F")
data_df_CYBZ=w.wsd("399006.SZ", "close", "2018-12-01", "2019-02-28", "Fill=Previous;PriceAdj=F")

data_df_QTJY=pd.DataFrame(data_df_QTJY.Data,columns=data_df_QTJY.Times,index=data_df_QTJY.Fields)
data_df_QTJY=data_df_QTJY.T
data_df_CYBZ=pd.DataFrame(data_df_CYBZ.Data,columns=data_df_CYBZ.Times,index=data_df_CYBZ.Fields)
data_df_CYBZ=data_df_CYBZ.T

data_df_QTJY['return']=np.log(data_df_QTJY['CLOSE']/data_df_QTJY['CLOSE'].shift(1))
data_df_QTJY['index_return']=np.log(data_df_CYBZ['CLOSE']/data_df_CYBZ['CLOSE'].shift(1))
data_df_QTJY=data_df_QTJY.dropna()

##樣本資料(Xi,Yi)，需要轉換成陣列(列表)形式
xi=np.array(data_df_QTJY['index_return'])
yi=np.array(data_df_QTJY['return'])
xi
##需要擬合的函式func 
def func(p,x):
    k,b=p
    return k*x+b
 
##偏差函式：x,y都是列表:這裡的x,y更上面的Xi,Yi中是一一對應的
def error(p,x,y):
    return func(p,x)-y
 
#k,b的初始值，可以任意設定,經過幾次試驗，發現p0的值會影響cost的值：Para[1]
p0=[0.5,0.03]
 
#把error函式中除了p0以外的引數打包到args中(使用要求)
Para=leastsq(error,p0,args=(xi,yi))
Para
#讀取結果
k,b=Para[0]
print("k=",k,"b=",b)
 
#畫樣本點
plt.figure(figsize=(8,6)) ##指定影像比例： 8：6
plt.scatter(xi,yi,color="green",label="樣本資料",linewidth=2) 
 
#畫擬合直線
x=np.linspace(-0.1,0.1,100) ##在150-190直接畫100個連續點
y=k*x+b ##函式式
plt.plot(x,y,color="red",label="擬合直線",linewidth=2) 
plt.legend() #繪製圖例
plt.show()
def get_ar(k):
    # 計算視窗期 AR
    data_df_QTJY=w.wsd("300359.SZ", "close", "2019-04-01", "2019-04-15", "PriceAdj=F")
    data_df_CYBZ=w.wsd("399006.SZ", "close", "2019-04-01", "2019-04-15", "PriceAdj=F")
    data_df_QTJY=pd.DataFrame(data_df_QTJY.Data,columns=data_df_QTJY.Times,index=data_df_QTJY.Fields)
    data_df_QTJY=data_df_QTJY.T
    data_df_CYBZ=pd.DataFrame(data_df_CYBZ.Data,columns=data_df_CYBZ.Times,index=data_df_CYBZ.Fields)
    data_df_CYBZ=data_df_CYBZ.T
    data_df_QTJY['return']=np.log(data_df_QTJY['CLOSE']/data_df_QTJY['CLOSE'].shift(1))
    data_df_QTJY['index_return']=np.log(data_df_CYBZ['CLOSE']/data_df_CYBZ['CLOSE'].shift(1))
    data_df_QTJY=data_df_QTJY.dropna()
    data_df_QTJY['moni_return']=(data_df_QTJY['index_return']-0.03/365)*k
    data_df_QTJY['ar']=data_df_QTJY['return']-data_df_QTJY['moni_return']
    
    return data_df_QTJY

k=0.9595205659434853
df=get_ar(k=k)    
def get_ar(k):
    # 計算視窗期 AR
    data_df_QTJY=w.wsd("300359.SZ", "close", "2019-04-01", "2019-04-15", "PriceAdj=F")
    data_df_CYBZ=w.wsd("399006.SZ", "close", "2019-04-01", "2019-04-15", "PriceAdj=F")
    data_df_QTJY=pd.DataFrame(data_df_QTJY.Data,columns=data_df_QTJY.Times,index=data_df_QTJY.Fields)
    data_df_QTJY=data_df_QTJY.T
    data_df_CYBZ=pd.DataFrame(data_df_CYBZ.Data,columns=data_df_CYBZ.Times,index=data_df_CYBZ.Fields)
    data_df_CYBZ=data_df_CYBZ.T
    data_df_QTJY['return']=np.log(data_df_QTJY['CLOSE']/data_df_QTJY['CLOSE'].shift(1))
    data_df_QTJY['index_return']=np.log(data_df_CYBZ['CLOSE']/data_df_CYBZ['CLOSE'].shift(1))
    data_df_QTJY=data_df_QTJY.dropna()
    data_df_QTJY['moni_return']=(data_df_QTJY['index_return']-0.03/365)*k
    data_df_QTJY['ar']=data_df_QTJY['return']-data_df_QTJY['moni_return']
    
    return data_df_QTJY

k=0.9595205659434853
df=get_ar(k=k)    

# t檢驗
from scipy import stats
result=stats.ttest_1samp(df['ar'],0)
car=[0.0986 ,0.0933 ,0.0512 ,0.0344 ,0.1287 ,0.0957 ,0.0837 ,0.1155 ,0.0916]
print('car顯著性檢驗',stats.ttest_1samp(car,1))
print('單樣本檢驗結果',result)
result2=stats.ttest_ind(df['return'],df['moni_return'],equal_var=False)
print('雙樣本檢驗結果',result2)
result3=stats.levene(df['return'],df['moni_return'])
print('原假設方差不同',result3)
x=[1,2,3,4,5,6,7,8,9]
plt.plot(x,df['return'],color='blue',label='actual return')
plt.plot(x,df['moni_return'],color='black',label='semulated_return')
plt.plot(x,np.zeros(9),color='red')
plt.legend() #繪製圖例
plt.show()
x=[1,2,3,4,5,6,7,8,9]
plt.plot(x,df['ar'],color='blue',label='ar')
car=[0.0986 ,0.0933 ,0.0512 ,0.0344 ,0.1287 ,0.0957 ,0.0837 ,0.1155 ,0.0916]
plt.plot(x,car,color='black',label='car')
plt.plot(x,np.zeros(9),color='red')
plt.legend() #繪製圖例
plt.show()

相關文章