In [1]:
%matplotlib inline
import numpy as np
from matplotlib.pyplot import plot
from matplotlib.pyplot import show
In [2]:
# 首先读入两只股票的收盘价，并计算收益率
bhp_cp = np.loadtxt('BHP.csv', delimiter=',', usecols=(6,), unpack=True)
vale_cp = np.loadtxt('VALE.csv', delimiter=',', usecols=(6,), unpack=True)
bhp_returns = np.diff(bhp_cp) / bhp_cp[:-1]
vale_returns = np.diff(vale_cp) / vale_cp[:-1]
协方差描述的是两个变量共同变化的趋势，其实就是归一化前的相关系数。 使用cov函数计算股票收益率的协方差矩阵
In [3]:
covariance = np.cov(bhp_returns, vale_returns)
print 'Covariance:\n', covariance
In [4]:
# 查看协方差矩阵对角线的元素
print 'Covariance diagonal:\n', covariance.diagonal()
# 计算矩阵的迹，即对角线之和
print 'Covariance trace:\n', covariance.trace()
# 计算相关系数，相关系数是协方差除以各自标准差的乘积
print 'Correlation coefficient:\n', covariance / (bhp_returns.std() * vale_returns.std())
用相关系数来度量两只股票的相关程度。相关系数的取值范围在-1到1之间，一组数据域自身的相关系数为1.使用corrcoef函数计算相关系数。
In [5]:
# 使用corrcoef计算更加精确
print 'Correlation coefficient:\n', np.corrcoef(bhp_returns, vale_returns)
相关系数矩阵关于对角线对称，BHP与VALE的相关系数等于VALE和BHP的相关系数。看起来0.68的相关系数表示他们的相关程度似乎不是很强。
In [6]:
difference = bhp_cp - vale_cp
avg = np.mean(difference)
dev = np.std(difference)
# 检查最后一次收盘价是否在同步状态
print "Out of sync : ", np.abs(difference[-1] - avg) > 2*dev
这说明，最后一次收盘价不再同步状态，我们暂时不能进行交易
In [7]:
# 绘制收益率曲线
t = np.arange(len(bhp_returns))
plot(t, bhp_returns, lw=1)
plot(t, vale_returns, lw=2)
show()
In [8]:
# 用三次多项式去拟合两只股票收盘价的差价
t = np.arange(len(bhp_cp))
poly = np.polyfit(t, bhp_cp-vale_cp, 3)
print "Polynomial fit\n", poly
In [9]:
# 用刚才得到的多项式对象，推断下一个值
print "Next value: ", np.polyval(poly, t[-1]+1)
理想情况下，BHP和VALE股票收盘价的差价越小越好。在极限情况下，差值可以在某个点为0。用roots函数找到拟合多项式函数在什么时候达到0。
In [10]:
print "Roots: ", np.roots(poly)
In [11]:
# 极值位于导数为0的点
der = np.polyder(poly)
print "Dervative:\n", der
# 得到多项式导函数的系数
In [12]:
# 求出导函数的根，即找出原多项式函数的极值点
print "Extremas: ", np.roots(der)
In [13]:
# 通过argmax和argmin函数找到最大最小值点来检查结果
vals = np.polyval(poly, t)
print "Maximum index: ", np.argmax(vals)
print "Minimum index: ", np.argmin(vals)
In [14]:
plot(t, bhp_cp-vale_cp)
plot(t, vals)
show()
我们需要在成交量前面乘上一个有收盘价变化决定的正负号。
In [15]:
cp, volume = np.loadtxt('BHP.csv', delimiter=',', usecols=(6,7), unpack=True)
change = np.diff(cp)
print "Change:", change
使用NumPy的sign函数返回每个元素的正负号。
In [16]:
signs = np.sign(change)
print "Signs:\n", signs
使用Numpy的piecewise函数获取数组元素的正负。piecewise(分段)，可以根据给定取值，得到分段。
In [17]:
pieces = np.piecewise(change, [change<0, change>0], [-1,1])
print "Pieces:\n", pieces
In [18]:
# 检查两次输出是否一致
print "Arrays equal?", np.array_equal(signs, pieces)
In [19]:
# OBV值的计算依赖于前一日的收盘价
print "On balance volume: \n", volume[1:]*signs
使用vectorize函数可以减少你的程序中使用循环的次数，NumPy中的vectorize函数相当于Python中的map函数。我们用它来计算单个交易日的利润。
In [20]:
# 读入数据
# op is opening price,hp is the highest price
# lp is the lowest price, cp is closing price
op, hp, lp, cp = np.loadtxt('BHP.csv', delimiter=',', usecols=(3,4,5,6), unpack=True)
我们尝试以比开盘价稍低一点的价格买入股票。如果这个价格不在当日的股价范围内，则尝试买入失败，没有获利，也没有亏损，我们返回0。否则，我们将以当日收盘价卖出，所获得的利润即买入卖出的差价。
In [24]:
def calc_profit(op, high, low, close):
# 以开盘价买入,这里不考虑买入多少股
buy = op
if low < buy < high:
return (close-buy) / buy
else:
return 0
In [25]:
# 矢量化一个函数，这样可以避免使用循环
func = np.vectorize(calc_profit)
profits = func(op, hp, lp, cp)
print 'Profits:\n', profits
In [23]:
# 我们选择非零利润的交易日并计算平均值
real_trades = profits[profits != 0]
print 'Number of trades:\n', len(real_trades), round(100.0 * len(real_trades)/len(cp), 2),"%"
print "Average profit/loss % :", round(np.mean(real_trades) * 100, 2)
In [ ]: