精通-Python-金融第二版（一）
精通 Python 金融第二版（一）
原文：
zh.annas-archive.org/md5/8b046e39ce2c1a10ac13fd89834aaadc译者：飞龙
前言
本书的第二版Mastering Python for Finance将指导您使用下一代方法在金融行业中进行复杂的金融计算。您将通过利用公开可用的工具来掌握 Python 生态系统，成功进行研究和建模，并学习如何使用高级示例来管理风险。
您将首先设置一个 Jupyter 笔记本，以实现本书中的任务。您将学习如何使用流行的库（如 TensorFlow、Keras、NumPy、SciPy、scikit-learn 等）做出高效而强大的数据驱动金融决策。您还将学习如何通过掌握股票、期权、利率及其衍生品以及使用计算方法进行风险分析等概念来构建金融应用程序。有了这些基础，您将学习如何对时间序列数据进行统计分析，并了解如何利用高频数据来设计交易策略，从而构建算法交易*台。您将学习通过实施事件驱动的回测系统来验证您的交易策略，并衡量其性能。最后，您将探索在金融领域应用的机器学习和深度学习技术。
通过本书，您将学会如何将 Python 应用于金融行业中的不同范式，并进行高效的数据分析。
这本书是为谁准备的
如果您是金融或数据分析师，或者是金融行业的软件开发人员，有兴趣使用高级 Python 技术进行量化方法，那么这本书就是您需要的！如果您想要使用智能机器学习技术扩展现有金融应用程序的功能，您也会发现这本书很有用。
充分利用本书
需要有 Python 的先验经验。
下载示例代码文件
您可以从www.packt.com的帐户中下载本书的示例代码文件。如果您在其他地方购买了这本书，您可以访问www.packt.com/support并注册，文件将直接发送到您的邮箱。
您可以按照以下步骤下载代码文件：
-
登录或注册www.packt.com。
-
选择“支持”选项卡。
-
单击“代码下载和勘误”。
-
在搜索框中输入书名，然后按照屏幕上的说明操作。
下载文件后，请确保使用最新版本的软件解压或提取文件夹：
-
WinRAR/7-Zip for Windows
-
Zipeg/iZip/UnRarX for Mac
-
7-Zip/PeaZip for Linux
该书的代码包也托管在 GitHub 上，网址为github.com/PacktPublishing/Mastering-Python-for-Finance-Second-Edition。如果代码有更新，将在现有的 GitHub 存储库上进行更新。
我们还有其他代码包，来自我们丰富的图书和视频目录，可在github.com/PacktPublishing/上找到。去看看吧！
下载彩色图像
我们还提供了一个 PDF 文件，其中包含本书中使用的屏幕截图/图表的彩色图像。您可以在这里下载：www.packtpub.com/sites/default/files/downloads/9781789346466_ColorImages.pdf。
使用的约定
本书中使用了许多文本约定。
CodeInText：表示文本中的代码词、数据库表名、文件夹名、文件名、文件扩展名、路径名、虚拟 URL、用户输入和 Twitter 句柄。例如：“默认情况下，pandas 的.plot()命令使用matplotlib库来显示图形。”
代码块设置如下：
In [ ]:
%matplotlib inline
import quandl
quandl.ApiConfig.api_key = QUANDL_API_KEY
df = quandl.get('EURONEXT/ABN.4')
daily_changes = df.pct_change(periods=1)
daily_changes.plot();
当我们希望引起您对代码块的特定部分的注意时，相关行或项目将以粗体显示：
2015-02-26 TICK WIKI/AAPL open: 128.785 close: 130.415
2015-02-26 FILLED BUY 1 WIKI/AAPL at 128.785
2015-02-26 POSITION value:-128.785 upnl:1.630 rpnl:0.000
2015-02-27 TICK WIKI/AAPL open: 130.0 close: 128.46
任何命令行输入或输出都是按照以下格式编写的：
$ cd my_project_folder
$ virtualenv my_env
粗体：表示一个新术语，一个重要词或屏幕上看到的词。例如，菜单或对话框中的单词会以这种方式出现在文本中。这是一个例子：“要启动你的第一个笔记本，选择新建，然后Python 3。”
警告或重要提示会显示为这样。提示和技巧会显示为这样。
第一部分：开始使用 Python
本节将帮助我们在准备在本书中运行代码示例之前，在我们的机器上设置 Python。
本节将只包含一个章节：
- 第一章，使用 Python 进行财务分析概述
第一章：使用 Python 进行金融分析概述
自从我之前的书《精通 Python 金融》出版以来，Python 本身和许多第三方库都有了重大升级。许多工具和功能已被弃用，取而代之的是新的工具和功能。本章将指导您如何获取最新的可用工具，并准备本书其余部分将使用的环境。
本书中涵盖的大部分数据集将使用 Quandl。Quandl 是一个提供金融、经济和替代数据的*台。这些数据来源于各种数据发布者，包括联合国、世界银行、中央银行、交易所、投资研究公司，甚至 Quandl 社区的成员。通过 Python Quandl 模块，您可以轻松下载数据集并进行金融分析，以得出有用的见解。
我们将使用pandas模块探索时间序列数据操作。pandas中的两个主要数据结构是 Series 对象和 DataFrame 对象。它们可以一起用于绘制图表和可视化复杂信息。本章将涵盖金融时间序列计算和分析的常见方法。
本章的目的是为您设置工作环境，使用本书中将使用的库。多年来，像任何软件包一样，pandas模块已经发生了巨大的演变，许多重大变化。多年前编写的代码与旧版本的pandas接口将不再起作用，因为许多方法已被弃用。本书中使用的pandas版本是 0.23。本书中编写的代码符合这个版本的pandas。
在本章中，我们将涵盖以下内容：
-
为您的环境设置 Python、Jupyter、Quandl 和其他库
-
从 Quandl 下载数据集并绘制您的第一个图表
-
绘制最后价格、成交量和蜡烛图
-
计算和绘制每日百分比和累积收益
-
绘制波动率、直方图和 Q-Q 图
-
可视化相关性并生成相关矩阵
-
可视化简单移动*均线和指数移动*均线
获取 Python
在撰写本文时，最新的 Python 版本是 3.7.0。您可以从官方 Python 网站www.python.org/downloads/下载 Windows、macOS X、Linux/UNIX 和其他操作系统的最新版本。按照安装说明在您的操作系统上安装基本的 Python 解释器。
安装过程应该将 Python 添加到您的环境路径中。要检查已安装的 Python 版本，请在 macOS X/Linux 上的终端中输入以下命令，或者在 Windows 上的命令提示符中输入以下命令：
$ python --version
Python 3.7.0
为了方便安装 Python 库，考虑使用 Anaconda（www.anaconda.com/download/）、Miniconda（conda.io/miniconda.html）或 Enthought Canopy（www.enthought.com/product/enthought-python-distribution/）等一体化 Python 发行版。然而，高级用户可能更喜欢控制哪些库与他们的基本 Python 解释器一起安装。
准备虚拟环境
此时，建议设置 Python 虚拟环境。虚拟环境允许您管理特定项目所需的单独包安装，隔离其他环境中安装的包。
要在终端窗口中安装虚拟环境包，请输入以下命令：
$ pip install virtualenv
在某些系统上，Python 3 可能使用不同的pip可执行文件，并且可能需要通过替代的pip命令进行安装；例如：$ pip3 install virtualenv。
要创建虚拟环境，请转到项目目录并运行virtualenv。例如，如果您的项目文件夹的名称是my_project_folder，请输入以下内容：
$ cd my_project_folder
$ virtualenv my_venv
virtualenv my_venv将在当前工作目录中创建一个文件夹，其中包含您之前安装的基本 Python 解释器的 Python 可执行文件，以及pip库的副本，您可以使用它来安装其他软件包。
在使用新的虚拟环境之前，需要激活它。在 macOS X 或 Linux 终端中，输入以下命令：
$ source my_venv/bin/activate
在 Windows 上，激活命令如下：
$ my_project_folder\my_venv\Scripts\activate
当前虚拟环境的名称现在将显示在提示的左侧（例如，(my_venv) current_folder$），以让您知道所选的 Python 环境已激活。从同一终端窗口进行的软件包安装将放在my_venv文件夹中，与全局 Python 解释器隔离开来。
虚拟环境可以帮助防止冲突，如果您有多个应用程序使用相同模块但来自不同版本。这一步（创建虚拟环境）完全是可选的，因为您仍然可以使用默认的基本解释器来安装软件包。
运行 Jupyter Notebook
Jupyter Notebook 是一个基于浏览器的交互式计算环境，用于创建、执行和可视化各种编程语言的交互式数据。它以前被称为IPython Notebook。IPython 仍然存在作为 Python shell 和 Jupyter 的内核。Jupyter 是一个开源软件，所有人都可以免费使用和学习各种主题，从基本编程到高级统计学或量子力学。
要安装 Jupyter，在终端窗口中输入以下命令：
$ pip install jupyter
安装后，使用以下命令启动 Jupyter：
$ jupyter notebook
...
Copy/paste this URL into your browser when you connect for the first time, to login with a token:
http://127.0.0.1:8888/?token=27a16ee4d6042a53f6e31161449efcf7e71418f23e17549d
观察您的终端窗口。当 Jupyter 启动时，控制台将提供有关其运行状态的信息。您还应该看到一个 URL。将该 URL 复制到 Web 浏览器中，即可进入 Jupyter 计算界面。
由于 Jupyter 在您发出前面的命令的目录中启动，Jupyter 将列出工作目录中所有保存的笔记本。如果这是您在该目录中工作的第一次，列表将为空。
要启动您的第一个笔记本，请选择 New，然后选择 Python 3。一个新的 Jupyter Notebook 将在新窗口中打开。今后，本书中的大多数计算将在 Jupyter 中进行。
Python Enhancement Proposal
Python 编程语言中的任何设计考虑都被记录为Python Enhancement Proposal（PEP）。已经编写了数百个 PEP，但您可能应该熟悉的是PEP 8，这是 Python 开发人员编写更好、更可读代码的风格指南。PEP 的官方存储库是github.com/python/peps。
什么是 PEP？
PEP 是一系列编号的设计文档，描述与 Python 相关的特性、过程或环境。每个 PEP 都在一个文本文件中进行精心维护，包含特定特性的技术规范及其存在的原因。例如，PEP 0 用作所有 PEP 的索引，而 PEP 1 提供了 PEP 的目的和指南。作为软件开发人员，我们经常阅读代码而不是编写代码。为了创建清晰、简洁和可读的代码，我们应该始终使用编码约定作为风格指南。PEP 8 是一组编写 Python 代码的风格指南。您可以在www.python.org/dev/peps/pep-0008/上了解更多关于 PEP 8 的信息。
Python 之禅
PEP 20 体现了 Python 之禅，这是一组指导 Python 编程语言设计的 20 个软件原则。要显示这个彩蛋，在 Python shell 中输入以下命令：
>> import this
The Zen of Python, by Tim Peters Beautiful is better than ugly. Explicit is better than implicit.
Simple is better than complex.
Complex is better than complicated.
Flat is better than nested.
Sparse is better than dense.
Readability counts.
Special cases aren't special enough to break the rules.
Although practicality beats purity.
Errors should never pass silently.
Unless explicitly silenced.
In the face of ambiguity, refuse the temptation to guess.
There should be one-- and preferably only one --obvious way to do it.
Although that way may not be obvious at first unless you're Dutch.
Now is better than never.
Although never is often better than *right* now.
If the implementation is hard to explain, it's a bad idea.
If the implementation is easy to explain, it may be a good idea.
Namespaces are one honking great idea -- let's do more of those!
只显示了 20 条格言中的 19 条。你能猜出最后一条是什么吗？我留给你的想象！
Quandl 简介
Quandl 是一个提供金融、经济和替代数据的*台。这些数据来源于各种数据发布者，包括联合国、世界银行、中央银行、交易所和投资研究公司。
使用 Python Quandl 模块，您可以轻松地将金融数据集导入 Python。Quandl 提供免费数据集，其中一些是样本。访问高级数据产品需要付费。
为您的环境设置 Quandl
Quandl包需要最新版本的 NumPy 和pandas。此外，我们将在本章的其余部分需要matplotlib。
要安装这些包，请在终端窗口中输入以下代码：
$ pip install quandl numpy pandas matplotlib
多年来，pandas库发生了许多变化。为旧版本的pandas编写的代码可能无法与最新版本一起使用，因为有许多已弃用的内容。我们将使用的pandas版本是 0.23。要检查您正在使用的pandas版本，请在 Python shell 中键入以下命令：
>>> import pandas
>>> pandas.__version__'0.23.3'
使用 Quandl 请求数据时需要一个 API（应用程序编程接口）密钥。
如果您没有 Quandl 账户，请按以下步骤操作：
- 打开浏览器，在地址栏中输入
www.quandl.com。这将显示以下页面：
-
选择注册并按照说明创建一个免费账户。成功注册后，您将会看到您的 API 密钥。
-
复制此密钥并将其安全地保存在其他地方，因为您以后会需要它。否则，您可以在您的账户设置中再次检索此密钥。
-
请记住检查您的电子邮件收件箱，查看欢迎消息并验证您的 Quandl 账户，因为继续使用 API 密钥需要验证和有效的 Quandl 账户。
匿名用户每 10 分钟最多可以调用 20 次，每天最多可以调用 50 次。经过身份验证的免费用户每 10 秒最多可以调用 300 次，每 10 分钟最多可以调用 2,000 次，每天最多可以调用 50,000 次。
绘制时间序列图表
在图表上可视化时间序列数据是一种简单而有效的分析技术，通过它我们可以推断出某些假设。本节将指导您完成从 Quandl 下载股价数据集并在价格和成交量图上绘制的过程。我们还将介绍绘制蜡烛图表，这将为我们提供比线图更多的信息。
从 Quandl 检索数据
从 Quandl 中获取数据到 Python 是相当简单的。假设我们对来自 Euronext 股票交易所的 ABN Amro Group 感兴趣。在 Quandl 中的股票代码是EURONEXT/ABN。在 Jupyter 笔记本单元格中，运行以下命令：
In [ ]:
import quandl
# Replace with your own Quandl API key
QUANDL_API_KEY = 'BCzkk3NDWt7H9yjzx-DY'
quandl.ApiConfig.api_key = QUANDL_API_KEY
df = quandl.get('EURONEXT/ABN')
将 Quandl API 密钥存储在常量变量中是一个好习惯。这样，如果您的 API 密钥发生变化，您只需要在一个地方更新它！
导入quandl包后，我们将 Quandl API 密钥存储在常量变量QUANDL_API_KEY中，这将在本章的其余部分中重复使用。这个常量值用于设置 Quandl 模块的 API 密钥，每次导入quandl包时只需要执行一次。接下来的一行调用quandl.get()方法，将 ABN 数据集从 Quandl 直接下载到我们的df变量中。请注意，EURONEXT是数据提供者 Euronext 股票交易所的缩写。
默认情况下，Quandl 将数据集检索到pandas DataFrame 中。我们可以按以下方式检查 DataFrame 的头和尾：
In [ ]:
df.head()
Out[ ]:
Open High Low Last Volume Turnover
Date
2015-11-20 18.18 18.43 18.000 18.35 38392898.0 7.003281e+08
2015-11-23 18.45 18.70 18.215 18.61 3352514.0 6.186446e+07
2015-11-24 18.70 18.80 18.370 18.80 4871901.0 8.994087e+07
2015-11-25 18.85 19.50 18.770 19.45 4802607.0 9.153862e+07
2015-11-26 19.48 19.67 19.410 19.43 1648481.0 3.220713e+07
In [ ]:
df.tail()
Out[ ]:
Open High Low Last Volume Turnover
Date
2018-08-06 23.50 23.59 23.29 23.34 1126371.0 2.634333e+07
2018-08-07 23.59 23.60 23.31 23.33 1785613.0 4.177652e+07
2018-08-08 24.00 24.39 23.83 24.14 4165320.0 1.007085e+08
2018-08-09 24.40 24.46 24.16 24.37 2422470.0 5.895752e+07
2018-08-10 23.70 23.94 23.28 23.51 3951850.0 9.336493e+07
默认情况下，head()和tail()命令将分别显示 DataFrame 的前五行和后五行。您可以通过在参数中传递一个数字来定义要显示的行数。例如，head(100)将显示 DataFrame 的前 100 行。
对于get()方法没有设置任何额外参数，将检索整个时间序列数据集，从上一个工作日一直回溯到 2015 年 11 月，每天一次。
要可视化这个 DataFrame，我们可以使用plot()命令绘制图表：
In [ ]:
%matplotlib inline
import matplotlib.pyplot as plt
df.plot();
最后一个命令输出一个简单的图表：
pandas的plot()方法返回一个 Axes 对象。这个对象的字符串表示形式与plot()命令一起打印在控制台上。要抑制这些信息，可以在最后一个语句的末尾添加一个分号(😉。或者，可以在单元格的底部添加一个pass语句。另外，将绘图函数分配给一个变量也会抑制输出。
默认情况下，pandas中的plot()命令使用matplotlib库来显示图表。如果出现错误，请检查是否安装了该库，并且至少调用了%matplotlib inline。您可以自定义图表的外观和感觉。有关pandas DataFrame 中plot命令的更多信息，请参阅pandas文档pandas.pydata.org/pandas-docs/stable/generated/pandas.DataFrame.plot.html。
绘制价格和成交量图表
当没有参数提供给plot()命令时，将使用目标 DataFrame 的所有列绘制一条线图，放在同一张图上。这会产生一个混乱的视图，无法提供太多信息。为了有效地从这些数据中提取见解，我们可以绘制一个股票的财务图表，显示每日收盘价与交易量的关系。为了实现这一点，输入以下命令：
In [ ]:
prices = df['Last']
volumes = df['Volume']
上述命令将我们感兴趣的数据分别存储到closing_prices和volumes变量中。我们可以使用head()和tail()命令查看生成的pandas Series 数据类型的前几行和最后几行：
In [ ]:
prices.head()
Out[ ]:
Date
2015-11-20 18.35
2015-11-23 18.61
2015-11-24 18.80
2015-11-25 19.45
2015-11-26 19.43
Name: Last, dtype: float64
In [ ]:
volumes.tail()
Out[ ]:
Date
2018-08-03 1252024.0
2018-08-06 1126371.0
2018-08-07 1785613.0
2018-08-08 4165320.0
2018-08-09 2422470.0
Name: Volume, dtype: float64
要找出特定变量的类型，使用type()命令。例如，type(volumes)会产生pandas.core.series.Series，告诉我们volumes变量实际上是一个pandas Series 数据类型对象。
注意，数据可从 2018 年一直回溯到 2015 年。现在我们可以绘制价格和成交量图表：
In [ ]:
# The top plot consisting of daily closing prices
top = plt.subplot2grid((4, 4), (0, 0), rowspan=3, colspan=4)
top.plot(prices.index, prices, label='Last')
plt.title('ABN Last Price from 2015 - 2018')
plt.legend(loc=2)
# The bottom plot consisting of daily trading volume
bottom = plt.subplot2grid((4, 4), (3,0), rowspan=1, colspan=4)
bottom.bar(volumes.index, volumes)
plt.title('ABN Daily Trading Volume')
plt.gcf().set_size_inches(12, 8)
plt.subplots_adjust(hspace=0.75)
这将产生以下图表：
在第一行中，subplot2grid命令的第一个参数(4,4)将整个图表分成 4x4 的网格。第二个参数(0,0)指定给定的图表将锚定在图表的左上角。关键字参数rowspan=3表示图表将占据网格中可用的 4 行中的 3 行，实际上占据了图表的 75%高度。关键字参数colspan=4表示图表将占据网格的所有 4 列，使用了所有可用的宽度。该命令返回一个matplotlib轴对象，我们将使用它来绘制图表的上部分。
在第二行，plot()命令呈现了上部图表，x轴上是日期和时间值，y轴上是价格。在接下来的两行中，我们指定了当前图表的标题，以及放在左上角的时间序列数据的图例。
接下来，我们执行相同的操作，在底部图表上呈现每日交易量，指定一个 1 行 4 列的网格空间，锚定在图表的左下角。
在legend()命令中，loc关键字接受一个整数值作为图例的位置代码。值为2表示左上角位置。有关位置代码的表格，请参阅matplotlib的图例文档matplotlib.org/api/legend_api.html?highlight=legend#module-matplotlib.legend。
为了使我们的图形看起来更大，我们调用set_size_inches()命令将图形设置为宽 9 英寸、高 6 英寸，从而产生一个长方形的图形。前面的gcf()命令简单地表示获取当前图形。最后，我们调用带有hspace参数的subplots_adjust()命令，以在顶部和底部子图之间添加一小段高度。
subplots_adjust()命令调整了子图布局。可接受的参数有left、right、bottom、top、wspace和hspace。有关这些参数的更多信息，请参阅matplotlib.org/api/_as_gen/matplotlib.pyplot.subplots_adjust.html中的matplotlib文档。
绘制蜡烛图
蜡烛图是另一种流行的金融图表类型，它显示了比单一价格更多的信息。蜡烛图代表了特定时间点的每个标记，其中包含四个重要的信息：开盘价、最高价、最低价和收盘价。
matplotlib.finance模块已被弃用。相反，我们可以使用另一个包mpl_finance，其中包含了提取的代码。要安装此包，在您的终端窗口中，输入以下命令：
$ pip install mpl-finance
为了更仔细地可视化蜡烛图，我们将使用 ABN 数据集的子集。在下面的示例中，我们从 Quandl 查询 2018 年 7 月的每日价格作为我们的数据集，并绘制一个蜡烛图，如下所示：
In [ ]:
%matplotlib inline
import quandl
from mpl_finance import candlestick_ohlc
import matplotlib.dates as mdates
import matplotlib.pyplot as plt
quandl.ApiConfig.api_key = QUANDL_API_KEY
df_subset = quandl.get('EURONEXT/ABN',
start_date='2018-07-01',
end_date='2018-07-31')
df_subset['Date'] = df_subset.index.map(mdates.date2num)
df_ohlc = df_subset[['Date','Open', 'High', 'Low', 'Last']]
figure, ax = plt.subplots(figsize = (8,4))
formatter = mdates.DateFormatter('%Y-%m-%d')
ax.xaxis.set_major_formatter(formatter)
candlestick_ohlc(ax,
df_ohlc.values,
width=0.8,
colorup='green',
colordown='red')
plt.show()
这将产生一个蜡烛图，如下截图所示：
您可以在quandl.get()命令中指定start_date和end_date参数，以检索所选日期范围的数据集。
从 Quandl 检索的价格被放置在一个名为df_dataset的变量中。由于matplotlib的plot函数需要自己的格式，mdates.date2num命令将包含日期和时间的索引值转换，并将它们放在一个名为Date的新列中。
蜡烛图的日期、开盘价、最高价、最低价和收盘价数据列被明确提取为df_ohlc变量中的 DataFrame。plt.subplots()创建了一个宽 8 英寸、高 4 英寸的图形。x轴上的标签被格式化为人类可读的格式。
我们的数据现在已准备好通过调用candlestick_ohlc()命令作为蜡烛图来绘制，蜡烛图的宽度为 0.8（或全天宽度的 80%）。收盘价高于开盘价的上涨标记以绿色表示，而收盘价低于开盘价的下跌标记以红色表示。最后，我们添加了plt.show()命令来显示蜡烛图。
在时间序列数据上执行金融分析
在本节中，我们将可视化金融分析中使用的时间序列数据的一些统计属性。
绘制收益
安全性表现的经典指标之一是其在先前时期的收益。在pandas中计算收益的一种简单方法是pct_change，它计算了每行在 DataFrame 中的前一行的百分比变化。
在下面的示例中，我们使用 ABN 股票数据绘制了每日百分比收益的简单图表：
In [ ]:
%matplotlib inline
import quandl
quandl.ApiConfig.api_key = QUANDL_API_KEY
df = quandl.get('EURONEXT/ABN.4')
daily_changes = df.pct_change(periods=1)
daily_changes.plot();
每日百分比收益的折线图如下所示：
在quandl.get()方法中，我们使用后缀符号.4来指定仅检索数据集的第四列，其中包含最后的价格。在调用pct_change时，period参数指定了要移动以形成百分比变化的周期数，默认为1。
我们可以使用column_index参数和列的索引来代替使用后缀符号来指定要下载的数据集的列。例如，quandl.get('EURONEXT/ABN.4')与调用quandl.get('EURONEXT/ABN', column_index=4)是相同的。
绘制累积收益
为了了解我们的投资组合的表现，我们可以在一段时间内对其收益进行求和。pandas的cumsum方法返回 DataFrame 的累积和。
在下面的例子中，我们绘制了之前计算的 ABN 的daily_changes的累积和：
In [ ]:
df_cumsum = daily_changes.cumsum()
df_cumsum.plot();
这给我们以下输出图表：
绘制直方图
直方图告诉我们数据的分布情况。在这个例子中，我们对 ABN 的每日收益的分布情况感兴趣。我们在一个具有 50 个箱子大小的 DataFrame 上使用hist()方法：
In [ ]:
daily_changes.hist(bins=50, figsize=(8, 4));
直方图输出如下：
当pandas DataFrame 中有多个数据列时，hist()方法将自动在单独的图表上绘制每个直方图。
我们可以使用describe()方法来总结数据集分布的中心趋势、离散度和形状：
In [ ]:
daily_changes.describe()
Out[ ]:
Last
count 692.000000
mean 0.000499
std 0.016701
min -0.125527
25% -0.007992
50% 0.000584
75% 0.008777
max 0.059123
从直方图中可以看出，收益倾向于围绕着 0.0 的均值分布，或者确切地说是0.000499。除了这个微小的右偏移，数据看起来相当对称和正态分布。标准偏差为0.016701。百分位数告诉我们，25%的点在-0.007992以下，50%在0.000584以下，75%在0.008777以下。
绘制波动率
分析收益分布的一种方法是测量其标准偏差。标准偏差是均值周围离散度的度量。过去收益的高标准偏差值表示股价波动的历史波动性较高。
pandas的rolling()方法帮助我们在一段时间内可视化特定的时间序列操作。为了计算我们计算的 ABN 数据集的收益百分比的标准偏差，我们使用std()方法，它返回一个 DataFrame 或 Series 对象，可以用来绘制图表。下面的例子说明了这一点：
In [ ]:
df_filled = df.asfreq('D', method='ffill')
df_returns = df_filled.pct_change()
df_std = df_returns.rolling(window=30, min_periods=30).std()
df_std.plot();
这给我们以下波动率图：
我们原始的时间序列数据集不包括周末和公共假期，在使用rolling()方法时必须考虑这一点。df.asfreq()命令将时间序列数据重新索引为每日频率，在缺失的索引位置创建新的索引。method参数的值为ffill，指定我们在重新索引时将最后一个有效观察结果向前传播，以替代缺失值。
在rolling()命令中，我们指定了window参数的值为 30，这是用于计算统计量的观察次数。换句话说，每个期间的标准偏差是用样本量 30 来计算的。由于前 30 行没有足够的样本量来计算标准偏差，我们可以通过将min_periods指定为30来排除这些行。
选择的值 30 接*月度收益的标准偏差。请注意，选择更宽的窗口期代表着被测量的数据量较少。
一个分位数-分位数图
Q-Q（分位数-分位数）图是一个概率分布图，其中两个分布的分位数相互绘制。如果分布是线性相关的，Q-Q 图中的点将位于一条直线上。与直方图相比，Q-Q 图帮助我们可视化偏离正态分布线的点，以及过度峰度的正负偏差。
scipy.stats的probplot()帮助我们计算并显示概率图的分位数。数据的最佳拟合线也被绘制出来。在下面的例子中，我们使用 ABN 股票数据集的最后价格，并计算每日百分比变化以绘制 Q-Q 图：
In [ ]:
%matplotlib inline
import quandl
from scipy import stats
from scipy.stats import probplot
quandl.ApiConfig.api_key = QUANDL_API_KEY
df = quandl.get('EURONEXT/ABN.4')
daily_changes = df.pct_change(periods=1).dropna()
figure = plt.figure(figsize=(8,4))
ax = figure.add_subplot(111)
stats.probplot(daily_changes['Last'], dist='norm', plot=ax)
plt.show();
这给我们以下的 Q-Q 图：
当所有点恰好落在红线上时，数据的分布意味着与正态分布完全对应。我们的大部分数据在分位数-2 和+2 之间几乎完全相关。在这个范围之外，分布的相关性开始有所不同，在尾部有更多的负偏斜。
下载多个时间序列数据
我们将单个 Quandl 代码作为字符串对象传递给quandl.get()命令的第一个参数，以下载单个数据集。要下载多个数据集，我们可以传递一个 Quandl 代码的列表。
在下面的例子中，我们对三家银行股票的价格感兴趣——ABN Amro、Banco Santander 和 Kas Bank。2016 年至 2017 年的两年价格存储在df变量中，只下载了最后价格：
In [ ]:
%matplotlib inline
import quandl
quandl.ApiConfig.api_key = QUANDL_API_KEY
df = quandl.get(['EURONEXT/ABN.4',
'EURONEXT/SANTA.4',
'EURONEXT/KA.4'],
collapse='monthly',
start_date='2016-01-01',
end_date='2017-12-31')
df.plot();
生成了以下图表：
默认情况下，quandl.get()返回每日价格。我们还可以指定数据集下载的其他类型频率。在这个例子中，我们指定collapse='monthly'来下载月度价格。
显示相关矩阵
相关性是两个变量之间线性关系有多密切的统计关联。我们可以对两个时间序列数据集的收益进行相关性计算，得到一个介于-1 和 1 之间的值。相关值为 0 表示两个时间序列的收益之间没有关系。接* 1 的高相关值表示两个时间序列数据的收益倾向于一起变动。接*-1 的低值表示收益倾向于相互反向变动。
在pandas中，corr()方法计算其提供的 DataFrame 中列之间的相关性，并将这些值输出为矩阵。在前面的例子中，我们在 DataFrame df中有三个可用的数据集。要输出收益的相关矩阵，运行以下命令：
In [ ]:
df.pct_change().corr()
Out[ ]:
EURONEXT/ABN - Last ... EURONEXT/KA - Last
EURONEXT/ABN - Last 1.000000 ... 0.096238
EURONEXT/SANTA - Last 0.809824 ... 0.058095
EURONEXT/KA - Last 0.096238 ... 1.000000
从相关矩阵输出中，我们可以推断出 ABN Amro 和 Banco Santander 股票在 2016 年至 2017 年的两年时间内高度相关，相关值为0.809824。
默认情况下，corr()命令使用 Pearson 相关系数来计算成对相关性。这相当于调用corr(method='pearson')。其他有效值是kendall和spearman，分别用于 Kendall Tau 和 Spearman 秩相关系数。
绘制相关性
也可以使用rolling()命令来可视化相关性。我们将使用 2016 年至 2017 年从 Quandl 获取的 ABN 和 SANTA 的每日最后价格。这两个数据集被下载到 DataFrame df中，并且其滚动相关性如下所示：
In [ ]:
%matplotlib inline
import quandl
quandl.ApiConfig.api_key = QUANDL_API_KEY
df = quandl.get(['EURONEXT/ABN.4', 'EURONEXT/SANTA.4'],
start_date='2016-01-01',
end_date='2017-12-31')
df_filled = df.asfreq('D', method='ffill')
daily_changes= df_filled.pct_change()
abn_returns = daily_changes['EURONEXT/ABN - Last']
santa_returns = daily_changes['EURONEXT/SANTA - Last']
window = int(len(df_filled.index)/2)
df_corrs = abn_returns\
.rolling(window=window, min_periods=window)\
.corr(other=santa_returns)
.dropna()
df_corrs.plot(figsize=(12,8));
以下是相关性图的截图：
df_filled变量包含一个 DataFrame，其索引以每日频率重新索引，并且准备好进行rolling()命令的缺失值前向填充。DataFrame daily_changes存储每日百分比收益，并且其列被提取为一个单独的 Series 对象，分别为abn_returns和santa_returns。window变量存储了两年数据集中每年的*均天数。这个变量被提供给rolling()命令的参数。参数window表示我们将执行一年的滚动相关性。参数min_periods表示当只有完整样本大小用于计算时才会计算相关性。在这种情况下，在df_corrs数据集中的第一年没有相关性值。最后，plot()命令显示了 2017 年全年每日收益的一年滚动相关性图表。
简单移动*均线
用于时间序列数据分析的常见技术指标是移动*均线。mean()方法可用于计算rolling()命令中给定窗口的值的*均值。例如，5 天的简单移动*均线（SMA）是最*五个交易日的价格的*均值，每天在一段时间内计算一次。同样，我们也可以计算一个更长期的 30 天简单移动*均线。这两个移动*均线可以一起使用以生成交叉信号。
在下面的示例中，我们下载 ABN 的每日收盘价，计算短期和长期 SMA，并在单个图表上可视化它们：
In [ ]:
%matplotlib inline
import quandl
import pandas as pd
quandl.ApiConfig.api_key = QUANDL_API_KEY
df = quandl.get('EURONEXT/ABN.4')
df_filled = df.asfreq('D', method='ffill')
df_last = df['Last']
series_short = df_last.rolling(window=5, min_periods=5).mean()
series_long = df_last.rolling(window=30, min_periods=30).mean()
df_sma = pd.DataFrame(columns=['short', 'long'])
df_sma['short'] = series_short
df_sma['long'] = series_long
df_sma.plot(figsize=(12, 8));
这产生了以下的图表：
我们使用 5 天的*均值作为短期 SMA，30 天作为长期 SMA。min_periods参数用于排除前几行，这些行没有足够的样本大小来计算 SMA。df_sma变量是一个新创建的pandas DataFrame，用于存储 SMA 计算。然后我们绘制一个 12 英寸乘 8 英寸的图表。从图表中，我们可以看到许多点，短期 SMA 与长期 SMA 相交。图表分析师使用交叉点来识别趋势并生成信号。5 和 10 的窗口期纯粹是建议值；您可以调整这些值以找到适合自己的解释。
指数移动*均线
在计算移动*均线时的另一种方法是指数移动*均线（EMA）。请记住，简单移动*均线在窗口期内为价格分配相等的权重。然而，在 EMA 中，最*的价格被分配比旧价格更高的权重。这种权重是以指数方式分配的。
pandas DataFrame 的ewm()方法提供了指数加权函数。span参数指定了衰减行为的窗口期。使用相同的 ABN 数据集绘制 EMA 如下：
In [ ]:
%matplotlib inline
import quandl
import pandas as pd
quandl.ApiConfig.api_key = QUANDL_API_KEY
df = quandl.get('EURONEXT/ABN.4')
df_filled = df.asfreq('D', method='ffill')
df_last = df['Last']
series_short = df_last.ewm(span=5).mean()
series_long = df_last.ewm(span=30).mean()
df_sma = pd.DataFrame(columns=['short', 'long'])
df_sma['short'] = series_short
df_sma['long'] = series_long
df_sma.plot(figsize=(12, 8));
这产生了以下的图表：
SMA 和 EMA 的图表模式基本相同。由于 EMA 对最*的数据赋予的权重比旧数据更高，因此它们对价格变动的反应比 SMA 更快。
除了不同的窗口期，您还可以尝试使用 SMA 和 EMA 价格的组合来得出更多见解！
摘要
在本章中，我们使用 Python 3.7 建立了我们的工作环境，并使用虚拟环境包来管理单独的包安装。pip命令是一个方便的 Python 包管理器，可以轻松下载和安装 Python 模块，包括 Jupyter、Quandl 和pandas。Jupyter 是一个基于浏览器的交互式计算环境，用于执行 Python 代码和可视化数据。有了 Quandl 账户，我们可以轻松获取高质量的时间序列数据集。这些数据来源于各种数据发布者。数据集直接下载到一个pandas DataFrame 对象中，使我们能够执行金融分析，如绘制每日百分比收益、直方图、Q-Q 图、相关性、简单移动*均线和指数移动*均线。
第二部分：金融概念
本节涵盖了金融工程实践者讨论的金融概念和数学模型。
在本节中，我们将介绍以下章节：
-
第二章，线性在金融中的重要性
-
第三章，金融中的非线性
-
第四章，定价期权的数值方法
-
第五章，建模利率和衍生品
-
第六章，时间序列数据的统计分析
第二章：金融中线性的重要性
非线性动力学在我们的世界中起着至关重要的作用。由于更容易研究和更容易建模的能力，线性模型经常在经济学中使用。在金融领域，线性模型被广泛用于帮助定价证券和执行最优投资组合分配，以及其他有用的事情。金融建模中线性的一个重要方面是它保证问题在全局最优解处终止。
为了进行预测和预测，回归分析在统计学领域被广泛使用，以估计变量之间的关系。由于 Python 具有丰富的数学库是其最大的优势之一，因此 Python 经常被用作科学脚本语言来帮助解决这些问题。像 SciPy 和 NumPy 这样的模块包含各种线性回归函数，供数据科学家使用。
在传统的投资组合管理中，资产配置遵循线性模式，投资者有个人的投资风格。我们可以将投资组合分配问题描述为一个线性方程组，包含等式或不等式。然后，这些线性系统可以以矩阵形式表示为Ax=B，其中A是已知的系数值，B是观察到的结果，x是我们想要找出的值的向量。往往，x包含最大化效用的最优证券权重。使用矩阵代数，我们可以使用直接或间接方法高效地解出x。
在本章中，我们将涵盖以下主题：
-
检查资本资产定价模型和证券市场线
-
使用回归解决证券市场线问题
-
检查 APT 模型并执行多元线性回归
-
理解投资组合中的线性优化
-
使用 Pulp 软件包执行线性优化
-
理解线性规划的结果
-
整数规划简介
-
使用二进制条件实现线性整数规划模型
-
使用矩阵线性代数解等式的线性方程组
-
使用 LU、Cholesky 和 QR 分解直接解线性方程组
-
使用 Jacobi 和 Gauss-Seidel 方法间接解线性方程组
资本资产定价模型和证券市场线
很多金融文献都专门讨论了资本资产定价模型（CAPM）。在本节中，我们将探讨突出金融中线性的重要性的关键概念。
在著名的 CAPM 中，描述了证券的风险和回报率之间的关系如下：
对于证券i，其回报被定义为R[i]，其 beta 被定义为β[i]。CAPM 将证券的回报定义为无风险利率R[f]和其 beta 与风险溢价的乘积之和。风险溢价可以被视为市场投资组合的超额回报，不包括无风险利率。以下是 CAPM 的可视化表示：
Beta 是股票系统风险的度量 - 无法分散的风险。实质上，它描述了股票回报与市场波动的敏感性。例如，beta 为零的股票无论市场走向如何都不会产生超额回报。它只能以无风险利率增长。beta 为 1 的股票表示该股票与市场完全同步。
beta 是通过将股票与市场回报的协方差除以市场回报的方差来数学推导的。
CAPM 模型衡量了投资组合篮子中每支股票的风险和回报之间的关系。通过概述这种关系的总和，我们可以得到在每个投资组合回报水*下产生最低投资组合风险的风险证券的组合或权重。希望获得特定回报的投资者将拥有一个最佳投资组合的组合，以提供可能的最低风险。最佳投资组合的组合位于一条称为有效边界的线上。
在有效边界上，存在一个切线点，表示最佳的可用最优投资组合，并以可能的最低风险换取最高的回报率。切线点处的最佳投资组合被称为市场投资组合。
从市场投资组合到无风险利率之间存在一条直线。这条线被称为资本市场线（CML）。CML 可以被认为是所有最优投资组合中最高夏普比率的夏普比率。夏普比率是一个风险调整后的绩效指标，定义为投资组合超额回报与标准差风险单位的比率。投资者特别感兴趣持有沿着 CML 线的资产组合。以下图表说明了有效边界、市场投资组合和 CML：
CAPM 研究中另一个有趣的概念是证券市场线（SML）。SML 绘制了资产的预期回报与其贝塔值的关系。对于贝塔值为 1 的证券，其回报完全匹配市场回报。任何定价高于 SML 的证券被认为是被低估的，因为投资者期望在相同风险下获得更高的回报。相反，任何定价低于 SML 的证券被认为是被高估的。
假设我们有兴趣找到证券的贝塔值β[i]。我们可以对公司的股票回报R[i]与市场回报R[M]进行回归分析，同时加上一个截距α，形成R[i]=α+βR[M]的方程。
考虑以下一组在五个时间段内测得的股票回报和市场回报数据：
| 时间段 | 股票回报 | 市场回报 |
|---|---|---|
| 1 | 0.065 | 0.055 |
| 2 | 0.0265 | -0.09 |
| 3 | -0.0593 | -0.041 |
| 4 | -0.001 | 0.045 |
| 5 | 0.0346 | 0.022 |
使用 SciPy 的stats模块，我们将对 CAPM 模型进行最小二乘回归，并通过在 Python 中运行以下代码来得出α和β[i]的值：
In [ ]:
"""
Linear regression with SciPy
"""
from scipy import stats
stock_returns = [0.065, 0.0265, -0.0593, -0.001, 0.0346]
mkt_returns = [0.055, -0.09, -0.041, 0.045, 0.022]
beta, alpha, r_value, p_value, std_err = \
stats.linregress(stock_returns, mkt_returns)
scipty.stats.linregress函数返回五个值：回归线的斜率、回归线的截距、相关系数、零斜率假设的假设检验的 p 值，以及估计的标准误差。我们有兴趣通过打印beta和alpha的值来找到线的斜率和截距，分别为：
In [ ]:
print(beta, alpha)
Out[ ]:
0.5077431878770808 -0.008481900352462384
股票的贝塔值为 0.5077，α几乎为零。
描述 SML 的方程可以写成：
术语E[R[M]]−R[f]是市场风险溢价，E[R[M]]是市场投资组合的预期回报。R[f]是无风险利率的回报，E[R[i]]是资产i的预期回报，β[i]是资产的贝塔值。
假设无风险利率为 5%，市场风险溢价为 8.5%。股票的预期回报率是多少？根据 CAPM，贝塔值为 0.5077 的股票将有 0.5077×8.5%的风险溢价，即 4.3%。无风险利率为 5%，因此股票的预期回报率为 9.3%。
如果在同一时间段内观察到证券的回报率高于预期的股票回报（例如，10.5%），则可以说该证券被低估，因为投资者可以期望在承担相同风险的情况下获得更高的回报。
相反，如果观察到证券的回报率低于 SML 所暗示的预期回报率（例如，7%），则可以说该证券被高估。投资者在承担相同风险的情况下获得了降低的回报。
套利定价理论模型
CAPM 存在一些局限性，比如使用均值-方差框架和事实上回报被一个风险因素（市场风险因素）捕获。在一个分散投资组合中，各种股票的非系统风险相互抵消，基本上被消除了。
套利定价理论（APT）模型被提出来解决这些缺点，并提供了一种除了均值和方差之外确定资产价格的一般方法。
APT 模型假设证券回报是根据多因素模型生成的，这些模型由几个系统风险因素的线性组合组成。这些因素可能是通货膨胀率、GDP 增长率、实际利率或股利。
根据 APT 模型的均衡资产定价方程如下：
在这里，E[R[i]]是第i个证券的预期回报率，α[i]是如果所有因素都可以忽略时第i个股票的预期回报，β[i,j]是第i个资产对第j个因素的敏感性，F[j]是影响第i个证券回报的第j个因素的值。
由于我们的目标是找到α[i]和β的所有值，我们将在 APT 模型上执行多元线性回归。
因子模型的多元线性回归
许多 Python 包（如 SciPy）都带有几种回归函数的变体。特别是，statsmodels包是 SciPy 的补充，具有描述性统计信息和统计模型的估计。Statsmodels 的官方页面是www.statsmodels.org。
如果您的 Python 环境中尚未安装 Statsmodels，请运行以下命令进行安装：
$ pip install -U statsmodels
如果您已经安装了一个现有的包，-U开关告诉pip将选定的包升级到最新可用版本。
在这个例子中，我们将使用statsmodels模块的ols函数执行普通最小二乘回归，并查看其摘要。
假设您已经实现了一个包含七个因素的 APT 模型，返回Y的值。考虑在九个时间段t[1]到t[9]内收集的以下数据集。X[1]到X[7]是在每个时间段观察到的自变量。因此，回归问题的结构如下：
可以使用以下代码对X和Y的值进行简单的普通最小二乘回归：
In [ ]:
"""
Least squares regression with statsmodels
"""
import numpy as np
import statsmodels.api as sm
# Generate some sample data
num_periods = 9
all_values = np.array([np.random.random(8) \
for i in range(num_periods)])
# Filter the data
y_values = all_values[:, 0] # First column values as Y
x_values = all_values[:, 1:] # All other values as X
x_values = sm.add_constant(x_values) # Include the intercept
results = sm.OLS(y_values, x_values).fit() # Regress and fit the model
让我们查看回归的详细统计信息：
In [ ]:
print(results.summary())
OLS 回归结果将输出一个相当长的统计信息表。然而，我们感兴趣的是一个特定部分，它给出了我们 APT 模型的系数：
===================================================================
coef std err t P>|t| [0.025
-------------------------------------------------------------------
const 0.7229 0.330 2.191 0.273 -3.469
x1 0.4195 0.238 1.766 0.328 -2.599
x2 0.4930 0.176 2.807 0.218 -1.739
x3 0.1495 0.102 1.473 0.380 -1.140
x4 -0.1622 0.191 -0.847 0.552 -2.594
x5 -0.6123 0.172 -3.561 0.174 -2.797
x6 -0.2414 0.161 -1.499 0.375 -2.288
x7 -0.5079 0.200 -2.534 0.239 -3.054
coef列给出了我们回归的系数值，包括c常数和X[1]到X[7]。同样，我们可以使用params属性来显示这些感兴趣的系数：
In [ ]:
print(results.params)
Out[ ]:
[ 0.72286627 0.41950411 0.49300959 0.14951292 -0.16218313 -0.61228465 -0.24143028 -0.50786377]
两个函数调用以相同的顺序产生了 APT 模型的相同系数值。
线性优化
在 CAPM 和 APT 定价理论中，我们假设模型是线性的，并使用 Python 中的回归来解决预期的证券价格。
随着我们投资组合中证券数量的增加，也会引入一定的限制。投资组合经理在追求投资者规定的某些目标时会受到这些规则的约束。
线性优化有助于克服投资组合分配的问题。优化侧重于最小化或最大化目标函数的值。一些例子包括最大化回报和最小化波动性。这些目标通常受到某些规定的约束，例如不允许空头交易规则，或者对要投资的证券数量的限制。
不幸的是，在 Python 中，没有一个官方的包支持这个解决方案。但是，有第三方包可用，其中包含线性规划的单纯形算法的实现。为了演示目的，我们将使用 Pulp，一个开源线性规划建模器，来帮助我们解决这个特定的线性规划问题。
获取 Pulp
您可以从github.com/coin-or/pulp获取 Pulp。该项目页面包含了一份全面的文档列表，以帮助您开始优化过程。
您还可以使用pip包管理器获取 Pulp 包：
$ pip install pulp
线性规划的最大化示例
假设我们有兴趣投资两种证券X和Y。我们想要找出每三单位X和两单位Y的实际投资单位数，使得总投资单位数最大化，如果可能的话。然而，我们的投资策略有一定的限制：
-
对于每 2 单位X和 1 单位Y的投资，总量不得超过 100
-
对于每单位X和Y的投资，总量不得超过 80
-
允许投资X的总量不得超过 40
-
不允许对证券进行空头交易
最大化问题可以用数学表示如下：
受限于：
通过在x和y图上绘制约束条件，可以看到一组可行解，由阴影区域给出：
该问题可以用pulp包在 Python 中进行翻译，如下所示：
In [ ]:
"""
A simple linear optimization problem with 2 variables
"""
import pulp
x = pulp.LpVariable('x', lowBound=0)
y = pulp.LpVariable('y', lowBound=0)
problem = pulp.LpProblem(
'A simple maximization objective',
pulp.LpMaximize)
problem += 3*x + 2*y, 'The objective function'
problem += 2*x + y <= 100, '1st constraint'
problem += x + y <= 80, '2nd constraint'
problem += x <= 40, '3rd constraint'
problem.solve()
LpVariable函数声明要解决的变量。LpProblem函数用问题的文本描述和优化类型初始化问题，本例中是最大化方法。+=操作允许添加任意数量的约束，以及文本描述。最后，调用.solve()方法开始执行线性优化。要显示优化器解决的值，使用.variables()方法循环遍历每个变量并打印出其varValue。
当代码运行时生成以下输出：
In [ ]:
print("Maximization Results:")
for variable in problem.variables():
print(variable.name, '=', variable.varValue)
Out[ ]:
Maximization Results:
x = 20.0
y = 60.0
结果显示，在满足给定的一组约束条件的情况下，当x的值为 20，y的值为 60 时，可以获得最大值 180。
线性规划的结果
线性优化有三种结果，如下：
-
线性规划的局部最优解是一个可行解，其目标函数值比其附*的所有其他可行解更接*。它可能是也可能不是全局最优解，即优于每个可行解的解。
-
如果找不到解决方案，线性规划是不可行的。
-
如果最优解是无界的或无限的，线性规划是无界的。
整数规划
在我们之前调查的简单优化问题中，线性规划的最大化示例，变量被允许是连续的或分数的。如果使用分数值或结果不现实怎么办？这个问题被称为线性整数规划问题，其中所有变量都受限于整数。整数变量的一个特殊情况是二进制变量，可以是 0 或 1。在给定一组选择时，二进制变量在模型决策时特别有用。
整数规划模型经常用于运筹学中来模拟现实工作问题。通常情况下，将非线性问题陈述为线性或甚至二进制的问题需要更多的艺术而不是科学。
整数规划的最小化示例
假设我们必须从三家经销商那里购买 150 份某种场外奇特证券。经销商X报价每份合同 500 美元加上 4000 美元的手续费，无论卖出的合同数量如何。经销商Y每份合同收费 450 美元，加上 2000 美元的交易费。经销商Z每份合同收费 450 美元，加上 6000 美元的费用。经销商X最多销售 100 份合同，经销商Y最多销售 90 份，经销商Z最多销售 70 份。从任何经销商那里交易的最低交易量是 30 份合同。我们应该如何最小化购买 150 份合同的成本？
使用pulp包，让我们设置所需的变量：
In [ ]:
"""
An example of implementing an integer
programming model with binary conditions
"""
import pulp
dealers = ['X', 'Y', 'Z']
variable_costs = {'X': 500, 'Y': 350, 'Z': 450}
fixed_costs = {'X': 4000, 'Y': 2000, 'Z': 6000}
# Define PuLP variables to solve
quantities = pulp.LpVariable.dicts('quantity',
dealers,
lowBound=0,
cat=pulp.LpInteger)
is_orders = pulp.LpVariable.dicts('orders',
dealers,
cat=pulp.LpBinary)
dealers变量只是包含用于稍后引用列表和字典的字典标识符的字典。variable_costs和fixed_costs变量是包含每个经销商收取的相应合同成本和费用的字典对象。Pulp 求解器解决了quantities和is_orders的值，这些值由LpVariable函数定义。dicts()方法告诉 Pulp 将分配的变量视为字典对象，使用dealers变量进行引用。请注意，quantities变量具有一个下限（0），防止我们在任何证券中进入空头头寸。is_orders值被视为二进制对象，指示我们是否应该与任何经销商进行交易。
对建模这个整数规划问题的最佳方法是什么？乍一看，通过应用这个方程似乎相当简单：
其中以下内容为真：
该方程简单地陈述了我们希望最小化总成本，并使用二进制变量isOrder[i]来确定是否考虑从特定经销商购买的成本。
让我们在 Python 中实现这个模型：
In [ ]:
"""
This is an example of implementing an integer programming model
with binary variables the wrong way.
"""
# Initialize the model with constraints
model = pulp.LpProblem('A cost minimization problem',
pulp.LpMinimize)
model += sum([(variable_costs[i] * \
quantities[i] + \
fixed_costs[i])*is_orders[i] \
for i in dealers]), 'Minimize portfolio cost'
model += sum([quantities[i] for i in dealers]) == 150\
, 'Total contracts required'
model += 30 <= quantities['X'] <= 100\
, 'Boundary of total volume of X'
model += 30 <= quantities['Y'] <= 90\
, 'Boundary of total volume of Y'
model += 30 <= quantities['Z'] <= 70\
, 'Boundary of total volume of Z'
model.solve() # You will get an error running this code!
当我们运行求解器时会发生什么？看一下：
Out[ ]:
TypeError: Non-constant expressions cannot be multiplied
事实证明，我们试图对两个未知变量quantities和is_order进行乘法运算，无意中导致我们执行了非线性规划。这就是在执行整数规划时遇到的陷阱。
我们应该如何解决这个问题？我们可以考虑使用二进制变量，如下一节所示。
具有二进制条件的整数规划
制定最小化目标的另一种方法是将所有未知变量线性排列，使它们是可加的：
与先前的目标方程相比，我们将获得相同的固定成本值。但是，未知变量quantity[i]仍然在方程的第一项中。因此，需要将quantity[i]变量作为isOrder[i]的函数来求解，约束如下所述：
让我们在 Python 中应用这些公式：
In [ ]:
"""
This is an example of implementing an
IP model with binary variables the correct way.
"""
# Initialize the model with constraints
model = pulp.LpProblem('A cost minimization problem',
pulp.LpMinimize)
model += sum(
[variable_costs[i]*quantities[i] + \
fixed_costs[i]*is_orders[i] for i in dealers])\
, 'Minimize portfolio cost'
model += sum([quantities[i] for i in dealers]) == 150\
, 'Total contracts required'
model += is_orders['X']*30 <= quantities['X'] <= \
is_orders['X']*100, 'Boundary of total volume of X'
model += is_orders['Y']*30 <= quantities['Y'] <= \
is_orders['Y']*90, 'Boundary of total volume of Y'
model += is_orders['Z']*30 <= quantities['Z'] <= \
is_orders['Z']*70, 'Boundary of total volume of Z'
model.solve()
当我们尝试运行求解器时会发生什么？让我们看看：
In [ ]:
print('Minimization Results:')
for variable in model.variables():
print(variable, '=', variable.varValue)
print('Total cost:', pulp.value(model.objective))
Out[ ]:
Minimization Results:
orders_X = 0.0
orders_Y = 1.0
orders_Z = 1.0
quantity_X = 0.0
quantity_Y = 90.0
quantity_Z = 60.0
Total cost: 66500.0
输出告诉我们，从经销商Y购买 90 份合同和从经销商Z购买 60 份合同可以以最低成本 66,500 美元满足所有其他约束。
正如我们所看到的，需要在整数规划模型的设计中进行仔细规划，以便得出准确的解决方案，使其在决策中有用。
使用矩阵解决线性方程
在前一节中，我们看到了如何解决带有不等式约束的线性方程组。如果一组系统线性方程有确定性约束，我们可以将问题表示为矩阵，并应用矩阵代数。矩阵方法以紧凑的方式表示多个线性方程，同时使用现有的矩阵库函数。
假设我们想要建立一个包含三种证券a、b和c的投资组合。投资组合的分配必须满足一定的约束条件：必须持有证券a的多头头寸 6 单位。对于每两单位的证券a，必须投资一单位的证券b和一单位的证券c，净头寸必须是多头四单位。对于每一单位的证券a，必须投资三单位的证券b和两单位的证券c，净头寸必须是多头五单位。
要找出要投资的证券数量，我们可以用数学方式表述问题，如下：
在所有系数可见的情况下，方程如下：
让我们把方程的系数表示成矩阵形式：
线性方程现在可以陈述如下：
要解出包含要投资的证券数量的x向量，需要取矩阵A的逆，方程写为：
使用 NumPy 数组，A和B矩阵分配如下：
In [ ]:
"""
Linear algebra with NumPy matrices
"""
import numpy as np
A = np.array([[2, 1, 1],[1, 3, 2],[1, 0, 0]])
B = np.array([4, 5, 6])
我们可以使用 NumPy 的linalg.solve函数来解决一组线性标量方程：
In [ ]:
print(np.linalg.solve(A, B))
Out[ ]:
[ 6\. 15\. -23.]
投资组合需要持有 6 单位的证券a的多头头寸，15 单位的证券b，和 23 单位的证券c的空头头寸。
在投资组合管理中，我们可以使用矩阵方程系统来解决给定一组约束条件的证券的最佳权重分配。随着投资组合中证券数量的增加，A矩阵的大小增加，计算A的矩阵求逆变得计算成本高昂。因此，人们可以考虑使用 Cholesky 分解、LU 分解、QR 分解、雅各比方法或高斯-赛德尔方法等方法，将A矩阵分解为更简单的矩阵进行因式分解。
LU 分解
LU 分解，又称下三角-上三角分解，是解决方阵线性方程组的方法之一。顾名思义，LU 分解将矩阵A分解为两个矩阵的乘积：一个下三角矩阵L和一个上三角矩阵U。分解可以表示如下：
在这里，我们可以看到a=l[11]u[11]，b=l[11]u[12]，依此类推。下三角矩阵是一个矩阵，它在其下三角中包含值，其余的上三角中填充了零。上三角矩阵的情况相反。
LU 分解方法相对于 Cholesky 分解方法的明显优势在于它适用于任何方阵。后者只适用于对称和正定矩阵。
回想一下前面的例子，使用矩阵解线性方程，一个 3 x 3 的A矩阵。这次，我们将使用 SciPy 模块的linalg包来执行 LU 分解，使用以下代码：
In [ ]:
"""
LU decomposition with SciPy
"""
import numpy as np
import scipy.linalg as linalg
# Define A and B
A = np.array([
[2., 1., 1.],
[1., 3., 2.],
[1., 0., 0.]])
B = np.array([4., 5., 6.])
# Perform LU decomposition
LU = linalg.lu_factor(A)
x = linalg.lu_solve(LU, B)
要查看x的值，请执行以下命令：
In [ ]:
print(x)
Out[ ]:
[ 6\. 15\. -23.]
我们得到了a、b和c的值分别为6、15和-23。
请注意，我们在这里使用了scipy.linalg的lu_factor()方法，它给出了A矩阵的置换 LU 分解作为LU变量。我们使用了lu_solve()方法，它接受置换的 LU 分解和B向量来解方程组。
我们可以使用scipy.linalg的lu()方法显示A矩阵的 LU 分解。lu()方法返回三个变量——置换矩阵P，下三角矩阵L和上三角矩阵U——分别返回：
In [ ]:
import scipy
P, L, U = scipy.linalg.lu(A)
print('P=\n', P)
print('L=\n', L)
print('U=\n', U)
当我们打印出这些变量时，我们可以得出 LU 分解和A矩阵之间的关系如下：
LU 分解可以看作是在两个更简单的矩阵上执行的高斯消元的矩阵形式：上三角矩阵和下三角矩阵。
Cholesky 分解
Cholesky 分解是解线性方程组的另一种方法。它可以比 LU 分解快得多，并且使用的内存要少得多，因为它利用了对称矩阵的性质。然而，被分解的矩阵必须是 Hermitian（或者是实对称的并且是方阵）和正定的。这意味着A矩阵被分解为A=LLT*，其中*L*是一个下三角矩阵，对角线上有实数和正数，*LT是L的共轭转置。
让我们考虑另一个线性方程组的例子，其中A矩阵既是 Hermitian 又是正定的。同样，方程的形式是Ax=B，其中A和B取以下值：
让我们将这些矩阵表示为 NumPy 数组：
In [ ]:
"""
Cholesky decomposition with NumPy
"""
import numpy as np
A = np.array([
[10., -1., 2., 0.],
[-1., 11., -1., 3.],
[2., -1., 10., -1.],
[0., 3., -1., 8.]])
B = np.array([6., 25., -11., 15.])
L = np.linalg.cholesky(A)
numpy.linalg的cholesky()函数将计算A矩阵的下三角因子。让我们查看下三角矩阵：
In [ ]:
print(L)
Out[ ]:
[[ 3.16227766 0\. 0\. 0\. ]
[-0.31622777 3.3015148 0\. 0\. ]
[ 0.63245553 -0.24231301 3.08889696 0\. ]
[ 0\. 0.9086738 -0.25245792 2.6665665 ]]
为了验证 Cholesky 分解的结果是否正确，我们可以使用 Cholesky 分解的定义，将L乘以它的共轭转置，这将使我们回到A矩阵的值：
In [ ]:
print(np.dot(L, L.T.conj())) # A=L.L*
Out [ ]:
[[10\. -1\. 2\. 0.]
[-1\. 11\. -1\. 3.]
[ 2\. -1\. 10\. -1.]
[ 0\. 3\. -1\. 8.]]
在解出x之前，我们需要将L^Tx解为y。让我们使用numpy.linalg的solve()方法：
In [ ]:
y = np.linalg.solve(L, B) # L.L*.x=B; When L*.x=y, then L.y=B
要解出x，我们需要再次使用L的共轭转置和y来解：
In [ ]:
x = np.linalg.solve(L.T.conj(), y) # x=L*'.y
让我们打印出x的结果：
In [ ]:
print(x)
Out[ ]:
[ 1\. 2\. -1\. 1.]
输出给出了我们的a、b、c和d的x的值。
为了证明 Cholesky 分解给出了正确的值，我们可以通过将A矩阵乘以x的转置来验证答案，从而得到B的值：
In [ ] :
print(np.mat(A) * np.mat(x).T) # B=Ax
Out[ ]:
[[ 6.]
[ 25.]
[-11.]
[ 15.]]
这表明通过 Cholesky 分解得到的x的值将与B给出的相同。
QR 分解
QR 分解，也称为QR 分解，是使用矩阵解线性方程的另一种方法，非常类似于 LU 分解。要解的方程是Ax=B的形式，其中矩阵A=QR。然而，在这种情况下，A是正交矩阵Q和上三角矩阵R的乘积。QR 算法通常用于解线性最小二乘问题。
正交矩阵具有以下特性：
-
它是一个方阵。
-
将正交矩阵乘以其转置返回单位矩阵：
- 正交矩阵的逆等于其转置：
单位矩阵也是一个方阵，其主对角线包含 1，其他位置包含 0。
现在问题Ax=B可以重新表述如下：
使用 LU 分解示例中的相同变量，我们将使用scipy.linalg的qr()方法来计算我们的Q和R的值，并让y变量代表我们的BQ^T的值，代码如下：
In [ ]:
"""
QR decomposition with scipy
"""
import numpy as np
import scipy.linalg as linalg
A = np.array([
[2., 1., 1.],
[1., 3., 2.],
[1., 0., 0]])
B = np.array([4., 5., 6.])
Q, R = scipy.linalg.qr(A) # QR decomposition
y = np.dot(Q.T, B) # Let y=Q'.B
x = scipy.linalg.solve(R, y) # Solve Rx=y
注意Q.T只是Q的转置，也就是Q的逆：
In [ ]:
print(x)
Out[ ]:
[ 6\. 15\. -23.]
我们得到了与 LU 分解示例中相同的答案。
使用其他矩阵代数方法求解
到目前为止，我们已经看过了使用矩阵求逆、LU 分解、Cholesky 分解和 QR 分解来解线性方程组。如果A矩阵中的财务数据规模很大，可以通过多种方案进行分解，以便使用矩阵代数更快地收敛。量化投资组合分析师应该熟悉这些方法。
在某些情况下，我们寻找的解可能不会收敛。因此，您可能需要考虑使用迭代方法。解决线性方程组的常见迭代方法包括雅各比方法、高斯-赛德尔方法和 SOR 方法。我们将简要介绍实现雅各比和高斯-赛德尔方法的示例。
雅各比方法
雅各比方法沿着其对角线元素迭代地解决线性方程组。当解收敛时，迭代过程终止。同样，要解决的方程式是Ax=B的形式，其中矩阵A可以分解为两个相同大小的矩阵，使得A=D+R。矩阵 D 仅包含 A 的对角分量，另一个矩阵 R 包含其余分量。让我们看一个 4 x 4 的A矩阵的例子：
然后迭代地获得解如下：
与高斯-赛德尔方法相反，在雅各比方法中，需要在每次迭代中使用x[n]的值来计算x[n+1]，并且不能被覆盖。这将占用两倍的存储空间。然而，每个元素的计算可以并行进行，这对于更快的计算是有用的。
如果A矩阵是严格不可约对角占优的，这种方法保证收敛。严格不可约对角占优矩阵是指每一行的绝对对角元素大于其他项的绝对值之和。
在某些情况下，即使不满足这些条件，雅各比方法也可以收敛。Python 代码如下：
In [ ]:
"""
Solve Ax=B with the Jacobi method
"""
import numpy as np
def jacobi(A, B, n, tol=1e-10):
# Initializes x with zeroes with same shape and type as B
x = np.zeros_like(B)
for iter_count in range(n):
x_new = np.zeros_like(x)
for i in range(A.shape[0]):
s1 = np.dot(A[i, :i], x[:i])
s2 = np.dot(A[i, i + 1:], x[i + 1:])
x_new[i] = (B[i] - s1 - s2) / A[i, i]
if np.allclose(x, x_new, tol):
break
x = x_new
return x
考虑 Cholesky 分解示例中的相同矩阵值。我们将在我们的jacobi函数中使用 25 次迭代来找到x的值：
In [ ] :
A = np.array([
[10., -1., 2., 0.],
[-1., 11., -1., 3.],
[2., -1., 10., -1.],
[0.0, 3., -1., 8.]])
B = np.array([6., 25., -11., 15.])
n = 25
初始化值后，我们现在可以调用函数并求解x：
In [ ]:
x = jacobi(A, B, n)
print('x', '=', x)
Out[ ]:
x = [ 1\. 2\. -1\. 1.]
我们求解了x的值，这与 Cholesky 分解的答案类似。
高斯-赛德尔方法
高斯-赛德尔方法与雅各比方法非常相似。这是使用迭代过程以Ax=B形式的方程解决线性方程组的另一种方法。在这里，A矩阵被分解为A=L+U，其中A矩阵是下三角矩阵L和上三角矩阵U的和。让我们看一个 4 x 4 A矩阵的例子：
然后通过迭代获得解决方案，如下所示：
使用下三角矩阵L，其中零填充上三角，可以在每次迭代中覆盖x[n]的元素，以计算x[n+1]。这样做的好处是使用雅各比方法时所需的存储空间减少了一半。
高斯-赛德尔方法的收敛速度主要取决于A矩阵的性质，特别是如果需要严格对角占优或对称正定的A矩阵。即使不满足这些条件，高斯-赛德尔方法也可能收敛。
高斯-赛德尔方法的 Python 实现如下：
In [ ]:
"""
Solve Ax=B with the Gauss-Seidel method
"""
import numpy as np
def gauss(A, B, n, tol=1e-10):
L = np.tril(A) # returns the lower triangular matrix of A
U = A-L # decompose A = L + U
L_inv = np.linalg.inv(L)
x = np.zeros_like(B)
for i in range(n):
Ux = np.dot(U, x)
x_new = np.dot(L_inv, B - Ux)
if np.allclose(x, x_new, tol):
break
x = x_new
return x
在这里，NumPy 的tril()方法返回下三角A矩阵，从中我们可以找到下三角U矩阵。将剩余的值迭代地插入x，将会得到以下解，其中由tol定义了一些容差。
让我们考虑雅各比方法和乔列斯基分解示例中的相同矩阵值。我们将在我们的gauss()函数中使用最多 100 次迭代来找到x的值，如下所示：
In [ ]:
A = np.array([
[10., -1., 2., 0.],
[-1., 11., -1., 3.],
[2., -1., 10., -1.],
[0.0, 3., -1., 8.]])
B = np.array([6., 25., -11., 15.])
n = 100
x = gauss(A, B, n)
让我们看看我们的x值是否与雅各比方法和乔列斯基分解中的值匹配：
In [ ]:
print('x', '=', x)
Out[ ]:
x = [ 1\. 2\. -1\. 1.]
我们解出了x的值，这些值与雅各比方法和乔列斯基分解的答案类似。
总结
在本章中，我们简要介绍了 CAPM 模型和 APT 模型在金融中的应用。在 CAPM 模型中，我们访问了 CML 的有效边界，以确定最佳投资组合和市场投资组合。然后，我们使用回归解决了 SML，这有助于我们确定资产是被低估还是被高估。在 APT 模型中，我们探讨了除使用均值方差框架之外，各种因素如何影响证券回报。我们进行了多元线性回归，以帮助我们确定导致证券价格估值的因素的系数。
在投资组合配置中，投资组合经理通常被投资者授权实现一组目标，同时遵循某些约束。我们可以使用线性规划来建模这个问题。使用 Pulp Python 包，我们可以定义一个最小化或最大化的目标函数，并为我们的问题添加不等式约束以解决未知变量。线性优化的三种结果可以是无界解、仅有一个解或根本没有解。
线性优化的另一种形式是整数规划，其中所有变量都受限于整数，而不是分数值。整数变量的特殊情况是二进制变量，它可以是 0 或 1，特别适用于在给定一组选择时建模决策。我们研究了一个包含二进制条件的简单整数规划模型，并看到了遇到陷阱有多容易。需要仔细规划整数规划模型的设计，以便它们在决策中有用。
投资组合分配问题也可以表示为一个具有相等性的线性方程组，可以使用矩阵形式的Ax=B来求解。为了找到x的值，我们使用了各种类型的A矩阵分解来求解A^(−1)B。矩阵分解方法有两种类型，直接和间接方法。直接方法在固定次数的迭代中执行矩阵代数运算，包括 LU 分解、Cholesky 分解和 QR 分解方法。间接或迭代方法通过迭代计算x的下一个值，直到达到一定的精度容差。这种方法特别适用于计算大型矩阵，但也面临着解不收敛的风险。我们使用的间接方法有雅各比方法和高斯-赛德尔方法。
在下一章中，我们将讨论金融中的非线性建模。
第三章：金融中的非线性
*年来，对经济和金融理论中的非线性现象的研究越来越受到关注。由于非线性串行依赖在许多金融时间序列的回报中起着重要作用，这使得证券估值和定价非常重要，从而导致对金融产品的非线性建模研究增加。
金融业的从业者使用非线性模型来预测波动性、价格衍生品，并计算风险价值（VAR）。与线性模型不同，线性代数用于寻找解决方案，非线性模型不一定推断出全局最优解。通常采用数值根查找方法来收敛到最*的局部最优解，即根。
在本章中，我们将讨论以下主题：
-
非线性建模
-
非线性模型的例子
-
根查找算法
-
根查找中的 SciPy 实现
非线性建模
尽管线性关系旨在以最简单的方式解释观察到的现象，但许多复杂的物理现象无法用这样的模型来解释。非线性关系定义如下：
尽管非线性关系可能很复杂，但为了充分理解和建模它们，我们将看一下在金融和时间序列模型的背景下应用的例子。
非线性模型的例子
许多非线性模型已被提出用于学术和应用研究，以解释线性模型无法解释的经济和金融数据的某些方面。金融领域的非线性文献实在太广泛和深刻，无法在本书中得到充分解释。在本节中，我们将简要讨论一些非线性模型的例子，这些模型可能在实际应用中遇到：隐含波动率模型、马尔可夫转换模型、阈值模型和*滑转换模型。
隐含波动率模型
也许最广泛研究的期权定价模型之一是 Black-Scholes-Merton 模型，或简称 Black-Scholes 模型。看涨期权是在特定价格和时间购买基础证券的权利，而不是义务。看跌期权是在特定价格和时间出售基础证券的权利，而不是义务。Black-Scholes 模型有助于确定期权的公*价格，假设基础证券的回报服从正态分布（N(.)）或资产价格呈对数正态分布。
该公式假定以下变量：行权价（K）、到期时间（T）、无风险利率（r）、基础收益的波动率（σ）、基础资产的当前价格（S）和其收益（q）。看涨期权的数学公式，C(S,t)，表示如下：
在这里：
通过市场力量，期权的价格可能偏离从 Black-Scholes 公式推导出的价格。特别是，实现波动性（即从历史市场价格观察到的基础收益的波动性）可能与由 Black-Scholes 模型隐含的波动性值不同，这由σ表示。
回想一下在第二章中讨论的资本资产定价模型（CAPM），即金融中的线性重要性。一般来说，具有更高回报的证券表现出更高的风险，这表明回报的波动性或标准差。
由于波动性在证券定价中非常重要，因此已经提出了许多波动性模型进行研究。其中一种模型是期权价格的隐含波动率建模。
假设我们绘制了给定特定到期日的 Black-Scholes 公式给出的股票期权的隐含波动率值。一般来说，我们得到一个常被称为波动率微笑的曲线，因为它的形状：
隐含波动率通常在深度实值和虚值期权上最高，受到大量投机驱动，而在*值期权上最低。
期权的特征解释如下：
-
实值期权（ITM）：当认购期权的行权价低于标的资产的市场价格时，被视为实值期权。当认沽期权的行权价高于标的资产的市场价格时，被视为实值期权。实值期权在行使时具有内在价值。
-
虚值期权（OTM）：当认购期权的行权价高于标的资产的市场价格时，被视为虚值期权。当认沽期权的行权价低于标的资产的市场价格时，被视为虚值期权。虚值期权在行使时没有内在价值，但可能仍具有时间价值。
-
*值期权（ATM）：当期权的行权价与标的资产的市场价格相同时，被视为*值期权。*值期权在行使时没有内在价值，但可能仍具有时间价值。
从前述波动率曲线中，隐含波动率建模的一个目标是寻找可能的最低隐含波动率值，或者换句话说，找到根。一旦找到，就可以推断出特定到期日的*值期权的理论价格，并与市场价格进行比较，以寻找潜在的机会，比如研究接**值期权或远虚值期权。然而，由于曲线是非线性的，线性代数无法充分解决根的问题。我们将在下一节根查找算法中看一些根查找方法。
马尔可夫转换模型
为了对经济和金融时间序列中的非线性行为进行建模，可以使用马尔可夫转换模型来描述不同世界或状态下的时间序列。这些状态的例子可能是一个波动状态，就像在 2008 年全球经济衰退中看到的，或者是稳步复苏经济的增长状态。能够在这些结构之间转换的能力让模型能够捕捉复杂的动态模式。
股票价格的马尔可夫性质意味着只有当前价值对于预测未来是相关的。过去的股价波动对于当前的出现方式是无关紧要的。
让我们以m=2个状态的马尔可夫转换模型为例：
ϵ[t]是一个独立同分布（i.i.d）白噪声。白噪声是一个均值为零的正态分布随机过程。同样的模型可以用虚拟变量表示：
马尔可夫转换模型的应用包括代表实际 GDP 增长率和通货膨胀率动态。这些模型反过来推动利率衍生品的估值模型。从前一状态i转换到当前状态j的概率可以写成如下形式：
阈自回归模型
一个流行的非线性时间序列模型类别是阈自回归（TAR）模型，它看起来与马尔可夫转换模型非常相似。使用回归方法，简单的 AR 模型可以说是最流行的模型来解释非线性行为。阈值模型中的状态是由其自身时间序列的过去值d相对于阈值c确定的。
以下是一个自激励 TAR（SETAR）模型的例子。SETAR 模型是自激励的，因为在不同制度之间的切换取决于其自身时间序列的过去值：
使用虚拟变量，SETAR 模型也可以表示如下：
TAR 模型的使用可能会导致状态之间出现急剧的转变，这由阈值变量c控制。
*滑转换模型
阈值模型中的突然制度变化似乎与现实世界的动态不符。通过引入一个*滑变化的连续函数，可以克服这个问题，从一个制度*滑地过渡到另一个制度。SETAR 模型成为一个逻辑*滑转换阈值自回归（LSTAR）模型，其中使用逻辑函数G(y[t−1];γ,c)：
SETAR 模型现在变成了 LSTAR 模型，如下方程所示：
参数γ控制从一个制度到另一个制度的*滑过渡。对于γ的大值，过渡是最快的，因为y[t−d]接*阈值变量c。当γ=0 时，LSTAR 模型等同于简单的AR(1)单制度模型。
根查找算法
在前面的部分，我们讨论了一些用于研究经济和金融时间序列的非线性模型。从连续时间给定的模型数据，因此意图是搜索可能推断有价值信息的极值点。使用数值方法，如根查找算法，可以帮助我们找到连续函数f的根，使得f(x)=0，这可能是函数的极大值或极小值。一般来说，一个方程可能包含多个根，也可能根本没有根。
在非线性模型上使用根查找方法的一个例子是前面讨论的 Black-Scholes 隐含波动率建模，在隐含波动率模型部分。期权交易员有兴趣根据 Black-Scholes 模型推导隐含价格，并将其与市场价格进行比较。在下一章，期权定价的数值方法，我们将看到如何将根查找方法与数值期权定价程序结合起来，以根据特定期权的市场价格创建一个隐含波动率模型。
根查找方法使用迭代程序，需要一个起始点或根的估计。根的估计可能会收敛于一个解，收敛于一个不需要的根，或者根本找不到解决方案。因此，找到根的良好*似是至关重要的。
并非每个非线性函数都可以使用根查找方法解决。下图显示了一个连续函数的例子，
，在这种情况下，根查找方法可能无法找到解。在x=0和x=2处，y值在-20 到 20 的范围内存在不连续性：
并没有固定的规则来定义良好的*似。建议在开始根查找迭代程序之前，先确定下界和上界的搜索范围。我们当然不希望在错误的方向上无休止地搜索我们的根。
增量搜索
解决非线性函数的一种粗糙方法是进行增量搜索。使用任意起始点a，我们可以获得每个dx增量的f(a)值。我们假设f(a+dx)，f(a+2dx)，f(a+3dx)…的值与它们的符号指示的方向相同。一旦符号改变，解决方案被认为已找到。否则，当迭代搜索越过边界点b时，迭代搜索终止。
迭代的根查找方法的图示示例如下图所示：
可以从以下 Python 代码中看到一个例子：
In [ ]:
"""
An incremental search algorithm
"""
import numpy as np
def incremental_search(func, a, b, dx):
"""
:param func: The function to solve
:param a: The left boundary x-axis value
:param b: The right boundary x-axis value
:param dx: The incremental value in searching
:return:
The x-axis value of the root,
number of iterations used
"""
fa = func(a)
c = a + dx
fc = func(c)
n = 1
while np.sign(fa) == np.sign(fc):
if a >= b:
return a - dx, n
a = c
fa = fc
c = a + dx
fc = func(c)
n += 1
if fa == 0:
return a, n
elif fc == 0:
return c, n
else:
return (a + c)/2., n
在每次迭代过程中，a将被c替换，并且在下一次比较之前，c将被dx递增。如果找到了根，那么它可能位于a和c之间，两者都包括在内。如果解决方案不在任何一个点上，我们将简单地返回这两个点的*均值作为最佳估计。变量n跟踪经历了寻找根的过程的迭代次数。
我们将使用具有解析解
的方程来演示和测量我们的根查找器，其中x被限制在-5 和 5 之间。给出了一个小的dx值 0.001，它也充当精度工具。较小的dx值产生更好的精度，但也需要更多的搜索迭代：
In [ ]:
# The keyword 'lambda' creates an anonymous function
# with input argument x
y = lambda x: x**3 + 2.*x**2 - 5.
root, iterations = incremental_search (y, -5., 5., 0.001)
print("Root is:", root)
print("Iterations:", iterations)
Out[ ]:
Root is: 1.2414999999999783
Iterations: 6242
增量搜索根查找方法是根查找算法基本行为的基本演示。当由dx定义时，精度最佳，并且在最坏的情况下需要极长的计算时间。要求的精度越高，解决方案收敛所需的时间就越长。出于实际原因，这种方法是所有根查找算法中最不受欢迎的，我们将研究替代方法来找到我们方程的根，以获得更好的性能。
二分法
二分法被认为是最简单的一维根查找算法。一般的兴趣是找到连续函数f的值x，使得f(x)=0。
假设我们知道区间a和b的两个点，其中a<b，并且连续函数上有f(a)<0和f(b)>0，则取该区间的中点作为c，其中
；然后二分法计算该值f(c)。
让我们用以下图示来说明沿着非线性函数设置点的情况：
由于f(a)的值为负，f(b)的值为正，二分法假设根x位于a和b之间，并给出f(x)=0。
如果f(c)=0或者非常接*零，通过预定的误差容限值，就宣布找到了一个根。如果f(c)<0，我们可以得出结论，根存在于c和b的区间，或者a和c的区间。
在下一次评估中，c将相应地替换为a或b。随着新区间缩短，二分法重复相同的评估，以确定c的下一个值。这个过程继续，缩小ab的宽度，直到根被认为找到。
使用二分法的最大优势是保证在给定预定的误差容限水*和允许的最大迭代次数下收敛到根的*似值。应该注意的是，二分法不需要未知函数的导数知识。在某些连续函数中，导数可能是复杂的，甚至不可能计算。这使得二分法对于处理不*滑函数非常有价值。
由于二分法不需要来自连续函数的导数信息，其主要缺点是在迭代评估中花费更多的计算时间。此外，由于二分法的搜索边界位于a和b区间内，因此需要一个良好的*似值来确保根落在此范围内。否则，可能会得到不正确的解，甚至根本没有解。使用较大的a和b值可能会消耗更多的计算时间。
二分法被认为是稳定的，无需使用初始猜测值即可收敛。通常，它与其他方法结合使用，例如更快的牛顿法，以便快速收敛并获得精度。
二分法的 Python 代码如下。将其保存为bisection.py：
In [ ]:
"""
The bisection method
"""
def bisection(func, a, b, tol=0.1, maxiter=10):
"""
:param func: The function to solve
:param a: The x-axis value where f(a)<0
:param b: The x-axis value where f(b)>0
:param tol: The precision of the solution
:param maxiter: Maximum number of iterations
:return:
The x-axis value of the root,
number of iterations used
"""
c = (a+b)*0.5 # Declare c as the midpoint ab
n = 1 # Start with 1 iteration
while n <= maxiter:
c = (a+b)*0.5
if func(c) == 0 or abs(a-b)*0.5 < tol:
# Root is found or is very close
return c, n
n += 1
if func(c) < 0:
a = c
else:
b = c
return c, n
In [ ]:
y = lambda x: x**3 + 2.*x**2 - 5
root, iterations = bisection(y, -5, 5, 0.00001, 100)
print("Root is:", root)
print("Iterations:", iterations)
Out[ ]:
Root is: 1.241903305053711
Iterations: 20
再次，我们将匿名的lambda函数绑定到y变量，带有输入参数x，并尝试解决
方程，与之前一样，在-5 到 5 之间的区间内，精度为 0.00001，最大迭代次数为 100。
正如我们所看到的，与增量搜索方法相比，二分法给出了更好的精度，迭代次数更少。
牛顿法
牛顿法，也称为牛顿-拉弗森法，使用迭代过程来求解根，利用函数的导数信息。导数被视为一个要解决的线性问题。函数的一阶导数f′代表切线。给定x的下一个值的*似值，记为x[1]，如下所示：
在这里，切线与x轴相交于x[1]，产生y=0。这也表示关于x[1]的一阶泰勒展开，使得新点
解决以下方程：
重复这个过程，x取值为x[1]，直到达到最大迭代次数，或者x[1]和x之间的绝对差在可接受的精度水*内。
需要一个初始猜测值来计算f(x)和f'(x)的值。收敛速度是二次的，被认为是以极高的精度获得解决方案的非常快速的速度。
牛顿法的缺点是它不能保证全局收敛到解决方案。当函数包含多个根或算法到达局部极值且无法计算下一步时，就会出现这种情况。由于该方法需要知道其输入函数的导数，因此需要输入函数可微。然而，在某些情况下，函数的导数是无法知道的，或者在数学上很容易计算。
牛顿法的图形表示如下截图所示。x[0]是初始x值。评估f(x[0])的导数，这是一个切线，穿过x轴在x[1]处。迭代重复，评估点x[1]，x[2]，x[3]等处的导数：
Python 中牛顿法的实现如下：
In [ ]:
"""
The Newton-Raphson method
"""
def newton(func, df, x, tol=0.001, maxiter=100):
"""
:param func: The function to solve
:param df: The derivative function of f
:param x: Initial guess value of x
:param tol: The precision of the solution
:param maxiter: Maximum number of iterations
:return:
The x-axis value of the root,
number of iterations used
"""
n = 1
while n <= maxiter:
x1 = x - func(x)/df(x)
if abs(x1 - x) < tol: # Root is very close
return x1, n
x = x1
n += 1
return None, n
我们将使用二分法示例中使用的相同函数，并查看牛顿法的结果：
In [ ]:
y = lambda x: x**3 + 2*x**2 - 5
dy = lambda x: 3*x**2 + 4*x
root, iterations = newton(y, dy, 5.0, 0.00001, 100)
print("Root is:", root)
print("Iterations:", iterations)
Out [ ]:
Root is: 1.241896563034502
Iterations: 7
注意除零异常！在 Python 2 中，使用值如 5.0，而不是 5，让 Python 将变量识别为浮点数，避免了将变量视为整数进行计算的问题，并给出了更好的精度。
使用牛顿法，我们获得了一个非常接*的解，迭代次数比二分法少。
割线法
使用割线法来寻找根。割线是一条直线，它与曲线的两个点相交。在割线法中，画一条直线连接连续函数上的两个点，使其延伸并与x轴相交。这种方法可以被视为拟牛顿法。通过连续画出这样的割线，可以逼*函数的根。
割线法在以下截图中以图形方式表示。需要找到两个x轴值的初始猜测，a和b，以找到f(a)和f(b)。从f(b)到f(a)画一条割线y，并在x轴上的点c处相交，使得：
因此，c的解如下：
在下一次迭代中，a和b将分别取b和c的值。该方法重复自身，为x轴值a和b、b和c、c和d等画出割线。当达到最大迭代次数或b和c之间的差异达到预先指定的容限水*时，解决方案终止，如下图所示：
割线法的收敛速度被认为是超线性的。其割线法比二分法收敛速度快得多，但比牛顿法慢。在牛顿法中，浮点运算的数量在每次迭代中占用的时间是割线法的两倍，因为需要计算函数和导数。由于割线法只需要在每一步计算其函数，因此在绝对时间上可以认为更快。
割线法的初始猜测值必须接*根，否则无法保证收敛到解。
割线法的 Python 代码如下所示：
In [ ]:
"""
The secant root-finding method
"""
def secant(func, a, b, tol=0.001, maxiter=100):
"""
:param func: The function to solve
:param a: Initial x-axis guess value
:param b: Initial x-axis guess value, where b>a
:param tol: The precision of the solution
:param maxiter: Maximum number of iterations
:return:
The x-axis value of the root,
number of iterations used
"""
n = 1
while n <= maxiter:
c = b - func(b)*((b-a)/(func(b)-func(a)))
if abs(c-b) < tol:
return c, n
a = b
b = c
n += 1
return None, n
再次重用相同的非线性函数，并返回割线法的结果：
In [ ]:
y = lambda x: x**3 + 2.*x**2 - 5.
root, iterations = secant(y, -5.0, 5.0, 0.00001, 100)
print("Root is:", root)
print("Iterations:", iterations)
Out[ ]:
Root is: 1.2418965622558549
Iterations: 14
尽管所有先前的根查找方法都给出了非常接*的解，割线法与二分法相比，迭代次数更少，但比牛顿法多。
组合根查找方法
完全可以使用前面提到的根查找方法的组合来编写自己的根查找算法。例如，可以使用以下实现：
-
使用更快的割线法将问题收敛到预先指定的误差容限值或最大迭代次数
-
一旦达到预先指定的容限水*，就切换到使用二分法，通过每次迭代将搜索区间减半，直到找到根
布伦特法或Wijngaarden-Dekker-Brent方法结合了二分根查找方法、割线法和反向二次插值。该算法尝试在可能的情况下使用割线法或反向二次插值，并在必要时使用二分法。布伦特法也可以在 SciPy 的scipy.optimize.brentq函数中找到。
根查找中的 SciPy 实现
在开始编写根查找算法来解决非线性甚至线性问题之前，先看看scipy.optimize方法的文档。SciPy 包含一系列科学计算函数，作为 Python 的扩展。这些开源算法很可能适合您的应用程序。
查找标量函数的根
scipy.optimize模块中可以找到一些根查找函数，包括bisect、newton、brentq和ridder。让我们使用 SciPy 的实现来设置我们在增量搜索部分中讨论过的示例：
In [ ]:
"""
Documentation at
http://docs.scipy.org/doc/scipy/reference/optimize.html
"""
import scipy.optimize as optimize
y = lambda x: x**3 + 2.*x**2 - 5.
dy = lambda x: 3.*x**2 + 4.*x
# Call method: bisect(f, a, b[, args, xtol, rtol, maxiter, ...])
print("Bisection method:", optimize.bisect(y, -5., 5., xtol=0.00001))
# Call method: newton(func, x0[, fprime, args, tol, ...])
print("Newton's method:", optimize.newton(y, 5., fprime=dy))
# When fprime=None, then the secant method is used.
print("Secant method:", optimize.newton(y, 5.))
# Call method: brentq(f, a, b[, args, xtol, rtol, maxiter, ...])
print("Brent's method:", optimize.brentq(y, -5., 5.))
当运行上述代码时，将生成以下输出：
Out[ ]:
Bisection method: 1.241903305053711
Newton's method: 1.2418965630344798
Secant method: 1.2418965630344803
Brent's method: 1.241896563034559
我们可以看到，SciPy 的实现给出了与我们推导的答案非常相似的答案。
值得注意的是，SciPy 对每个实现都有一组明确定义的条件。例如，在文档中二分法例程的函数调用如下所示：
scipy.optimize.bisect(f, a, b, args=(), xtol=1e-12, rtol=4.4408920985006262e-16, maxiter=100, full_output=False, disp=True)
该函数将严格评估函数f以返回函数的零点。f(a)和f(b)不能具有相同的符号。在某些情况下，很难满足这些约束条件。例如，在解决非线性隐含波动率模型时，波动率值不能为负。在活跃市场中，如果不修改基础实现，几乎不可能找到波动率函数的根或零点。在这种情况下，实现我们自己的根查找方法也许可以让我们更加掌控我们的应用程序应该如何执行。
一般非线性求解器
scipy.optimize模块还包含了我们可以利用的多维通用求解器。root和fsolve函数是一些具有以下函数属性的例子：
-
root(fun, x0[, args, method, jac, tol, ...])：这找到向量函数的根。 -
fsolve(func, x0[, args, fprime, ...])：这找到函数的根。
输出以字典对象的形式返回。使用我们的示例作为这些函数的输入，我们将得到以下输出：
In [ ]:
import scipy.optimize as optimize
y = lambda x: x**3 + 2.*x**2 - 5.
dy = lambda x: 3.*x**2 + 4.*x
print(optimize.fsolve(y, 5., fprime=dy))
Out[ ]:
[1.24189656]
In [ ]:
print(optimize.root(y, 5.))
Out[ ]:
fjac: array([[-1.]])
fun: array([3.55271368e-15])
message: 'The solution converged.'
nfev: 12
qtf: array([-3.73605502e-09])
r: array([-9.59451815])
status: 1
success: True
x: array([1.24189656])
使用初始猜测值5，我们的解收敛到了根1.24189656，这与我们迄今为止得到的答案非常接*。当我们选择图表另一侧的值时会发生什么？让我们使用初始猜测值-5：
In [ ]:
print(optimize.fsolve(y, -5., fprime=dy))
Out[ ]:
[-1.33306553]
c:\python37\lib\site-packages\scipy\optimize\minpack.py:163: RuntimeWarning: The iteration is not making good progress, as measured by the
improvement from the last ten iterations.
warnings.warn(msg, RuntimeWarning)
In [ ]:
print(optimize.root(y, -5.))
Out[ ]:
fjac: array([[-1.]])
fun: array([-3.81481496])
message: 'The iteration is not making good progress, as measured by the \n improvement from the last ten iterations.'
nfev: 28
qtf: array([3.81481521])
r: array([-0.00461503])
status: 5
success: False
x: array([-1.33306551])
从显示输出中可以看出，算法没有收敛，并返回了一个与我们先前答案略有不同的根。如果我们在图表上看方程，会发现曲线上有许多点非常接*根。需要一个根查找器来获得所需的精度水*，而求解器则试图以最快的时间解决最*的答案。
总结
在本章中，我们简要讨论了经济和金融中非线性的持久性。我们看了一些在金融中常用的非线性模型，用来解释线性模型无法解释的数据的某些方面：Black-Scholes 隐含波动率模型、Markov 转换模型、阈值模型和*滑转换模型。
在 Black-Scholes 隐含波动率建模中，我们讨论了波动率微笑，它由通过 Black-Scholes 模型从特定到期日的看涨或看跌期权的市场价格推导出的隐含波动率组成。您可能会对寻找可能的最低隐含波动率值感兴趣，这对于推断理论价格并与潜在机会的市场价格进行比较可能是有用的。然而，由于曲线是非线性的，线性代数无法充分解决最优点的问题。为此，我们将需要使用根查找方法。
根查找方法试图找到函数的根或零点。我们讨论了常见的根查找方法，如二分法、牛顿法和割线法。使用根查找算法的组合可能有助于更快地找到复杂函数的根。Brent 方法就是一个例子。
我们探讨了scipy.optimize模块中包含的这些根查找方法的功能，尽管有约束条件。其中一个约束条件要求两个边界输入值被评估为一个负值和一个正值的对，以便解收敛成功。在隐含波动率建模中，这种评估几乎是不可能的，因为波动率没有负值。实现我们自己的根查找方法也许可以让我们更加掌控我们的应用程序应该如何执行。
使用通用求解器是另一种寻找根的方法。它们可能会更快地收敛到我们的解，但这种收敛并不由初始给定的值保证。
非线性建模和优化本质上是一项复杂的任务，没有通用的解决方案或确定的方法来得出结论。本章旨在介绍金融领域的非线性研究。
在下一章中，我们将介绍常用于期权定价的数值方法。通过将数值程序与寻根算法配对，我们将学习如何利用股票期权的市场价格构建隐含波动率模型。
