基于 Python 实现股票配对交易统计套利

作者：老余捞鱼

原创不易，转载请标明出处及原作者。

写在前面的话：本文主要介绍使用 Python 进行配对交易统计套利的策略，通过识别协整的股票对，并利用统计方法来预测价格回归均值，从而实现市场中性交易。

一、什么是配对交易统计套利

统计套利（Statistical Arbitrage，简称Stat Arb）是金融领域中一种利用数学模型和统计方法来识别并利用市场价格差异的高级交易策略。这种方法具有高度的定量性和分析性，通常涉及广泛分散的证券投资组合。

其核心特点一是依赖统计分析：统计套利依赖于大量的历史数据和复杂的统计分析方法，如协整分析、相关性分析等，来识别价格之间的关系和潜在的套利机会。二是自动化交易：由于需要处理大量数据并及时捕捉市场机会，统计套利策略往往通过计算机算法实现自动化或半自动化交易。三是低风险高收益：理论上，统计套利策略能够在较低风险的情况下获得较为稳定的收益，因为它依赖于价格回归均值的规律。

简单来说，配对交易统计套利是一种市场中性策略，利用统计技术来识别历史上表现出很强相关性的两只股票。当这两只股票的价格关系出现偏离时，我们预计其将回归均值。为此，我们买入表现不佳的股票，同时卖出表现出色的股票。如果我们的均值回归假设成立，价格应该会回归到长期平均水平，从而实现盈利交易。但是，如果价格背离不是暂时的，而是由结构性因素造成的，那么就有很大的亏损风险。

二、初始步骤

第一步是确定股票范围，并找出具有高度相关性的股票对。关键在于，这种相关性必须基于经济关系，例如在同一行业内运营的公司；否则，这种相关性可能是虚假的。在本次分析中，我选择了 NSE-100 指数中所有被归类为“金融服务”公司的成分股。这一选择为我们提供了一份包含 25 家公司的初始名单。然而，在剔除了每日定价数据不足 10 年的公司后，我们最终得到了 15 只股票。接下来，我们获取了这 15 只股票的每日收盘价，并将数据划分为训练集和测试集。这种划分方式确保了我们对协整对的选择是基于训练数据集，而回溯测试则使用样本外测试数据集。第一步，我们使用皮尔逊相关系数来初步了解这些股票之间的关系，然后使用 statsmodels.tsa.stattools 中的 coint 函数来识别协整对。coint 函数会返回每对股票的协整检验 p 值。这些 p 值被存储在一个数组中，并以热图的形式呈现。p 值小于 0.05 表示我们可以拒绝零假设，这意味着两个不同符号的时间序列可能存在协整关系。代码如下：

# make the necessary imports
import pandas as pd
import numpy as np
import matplotlib.pyplot as plt
import seaborn as sns
import datetime
idx = pd.IndexSlice
import statsmodels.api as smfrom statsmodels.regression.linear_model import OLS
from statsmodels.tsa.stattools import adfuller
from statsmodels.tsa.stattools import coint
from sklearn.model_selection import train_test_split%matplotlib inline


%config InlineBackend.figure_format = ‘retina’# read the matadata csv


nifty_meta = pd.read_csv('data/nifty_meta.csv')
nifty_meta.head(2)# get the ticker list with industry is equal to FINANCIAL SERVICES
tickers = list(nifty_meta[nifty_meta.Industry=='FINANCIAL SERVICES'].Symbol)

print(tickers)
print(len(tickers))# start and end dates for backtesting

fromdate = datetime.datetime(2010, 1, 1)
todate = datetime.datetime(2020, 6, 15)# read back the pricing data
prices = pd.read_csv('data/prices.csv', index_col=['ticker','date'], parse_dates=True)
prices.head(2)# remove tickers where we have less than 10 years of data.
min_obs = 2520
nobs = prices.groupby(level='ticker').size()
keep = nobs[nobs>min_obs].indexprices = prices.loc[idx[keep,:], :]
prices.info()# final tickers list
TICKERS = list(prices.index.get_level_values('ticker').unique())

print(len(TICKERS))
print(TICKERS)# unstack and take close price

close = prices.unstack('ticker')['close'].sort_index()
close = close.dropna()
close.head(2)# train test split 
train_close, test_close = train_test_split(close, test_size=0.5, shuffle=False)# quick view of head and tail of train set
train_close.head(2).append(train_close.tail(2))# Pearson correlation to get the basic idea about the relationship
fig, ax = plt.subplots(figsize=(10,7))
sns.heatmap(train_close.pct_change().corr(method ='pearson'), ax=ax, cmap='coolwarm', annot=True, fmt=".2f") #spearman
ax.set_title('Assets Correlation Matrix')
plt.savefig('images/chart1', dpi=300)# function to find cointegrated pairs

def find_cointegrated_pairs(data):
    n = data.shape[1]
    pvalue_matrix = np.ones((n, n))
    keys = data.keys()
    pairs = []
    for i in range(n):
        for j in range(i+1, n):
            result = coint(data[keys[i]], data[keys[j]])
            pvalue_matrix[i, j] = result[1]
            if result[1] < 0.05:
                pairs.append((keys[i], keys[j]))
    return pvalue_matrix, pairs# calculate p-values and plot as a heatmap

pvalues, pairs = find_cointegrated_pairs(train_close)
print(pairs)
fig, ax = plt.subplots(figsize=(10,7))
sns.heatmap(pvalues, xticklabels = train_close.columns,
                yticklabels = train_close.columns, cmap = 'RdYlGn_r', annot = True, fmt=".2f",
                mask = (pvalues >= 0.99))

ax.set_title('Assets Cointregration Matrix p-values Between Pairs')
plt.tight_layout()
plt.savefig('images/chart2', dpi=300)

代码解释：

1. 我们首先导入所需的 Python 库。
2. 接下来，我们加载 nifty_meta.csv 文件，创建行业分类为金融服务的股票列表，并定义开始和结束日期。
3. 然后，我们从 prices.csv 中读取每日定价数据，并排除数据不足 10 年的股票。
4. 定价数据帧经过重塑，只包含每日收盘价，然后分为训练数据集和测试数据集，各占数据的 50%。
5. 我们绘制了每日收益率的皮尔逊相关性，以初步了解股票之间的关系。
6. 定义了一个名为 find_cointegrated_pairs 的函数，用于识别协整对及其相关的 p 值。
7. 最后，我们通过绘制热图将 p 值可视化，如下图。

皮尔逊相关系数（Pearson correlation coefficient）从 +1 到 -1 不等，是两个变量之间关系的线性指标。+1表示强正相关，0表示无相关，-1表示强负相关。从上面的热图可以看出，有几对变量呈现出很强的正相关性，许多对的 p 值都低于 0.05。这表明，对于这些配对它们可能存在协整关系。

三、对所选数据对进行静态检验

在这一阶段，我们的策略有几个潜在的交易对，它们的 p 值都在 0.05 以下。选择正确的交易对至关重要，因为如果价格走势完全一致，策略就不会有效。要使策略盈利，价格必须出现背离，然后回归均值。

接下来，我们来深入分析 BANKBARODA 和 SBIN 这两只股票，并通过增强 Dickey-Fuller 检验来进一步考察它们之间价差的平稳性。价差必须是平稳的。当一个时间序列的参数（如均值和方差）在一段时间内保持恒定，且不存在单位根时，该时间序列就被视为平稳的。我们首先运用 OLS 回归法来计算这两只股票之间的套期保值比率。然后，依据套期保值比率，我们会计算出差价，并进行增强 Dickey-Fuller 检验。

# final pair to test strategy
asset1 = ‘BANKBARODA’
asset2 = ‘SBIN’

# create a train dataframe of 2 assets
train = pd.DataFrame()
train['asset1'] = train_close[asset1]
train['asset2'] = train_close[asset2]

# visualize closing prices
ax = train[['asset1','asset2']].plot(figsize=(12, 6), title = 'Daily Closing Prices for {} and {}'.format(asset1,asset2))
ax.set_ylabel("Closing Price")
ax.grid(True);
plt.savefig('images/chart3', dpi=300)

# run OLS regression
model=sm.OLS(train.asset2, train.asset1).fit()

# print regression summary results
plt.rc('figure', figsize=(12, 7))
plt.text(0.01, 0.05, str(model.summary()), {'fontsize': 16}, fontproperties = 'monospace')
plt.axis('off')
plt.tight_layout()
plt.subplots_adjust(left=0.2, right=0.8, top=0.7, bottom=0.1)
plt.savefig('images/chart4', dpi=300);print('Hedge Ratio = ', model.params[0])

# calculate spread
spread = train.asset2 - model.params[0] * train.asset1

# Plot the spread
ax = spread.plot(figsize=(12, 6), title = "Pair's Spread")
ax.set_ylabel("Spread")
ax.grid(True);

# conduct Augmented Dickey-Fuller test
adf = adfuller(spread, maxlag = 1)
print('Critical Value = ', adf[0])

# probablity critical values
print(adf[4])

代码解释：

1. 指定 BANKBARODA 为 asset1 （资产 1），SBIN 为asset2（资产 2）
2. 从训练数据集中创建一个包含这两只股票收盘价的数据帧，并将其可视化。
3. 进行 OLS 回归，得出斜率系数，作为我们的对冲比率。
4. 计算价差并绘制可视化图。
5. 执行增强 Dickey-Fuller 检验，以评估价差的静态性，并检查是否存在单位根。

从上图（所选股票对的每日收盘价）可以明显看出，这两只股票的收盘价走势趋于接近。下图是成果输出：

下图是对的传播（Pair’s Spread）

OLS 回归得出的高 R 平方值和接近零的 p 值，显示这两只股票之间存在紧密的相关性。价差似乎是平稳的，增强型 Dickey-Fuller 检验的临界值为-3.459，低于 1%显著性水平的-3.435。这使我们能够否定价差包含单位根的假设，从而确认价差本质上是稳定的。

四、利用 Z 值生成交易信号

到目前为止，我们已经使用训练数据集为我们的策略选择了股票配对。接下来，我们将依靠测试数据集生成交易信号，并使用样本外数据集进行回溯测试。我们将使用两只股票价格比率的 z 值来生成交易信号，并定义上下限阈值。Z 值表示当前价格与平均价格的偏离程度。如果 z 值为正且超过上下限，则表示股价高于平均值。因此，我们预计价格会下跌，促使我们做空（卖出）这只股票，做多（买入）另一只股票。

# calculate z-score
def zscore(series):
 return (series — series.mean()) / np.std(series)

# create a dataframe for trading signals
signals = pd.DataFrame()
signals['asset1'] = test_close[asset1] 
signals['asset2'] = test_close[asset2]
ratios = signals.asset1 / signals.asset2

# calculate z-score and define upper and lower thresholds
signals['z'] = zscore(ratios)
signals['z upper limit'] = np.mean(signals['z']) + np.std(signals['z'])
signals['z lower limit'] = np.mean(signals['z']) - np.std(signals['z'])

# create signal - short if z-score is greater than upper limit else long
signals['signals1'] = 0
signals['signals1'] = np.select([signals['z'] > \
                                 signals['z upper limit'], signals['z'] < signals['z lower limit']], [-1, 1], default=0)

# we take the first order difference to obtain portfolio position in that stock
signals['positions1'] = signals['signals1'].diff()
signals['signals2'] = -signals['signals1']
signals['positions2'] = signals['signals2'].diff()

# verify datafame head and tail
signals.head(3).append(signals.tail(3))

# visualize trading signals and position
fig=plt.figure(figsize=(14,6))
bx = fig.add_subplot(111)   
bx2 = bx.twinx()

#plot two different assets
l1, = bx.plot(signals['asset1'], c='#4abdac')
l2, = bx2.plot(signals['asset2'], c='#907163')u1, = bx.plot(signals['asset1'][signals['positions1'] == 1], lw=0, marker='^', markersize=8, c='g',alpha=0.7)

代码解释：

1. 定义一个名为 zscore 的函数来计算 z 分数。
2. 使用测试数据集中的收盘价，为两只股票创建名为信号的数据帧，并计算它们的价格比。
3. 计算该比率的 Z 值，并根据正负一个标准差确定上下限。
4. 根据这一逻辑添加信号列：如果 z 分数超过上阈值，则赋值 -1 （短信号）；如果 z 分数低于下阈值，则赋值 +1 （长信号）；否则，默认为 0 表示无信号。
5. 计算信号列的一阶差值，以确定股票仓位。数值 +1 表示多头头寸，-1 表示空头头寸，0 表示无头寸。
6. 第二个信号将与第一个信号相反：做多一只股票，同时做空另一只股票。同样，计算第二个信号的一阶差值，创建第二个仓位列。
7. 最后，将这两只股票的价格及其在投资组合中的多头和空头头寸可视化，交易信号和头寸如下图所示。

五、投资组合损益的计算

我们将从10万美元的初始投资开始，并根据这笔资金确定每只股票的最大买入股数。在每个交易日，第一只股票的盈亏总额将根据该股票的持股总价值和分配给它的现金余额计算得出。同样，第二只股票的盈亏将根据其持有股票的总价值和该股票的现金头寸计算得出。值得注意的是，我们保持市场中性立场，即同时以大致相等的资金量做多和做空。为了确定总体盈亏，我们将这两个值相加。根据股票 1 和股票 2 的仓位，我们将计算它们各自的每日收益。此外，我们还将在最终的投资组合数据框中加入一列 z-分数，其中包含上下限指标，以便直观显示。

# initial capital to calculate the actual pnl
initial_capital = 100000

# shares to buy for each position
positions1 = initial_capital// max(signals['asset1'])
positions2 = initial_capital// max(signals['asset2'])

# since there are two assets, we calculate each asset Pnl 
# separately and in the end we aggregate them into one portfolio
portfolio = pd.DataFrame()
portfolio['asset1'] = signals['asset1']
portfolio['holdings1'] = signals['positions1'].cumsum() * signals['asset1'] * positions1
portfolio['cash1'] = initial_capital - (signals['positions1'] * signals['asset1'] * positions1).cumsum()
portfolio['total asset1'] = portfolio['holdings1'] + portfolio['cash1']
portfolio['return1'] = portfolio['total asset1'].pct_change()
portfolio['positions1'] = signals['positions1']

# pnl for the 2nd asset
portfolio['asset2'] = signals['asset2']
portfolio['holdings2'] = signals['positions2'].cumsum() * signals['asset2'] * positions2
portfolio['cash2'] = initial_capital - (signals['positions2'] * signals['asset2'] * positions2).cumsum()
portfolio['total asset2'] = portfolio['holdings2'] + portfolio['cash2']
portfolio['return2'] = portfolio['total asset2'].pct_change()
portfolio['positions2'] = signals['positions2']

# total pnl and z-score
portfolio['z'] = signals['z']
portfolio['total asset'] = portfolio['total asset1'] + portfolio['total asset2']
portfolio['z upper limit'] = signals['z upper limit']
portfolio['z lower limit'] = signals['z lower limit']
portfolio = portfolio.dropna()

# plot the asset value change of the portfolio and pnl along with z-score
fig = plt.figure(figsize=(14,6),)
ax = fig.add_subplot(111)
ax2 = ax.twinx()
l1, = ax.plot(portfolio['total asset'], c='g')
l2, = ax2.plot(portfolio['z'], c='black', alpha=0.3)
b = ax2.fill_between(portfolio.index,portfolio['z upper limit'],\
                portfolio['z lower limit'], \
                alpha=0.2,color='#ffb48f')
ax.set_ylabel('Asset Value')
ax2.set_ylabel('Z Statistics',rotation=270)
ax.yaxis.labelpad=15
ax2.yaxis.labelpad=15
ax.set_xlabel('Date')
ax.xaxis.labelpad=15
plt.title('Portfolio Performance with Profit and Loss')
plt.legend([l2,b,l1],['Z Statistics',
                      'Z Statistics +-1 Sigma',
                      'Total Portfolio Value'],loc='upper left');
plt.savefig('images/chart8', dpi=300);

# calculate CAGR
final_portfolio = portfolio['total asset'].iloc[-1]
delta = (portfolio.index[-1] - portfolio.index[0]).days
print('Number of days = ', delta)
YEAR_DAYS = 365
returns = (final_portfolio/initial_capital) ** (YEAR_DAYS/delta) - 1print('CAGR = {:.3f}%' .format(returns * 100))

代码解释：

1. 从10万美元的初始资金开始，确定每种股票的购买股数。
2. 然后，计算第一只股票的持仓量，将其持仓量的累计总和乘以股票价格和股票总数。
3. 从初始现金额中减去持有价值，计算现金余额。投资组合中股票的总仓位就是现金价值和持有价值的总和。
4. 根据股票总仓位确定每日收益。
5. 对第二只股票重复步骤 1 至 4，然后将两种资产的仓位合并，计算出投资组合的整体价值。
6. 包括 z 分数和上下阈值，以便可视化。
7. 将投资组合的表现与 Z 值及其上下临界值一起可视化。
8. 计算投资组合的复合年增长率 (CAGR)，投资组合损益情况见下图。

六、复合年增长率（CAGR）及观点总结

这一战略的复合年增长率（CAGR）为 16.5%，似乎还不错。不过，在得出任何明确结论之前，必须考虑几个重要因素。主要考虑因素包括：

1. 鉴于这是一项市场中性策略，很大程度上取决于我们执行卖空的能力，而这可能受到各种因素的限制，包括监管限制。
2. 该分析未考虑与交易相关的成本、市场滑点或与借入证券相关的费用。市场中性策略通常涉及高频率的交易。
3. 使用历史数据作为预测未来业绩的依据，始终存在固有的局限性。

重要的是要记住，只有在彻底评估了所有关键的绩效指标，包括其实用性以及扣除费用和开支后的回报后，才能做出实施战略的决定。

• 配对交易是一种市场中性策略：通过找到历史上表现出高度相关性的两只股票，当其价格关系出现偏离时，预计其将回归均值，实现盈利。
• 协整是配对交易的关键：通过计算相关性和协整性，可以识别出具有长期平衡关系的股票对，这对于成功的配对交易至关重要。
• 统计方法的应用：利用 Python 中的统计库，如 statsmodels，可以进行相关性分析、协整检验以及回归分析，为交易策略提供数据支持。
• 交易信号的生成：通过计算 z-score 并设定上下限阈值，可以生成交易信号，指导何时做多和做空。
• 投资组合的管理和评估：在实施策略时，需要考虑资金的分配、仓位的管理以及投资组合的整体表现，并通过 CAGR 来评估策略的长期收益能力。
• 实际应用中的考虑因素：在将策略应用于实际交易之前，必须考虑到实际交易中可能遇到的各种成本和限制，包括交易成本、市场滑点、借入证券的费用以及监管限制等。

感谢您阅读到最后，希望本文能给您带来新的收获。祝您投资成功！如果对文中的内容有任何疑问，请给我留言，必复。

---

本文内容仅仅是技术探讨和学习，并不构成任何投资建议。

转发请注明原作者和出处。