原文链接：http://tecdat.cn/?p=19469
本文将分析工业指数（DJIA）。工业指数（DIJA）是一个股市指数，表明30家大型上市公司的价值。工业指数（DIJA）的价值基于每个组成公司的每股股票价格之和。
获取数据
利用quantmod软件包中提供的getSymbols（）函数，我们可以获得2007年至2018年底的工业平均指数。
getSymbols("^DJI", from = "2007-01-01", to = "2019-01-01") dim(DJI) ## \[1\] 3020 6 class(DJI) ## \[1\] "xts" "zoo"
让我们看一下DJI xts对象，它提供了六个时间序列，我们可以看到。
head(DJI) ## DJI.Open DJI.High DJI.Low DJI.Close DJI.Volume DJI.Adjusted ## 2007-01-03 12459.54 12580.35 12404.82 12474.52 327200000 12474.52 ## 2007-01-04 12473.16 12510.41 12403.86 12480.69 259060000 12480.69 ## 2007-01-05 12480.05 12480.13 12365.41 12398.01 235220000 12398.01 ## 2007-01-08 12392.01 12445.92 12337.37 12423.49 223500000 12423.49 ## 2007-01-09 12424.77 12466.43 12369.17 12416.60 225190000 12416.60 ## 2007-01-10 12417.00 12451.61 12355.63 12442.16 226570000 12442.16 tail(DJI) ## DJI.Open DJI.High DJI.Low DJI.Close DJI.Volume DJI.Adjusted ## 2018-12-21 22871.74 23254.59 22396.34 22445.37 900510000 22445.37 ## 2018-12-24 22317.28 22339.87 21792.20 21792.20 308420000 21792.20 ## 2018-12-26 21857.73 22878.92 21712.53 22878.45 433080000 22878.45 ## 2018-12-27 22629.06 23138.89 22267.42 23138.82 407940000 23138.82 ## 2018-12-28 23213.61 23381.88 22981.33 23062.40 336510000 23062.40 ## 2018-12-31 23153.94 23333.18 23118.30 23327.46 288830000 23327.46
更准确地说，我们有可用的OHLC（开盘，高，低，收盘）指数值，调整后的收盘价和交易量。在这里，我们可以看到生成的相应图表。
我们在此分析调整后的收盘价。
DJI\[,"DJI.Adjusted"\]
简单对数收益率
简单的收益定义为：
对数收益率定义为：
我们计算对数收益率。
CalculateReturns(dj_close, method = "log")
让我们看看。
head(dj_ret) ## DJI.Adjusted ## 2007-01-04 0.0004945580 ## 2007-01-05 -0.0066467273 ## 2007-01-08 0.0020530973 ## 2007-01-09 -0.0005547987 ## 2007-01-10 0.0020564627 ## 2007-01-11 0.0058356461 tail(dj_ret) ## DJI.Adjusted ## 2018-12-21 -0.018286825 ## 2018-12-24 -0.029532247 ## 2018-12-26 0.048643314 ## 2018-12-27 0.011316355 ## 2018-12-28 -0.003308137 ## 2018-12-31 0.011427645
给出了下面的图。
可以看到波动率的急剧上升和下降。第3部分将对此进行深入验证。
辅助函数
我们需要一些辅助函数来简化一些基本的数据转换，摘要和绘图。
1.从xts转换为带有year and value列的数据框。这样就可以进行年度总结和绘制。
df\_t <- data.frame(year = factor(year(index(data\_xts))), value = coredata(data_xts)) colnames(df_t) <- c( "year", "value")
2.摘要统计信息，用于存储为数据框列的数据。
rownames(basicStats(rnorm(10,0,1))) # 基本统计数据输出行名称 with(dataset, tapply(value, year, basicStats))
3.返回关联的列名。
colnames(basicstats\[r, which(basicstats\[r,\] > threshold), drop = FALSE\])
4.基于年的面板箱线图。
p <- ggplot(data = data, aes(x = year, y = value)) + theme\_bw() + theme(legend.position = "none") + geom\_boxplot(fill = "blue")
5.密度图，以年份为基准。
p <- ggplot(data = data, aes(x = value)) + geom_density(fill = "lightblue") p <- p + facet_wrap(. ~ year)
6.基于年份的QQ图。
p <- ggplot(data = dataset, aes(sample = value)) + stat\_qq(colour = "blue") + stat\_qq_line() p <- p + facet_wrap(. ~ year)
7. Shapiro检验
pvalue <- function (v) { shapiro.test(v)$p.value }
每日对数收益率探索性分析
我们将原始的时间序列转换为具有年和值列的数据框。这样可以按年简化绘图和摘要。
head(ret_df) ## year value ## 1 2007 0.0004945580 ## 2 2007 -0.0066467273 ## 3 2007 0.0020530973 ## 4 2007 -0.0005547987 ## 5 2007 0.0020564627 ## 6 2007 0.0058356461 tail(ret_df) ## year value ## 3014 2018 -0.018286825 ## 3015 2018 -0.029532247 ## 3016 2018 0.048643314 ## 3017 2018 0.011316355 ## 3018 2018 -0.003308137 ## 3019 2018 0.011427645
基本统计摘要
给出了基本统计摘要。
## 2007 2008 2009 2010 2011 ## nobs 250.000000 253.000000 252.000000 252.000000 252.000000 ## NAs 0.000000 0.000000 0.000000 0.000000 0.000000 ## Minimum -0.033488 -0.082005 -0.047286 -0.036700 -0.057061 ## Maximum 0.025223 0.105083 0.066116 0.038247 0.041533 ## 1. Quartile -0.003802 -0.012993 -0.006897 -0.003853 -0.006193 ## 3. Quartile 0.005230 0.007843 0.008248 0.004457 0.006531 ## Mean 0.000246 -0.001633 0.000684 0.000415 0.000214 ## Median 0.001098 -0.000890 0.001082 0.000681 0.000941 ## Sum 0.061427 -0.413050 0.172434 0.104565 0.053810 ## SE Mean 0.000582 0.001497 0.000960 0.000641 0.000837 ## LCL Mean -0.000900 -0.004580 -0.001207 -0.000848 -0.001434 ## UCL Mean 0.001391 0.001315 0.002575 0.001678 0.001861 ## Variance 0.000085 0.000567 0.000232 0.000104 0.000176 ## Stdev 0.009197 0.023808 0.015242 0.010182 0.013283 ## Skewness -0.613828 0.224042 0.070840 -0.174816 -0.526083 ## Kurtosis 1.525069 3.670796 2.074240 2.055407 2.453822 ## 2012 2013 2014 2015 2016 ## nobs 250.000000 252.000000 252.000000 252.000000 252.000000 ## NAs 0.000000 0.000000 0.000000 0.000000 0.000000 ## Minimum -0.023910 -0.023695 -0.020988 -0.036402 -0.034473 ## Maximum 0.023376 0.023263 0.023982 0.038755 0.024384 ## 1. Quartile -0.003896 -0.002812 -0.002621 -0.005283 -0.002845 ## 3. Quartile 0.004924 0.004750 0.004230 0.005801 0.004311 ## Mean 0.000280 0.000933 0.000288 -0.000090 0.000500 ## Median -0.000122 0.001158 0.000728 -0.000211 0.000738 ## Sum 0.070054 0.235068 0.072498 -0.022586 0.125884 ## SE Mean 0.000470 0.000403 0.000432 0.000613 0.000501 ## LCL Mean -0.000645 0.000139 -0.000564 -0.001298 -0.000487 ## UCL Mean 0.001206 0.001727 0.001139 0.001118 0.001486 ## Variance 0.000055 0.000041 0.000047 0.000095 0.000063 ## Stdev 0.007429 0.006399 0.006861 0.009738 0.007951 ## Skewness 0.027235 -0.199407 -0.332766 -0.127788 -0.449311 ## Kurtosis 0.842890 1.275821 1.073234 1.394268 2.079671 ## 2017 2018 ## nobs 251.000000 251.000000 ## NAs 0.000000 0.000000 ## Minimum -0.017930 -0.047143 ## Maximum 0.014468 0.048643 ## 1. Quartile -0.001404 -0.005017 ## 3. Quartile 0.003054 0.005895 ## Mean 0.000892 -0.000231 ## Median 0.000655 0.000695 ## Sum 0.223790 -0.057950 ## SE Mean 0.000263 0.000714 ## LCL Mean 0.000373 -0.001637 ## UCL Mean 0.001410 0.001175 ## Variance 0.000017 0.000128 ## Stdev 0.004172 0.011313 ## Skewness -0.189808 -0.522618 ## Kurtosis 2.244076 2.802996
在下文中，我们对上述一些相关指标进行了具体评论。
平均值
每日对数收益率具有正平均值的年份是：
filter_stats(stats, "Mean", 0) ## \[1\] "2007" "2009" "2010" "2011" "2012" "2013" "2014" "2016" "2017"
按升序排列。
## 2008 2018 2015 2011 2007 2012 2014 ## Mean -0.001633 -0.000231 -9e-05 0.000214 0.000246 0.00028 0.000288 ## 2010 2016 2009 2017 2013 ## Mean 0.000415 5e-04 0.000684 0.000892 0.000933
中位数
正中位数是：
filter\_stats(dj\_stats, "Median", 0) ## \[1\] "2007" "2009" "2010" "2011" "2013" "2014" "2016" "2017" "2018"
以升序排列。
## 2008 2015 2012 2017 2010 2018 2014 ## Median -0.00089 -0.000211 -0.000122 0.000655 0.000681 0.000695 0.000728 ## 2016 2011 2009 2007 2013 ## Median 0.000738 0.000941 0.001082 0.001098 0.001158
偏度
偏度（Skewness）可以用来度量随机变量概率分布的不对称性。
公式： 
其中 
是均值， 
是标准差。
几何意义：
偏度的取值范围为(-∞,+∞)
当偏度<0时，概率分布图左偏（也叫负偏分布，其偏度<0）。
当偏度=0时，表示数据相对均匀的分布在平均值两侧，不一定是绝对的对称分布。
当偏度>0时，概率分布图右偏（也叫正偏分布，其偏度>0）。
例如上图中，左图形状左偏，右图形状右偏。
每日对数收益出现正偏的年份是：
## \[1\] "2008" "2009" "2012"
按升序返回对数偏度。
stats\["Skewness",order(stats\["Skewness", ## 2007 2011 2018 2016 2014 2013 ## Skewness -0.613828 -0.526083 -0.522618 -0.449311 -0.332766 -0.199407 ## 2017 2010 2015 2012 2009 2008 ## Skewness -0.189808 -0.174816 -0.127788 0.027235 0.07084 0.224042
峰度
峰度（Kurtosis）可以用来度量随机变量概率分布的陡峭程度。
公式：
其中 
是均值， 
是标准差。
几何意义：
峰度的取值范围为[1,+∞)，完全服从正态分布的数据的峰度值为 3，峰度值越大，概率分布图越高尖，峰度值越小，越矮胖。
例如上图中，左图是标准正太分布，峰度=3，右图的峰度=4，可以看到右图比左图更高尖。
通常我们将峰度值减去3，也被称为超值峰度（Excess Kurtosis），这样正态分布的峰度值等于0，当峰度值>0，则表示该数据分布与正态分布相比较为高尖，当峰度值<0，则表示该数据分布与正态分布相比较为矮胖。
点击标题查阅往期内容
R语言风险价值：ARIMA，GARCH，Delta-normal法滚动估计VaR（Value at Risk）和回测分析股票数据
01
02
03
04
每日对数收益出现超值峰度的年份是：
## \[1\] "2007" "2008" "2009" "2010" "2011" "2012" "2013" "2014" "2015" "2016" ## \[11\] "2017" "2018"
按升序返回超值峰度。
## 2012 2014 2013 2015 2007 2010 2009 ## Kurtosis 0.84289 1.073234 1.275821 1.394268 1.525069 2.055407 2.07424 ## 2016 2017 2011 2018 2008 ## Kurtosis 2.079671 2.244076 2.453822 2.802996 3.670796
2018年的峰度最接近2008年。
箱形图
我们可以看到2008年出现了最极端的值。从2009年开始，除了2011年和2015年以外，其他所有值的范围都变窄了。但是，与2017年和2018年相比，产生极端值的趋势明显改善。
密度图
densityplot(ret_df)
2007年具有显着的负偏。2008年的特点是平坦。2017年的峰值与2018年的平坦度和左偏一致。
R语言股票市场指数：ARMA-GARCH模型和对数收益率数据探索性分析（中）：https://developer.aliyun.com/article/1491668
