书名:量化交易学习指南——基于R语言
ISBN:978-7-115-49874-8
本书由人民邮电出版社发行数字版。版权所有,侵权必究。
您购买的人民邮电出版社电子书仅供您个人使用,未经授权,不得以任何方式复制和传播本书内容。
我们愿意相信读者具有这样的良知和觉悟,与我们共同保护知识产权。
如果购买者有侵权行为,我们可能对该用户实施包括但不限于关闭该帐号等维权措施,并可能追究法律责任。
著 [印度] 帕勒姆•吉特(Param Jeet)
[印度] 普拉桑特•瓦次(Prashant Vats)
译 曾永艺 许健男
责任编辑 胡俊英
人民邮电出版社出版发行 北京市丰台区成寿寺路11号
邮编 100164 电子邮件 315@ptpress.com.cn
读者服务热线:(010)81055410
反盗版热线:(010)81055315
Copyright ©2017 Packt Publishing. First published in the English language under the title Learning Quantitative Finance with R.
All rights reserved.
本书由英国Packt Publishing公司授权人民邮电出版社出版。未经出版者书面许可,对本书的任何部分不得以任何方式或任何手段复制和传播。
版权所有,侵权必究。
R既是统计、挖掘、计算、分析、制图等方面的优秀工具,也是一个强大的开发与应用平台。在大数据时代,任何与数据相关的难题,都可以借助R语言来解决。而金融领域正是与数据密切相关的行业,可以通过R这一工具来实现量化金融建模与量化交易。
本书包括9章内容,书中包含诸多真实的金融案例,旨在通过循序渐进的讲解帮助读者了解R,并学会在量化金融与量化交易中使用R。本书还介绍了有关统计建模、计量分析与小波分析、时间序列建模、算法交易、基于机器学习的交易、风险管理、最优化、衍生品定价等重要内容。
本书适合对R及其应用感兴趣的读者阅读,尤其适合想要在量化交易中使用R的读者学习。本书并不要求读者具备R编程的知识,但希望读者对数学分析有一些了解。
Param Jeet从印度理工学院马德拉斯分校(IITM)获得数学博士学位。Param Jeet博士在多个国际刊物上发表多篇数学方面的研究论文。过去几年,Param Jeet博士进入数据分析行业,并以数据科学家的身份为多家顶级的跨国公司工作或提供咨询服务。
感谢我的父母亲S. Dhayan Singh和Jeet Kaur,在我生命中的每个阶段给我提供一贯的支持;感谢我的妻子Manpreet Kaur,总是把我摆在第一位,全力支持并鼓励我撰写本书;感谢我的孩子Kavan Singh,他天真无邪的笑容让我倍加珍爱家人并努力工作。还要感谢我的博士论文指导老师Satyajit Roy教授这么多年来给我指导。感谢同事和朋友,如果没有他们的帮助,这本书绝不可能完成。我希望能将我的知识分享给那些迫切想要用R学习量化金融的读者。
Prashant Vats从印度理工学院获得数学硕士学位。Prashant进入数据分析行业已有10年时间,作为数据科学家供职于多家顶级的跨国公司,提供多领域的咨询服务。
感谢我的父母亲Devendra K. Singh和Sushila Sinha,允许我追逐梦想并在我的职业生涯中给我以支持;感谢我的妻子Namrata,在人生的每个阶段总陪在我身边,支持我撰写本书;感谢我的孩子Aahan Vats,他的微笑总是给我以灵感。还要感谢我的导师、同事和朋友,如果没有他们的帮助,这本书不可能完成。我希望能将自己的知识分享给更多有需要的人。
Manuel Amunategui是一位应用数据科学家,为多个行业(包括医疗、金融和销售)开发企业级的预测分析解决方案。在此之前,他作为量化开发人员在华尔街为一家大规模的股票和期权做市商工作6年,并曾作为软件开发人员在微软公司工作4年。他拥有美国西北大学预测分析专业的硕士学位以及国际训练学院国际管理专业的硕士学位。目前他是SpringML公司的数据科学副总裁。这是一家提供先进的预测性CRM分析建议、仪表盘和自动报告工具的创业公司。SpringML公司的客户包括谷歌云平台、雪佛兰、雅马哈、Tesoro和Salesforce。他还是数据科学的倡导者、博主/视频博主(amunategui.github.io)以及Udemy.com和O’Reilly Media的培训师。
本书借助统计语言R来阐释量化金融领域的实用例子。写作本书的目的在于向读者介绍如何利用R来学习量化金融。本书涉及从基础到高级的众多话题。具体而言,我们会介绍统计学、时间序列分析、小波分析以及它们在算法交易中的应用。在本书中,我们还会尽己所能地解释机器学习、风险管理、最优化、期权定价等领域知识的应用。
第1章“R语言入门”介绍R语言相关的基础知识,逐一说明R和R程序包的安装、数据类型、数据框、循环等知识。本章还说明如何编写和调用函数以及如何从不同格式的文档导入数据到R。本章的目的在于向读者提供关于R语言的基础性说明。
第2章“统计建模”介绍如何进行探索性分析(如常用的概率分布、相关系数、趋中性度量、异常值侦测等)以更好地理解数据。本章还介绍数据抽样和数据标准化/归一化的技术,这有助于我们更好地为分析步骤准备数据。在本章的最后,我们还说明如何进行假设检验和参数估计。
第3章“计量分析与小波分析”解释简单线性回归模型和多元线性回归模型,它们是构成很多统计分析技术的支柱。本章还介绍ANOVA和特征选择的相关知识。最后,我们利用小波分析构建若干模型。
在第4章“时间序列建模”中,我们借助例子展示如何使用ts()
、zoo()
和xts()
等函数将数据转换为时间序列数据,这是时间序列预测建模的基础性工作。接着,我们讨论不同的预测技术,如AR、ARIMA、GARCH和VGARCH等,并使用例子说明它们在R中如何实现。
第5章“算法交易”包含多个算法交易领域的生动示例,涵盖动量交易以及使用多种方法的配对交易等。本章还涉及CAPM、多因子模型、投资组合构建等方面的知识。
第6章“基于机器学习的交易”展示如何使用资本市场的数据构建机器学习算法模型,具体包括监督型学习算法和非监督型学习算法。
第7章“风险管理”讨论度量市场风险和投资组合风险的技术,并解释计算VAR的常见方法。本章还给出银行业内度量信用风险的最佳实践示例。
第8章“最优化”展示金融领域应用最优化技术的不同例子,如动态再平衡、前行测试、网格测试、基因算法等。
第9章“衍生品定价”介绍如何用R进行衍生品定价。本章涵盖普通期权定价、奇异期权定价、债券定价、信用利差和信用违约互换等知识。本章相对复杂,要求读者具备关于衍生品的基础知识。
首先,确保电脑上已经安装好R软件。书中的例子全部都用R来实现,相关代码可在R编辑器或控制台上直接执行。R是开源软件,可免费从官网下载并安装到任何类型的操作系统中。在该网站上可找到安装指南。一旦安装好R软件后,你就可以从第1章开始学习。每个章节都会说明本章所需的R程序包以及如何把它们安装并加载到R的工作环境中。
写作本书的目的在于向目标读者对学习R语言及其应用感兴趣的读者介绍知识。书中所用例子来自金融领域。这些例子有的比较基础,有的相对复杂,所用到R代码的复杂性也有所不同。本书并不要求读者已经掌握R语言编程的知识,但读者要具备一些数学分析的概念。即使读者已经能够熟练使用R,本书中用到的众多来自数据分析领域(特别是资本市场)的生动例子对这些读者肯定也会大有裨益。
在本书中,你会发现多种用来区分不同类型信息的文本格式。以下是一些文本格式的示例以及相应说明。
正文中的代码、数据表名、文件夹名、文件名、文件扩展名、路径、URL、用户输入、推特标签等使用如下格式:“我们会多次用到quantmod程序包。”
代码块示例如下:
getSymbols("^DJI", src = "yahoo")
dji <- DJI[, "DJI.Close"]
当我们需要你关注代码块的某个特定部分时,相应的代码行会加粗显示:
corr<-rollapply (data, 252, correlation, by.column = FALSE)
为了方便读者直接复制并执行相应代码,本书中输出的代码行前没有R控制台的命令提示符>,长表达式的换行也不注明+,而执行代码的输出结果前有两个注释符号,即##。
新术语和关键词以粗体表示,例如以读者在屏幕上、菜单栏或对话框中看到的关键词在文中格式如下:“点击Next按钮进入下一界面。”
本书由异步社区出品,社区(https://www.epubit.com/)为您提供相关资源和后续服务。
本书提供配套源代码,要获得该配套资源,请在异步社区本书页面中单击,跳转到下载界面,按提示进行操作即可。注意:为保证购书读者的权益,该操作会给出相关提示,要求输入提取码进行验证。
如果您是教师,希望获得教学配套资源,请在社区本书页面中直接联系本书的责任编辑。
作者和编辑尽最大努力来确保书中内容的准确性,但难免会存在疏漏。欢迎您将发现的问题反馈给我们,帮助我们提升图书的质量。
当您发现错误时,请登录异步社区,按书名搜索,进入本书页面,单击“提交勘误”,输入勘误信息,单击“提交”按钮即可。本书的作者和编辑会对您提交的勘误进行审核,确认并接受后,您将获赠异步社区的100积分。积分可用于在异步社区兑换优惠券、样书或奖品。
我们的联系邮箱是contact@epubit.com.cn。
如果您对本书有任何疑问或建议,请您发邮件给我们,并请在邮件标题中注明本书书名,以便我们更高效地做出反馈。
如果您有兴趣出版图书、录制教学视频,或者参与图书翻译、技术审校等工作,可以发邮件给我们;有意出版图书的作者也可以到异步社区在线提交投稿(直接访问www.epubit.com/selfpublish/submission即可)。
如果您是学校、培训机构或企业,想批量购买本书或异步社区出版的其他图书,也可以发邮件给我们。
如果您在网上发现有针对异步社区出品图书的各种形式的盗版行为,包括对图书全部或部分内容的非授权传播,请您将怀疑有侵权行为的链接发邮件给我们。您的这一举动是对作者权益的保护,也是我们持续为您提供有价值的内容的动力之源。
“异步社区”是人民邮电出版社旗下IT专业图书社区,致力于出版精品IT技术图书和相关学习产品,为作译者提供优质出版服务。异步社区创办于2015年8月,提供大量精品IT技术图书和电子书,以及高品质技术文章和视频课程。更多详情请访问异步社区官网https://www.epubit.com。
“异步图书”是由异步社区编辑团队策划出版的精品IT专业图书的品牌,依托于人民邮电出版社近30年的计算机图书出版积累和专业编辑团队,相关图书在封面上印有异步图书的LOGO。异步图书的出版领域包括软件开发、大数据、AI、测试、前端、网络技术等。
异步社区
微信服务号
在本章,我们将讨论R语言的基础概念,这是学习本书其他章节的基础,但我们并不会深入讨论关于R的每个概念。本章的目标读者是那些并不了解R语言、想要在量化金融领域寻求职业机会或者应用R语言进行量化金融分析的初学者。本章会教你开始用R编写程序的基础知识,但假如你想要编写更加复杂的程序,则需要查看其他相关书籍。
本章包括以下内容。
许多统计软件都可以用来解决量化金融问题,但R本质上不是一个统计软件包,而是一种灵活且强大的语言,我们可用R来完成高质量的分析工作。
想要用R,人们无需成为专业编程人员或计算机专家。了解基础的编程知识对学习R当然有所助益,但这并不是先决条件。
R的绘图功能少有匹敌,这是R语言的一个突出的优点。
R语言的另一个优点在于它有着众多的扩展程序包。对于任何一个统计概念,在R中很有可能已经存在对应的程序包。R有许多专门用来进行统计和量化金融分析的工具。
R具备扩展性,提供很多功能来帮助量化金融领域的开发者编写解决分析性问题的工具或方法。当学术界有新的研究发表时,可能就会出现新的相关程序包。R的开发社区非常活跃且易于接近,人们热衷于为各种新概念开发出新的程序包。这都使得R和量化金融领域出现的前沿概念能保持同步。
R生来就是为了处理数据,但在它诞生之时,并不存在所谓的大数据。处理大数据带来的额外挑战包括数据多样化(如文本数据、测量数据等)、数据安全、存储器、CPU I/O RSC要求、多机处理等。R中用来应对大数据挑战的相关技术包括映射-归约(map-reducing)、内存中处理、流数据处理、降频采样(down sampling)、分块处理(chunking)等。
总之,R是免费软件,有着优秀的数据处理和制图功能,在网上存在海量的关于R程序包的帮助文档或技术文档。因此,R是高效且容易学习的工具。对于那些想要在量化金融领域寻求职业机会的人们来说,学习R是时代的需求。
在本节中,我们将讨论如何为不同的操作系统(Windows、Linux和Mac)下载并安装R软件。
打开网页浏览器,登录R官网,你可以基于不同的操作系统下载相应版本的R。
对于Windows版本,依次点击Download R for Windows链接、base链接和Download R 3.#.# for Windows链接,从而给你的Windows操作系统下载相应版本的R。双击打开下载的安装包,选择你偏好的语言选项,接下来的安装步骤依次如下。
1.安装向导。
2.许可协议。
3.选定你打算安装到的目标文件夹。
4.选择安装组件。基于你的系统设置相应选项,如果你不知道系统设置,则选择全部选项。
5.如果你想要自定义安装,则点击选项(option)按钮。
6.根据你的需要选择R的运行选项和桌面快捷方式选项。
至此完成Windows环境下R软件的下载和安装。
类似地,你也可以在Linux和Mac操作环境下下载相应的安装程序,然后循序完成不同的安装选项设置,直至完成安装。
程序包(package)由R函数、编译代码和示例数据集等组成,它们的存储目录被称作程序库(library)。R会默认安装一组程序包,如果你要使用其他程序包则必须自行添加。
我们使用如下代码来检测系统中存在哪些R程序包:
.libPaths()
用上面的命令来获得或设定R程序库的路径,它会输出类似如下的结果:
## "C:/Program Files/R/R-3.3.1/library"
执行如下代码会列出所有可用的程序包:
library()
存在两种安装新程序包的方法:直接从CRAN安装和手动安装。
CRAN是“综合性R存档网络(Comprehensive R Archive Network)”的首字母简称,它是由FTP服务器构成的全球网络,而在这些FTP服务器中存储着完全相同的各种版本的R代码和文档。
下面的代码用来直接从CRAN下载并安装R程序包(用户需要选定合适的镜像):
install.packages("包名")
举例说明,如果你要安装ggplot2
包和forecast
包,则代码如下:
install.packages("ggplot2")
install.packages("forecast")
手动下载所需的程序包,然后用install.packages()
将程序压缩包解压存入系统中的指定位置(如"/DATA/RPACKAGES/")。譬如,我们想要安装ggplot2
包,则使用如下代码可完成安装并将其加载到R的当前工作环境中(其他程序包也可用类似方法安装):
install.packages("ggplot2", lib = "/data/Rpackages/")
library(ggplot2, lib.loc = "/data/Rpackages/")
在编程语言中,人们需要用不同的变量来存储不同的信息。变量(variable)是用来存储值的内存位置。创建变量就是在内存中预订存储空间。用户可能需要存储不同类型的数据,如字符、浮点数、布尔值等。操作系统根据数据类型分配内存并决定在预订的内存位置存入哪些数据。
你在R中碰到的一切都被称作对象(object)。
R有5种基本的数据对象类型,又称原子对象(atomic object),而其他数据对象都基于原子对象。接下来我们将会给出关于基本对象的一些示例并验证它的类(class)。
我们将字符取值赋予一个变量并验证它的类,代码及结果如下:
a <- "hello"
class(a)
## [1] "character"
我们将数值赋予一个变量并验证它的类,代码及结果如下:
b <- 2.5
class(b)
## [1] "numeric"
我们将整数赋予一个变量并验证它的类,代码及结果如下:
c <- 6L
class(c)
## [1] "integer"
我们将复数赋予一个变量并验证它的类,代码及结果如下:
d <- 1 + 2i
class(d)
## [1] "complex"
我们将逻辑取值赋予一个变量并验证它的类,代码及结果如下:
e <- TRUE
class(e)
## [1] "logical"
R中基本的数据类型被称作向量(vector),它由相同类型的对象组成。向量不能同时包含两种不同类型的对象,比如同时包含字符对象和数值对象。但列表(list)是个例外,它可以同时包含多种类型的对象,比如同时包含字符对象、数值对象和列表对象。
接下来我们讨论R中常见的数据类型,并给出相应示例。
前面已经给出向量的定义。如果我们想构建一个包含多个元素的向量,我们可用c()
函数将元素合并成一个向量,如:
a <- "Quantitative"
b <- "Finance"
c(a, b)
## [1] "Quantitative" "Finance"
类似地,
Var <- c(1, 2, 3)
Var
## [1] 1 2 3
列表是一种包含多种类型对象(如向量,甚至列表)的R对象。让我们用如下代码构建一个列表并将其列印出来:
List1 <- list(c(4, 5, 6), "Hello", 24.5)
List1
## [[1]]
## [1] 4 5 6
##
## [[2]]
## [1] "Hello"
##
## [[3]]
## [1] 24.5
我们可根据需要提取列表中的元素,比如提取上述列表List1
中的第二个元素:
List1[2]
## [[1]]
## [1] "Hello"
我们也可用c()
函数将两个列表合并起来,如:
list1 <- list(5, 6, 7)
list2 <- list("a", "b", "c")
Combined_list <- c(list1, list2)
Combined_list
## [[1]]
## [1] 5
##
## [[2]]
## [1] 6
##
## [[3]]
## [1] 7
##
## [[4]]
## [1] "a"
##
## [[5]]
## [1] "b"
##
## [[6]]
## [1] "c"
矩阵(matrix)是具有两个维度的矩形数据集,它可通过将向量输入到matrix()
函数中得到。使用如下命令创建一个2行3列的矩阵并将其打印出来:
M <- matrix(c(1, 2, 3, 4, 5, 6), nrow = 2, ncol = 3)
M
## [,1] [,2] [,3]
## [1,] 1 3 5
## [2,] 2 4 6
矩阵只有两个维度,数组(array)则可有着任意的维度。array()
函数有个dim
参数,设定所需的维度。使用如下命令创建一个数组并将其打印出来:
a <- array(1:18, dim = c(3, 3, 2))
a
## , , 1
##
## [,1] [,2] [,3]
## [1,] 1 4 7
## [2,] 2 5 8
## [3,] 3 6 9
##
## , , 2
##
## [,1] [,2] [,3]
## [1,] 10 13 16
## [2,] 11 14 17
## [3,] 12 15 18
因子(factor)是由向量构建的R对象,它将向量及其不同元素取值所对应的标签一起存储起来。不管向量是否为数值型、字符型,还是布尔型,标签通常为字符形式。
因子用factor()
函数创建,nlevels()
函数可得到因子水平数,如:
a <- c(2, 3, 4, 2, 3)
fact <- factor(a)
fact
## [1] 2 3 4 2 3
## Levels: 2 3 4
nlevels(fact)
## [1] 3
数据框(dataframe)是表格形式的数据对象,由相同长度的向量组成,但每个向量可存储不同类型的数据(如数值型、字符型、布尔型等)。
数据框可用data.frame()
函数生成,如:
data <- data.frame(
Name = c("Alex", "John", "Bob"),
Age = c(18, 20, 23),
Gender = c("M", "M", "M")
)
data
## Name Age Gender
## 1 Alex 18 M
## 2 John 20 M
## 3 Bob 23 M
在R中,我们可以读入存储在R环境之外的文档,也可以将数据输出到可由操作系统直接读取的文档中。R可读取和存储不同格式的文档(如CSV、Excel、TXT等)。在本节中,我们将讨论如何读取和存储不同格式的文档。
读入或存储文档的第一步就是明确工作目录,代码如下:
getwd()
这会给出你所运行R的当前工作目录。如果这不是你想要的目标目录,则可用setwd()
函数进行设定。比如,如下代码会将目录"C:/Users"设定为R的工作目录:
setwd("C:/Users")
CSV格式文档是用逗号分隔数据取值的纯文本文档。让我们考虑包含如下股票市场数据内容的CSV文档:
Date,Open,High,Low,Close,Volume,Adj Close
14-10-2016,2139.679932,2149.189941,2132.97998,2132.97998,3228150000,2132.97998
13-10-2016,2130.26001,2138.189941,2114.719971,2132.550049,3580450000,2132.550049
12-10-2016,2137.669922,2145.360107,2132.77002,2139.179932,2977100000,2139.179932
11-10-2016,2161.350098,2161.560059,2128.840088,2136.72998,3438270000,2136.72998
10-10-2016,2160.389893,2169.600098,2160.389893,2163.659912,2916550000,2163.659912
CSV文档名为ch01_sample.csv
,存储在工作目录的子文件夹data
中。在此我们使用相对路径读入上述文档,代码如下:
data <- read.csv("data/ch01_sample.csv")
data
## Date Open High Low Close Volume Adj.Close
## 1 14-10-2016 2140 2149 2133 2133 3.228e+09 2133
## 2 13-10-2016 2130 2138 2115 2133 3.580e+09 2133
## 3 12-10-2016 2138 2145 2133 2139 2.977e+09 2139
## 4 11-10-2016 2161 2162 2129 2137 3.438e+09 2137
## 5 10-10-2016 2160 2170 2160 2164 2.917e+09 2164
read.csv()
函数默认生成数据框,这可通过如下代码进行确认:
is.data.frame(data)
## [1] TRUE
接下来,我们在R中将相应的数据分析函数应用于这个数据框上。一旦完成分析,用户也可用write.csv()
函数将数据输出到目标文件夹的文档中,代码如下(输出结果从略):
write.csv(data, "data/ch01_result.csv", row.names = FALSE)
# 重新读入并列印数据,看看数据输出是否正确
output <- read.csv("data/ch01_result.csv")
output
Excel文档是最常见的存储数据的文档,扩展名为.xls或.xlsx。xlsx
包可用来在R环境下读取或写入.xlsx文档。安装xlsx
包依赖于Java,所以系统中要安装Java。xlsx
包可用如下命令进行安装:
install.packages("xlsx")
当运行上述命令时,R会询问打算使用的CRAN镜像,用户应选定目标镜像。我们可通过如下命令来验证xlsx
包是否安装成功:
any(grepl("xlsx", installed.packages()))
## [1] TRUE
如果显示结果为[1] TRUE
,则说明该包成功安装。
接着,我们运行如下代码载入xlsx
包:
library("xlsx")
现在让我们在Excel软件中将之前的ch01_sample.csv
文档存储为xlsx格式的文档ch01_sample.xlsx
,然后在R环境中用如下代码读入数据:
data <- read.xlsx("data/ch01_sample.xlsx", sheetIndex = 1)
data
## Date Open High Low Close Volume Adj.Close
## 1 14-10-2016 2140 2149 2133 2133 3.228e+09 2133
## 2 13-10-2016 2130 2138 2115 2133 3.580e+09 2133
## 3 12-10-2016 2138 2145 2133 2139 2.977e+09 2139
## 4 11-10-2016 2161 2162 2129 2137 3.438e+09 2137
## 5 10-10-2016 2160 2170 2160 2164 2.917e+09 2164
类似地,你也可用xlsx
包中的write.xlsx()
将R数据对象直接存储为.xlsx格式的文档(输出结果从略):
write.xlsx(data, "data/ch01_result.xlsx", row.names = FALSE)
# 重新读入并列印数据,看看数据输出是否正确
output <- read.xlsx("data/ch01_result.xlsx", sheetIndex = 1)
output
如今网络是数据的一个主要来源,我们可将网络上的数据直接读入到R环境中,例如:
# URL <- "http://ichart.finance.yahoo.com/table.csv?s=^GSPC"
# snp <- as.data.frame(read.csv(URL))
# head(snp)
# 译者注:书中给出的类似URL已无法使用,
# 作为替代,我们使用quantmod包从Yahoo下载数据
# install.packages("quantmod")
library(quantmod)
snp <- getSymbols(Symbols = "^GSPC", src = "yahoo",
from = "2015-10-15", to = "2016-10-15",
auto.assign = FALSE)
snp <- as.data.frame(snp)
colnames(snp) <- c("Open", "High", "Low", "Close", "Volume", "Adj.Close")
snp$Date <- rownames(snp)
snp <- snp[order(snp$Date, decreasing = TRUE), ]
snp <- snp[c("Date", "Open", "High", "Low", "Close", "Volume", "Adj.Close")]
rownames(snp) <- NULL
head(snp)
## Date Open High Low Close Volume Adj.Close
## 1 2016-10-14 2140 2149 2133 2133 3.228e+09 2133
## 2 2016-10-13 2130 2138 2115 2133 3.580e+09 2133
## 3 2016-10-12 2138 2145 2133 2139 2.977e+09 2139
## 4 2016-10-11 2161 2162 2129 2137 3.438e+09 2137
## 5 2016-10-10 2160 2170 2160 2164 2.917e+09 2164
## 6 2016-10-07 2164 2166 2145 2154 3.620e+09 2154
执行上述代码会从Yahoo Finance下载标普500指数("^GSPC
")的数据,处理为数据框格式,并用head()
函数列示出部分数据。类似地,我们也可通过执行下面的代码下载道琼斯工业指数("^DJI"
)的数据到R环境中:
dji <- getSymbols(
Symbols = "^DJI", src = "yahoo",
from = "2015-10-15", to = "2016-10-15",
auto.assign = FALSE
)
dji <- as.data.frame(dji)
colnames(dji) <- c("Open", "High", "Low", "Close", "Volume", "Adj.Close")
dji$Date <- rownames(dji)
dji <- dji[order(dji$Date, decreasing = TRUE), ]
dji <- dji[c("Date", "Open", "High", "Low", "Close", "Volume", "Adj.Close")]
rownames(dji) <- NULL
head(dji)
## Date Open High Low Close Volume Adj.Close
## 1 2016-10-14 18177 18261 18138 18138 87050000 18138
## 2 2016-10-13 18088 18138 17960 18099 83160000 18099
## 3 2016-10-12 18133 18194 18082 18144 72230000 18144
## 4 2016-10-11 18308 18312 18062 18129 88610000 18129
## 5 2016-10-10 18283 18400 18283 18329 72110000 18329
## 6 2016-10-07 18295 18320 18149 18240 82680000 18240
请注意在本书后面的示例中我们会经常用到snp
和dji
这两个关于股票市场指数的数据集。
关系型数据库以标准化的格式存储数据,为了进行统计分析,我们需要编写复杂的高级查询指令。但R能连接各种数据库(如MySQL、Oracle、SQL Server等),很容易提取数据表并将其转为数据框格式。一旦数据变为数据框格式,就可以利用R中内置函数和扩展R包对数据进行统计分析。
在本小节中,我们以MySQL数据库为例进行说明。R现有的程序包RMySQL
提供和MySQL数据库的连接。通过如下命令(安装并)加载R数据库接口程序包DBI
包和相应的RMySQL
包:
# install.packages(c("DBI", "RMySQL"))
library(DBI)
library(RMySQL)
安装并载入相应的程序包后,我们就可以创建一个和数据库建立连接的对象。它以用户名、密码、数据库名、本地主机名作为输入,示例如下:
mysqlconnection <- dbConnect( MySQL(), user = "...", password = "...",
dbname = "...", host = "..." )
当创建数据库连接之后,我们可执行下面的命令将该数据库中的表列示出来:
dbListTables(mysqlconnection)
我们可用dbSendQuery()
函数对数据库进行查询,查询结果用fetch()
函数返回给R并存储为数据框:
result <- dbSendQuery(mysqlconnection, "select * from <table name>")
result_df <- fetch(result)
result_df
我们还可以通过dbSendQuery()
函数给数据库发送指令,完成诸如过滤、更新行、插入数据、创建表、删除表等各种操作。
在本节中,我们先讨论如何编写基础的表达式,这是编写程序的核心要素。接着,我们还会讨论如何编写自定义函数。
R代码由一个或多个表达式(expression)组成,表达式是指执行特定任务的语句。比如,下面的表达式计算两个数之和:
4 + 5
## [1] 9
如果程序中存在多个表达式,它们会按照出现的顺序逐一执行。
现在,我们讨论几种基本的表达式。
最简单的表达式为常数表达式,其取值可为字符型或数值型。比如,100
是数值型常数表达式,"Hello World"是字符型常数表达式。
R语言支持标准的算术运算符,其中部分算术运算符如表1-1所示。
表1-1 R支持的部分算术运算符
运算符 |
运算 |
---|---|
+ |
加法 |
- |
减法 |
* |
乘法 |
/ |
除法 |
^ |
乘方 |
使用上述算术运算符,我们可生成算术表达式,示例如下(输出结果从略):
4 + 5
4 - 5
4 * 5
R遵循BODMAS规则,即“先计算括号内的以及先乘除后加减”的运算优先规则。在创建算术表达式时,人们可使用括号来避免歧义。
条件表达式对两个取值进行比较并返回逻辑值(TRUE
或FALSE
),R支持标准的比较运算符和条件连接运算符(表1-2)。
表1-2 R支持的比较运算符和条件连接运算符
运算符 |
运算 |
---|---|
== |
等于 |
> (>=) |
大于(大于等于) |
< (<=) |
小于(小于等于) |
!= |
不等于 |
&& |
逻辑和 |
|| |
逻辑或 |
! |
逻辑否 |
比如:10 > 5
,执行后返回TRUE
;5 > 10
,执行后返回FALSE
。
R中最常见也最有用的表达式是函数调用表达式。R中存在很多内置函数,用户也可以自定义函数。在本小节中,我们说明调用函数的基本结构。
函数调用表达式由函数名、括号以及括号内用逗号分隔的参数构成。参数是给函数执行特定任务时提供必要信息的表达式。在下文讨论如何构建自定义函数时,我们会给出一个函数调用表达式的例子。
R代码中包含关键词(keyword)和符号(symbol)。
符号是内存中存储对象的标签。当执行程序时,程序会从内存中取得存储对象的取值。R有很多已经预先定义好取值的符号,可在程序中直接调用。比如,letters
是会返回26个小写英文字母的向量。
通过赋值运算符<-
可将表达式的结果赋值给某个符号。比如,value <- 4 + 6
将加法运算的结果10赋值给符号value
并存储在内存中。
有些符号用来表示某些特殊的取值,这些符号无法被重新赋值,如下所示。
NA
:用来表示缺失值。Inf
:用来表示无穷。比如,1/0
就会返回Inf
。NaN
:用来表示无法定义的算术表达式的结果。比如,0/0
就会返回NaN
。NULL
:用来表示空值。TRUE
和FALSE
:这是布尔逻辑值,通常是比较变量取值的返回值。当编写R代码时,我们需要用很多符号来存储不同的信息。为了让代码更容易理解,我们要尽量给这些符号取有意义的名字。符号命名应该是不言自明的,太短的符号名会使代码更难以理解。
比如,我们可用DateOfBirth
或DOB
来存储生日的信息存储,但显然前者更好,让人更容易理解。
在本节中,我们将提供关于R中内置函数的若干例子,并编写一个完成特定任务的自定义函数。
函数(function)是为执行特定任务而放在一起的一组语句。R中内置有很多函数,用户也可以自定义函数。在R中,解释器根据要求将控制以及为完成任务而设置的函数参数交给函数对象。在完成任务之后,函数会将控制交还给解释器。
定义函数的语法如下:
function_name <- function(arg1, arg2, ...) {
function_body
}
具体解释如下所示。
这里提供两个内置函数及其返回值的示例:
mean(25:82)
## [1] 53.5
sum(41:68)
## [1] 1526
现在,让我们看看如何自定义函数。我们想要计算并打印出给定整数的平方序列,函数名为findingSqrFunc
,参数为value
:
findingSqrFunc <- function(value) {
for (i in 1:value) {
sqr <- i^2
print(sqr)
}
}
在执行上述代码之后,我们就可以调用该函数:
findingSqrFunc(4)
## [1] 1
## [1] 4
## [1] 9
## [1] 16
定义并调用一个无需参数的函数示例如下:
func_wo_args <- function() {
for (i in 1:3) {
print(i * 5)
}
}
func_wo_args()
## [1] 5
## [1] 10
## [1] 15
我们可按照定义函数时参数的顺序来给函数传递参数取值。但如果使用参数名,就可以不按顺序来传递参数的取值。下面给出定义并调用函数的示例:
func_w_args <- function(a, b, c) {
result <- a * b + c
print(result)
}
func_w_args(2, 3, 4)
## [1] 10
func_w_args(c = 4, b = 3, a = 4)
## [1] 16
在本节中,我们讨论运行R程序的不同方式。
遵循如下步骤在R代码窗口中运行程序。
1.打开R软件(双击桌面图标或从开始程序菜单打开)。
2.点击File菜单下的打开脚本(Open Script),然后定位到脚本文档所在的目录。
3.选定并打开你要运行的脚本文档,它会在R编辑器窗口中打开。
4.右键菜单->全选(Select All)(或使用快捷键Ctrl + A)。
5.右键菜单->运行代码行或选定代码(Run LineorSelection)(或使用快捷键Ctrl + R)。
6.程序的运行结果将出现在R控制台窗口中。
如果你想让R无须等待指令而直接执行多行代码,你可以使用source()
函数来加载整个脚本文档,步骤如下。
1.在R编辑器窗口中编写全部代码,然后存盘为目标文件夹下的脚本文档(如"D:/Rcode/first_prog.R"
)。代码示例如下:
a <- 5
print(a)
2.使用绝对路径加载R脚本文档,代码示例如下:
source("D:/Rcode/first_prog.R")
当然,我们也可先用setwd()
函数修改R工作目录再直接加载R脚本文档。加载脚本文档first_prog.R
会在R控制台窗口中输出5。
注释是程序的一部分,但在执行代码时解释器会忽略注释语句。注释用#
加以标注,如:
# 这是程序中的一句注释
循环(loop)是将需要重复的操作分组并安排执行顺序从而自动执行多步过程的指令。所有的编程语言都有内置的方式来实现重复执行指令或指令块。条件决策也是编程语言的重要构件,在R语言中这通过使用条件语句(conditional statement)if ... else
来实现。
本节先讨论if ... else
条件语句,接着再讨论循环。
让我们看看if
语句在R中是如何工作的。if
语句的通用语法如下:
if (表达式) {
语句
}
如果表达式成立,就会执行语句,否则就不执行语句。
表达式可以是逻辑或数值向量。当表达式为数值向量时,0被视作FALSE
,其他取值被视为TRUE
。
x <- 5
if (x > 0)
{
print("I am Positive")
}
执行上述代码,R将输出"I am Positive"
。
现在让我们看看if
和else
条件语句在R中是如何工作的。语法如下:
if (表达式) {
语句1
} else {
语句2
}
当if ()
内的表达式为FALSE
时,else
部分的语句2就会被执行,如:
x <- -5
if (x > 0)
{
print("I am Positive")
} else
{
print("I am Negative")
}
执行上述代码,R将输出"I am Negative"
。
通过计数器或索引,for
循环会重复执行预定的次数,每次循环后计数器或索引会自动增加。for
循环的语法如下:
for (val in sequence) {
语句
}
下面是个示例:
vals <- c(3, 6, 8, 9, 11, 16)
counter <- 0
for (val in vals) {
if (val %% 2 != 0) counter <- counter + 1
}
print(counter)
上述代码会计算向量vals
中奇数的个数,输出结果为3。
while
循环在每次循环开始前都会检验逻辑条件是否成立。while
循环的语法如下:
while (表达式) {
语句
}
在每次循环中,R都会对表达式进行求值,如果结果为TRUE
,则重复执行循环体中的语句。也就是说,while
循环会重复执行直至表达式的取值为FALSE
。示例如下:
helo <- "Hello"
counter <- 4
while (counter < 7) {
print(helo)
counter <- counter + 1
}
## [1] "Hello"
## [1] "Hello"
## [1] "Hello"
apply()
函数会对矩阵、向量、数组等数据对象的每行、每列或每个元素依次执行某项操作。例如我们可使用apply()
函数来计算矩阵每一行中各元素的和,代码及输出结果如下:
sample <- matrix(c(1:10), nrow = 5 , ncol = 2)
apply(sample, 1, sum)
## [1] 7 9 11 13 15
sapply()
函数会对数据集(如列表或向量)中的每个元素调用指定函数,示例代码如下:
sapply(1:5, function(x) x^3)
## [1] 1 8 27 64 125
上述代码计算得到自然数序列1到5的三次方。
循环控制语句(break
和next
)会变更循环执行的正常顺序,在此我们加以简要讨论。
break停止当前循环并将控制交给循环之后的下一个语句,例如:
helo <-c ("Hello")
counter <- 5
repeat {
print(helo)
counter <- counter + 1
if(counter > 8) {
break
}
}
由于break
语句的存在,上述代码只会打印"Hello"
4次,然后就会退出循环。在此,repeat
是另一类循环语句,它会重复执行直至循环体内触发停止循环的条件。
next
并不终止循环,而是跳过当前循环,直接进入下一个循环。示例如下:
vec <- c(2, 3, 4, 5, 6)
for (i in vec) {
if (i == 4) {
next
}
print(i)
}
## [1] 2
## [1] 3
## [1] 5
## [1] 6
在上述例子中,当循环至向量vec
的第3个元素4时,控制语句next
会跳过当前循环,进入下一个循环。因此,执行上述代码只会打印出向量vec
的元素2、3、5和6,而跳过元素4。
1.R的原子对象包括哪些?
2.R的向量指的是什么?
3.向量和列表的区别是什么?
4.数组和矩阵的区别是什么?
5.什么是数据框?它在R中的重要性体现在哪里?
6.在R中怎样读取和存储CSV和XLSX文档?
7.在R中怎样读取和存储股票市场数据?
8.解释用R连接关系型数据库的步骤。
9.什么是函数?它在R中的重要性体现在哪里?
10.什么是R的赋值运算符?
11.在R中如何调用函数?
12.在R中如何加载R脚本文档?
13.在R中for循环和while循环有哪些不同?
在本章中我们学到的内容概括如下。
if
和else
)。for
和while
)。本章讨论统计建模。统计建模是量化金融的驱动力,因此这是在R中学习量化金融的第一步。我们假定读者熟悉R语言基础编程知识并对相关的统计概念有一定了解。在本章中,我们不会详细讨论统计概念,而是展示如何用R进行统计建模。
本章包括以下内容。
概率分布(probability distributions)确定随机变量的取值是如何分布的。比如,抛硬币得到所有可能结果的集合是二项式分布。基于总体数据的大样本均值遵循正态分布,而正态分布是最常见和有用的分布。
统计学对这些概率分布的性质有很好的理解,可用来进行关于总体的推断。在本节中,我们会讨论几种最常见的概率分布并说明用来计算概率分布相关性质的R函数。
正态分布(normal distribution)是金融行业最广泛使用的概率分布。正态分布为钟型曲线,它的均值、中位数和众数都相同。正态分布表示为,其中:μ为样本的均值,
为样本的方差。如果均值等于0,方差等于1,则该正态分布就称作标准正态分布N(0,1)。
接下来讨论用来计算正态分布相关性质的主要函数。请注意,本章中全部计算用到的数据集为data/ch02_data.csv
,用如下代码将数据集导入到R中,存为变量Sampledata
,并列示出数据集最开始的5行。
Sampledata <- read.csv("data/ch02_data.csv")
head(Sampledata, n = 5)
## Date Open High Low Close Volume Adj.Close Return Flag Sentiments
## 1 12/14/2016 198.7 203.0 196.8 198.7 4144600 198.7 0.00 1 Good
## 2 12/13/2016 193.2 201.3 193.0 198.1 6816100 198.1 0.03 1 Bad
## 3 12/12/2016 192.8 194.4 191.2 192.4 615800 192.4 0.00 1 Good
## 4 12/9/2016 190.9 193.8 190.8 192.2 2719600 192.2 0.00 0 Bad
## 5 12/8/2016 192.1 192.5 189.5 192.3 3187300 192.3 0.00 0 Good
在样本数据中,Date
是数据收集的时间,Open
、High
、Low
和Close
分别为开盘价、最高价、最低价和收盘价。Volume
为成交量,Adj.Close
为复权后收盘价,Return
为根据今日和昨日复权后收盘价计算得到的收益率。Flag
和Sentiments
则是分析中用到的两个哑变量。
dnorm()
函数返回正态分布概率密度的取值,函数定义如下:
dnorm(x, mean, sd)
其中,x
为数值向量,mean
为正态分布的均值,sd
为正态分布的标准差。
执行如下代码将会生成显示正态分布概率密度曲线的散点图(图2-1):
y <- dnorm(
Sampledata$Return,
mean = mean(Sampledata$Return),
sd = sd(Sampledata$Return, na.rm = FALSE)
)
plot(Sampledata$Return, y)
图2-1 显示正态分布曲线高度的散点图
pnorm()
函数是正态分布的累积分布函数,给出随机变量的取值小于给定值的概率,函数如下:
pnorm(q, mean, sd)
运行如下代码:
pnorm(
.02,
mean = mean(Sampledata$Return),
sd = sd(Sampledata$Return, na.rm = FALSE)
)
## [1] 0.8402
代码给出的结果意味着日收益率Return
大约只有16%的概率会大于2%。
qnorm()
函数以给定概率值作为输入,返回累积概率等于给定概率值的随机数分位点的取值,函数定义如下:
qnorm(p, mean, sd)
其中,p
为给定概率值。
运行如下代码:
qnorm(
0.159837,
mean = mean(Sampledata$Return),
sd = sd(Sampledata$Return, na.rm = FALSE),
lower.tail = FALSE
)
## [1] 0.02
代码给出的结果意味着日收益率Return
大于等于2%的概率约为16%。
rnorm()
函数用来生成正态分布的随机变量,函数定义如下:
qnorm(n, mean, sd)
其中,n
为拟生成的随机数的个数。
运行如下代码将会生成5个满足正态分布的随机数,正态分布的均值和标准差分别等于样本数据中Return
变量的均值和标准差:
rnorm(
5,
mean = mean(Sampledata$Return),
sd = sd(Sampledata$Return, na.rm = FALSE)
)
## [1] -0.035468 0.013173 -0.014893 0.003836 -0.012155
在金融时间序列中,对数正态分布(lognormal distribution)往往发挥着比正态分布更为关键的作用。和2.1.1节“正态分布”一样,我们也会讨论计算对数正态分布相关性质的函数。
dlnorm()
函数用来计算对数正态分布的密度函数,其通用语法如下:
dlnorm(x, meanlog, sdlog)
我们使用如下代码计算样本数据交易量变量Volume
的密度函数图(图2-2):
y <- dlnorm(
Sampledata$Volume,
meanlog = mean(Sampledata$Volume),
sdlog = sd(Sampledata$Volume, na.rm = FALSE)
)
plot(Sampledata$Volume, y)
图2-2 对数正态分布密度函数的散点图
plnorm()
函数给出对数正态分布的累积分布函数,函数语法如下:
plnorm(q, meanlog, sdlog)
我们用如下代码计算出样本数据交易量变量Volume
所对应的累积分布函数(cdf)及其图示(图2-3):
y <- plnorm(
Sampledata$Volume,
meanlog = mean(Sampledata$Volume),
sdlog = sd(Sampledata$Volume, na.r = FALSE)
)
plot(Sampledata$Volume, y)
图2-3 对数正态累积分布函数的散点图
qlnorm()
函数用来计算对数正态分布的p
分位数,函数语法如下:
qlnorm(p, meanlog, sdlog)
rlnorm()
函数用来生成满足给定均值和标准差的对数正态分布的随机数,函数语法如下:
rlnorm((n, meanlog, sdlog)
泊松分布(Poisson distribution)是单位间隔时间内独立事件发生次数的概率分布。如果在单位间隔时间内事件发生的平均次数为λ,则单位间隔时间内发生x次事件的概率为:
其中,。
如果每分钟内平均有10只股票的收益率为正,那么我们可用如下代码计算出某特定分钟内收益率为正股票只数小于等于15只的概率(约等于95%):
ppois(15, lambda = 10)
## [1] 0.9513
类似地,我们也用如下代码算出超出的概率:
ppois(15, lambda = 10, lower = FALSE)
## [1] 0.04874
连续均匀分布(continuous uniform distribution)是从连续区间[a, b]中选出随机变量的概率分布,其密度函数如下:
其中,且当
或
时
。
让我们用如下代码生成10个在1到5之间均匀分布的随机数:
runif(10, min = 1, max = 5)
## [1] 3.401 2.883 3.635 4.535 4.443 1.979 4.760 2.943 4.303 3.273
多数常用的统计分布更加关注分布的中间区域,而忽视分布的尾部区域(包含着极值或异常值)。风险管理者面临的最大挑战就是开发出考虑到罕见和极端事件的风险模型。极值理论(Extreme value theory,EVT)试图对分布的尾部区域提供最佳的估计。
存在两类估计极值的模型,即用广义极值(Generalized Extreme Value, GEV)分布拟合的块极值模型(Block Maxima Models)和用广义帕累托分布(Generalized Pareto Distribution,GPD)拟合的超阈峰值模型(Peaks Over Threshold,POT),通常POT更为常用。本章给出一个POT例子,直接使用POT
程序包中的数据集ardieres
的子集。
为了找出尾部分布,首先需要确定阈值,代码如下:
# install.packages("POT")
library(POT)
data(ardieres)
abc <- ardieres[1:10000, ]
events <- clust(abc, u = 1.5, tim.cond = 8 / 365, clust.max = TRUE)
par(mfrow = c(2, 2))
mrlplot(events[, "obs"])
diplot(events)
tcplot(events[, "obs"], which = 1)
tcplot(events[, "obs"], which = 2)
由此可得到图2-4。
图2-4 EVT中对阈值选择的分析
根据对图2-4的分析,我们将阈值设置为7,接着就可以拟合GPD的参数,代码如下:
obs <- events[,"obs"]
ModelFit <- fitgpd(obs, thresh = 7, "pwmu")
ModelFit
## Estimator: PWMU
##
## Varying Threshold: FALSE
##
## Threshold Call: 7
## Number Above: 10
## Proportion Above: 0.077
##
## Estimates
## scale shape
## 3.671 -0.147
##
## Standard Error Type:
##
## Standard Errors
## scale shape
## 1.768 0.384
##
## Asymptotic Variance Covariance
## scale shape
## scale 3.126 -0.562
## shape -0.562 0.148
##
## Correlation
## scale shape
## scale 1.000 -0.827
## shape -0.827 1.000
##
## Optimization Information
## Convergence: NA
## Function Evaluations: NA
当构建金融模型时,我们可能面对非常大的数据集,这导致每次估计模型都要占用大量时间。在完成建模之后,对模型进行调整也会花费大量时间。在这种情况下,如果我们只用总体的一个随机样本数据进行建模就会更加容易且用时更少。在本节中,我们讨论如何对数据进行随机抽样和分层抽样,这在基于总体抽样数据的建模中起到关键作用。
随机抽样(random sampling)是指按照相同的概率抽取总体观测来构建样本,这可通过两种方式进行:不回置抽样和回置抽样。
不回置随机抽样的示例代码及其结果如下(注意:由于未设定随机数种子,每次运行结果都会不同):
RandomSample <- Sampledata[sample(1:nrow(Sampledata), 10, replace = FALSE), ]
RandomSample
## Date Open High Low Close Volume Adj.Close Return Flag Sentiments
## 10 12/1/2016 188.2 188.5 181.0 181.9 5112100 181.9 -0.04 0 Bad
## 3 12/12/2016 192.8 194.4 191.2 192.4 615800 192.4 0.00 1 Good
## 7 12/6/2016 185.5 186.6 182.7 185.9 3372000 185.9 -0.01 0 Good
## 45 10/12/2016 200.9 203.9 200.4 201.5 1970700 201.5 0.01 1 Good
## 33 10/28/2016 204.0 205.3 199.8 200.0 4255500 200.0 -0.02 0 Good
## 49 10/6/2016 202.5 204.2 200.2 201.0 4703400 201.0 -0.04 0 Good
## 1 12/14/2016 198.7 203.0 196.8 198.7 4144600 198.7 0.00 1 Good
## 17 11/21/2016 185.0 188.9 184.4 184.5 4344600 184.5 0.00 0 Good
## 39 10/20/2016 202.1 203.0 197.1 199.1 5072900 199.1 -0.02 0 Good
## 24 11/10/2016 191.1 191.6 180.4 185.4 6748500 185.4 -0.02 0 Bad
回置随机抽样的示例代码及其结果如下(注意:由于未设定随机数种子,每次运行结果都会不同):
RandomSample <- Sampledata[sample(1:nrow(Sampledata), 10, replace = TRUE), ]
RandomSample
## Date Open High Low Close Volume Adj.Close Return Flag Sentiments
## 13 11/28/2016 195.5 199.4 194.6 196.1 4487100 196.1 0.00 0 Good
## 8 12/5/2016 182.5 188.9 182.5 186.8 4066500 186.8 0.03 1 Bad
## 20 11/16/2016 182.6 184.7 181.2 183.9 3430400 183.9 0.00 1 Bad
## 34 10/27/2016 211.3 213.7 201.6 204.0 13066400 204.0 0.01 1 Bad
## 37 10/24/2016 201.0 203.9 200.2 202.8 2751600 202.8 0.01 1 Good
## 38 10/21/2016 198.6 201.6 197.4 200.1 2943400 200.1 0.00 1 Bad
## 42 10/17/2016 197.1 198.4 192.0 194.0 4554100 194.0 -0.01 0 Bad
## 6 12/7/2016 186.1 193.4 185.0 193.1 5441400 193.1 0.04 1 Bad
## 29 11/3/2016 189.0 191.5 187.0 187.4 2641400 187.4 0.00 0 Good
## 28 11/4/2016 189.0 193.5 186.0 190.6 5140900 190.6 0.02 1 Bad
在分层抽样(stratified sampling)中,我们将总体划分为不同组别(即“分层”)。然后从每个分层抽出一个随机样本。相比于简单随机抽样,分层抽样具备几个优点,如为了提高样本的代表性,用较小的样本得到更精确的估计。
首先,让我们看下根据变量Flag
和Sentiments
可将样本分为多少组别,代码及其结果如下:
table(Sampledata$Flag, Sampledata$Sentiments)
##
## Bad Good
## 0 13 14
## 1 12 11
接着,我们使用sampling
程序包中的strata()
函数从不同组别中选出随机样本(注意:由于未设定随机数种子,每次运行结果都会不同):
# install.packages("sampling")
library(sampling)
Stratsubset <- strata(
Sampledata,
c("Flag", "Sentiments"),
size = c(6, 5, 4, 3),
method = "srswor"
)
Stratsubset
## Flag Sentiments ID_unit Prob Stratum
## 15 1 Good 15 0.5455 1
## 23 1 Good 23 0.5455 1
## 27 1 Good 27 0.5455 1
## 37 1 Good 37 0.5455 1
## 45 1 Good 45 0.5455 1
## 47 1 Good 47 0.5455 1
## 2 1 Bad 2 0.4167 2
## 14 1 Bad 14 0.4167 2
## 28 1 Bad 28 0.4167 2
## 34 1 Bad 34 0.4167 2
## 50 1 Bad 50 0.4167 2
## 10 0 Bad 10 0.3077 3
## 30 0 Bad 30 0.3077 3
## 44 0 Bad 44 0.3077 3
## 48 0 Bad 48 0.3077 3
## 7 0 Good 7 0.2143 4
## 31 0 Good 31 0.2143 4
## 35 0 Good 35 0.2143 4
对于某个给定数据集,人们试图用数据的中间位置来概括数据,这就是所谓的趋中性度量(measure of central tendency)或汇总统计量(summary statistics)。存在几种趋中性的度量指标,如均值、中位数和众数。在不同的场合下,人们使用不同的趋中性度量,但以均值的使用最为广泛。接下来,我们给出R中计算趋中性度量的例子。
均值(mean)是样本的等权重平均数。比如,我们可用如下代码计算出Sampledata
数据集中Volume
变量的均值(即算术平均数):
mean(Sampledata$Volume)
中位数(median)是排序后向量中居中元素的取值,示例代码如下:
median(Sampledata$Volume)
众数(mode)是最大频率元素的取值。R中不存在计算众数的内置函数,但我们可以自行编写一个:
findmode <- function(x) {
uniqx <- unique(x)
uniqx[which.max(tabulate(match(x, uniqx)))]
}
findmode(Sampledata$Return)
执行上述代码会给出Sampledata
数据集中Return
变量的众数。
我们也可以用summary()
函数计算向量的基本统计量:
summary(Sampledata$Volume)
## Min. 1st Qu. Median Mean 3rd Qu. Max.
## 615800 3267625 4232100 4391172 5102300 13066400
这将同时给出数据的均值、中位数、最小值、最大值、第1四分位数、第3四分位数等统计量。
矩(moment)给出关于总体方差、偏度等性质的度量,示例代码如下:
# install.packages("e1071")
library(e1071)
moment(Sampledata$Volume, order = 3, center = TRUE)
通过设置moment()
函数的相关参数,上述代码计算Sampledata
数据集中Volume
变量的三阶(order = 3
)中心矩(center = TRUE
)。
峰度(kurtosis)衡量相对于正态分布来说,数据是否出现厚尾或薄尾分布(heavy-tailed or light-tailed)。如果数据集的峰度值较大,则多为厚尾分布(异常值较多);而峰度值较小,则多为薄尾分布(异常值较少)。计算得到的峰度值需要和正态分布的峰度值进行比较,并在此基础上进行解释。
我们用e1071
包中的kurtosis()
函数来计算Volume
变量峰度值:
kurtosis(Sampledata$Volume)
结果显示Volume
变量的峰度值为5.7771,这说明Volume
变量为尖峰态分布(leptokurtic)。
偏度(skewness)衡量分布的对称性。如果变量的均值小于中位数则说明分布是左偏的(left-skewed);而如果均值大于中位数则说明分布是右偏的(right-skewed)。
我们同样使用e1071
包中的skewness()
函数来计算Sampledata
数据中Volume
变量的偏度值:
skewness(Sampledata$Volume)
代码给出的结果是1.7237,这说明Volume
变量的分布是右偏的。
在量化金融领域,相关性(correlation)发挥着重要的作用,它是对两个金融变量之间线性关系的度量。相关性不仅可以用来确定不同金融变量之间的关系,而且还可用来预测金融工具价格的未来走势。
现在让我们用R来计算Sampledata
数据集不同类型价格数据之间的相关性,这可用来确定预测性金融模型所需变量的数量。我们先提取数据子集,然后再计算不同类型价格变量之间的相关系数:
x <- Sampledata[, 2:5]
cor(x, method = "pearson")
## Open High Low Close
## Open 1.0000 0.9621 0.9342 0.8786
## High 0.9621 1.0000 0.9527 0.9454
## Low 0.9342 0.9527 1.0000 0.9604
## Close 0.8786 0.9454 0.9604 1.0000
自相关性(autocorrelation)是数据序列与其自身的过去值或将来值之间的关系,它也被称作系列相关性或滞后相关性。自相关性在时间序列预测模型中发挥着重要作用。我们用acf()
函数来估计数据序列的自相关函数。
下面的代码给出数据序列与其不同滞后阶数取值之间的自相关图(图2-5):
acf(Sampledata$Volume)
图2-5 数据序列与其滞后值之间的自相关图
acf()
函数还有lag.max
、plot
等参数。
偏自相关性(partial autocorrelation)是在控制所有更低阶滞后项的影响后,时间序列的当前值与其滞后值之间的相关性,在时间序列模型中可用它来确定预测模型滞后变量的阶数。计算偏自相关性的示例代码如下(图2-6):
pacf(Sampledata$Volume)
图2-6 数据序列与其滞后值之间的偏自相关图
pacf()
函数也包含用来控制最大滞后阶数(lag.max
)和是否输出图形的参数(plot
)。
交叉相关性(cross-correlation)是对两组数据序列的相关性与位移间隔之间关系的度量。和自相关性、偏自相关性一样,交叉相关性在时间序列预测中也起到重要作用。我们用ccf()
函数计算交叉相关性,代码及输出如下(图2-7):
ccf(Sampledata$Volume, Sampledata$High, main = "ccf plot")
图2-7 两组数据系列之间的交叉相关图
假设检验(hypothesis testing)根据观测样本的检验统计量而拒绝或接受一个假设。我们不会深入分析假设检验的理论基础,而只讨论在R中如何实现各种类型的假设检验。
原假设为,其中
是总体均值的下边界。
假定投资者预期一只股票自上市以来日收益的均值大于$10,而30天样本日收益的平均值是$9.9,总体的标准差为$1.1,我们能否以0.05的显著性水平拒绝原假设吗?
让我们用如下代码计算检验统计量z
:
xbar <- 9.9 # 样本均值
mu0 <- 10 # 假设的总体均值
sig <- 1.1 # 总体的标准差
n <- 30 # 样本数
z <- (xbar - mu0) / (sig / sqrt(n)) # 检验统计量
z
## [1] -0.4979
接着,我们需要找出0.05显著性水平的临界值,计算代码如下:
alpha <- .05
z.alpha <- qnorm(1 - alpha)
-z.alpha
## [1] -1.645
样本的z检验统计量为-0.4979,大于0.05显著性水平要求的临界值-1.645,因此我们无法拒绝认为该股票日收益大于$10的原假设。
作为替代方法,我们也可以用pnorm()
函数计算检验统计量z
的p值,代码如下:
pnorm(z)
## [1] 0.3093
由于上面计算得到的p值大于0.05的显著性水平,我们同样无法拒绝原假设。
原假设为,其中的
是总体均值的上边界。
假定投资者预期一只股票自上市以来日收益的均值小于$5,而30天样本日收益的平均值是$5.1,总体的标准差为$0.25,我们能否以0.05的显著性水平拒绝原假设吗?
让我们用如下代码计算检验统计量z:
xbar <- 5.1 # 样本均值
mu0 <- 5 # 假设的总体均值
sig <- 0.25 # 总体的标准差
n <- 30 # 样本数
z <- (xbar - mu0) / (sig / sqrt(n)) # 检验统计量
z
## [1] 2.191
接着,我们需要找出0.05显著性水平的临界值,计算代码如下:
alpha <- .05
z.alpha <- qnorm(1 - alpha)
z.alpha
## [1] 1.645
0.05显著性水平所要求的临界值1.645小于样本的z检验统计量(为2.191),因此我们可以拒绝原假设。
我们也可以计算检验统计量z
的p值,代码如下:
pnorm(z, lower.tail = FALSE)
## [1] 0.01423
由于上面计算得到的p值小于0.05的显著性水平,我们同样可以拒绝原假设。
原假设为,其中的
是总体均值。
假定一只股票去年日收益均值为$2,而30天样本日收益的平均值是$1.5,总体的标准差为$0.1,我们能否以0.05的显著性水平拒绝认为两者均值不存在显著差异原假设吗?
让我们用如下代码计算检验统计量z
:
xbar <- 1.5 # 样本均值
mu0 <- 2 # 假设的总体均值
sig <- 0.1 # 总体的标准差
n <- 30 # 样本数
z <- (xbar - mu0) / (sig / sqrt(n)) # 检验统计量
z
## [1] -27.39
接着,让我们找出0.05显著性水平所对应的临界值,代码如下:
alpha <- .05
z.half.alpha <- qnorm(1 - alpha / 2)
c(-z.half.alpha, z.half.alpha)
## [1] -1.96 1.96
样本的z检验统计量(为-27.39)并未落在区间(-1.96, 1.96)之内,因此我们可以拒绝认为两者的均值不存在显著差异的原假设。
双尾p值计算如下:
2 * pnorm(z)
## [1] 4.012e-165
计算得到的p值小于0.05的显著性水平,因此我们可以拒绝原假设。
在前面3个例子中,我们都假定方差已知并使用正态分布进行假设检验。在接下来的例子中,我们不再给定方差,转而使用t分布进行假设检验。
原假设为,其中的
是总体均值的下边界。
假定投资者预期一只股票自上市以来日收益的均值大于$1,而30天样本日收益的平均值是$0.9,样本的标准差为$0.01,我们能否以0.05的显著性水平拒绝原假设吗?
在本例中,我们用如下代码计算检验统计量t
:
xbar <- 0.9 # 样本均值
mu0 <- 1 # 假设的总体均值
sig <- 0.1 # 样本的标准差
n <- 30 # 样本数
t <- (xbar - mu0) / (sig / sqrt(n)) # 检验统计量
t
## [1] -5.477
接着,让我们计算0.05显著性水平所对应的临界值,代码如下:
alpha <- .05
t.alpha <- qt(1 - alpha, df = n - 1)
-t.alpha
## [1] -1.699
样本的检验统计量为-5.477,小于0.05显著性水平要求的临界值-1.699,因此我们可以拒绝原假设。
作为替代方法,我们也可以计算检验统计量t
所对应的p值,代码如下:
pt(t, df = n - 1)
## [1] 3.37e-06
由于上面计算得到的p值小于0.05的显著性水平,我们同样得以拒绝原假设。
原假设为,其中的
是总体均值的上边界。
假定投资者预期一只股票自上市以来日收益的均值小于$3,而30天样本日收益的平均值是$3.1,样本的标准差为$0.2,我们能否以0.05的显著性水平拒绝原假设呢?
让我们用如下代码计算检验统计量t
:
xbar <- 3.1 # 样本均值
mu0 <- 3 # 假设的总体均值
sig <- 0.2 # 样本的标准差
n <- 30 # 样本数
t <- (xbar - mu0) / (sig / sqrt(n)) # 检验统计量
t
## [1] 2.739
现在让我们找出0.05显著性水平所对应检验统计量的临界值,计算代码如下:
alpha <- .05
t.alpha <- qt(1 - alpha, df = n - 1)
t.alpha
## [1] 1.699
临界值(为1.699)小于检验统计量(为2.739),因此我们拒绝原假设。
同样地,检验统计量所对应的p值为:
pt(t, df = n - 1, lower.tail = FALSE)
## [1] 0.005219
计算得到的p值小于0.05的显著性水平,同样得以拒绝原假设。
原假设为,其中的
是总体均值。
假定一只股票去年日收益均值为$2,而30天样本日收益的平均值是$1.9,样本的标准差为$0.1,我们能否以0.05的显著性水平拒绝认为两者均值不存在显著差异原假设吗?
让我们用如下代码计算检验统计量t
:
xbar <- 1.9 # 样本均值
mu0 <- 2 # 假设的总体均值
sig <- 0.1 # 样本的标准差
n <- 30 # 样本数
t <- (xbar - mu0) / (sig / sqrt(n)) # 检验统计量
t
## [1] -5.477
现在让我们找出临界区间,代码如下:
alpha <- .05
t.half.alpha <- qt(1 - alpha / 2, df = n - 1)
c(-t.half.alpha, t.half.alpha)
## [1] -2.045 2.045
t
检验统计量等于-5.477,并未落上面的临界区间内,因此我们得以拒绝原假设。
在本节中,我们将讨论几种估计参数的算法。
极大似然估计(Maximum Likelihood Estimation,MLE)是基于给定数据集估计模型参数的一种方法。
作为示例,让我们试着找出一个正态分布概率密度函数的参数估计值。
我们首先生成一组随机数系列,代码如下:
set.seed(100)
NO_values <- 100
Y <- rnorm(NO_values, mean = 5, sd = 1)
mean(Y)
## [1] 5.003
sd(Y)
## [1] 1.021
接着让我们定义一个计算对数似然率的函数:
LogL <- function(mu, sigma) {
A <- dnorm(Y, mu, sigma)
-sum(log(A))
}
接着应用stats4
程序包中的mle()
函数来估计参数均值和标准差:
library(stats4)
mle(LogL, start = list(mu = 2, sigma = 2)) # start参数设定mu和sigma的初始值
## Warning in dnorm(Y, mu, sigma): NaNs produced
## Warning in dnorm(Y, mu, sigma): NaNs produced
## Warning in dnorm(Y, mu, sigma): NaNs produced
## Warning in dnorm(Y, mu, sigma): NaNs produced
## Warning in dnorm(Y, mu, sigma): NaNs produced
## Warning in dnorm(Y, mu, sigma): NaNs produced
## Warning in dnorm(Y, mu, sigma): NaNs produced
## Warning in dnorm(Y, mu, sigma): NaNs produced
## Warning in dnorm(Y, mu, sigma): NaNs produced
## Warning in dnorm(Y, mu, sigma): NaNs produced
##
## Call:
## mle(minuslogl = LogL, start = list(mu = 2, sigma = 2))
##
## Coefficients:
## mu sigma
## 5.003 1.016
当mle()
函数对sigma
参数尝试取负值时,就会出现NaNs produced
的警告信息。这可通过给mle()
设定相应的参数来加以控制,如:
mle(
LogL,
start = list(mu = 2, sigma = 2),
method = "L-BFGS-B",
lower = c(-Inf, 0),
upper = c(Inf, Inf)
)
##
## Call:
## mle(minuslogl = LogL, start = list(mu = 2, sigma = 2), method = "L-BFGS-B",
## lower = c(-Inf, 0), upper = c(Inf, Inf))
##
## Coefficients:
## mu sigma
## 5.003 1.016
在线性回归模型(linear regression model)中,我们试图利用自变量/预测变量的取值 来预测因变量/反应变量的取值,计算得到通过给定散点的最佳拟合线(也称作回归线)。我们使用统计软件来估计回归线的系数。回归系数的截距项表示当预测变量取值为0时因变量的预测均值。此外,预测变量每增加一个单位,因变量就按估计系数的大小保持同步变动。现在,我们来估计以Sampledata
数据集中的Adj.Close
(调整后收盘价)为因变量,以Volume
(交易量)为自变量的线性回归方程。代码及回归估计结果如下:
Y <- Sampledata$Adj.Close
X <- Sampledata$Volume
fit <- lm(Y ~ X)
summary(fit)
##
## Call:
## lm(formula = Y ~ X)
##
## Residuals:
## Min 1Q Median 3Q Max
## -12.305 -5.163 -0.419 5.911 14.279
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 1.95e+02 2.43e+00 80.15 <2e-16 ***
## X -1.88e-07 5.06e-07 -0.37 0.71
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 6.93 on 48 degrees of freedom
## Multiple R-squared: 0.00287, Adjusted R-squared: -0.0179
## F-statistic: 0.138 on 1 and 48 DF, p-value: 0.712
我们也可以用类似的方法来估计其他回归模型(如多元或其他形式的回归模型)。
所谓异常值(outliers)是指远离其他数据的观测,它可能导致统计分析结果出现偏差,因此在分析中考虑异常值的存在并对其加以处理是非常重要的步骤。在R中存在多种侦测异常值的方法,本节讨论其中最常用的一种。
让我们为Sampledata
数据集中的Volume
变量生成一幅箱形图(boxplot),代码及示意图(图2-8)如下:
boxplot(Sampledata$Volume, main = "Volume", boxwex = 0.1)
图2-8 用来侦测异常值的箱形图
从图2-8中我们可以清楚地看到落在箱形图的虚线之外的两个异常值。
局部异常值因子(Local Outlier Factor, LOF)用来识别基于密度的局部异常值。在LOF分析中,将样本点的局部密度与其邻居的局部密度进行比较。如果与其邻居相比,一个样本点落在更加稀疏的区域,则该样本点被视为异常值。让我们用DMwR
程序包中的lofactor()
函数计算Sampledata
数据集的某些变量的LOF,代码及示意图(图2-9)如下:
# install.packages("DMwR")
library(DMwR)
Sampledata1 <- Sampledata[, 2:4]
outlier.scores <- lofactor(Sampledata1, k = 4)
plot(density(outlier.scores))
图2-9 用LOF法显示异常值
如果用户想要提取出排名前5的异常值,则可运行如下代码:
order(outlier.scores, decreasing = TRUE)[1:5]
## [1] 50 34 40 33 22
结果显示的是排名前5异常值所在的行号。
统计建模有时会用到很多变量,而这些变量的量纲可能各不相同,此时对变量进行标准化(standardization)变换至关重要。基于比较的目的,我们需要对这些变量进行标准化,使得它们的量纲大致相同。在R中可用scale()
函数来对变量进行标准化操作并计算z值。
scale()
函数的参数如下。
x
:数值对象。center
:如果为TRUE
,则从数值对象每列的取值中减去对应列的均值(忽略缺失值);如果为FALSE
,则不进行中心化操作。scale
:如果为TRUE
且center
参数也为TRUE
,则将中心化后每列的取值除以对应列的标准差,否则除以均方根;如果为FALSE
,则不更改量纲。如果我们想要中心化Sampledata
数据集中的Volume
变量,只需执行如下代码:
scale(Sampledata$Volume, center = TRUE, scale = FALSE)
如果我们想要标准化Sampledata
数据集中的Volume
变量,则需要执行如下代码:
scale(Sampledata$Volume, center = TRUE, scale = TRUE)
归一化(normalization)是指用变量的最小值和最大值来使得不同变量的量纲相同,计算公式如下:
如果我们想要归一化Sampledata
数据集中的Volume
变量,则可执行如下代码(输出结果从略):
normalized <- (Sampledata$Volume - min(Sampledata$Volume)) /
(max(Sampledata$Volume) - min(Sampledata$Volume))
normalized
1.在R中分别构建正态分布、泊松分布和均匀分布的例子。
2.在R中如何进行随机抽样和分层抽样?
3.趋中性的度量指标有哪些?在R中它们各自的实现函数是什么?
4.在R中如何计算峰度和偏度?
5.在R中如何完成已知总体方差和未知总体方差的假设检验?
6.在R中如何侦测异常值?
7.在R中如何对线性模型和MLE进行参数估计?
8.什么是变量的标准化和归一化变换?在R中如何实现?
在本章中,我们逐一讨论如下概念和问题:金融领域中最常用到的分布函数以及如何用R完成相应数值的计算;随机抽样和分层抽样;趋中性的度量指标;相关性以及时间序列模型中不同类型的相关性;已知方差和未知方差的假设检验(单尾和双尾);侦测异常值;参数估计以及使得变量的量纲大致可比的标准化/归一化变换。