前言

美赛已经结束4天了，一直忙于教资考试的准备，今天我终于抽空写了这篇C题思路复盘的博客。

题目大致要求

题目叫'Trading Strategies（交易策略）'，一共给了两个文件，分别是比特币和黄金价格随时间变化的CSV文件。大致要求可以分为以下四个步骤：

1. 基于截止至当日的价格情况建立模型，预测2021年9月10日原来的本金1000美元会变成多少钱？（Develop a model that gives the best daily trading strategy based only on price data up to that day. How much is the initial $1000 investment worth on 9/10/2021 using your model and strategy?）
2. 证明你的模型提供了最佳策略（Present evidence that your model provides the best strategy.）
3. 就是让你对手续费那两个参数做敏感性分析（Determine how sensitive the strategy is to transaction costs. How do transaction costs affect the strategy and results?）
4. 美赛的老一套，让你写一个什么备忘录，向投资者陈述你的战略、模型和结果（Communicate your strategy, model, and results to the trader in a memorandum of at most two pages.）

具体题目内容见美赛网站：2022 MCM Problem C (immchallenge.org)http://www.immchallenge.org/mcm/2022_MCM_Problem_C.pdf

分析题目 

对于问题1，思路无非是分为两个步骤，第一个步骤是预测出后几日的价格进而估计得到后几日的收益率，然后通过动态规划模型进行最优化的求解。具体思路如下：

首先先看数据集，你会发现数据集只有两列，时间一列和价格一列。这意味着知网或者其他数据库中查到的很多论文写的预测算法你都用不了了。为什么呢？因为人家写论文时候用的数据集都有其他的特征，而不是这里只有时间一列特征。所以那些什么SVM、贝叶斯网、向量自回归等要求多维特征的算法都不适用。查来查去，有推荐使用LSTM、神经网络模型还有xgboost算法，但我都是不是很熟悉，最终还是用了时间序列ARIMA模型来做的预测。为什么用ARIMA模型来做的预测，我给出以下几点原因：

• ARIMA模型属于单时间序列模型，符合题目数据集提供的数据以及题目给出的要求。
• ARIMA模型是一个很常见的时间序列模型，较为成熟也是我所熟悉的，我能拿它写不少东西。
• 我使用Python编程，而Python的statsmodels库中提供了ARIMA现成的API，不用另外写代码，而且调用API绘制图像的操作十分方便。

对于问题二呢，网上的思路都是说，改变什么参数，看看收益率会不会提高，如果不会的话就说明是最优的，但也没说清楚具体怎么操作。而我的想法不完全是这样的，我看了他们的思路突然想到了"梯度下降法"，那你说如果我用模型的两个参数（具体来讲一个是ARIMA模型的AR参数另一个是MA参数，作为梯度下降法改变模型的两个搜索方向，然后进行算法计算，理论上是完全可行的，于是乎我就这么干了。

对于问题三和问题四就不再多说了。但我也罗嗦一下，问敏感性分析，我除了做了问题三的手续费的敏感性分析外，我还针对ARIMA模型做了敏感度分析。问题四就总分总来写，总一段，分三段，也就是战略、模型和结果各一段写，最后写个总结，要注意的是英文备忘录是有格式的，得按照格式来，要有点仪式感，不然你就可能会被扣分。

还有就是别忘了一开始写数据清洗那一段，我用了基于时间序列的牛顿插值法进行插补，用箱型图法判断离群点。

接下来我贴出我的代码。

代码

import pandas as pd
from pylab import *
import statsmodels.api as sm
from datetime import datetime
from statsmodels.tsa.stattools import adfuller
from statsmodels.tsa.stattools import kpss
from statsmodels.stats.diagnostic import acorr_ljungbox


def adf_test(timeseries):
    print('Results of Dickey-Fuller Test:')
    dftest = adfuller(timeseries, autolag='AIC')
    dfoutput = pd.Series(dftest[0:4], index=['Test Statistic', 'p-value', '#Lags Used', 'Number of Observations Used'])
    for key, value in dftest[4].items():
       dfoutput['Critical Value (%s)'%key] = value
    print(dfoutput)


def kpss_test(timeseries):
    print('Results of KPSS Test:')
    kpsstest = kpss(timeseries, regression='c')
    kpss_output = pd.Series(kpsstest[0:3], index=['Test Statistic','p-value','Lags Used'])
    for key, value in kpsstest[3].items():
        kpss_output['Critical Value (%s)'%key] = value
    print(kpss_output)


dataGOLD = pd.read_csv("C:\\ProblemC\\LBMA-GOLD.csv")
dataMKPRU = pd.read_csv("C:\\ProblemC\\BCHAIN-MKPRU.csv")
dataGOLD['Date'] = pd.to_datetime(dataGOLD.Date)
dataMKPRU['Date'] = pd.to_datetime(dataMKPRU.Date)
dataGOLD = dataGOLD[dataGOLD['Date'] >= datetime.strptime('2016-9-11', "%Y-%m-%d")]
dataMKPRU = dataMKPRU[dataMKPRU['Date'] >= datetime.strptime('2016-9-11', "%Y-%m-%d")]
dataGOLD = dataGOLD[dataGOLD['Date'] <= datetime.strptime('2021-9-10', "%Y-%m-%d")]
dataMKPRU = dataMKPRU[dataMKPRU['Date'] <= datetime.strptime('2021-9-10', "%Y-%m-%d")]
dataGOLD = dataGOLD.sort_values('Date').reset_index(drop=True)
dataMKPRU = dataMKPRU.sort_values('Date').reset_index(drop=True)

adf_test(dataGOLD['Value'])
kpss_test(dataGOLD['Value'])

diffGOLD = dataGOLD['Value'].diff(1).dropna()
diffMKPRU = dataMKPRU['Value'].diff(1).dropna()

fig = plt.figure(figsize=(12, 8))
ax1 = fig.add_subplot(211)
sm.graphics.tsa.plot_acf(diffMKPRU, ax=ax1)
ax1.xaxis.set_ticks_position('top')
ax2 = fig.add_subplot(212)
sm.graphics.tsa.plot_pacf(diffMKPRU, ax=ax2)
ax2.xaxis.set_ticks_position('bottom')
fig.tight_layout()
plt.show()

fig = plt.figure(figsize=(12, 8))
ax11 = fig.add_subplot(211)
sm.graphics.tsa.plot_acf(diffGOLD, ax=ax11)
ax11.xaxis.set_ticks_position('top')
ax22 = fig.add_subplot(212)
sm.graphics.tsa.plot_pacf(diffGOLD, ax=ax22)
ax22.xaxis.set_ticks_position('bottom')
fig.tight_layout()
plt.show()

p_value = acorr_ljungbox(dataGOLD['Value'])
print(p_value)

results = sm.tsa.arma_order_select_ic(dataMKPRU['Value'], ic=['aic', 'bic'], trend='nc', max_ar=5, max_ma=5)
print('AIC', results.aic_min_order)
print('BIC', results.bic_min_order)

n_sample = dataGOLD['Value'].shape[0]
n_train = int(n_sample * 0.95) + 1
n_forecast = n_sample - n_train

ts_train = dataGOLD.iloc[:n_train]['Value']
ts_test = dataGOLD.iloc[n_train:]['Value']

arima = sm.tsa.SARIMAX(ts_train, order=(2, 1, 0))
model_results = arima.fit()
model_results.plot_diagnostics(figsize=(16, 12))
plt.show()
print(model_results.summary())

plt.title("Bitcoin price image after first-order difference")
plt.plot(dataMKPRU['Date'], dataMKPRU['Value'].diff(1))
plt.show()

plt.title("Gold price image after first-order difference")
plt.plot(dataGOLD['Date'], dataGOLD['Value'].diff(1))
plt.show()

代码绘图展示

我们使用matplotlib进行绘图，绘制了如下的图片进行展示：

箱型图法检测的两张图，如下：

价格随时间变化的趋势图两张：

一阶差分图两张： 

应用新策略后的趋势图一张：

 ARIMA模型预测图一张

ACF图和PACF图两张： 

各类图表：一阶差分图、正态分布直方图、QQ图等 

代码分析与第一问细节 

可以看到，statsmodels确实十分强大，各种API拿来就用，十分方便。

其中我们用ACF图和PACF图确定了参数的大致阶数，然后用BIC和AIC来定阶。代码只展示了黄金的ARIMA模型的阶数为（2，1，0），比特币模型的阶数只要改下变量就行了，就不在展示了，比特币模型的阶数最终定为（4，1，4）。

就写到这里吧，仅供参考和复盘。

最终获奖结果

更新于5月8日。

我们团队最终获得了M奖，祝贺所有人！！！

结语

 这是我们的论文链接，仅供参考！学术生涯不易，请投喂我们，论文只要1.9元，也不是很贵！

2022年美赛（MCM）C题M奖论文https://download.csdn.net/download/qq_41938259/85318021

END