第3章 预测股票市场收益

3.3 定义预测任务

图3-1通过蜡烛图来反映股票的交易趋势的变化。对应的代码为

T.ind <- function(quotes,tgt.margin=0.025,n.days=10) {  v <- apply(HLC(quotes),1,mean)  r <- matrix(NA,ncol=n.days,nrow=NROW(quotes))  ## The following statment is wrong in the book (page 109)!  for(x in 1:n.days) r[,x] <- Next(Delt(Cl(quotes),v,k=x),x)  x <- apply(r,1,function(x) sum(x[x > tgt.margin | x < -tgt.margin]))  if (is.xts(quotes)) xts(x,time(quotes)) else x}candleChart(last(GSPC,'3 months'),theme='white',TA=NULL)avgPrice <- function(p) apply(HLC(p),1,mean)addAvgPrice <- newTA(FUN=avgPrice,col=1,legend='AvgPrice')addT.ind <- newTA(FUN=T.ind,col='red',legend='tgtRet')addAvgPrice(on=1)addT.ind()

在R中运行上面的代码，得到如下的结果图。

3.3.2 预测的变量是什么

　　我们假设如果过去的某些行为p之后是另一个行为f,并且如果这一因果链经常发生，那假设这一现象未来会再次发生就是合理的。技术指标反映价格时间序列特征的数值汇总。在R的添加包中TTR中可以找到这些指标。这些技术指标可以获取价格序列的某些特征，例如价格是否波动太大，是否具有某种特定趋势等。目前为止，针对不同领域应用的技术指标仍然是一个研究课题。

　　变量选择方法通常分为两类：1） 变量过滤器；2）变量封装。前者的变量选择独立于变量选择之后的模型构建阶段。它一般使用变量的某些统计特性（如相关）来选择用于建模的数据集。变量封装方法在变量的选择过程中要包含后面构建的模型信息。它是一个循环选择的过程，在循环选择过程的每一步中，一组候选变量用于特定的模型病记录模型的相应的结果。给予选择记录的结果，用某些搜索操作确定一组新的候选输入变量自己，重复这个过程直到定义了最终变量自己的收敛准则满足。

　　下面的例子，给出一组输入变量，然后用一个方法来估计这些变量的重要性，然后基于他们的重要性来选择最适合问题的变量。

　　在下面的代码中，我们对一些TTR函数的输出进行了后续处理以获得一个单一值：

myATR <- function(x) ATR(HLC(x))[,'atr']mySMI <- function(x) SMI(HLC(x))[,'SMI']myADX <- function(x) ADX(HLC(x))[,'ADX']myAroon <- function(x) aroon(x[,c('High','Low')])$oscillatormyBB <- function(x) BBands(HLC(x))[,'pctB']myChaikinVol <- function(x) Delt(chaikinVolatility(x[,c("High","Low")]))[,1]myCLV <- function(x) EMA(CLV(HLC(x)))[,1]myEMV <- function(x) EMV(x[,c('High','Low')],x[,'Volume'])[,2]myMACD <- function(x) MACD(Cl(x))[,2]myMFI <- function(x) MFI(x[,c("High","Low","Close")], x[,"Volume"])mySAR <- function(x) SAR(x[,c('High','Close')]) [,1]myVolat <- function(x) volatility(OHLC(x),calc="garman")[,1]

　　上面的变量给出了预测指标T的未来值的初始预测变量集合。下面应用变量选择方法把这22个预测变量进行精简。随机森林也可以用于估计预测任务中变量的重要性。是通过计算每个变量被移除后随机森林的误差来计算机变量的重要性。

　　书中的案例，是把数据集分成2个独立的自己：1.一个数据集用于构建交易系统。2.用于最终的交易系统的测试。

训练集数据构建采用随机森林模型，代码如下：

T.ind <- function(quotes,tgt.margin=0.025,n.days=10) {  v <- apply(HLC(quotes),1,mean)  r <- matrix(NA,ncol=n.days,nrow=NROW(quotes))  ## The following statment is wrong in the book (page 109)!  for(x in 1:n.days) r[,x] <- Next(Delt(Cl(quotes),v,k=x),x)  x <- apply(r,1,function(x) sum(x[x > tgt.margin | x < -tgt.margin]))  if (is.xts(quotes)) xts(x,time(quotes)) else x}# 加载自定义的函数myATR <- function(x) ATR(HLC(x))[,'atr']mySMI <- function(x) SMI(HLC(x))[,'SMI']myADX <- function(x) ADX(HLC(x))[,'ADX']myAroon <- function(x) aroon(x[,c('High','Low')])$oscillatormyBB <- function(x) BBands(HLC(x))[,'pctB']myChaikinVol <- function(x) Delt(chaikinVolatility(x[,c("High","Low")]))[,1]myCLV <- function(x) EMA(CLV(HLC(x)))[,1]myEMV <- function(x) EMV(x[,c('High','Low')],x[,'Volume'])[,2]myMACD <- function(x) MACD(Cl(x))[,2]myMFI <- function(x) MFI(x[,c("High","Low","Close")], x[,"Volume"])mySAR <- function(x) SAR(x[,c('High','Close')]) [,1]myVolat <- function(x) volatility(OHLC(x),calc="garman")[,1]# 使用randomForest R添加包library(randomForest)# 根据几个初始化的技术指标建立随机森林模型data.model <- specifyModel(T.ind(GSPC) ~ Delt(Cl(GSPC),k=1:10) +        myATR(GSPC) + mySMI(GSPC) + myADX(GSPC) + myAroon(GSPC) +        myBB(GSPC)  + myChaikinVol(GSPC) + myCLV(GSPC) +        CMO(Cl(GSPC)) + EMA(Delt(Cl(GSPC))) + myEMV(GSPC) +        myVolat(GSPC)  + myMACD(GSPC) + myMFI(GSPC) + RSI(Cl(GSPC)) +       mySAR(GSPC) + runMean(Cl(GSPC)) + runSD(Cl(GSPC)))# 设置随机种子set.seed(1234)# 构建模型：1.指明模型，2.指明构建方法 3. 指明数据集 4.importance 随机森林开启变量估计rf <- buildModel(data.model,method='randomForest',             training.per=c(start(GSPC),index(GSPC["1999-12-31"])),             ntree=50, importance=T)ex.model <- specifyModel(T.ind(IBM) ~ Delt(Cl(IBM),k=1:3))data <- modelData(ex.model,data.window=c('2009-01-01','2009-08-10'))varImpPlot([e-mail],type=1)

　　现在我们需要确定一个界限值来选择重要性评分高的变量子集。纵观图3-2的结果，因为这是一个简单的用随机森林来选择变量的例子，所以我们将使用界限值10：

> imp <- importance([e-mail],type=1)> rownames(imp)[which(imp > 10)][3] "myATR.GSPC"   "myADX.GSPC"   "myEMV.GSPC"   "myVolat.GSPC" "myMACD.GSPC" [6] "myMFI.GSPC" 

　　函数importance()将得到每个变量的具体重要性分数。然后，我们用给定的界限值来筛选用于建模的变量子集。

> data.model <- specifyModel(T.ind(GSPC) ~ Delt(Cl(GSPC),k=1) + myATR(GSPC) + myADX(GSPC) +    myEMV(GSPC) + myVolat(GSPC)  + myMACD(GSPC) + mySAR(GSPC) + runMean(Cl(GSPC)) )Warning message:In sqrt(N/n * runSum(0.5 * log(OHLC[, 2]/OHLC[, 3])^2 - (2 * log(2) -  :  NaNs produced> set.seed(1234)> rf <- buildModel(data.model,method='randomForest',+              training.per=c(start(GSPC),index(GSPC["1999-12-31"])),+              ntree=50, importance=T)> varImpPlot([e-mail],type=1)

3.3.3 预测任务

　　在3.3.2节中，我们已经获得一个quantmod对象，该对象包含了用于建立模型的数据集。

　　目标变量为名义变量的预测问题称为“分类任务或者是判别任务”。分类任务和回归任务的主要区别是目标变量的类型。回归任务有一个数值型的目标变量，而分类任务的目标变量是取值为有限个的名义变量。

　　R的xts添加包不允许训练集数据的一列为名义变量，而其他列是数值型变量。为了克服这个问题，我们需要在xts之外完成整个建模的处理过程。我们会使用xts的功能来进行数据子集的选取和绘图。下面的代码为构造本节下面部分两个预测任务的预测模型所应用的数据结构。

# 设置训练数据集的时间窗口为1970/01/02~1999/12/31的十年数据Tdata.train <- as.data.frame(modelData(data.model,                       data.window=c('1970-01-02','1999-12-31')))# 设置测试数据集为近10年的数据Tdata.eval <- na.omit(as.data.frame(modelData(data.model,                       data.window=c('2000-01-01','2009-09-15'))))Tform <- as.formula('T.ind.GSPC ~ .')

3.3.4 模型评价准则

　　上一节描述的预测任务获得的模型可以给出未来市场方向的预测值。对于回归任务而言，这个预测值得一个数值，对分类问题而言，就是一个信号。

　　金融市场预测是罕见的事件驱动应用的一个例子。基于事件的预测任务通常由决策精确度指标和回溯精确度来衡量，这两个指标集中于所关注的事件的评估而不是常见事情的评估。决策精确度衡量模型给出的事件信号的正确百分比；而回溯精确度则是指模型给出的事件信号占事实存在的百分比。通过分类矩阵是通过比较模型的预测值和真实值。

3.4 预测模型

　　主要通过研究模型处理非线性回归问题的优劣来选择最佳模型，这就是我们案例中所要研究的问题。当然，许多其他的方法也可以用来解决这个问题。在研究股票收益这样一个领域中，任何一个更好的研究方法都必然需要经历更多的比较和更多的选择。

3.4.1 如何应用训练集数据来建模

　　复杂的事件序列问题通常有着不同的表现形式，比如序列由于波动大的阶段到相对比较平稳的阶段，或者是具有某种趋势倾向的序列。这些类型的序列通常称为是非平稳时间序列，在某些模型的基本假设条件下，非平稳序列将会带来更严重的问题。

　　模型需要能够适应当前数据的变化。

　　增量学习使用心得训练案例数据来改变原来的模型。这些方法通常称为增量学习，因为它更新模型以适应最新的信息而无需从头开始。目前还没有太多的建模技巧应用这种更新模型的方式。更一般的做法是使用新的数据来重新建立模型。
　　扩大和滑动窗口的方法都涉及到一个关键的决定：什么时候应该通过纳入更新的数据来改善或者调整模型？解决这个问题的方法：1.需要通过检查来估计模型的预测能力从什么时候开始降低如果我们管擦在某个时候，模型的预测性能突然降低，我们就可以认为在这时模型发生了变化。这种方法最主要的挑战在于给出模型性能变化的合理估计。第二种更简单一些的方法是在常规的事件的基础上更新模型，也就是说，每个w个测试用例，用更新的数据来建立新的模型，在本案例中，我们就采用这种方法。

　　总之，我们考虑每一个模型应用会用到以下三种方法：1所有的测试时间段都使用一个模型；2每隔w个测试用例更新数据增长窗口；3用每隔w天的护具滑动窗口。

3.4.2 建模工具

　　本节将会介绍建模的方法，将使用这些方法来完成我们的预测任务，并说明在R上如何使用这些方法。

3.4.2.1 人工神经网络

　　人工神经网络经常在金融预测中使用，因为运用人工神经网络可以处理高度非线性问题。R的添加包nnet可以实现前馈神经网络。这种类型的神经网络是最常用的，也是我们将要使用的神经网络。

　　人工神经网络由中相互联系的计算单元构成。每个神经元执行两次连续的计算：输入的线性组合；之后对前面的结果的非线性计算得到的输出值作为神经网络的下一个神经元的输入。每个神经元都有一个相关联的权重。要构建人工神经网络，先要建立网络体系结构，然后使用一种算法来计算出神经元之间的连接权重。

　　前馈神经网络按照层来组织神经元。第一层包含输入神经元，训练集的观测值通过这些输入神经元传递给网络。最后一层包含了任何情况下传递给神经网络输入神经元的神经网络预测值。在这两层之间，通常有一个或者多个隐藏层神经元。权重更新算法，比如反向传播法，试图获得能够优化某个误差标准的连接权重。也就是试图确保网络输出于提交给神经网络模型的训练集个案一致。这是通过在网络输入结点多次传入训练个案来进行迭代的一个过程，在网络输出结点获得预测值病计算出各自的预测误差之后，通过更新网络中的权重来减小模型的预测误差。这种迭代过程反复进行，知道满足一定的收敛准侧。

前馈神经网络既可以用于分类问题，也可以用于回归问题。

　　人工神经网络对预测问题中的变量的尺度敏感。这种情况下，我们有必要在运用神经网络之前，对数据进行转换。这样可以避免神经网络模型的性能受到变量尺度的影响。在我们的案例中，先进行数据的标准化处理，使所有变量均具有零均值和标准差为1.通过下面的公式我们可以很容易地对数据集中的每一列应用下列公式进行转换：

yi=xi−x¯σx

其中，x¯ 是原始变量X的均值，σx是变量X的标准差。

　　函数scale()可以用来进行上面的数据变换，在本书的R添加包中，还可以找到函数unscale()，这个函数可以进行标准化过程的逆运算，将转换后的数据转变为原来尺度的数据。
下面的例子说明了如何在R中获取和使用这种类型的神经网络。

> set.seed(1234)> library(nnet)> norm.data <- scale(Tdata.train)> nn <- nnet(Tform,norm.data[1:1000,],size=10,decay=0.01,maxit=1000,linout=T,trace=F)> norm.preds <- predict(nn,norm.data[1001:2000,])> preds <- unscale(norm.preds,norm.data)

　　注意：在默认的情况下，函数nnet()以在区间[−0.5,0.5]的随机值来设置结点之间链接的初始权重，这意味着连续两次运行有完全相同的参数的统一函数，得到的实际结果可能是不一样的。为了确保可以得到同样的结果，我们像以前一样，增加了函数set.seed()的调用，这样可以确保得到与书中一样的结果。nnet() 函数的前两个参数是所有的函数都要的参数。模型的函数形式通过公式化来具体化，训练集数据用于建立模型。这里的size用来指定隐藏层的结点个数。参数decay控制反向传播算法更新权重的更新率。最后，参数maxit控制权重收敛过程中所允许使用的最大迭代次数；参数linout=T告诉该函数是处理回归问题；参数trace=F用来避免一些和优化过程相关的结果被输出。很好奇优化的过程是什么，于是改成linout=T执行结果如下：

# weights:  101initial  value 3888.449360 iter  10 value 901.426850iter  20 value 852.958920iter  30 value 762.829104iter  40 value 676.779362iter  50 value 586.451259iter  60 value 542.854098iter  70 value 496.351118iter  80 value 440.983977iter  90 value 420.154494iter 100 value 409.918942iter 110 value 394.500011iter 120 value 379.637372iter 130 value 374.019078iter 140 value 368.425846iter 150 value 364.667234iter 160 value 358.498127iter 170 value 350.518014iter 180 value 346.602056iter 190 value 343.502375iter 200 value 340.811954iter 210 value 339.644652iter 220 value 338.983495iter 230 value 337.912043iter 240 value 336.640434iter 250 value 334.378456iter 260 value 333.175536iter 270 value 331.871760iter 280 value 331.458035iter 290 value 331.042287iter 300 value 330.810038iter 310 value 330.736837iter 320 value 330.722524iter 330 value 330.707416iter 340 value 330.570134iter 350 value 330.442832iter 360 value 330.416749iter 370 value 330.269002iter 380 value 329.977801iter 390 value 329.931536iter 400 value 329.917680iter 410 value 329.915387final  value 329.915345 converged

　　函数predict() 可以用来获得测试数据集神经网络预测值。取得这些预测值之后，我们使用R添加包提供的函数unscale()转换为原来尺度的数据。
下面将评估人工神经网络模型预测测试集信号的准确性。我们可以将数值型预测值转换成信号，然后使用3.3.4节中的统计方法。

> norm.preds <- predict(nn,norm.data[1001:2000,])> preds <- unscale(norm.preds,norm.data)> sigs.nn <- trading.signals(preds,0.1,-0.1)> true.sigs <- trading.signals(Tdata.train[1001:2000,'T.ind.GSPC'],0.1,-0.1)> sigs.PR(sigs.nn,true.sigs)    precision    recalls   0.1923077 0.2222222b   0.3009259 0.4140127s+b 0.2476415 0.3115727

　　函数sigs.PR()可以将我们关系的两类事件以及决策精度和回溯精度矩阵。这些值表面人工神经网络的预测性能不是很好。实际上，会得到相当低的预测精确度，回溯精确度值也不是很好。但是，后者比较差所导致的问题不是太严重，因为它基本上意味着失去机会而不是意味着有成本损失。另外一方面，预测精确度的较小值意味着该模型频繁给出错误信号。如果这些信号被用来交易，就可能导致严重的损失。

　　人工神经网络也可以用于分类问题。对于这些问题，在网络拓扑方面最重要的区别是一个单一的输出单位。我们将有和目标变量值（有时称为类变量）一样多的输出单位，这些输出单位的每一个都将产生各自类值的概率估计。这意味着，对于每一个测试个案，人工神经网络都可以产生一组概率值，每一个概率值相应于一个可能的类值。
　　使用nnet()函数来完成这个任务与使用该函数来解决回归问题很相似。使用训练集数据并使用下面的数据进行演示：

> set.seed(1234)> library(nnet)> signals <- trading.signals(Tdata.train[,'T.ind.GSPC'],0.1,-0.1)> norm.data <- data.frame(signals=signals,scale(Tdata.train[,-1]))> nn <- nnet(signals ~ .,norm.data[1:1000,],size=10,decay=0.01,maxit=1000,trace=F)> preds <- predict(nn,norm.data[1001:2000,],type='class')

　　这里的参数type="class"是用于获得测试集个案的类标签，而不是概率的估计值。上面的代码如之前一样，也没有将优化的过程输出。通过使用trace=T来输出优化的过程。在神经网络预测中，可以计算出模型的预测精确度和回溯精确度。代码如下所示：

> sigs.PR(preds,norm.data[1001:2000,1])    precision    recalls   0.3461538 0.2000000b   0.3380282 0.3057325s+b 0.2398374 0.1750742

　　上面的结果中，预测精确度和回溯精确度的值还比较低，但是都高于回归任务中相应的值。

3.4.2.2 支持向量机

　　支持向量机和人工神经网络一样，也是一种建模工具，可以处理回归和分类问题。给予其成功应用到多个领域和其强大的理论背景，支持向量机已经变到越来越多不同研究领域的关注。
　　支持向量机的基本思想是，将原始数据映射到一个新的高维空间中，在这个新的高维空间中，有可能应用线性模型来获得一个超平面来进行分离，例如在分类任务中，分离问题中的不同类别。将原始数据映射到这一新的空间是在所谓的核函数的帮助下进行的。支持向量机是作用在由核函数所引入的对称表示的线性机。
在新的对称表示下进行超平面分割，这是通过最大化不同类别之间个案的分割边际来进行的。

　　下面，我们将提供使用R中这类模型的简单例子。我们从使用添加包e1071中的函数进行回归任务开始，代码如下：

> library(e1071)> sv <- svm(Tform,Tdata.train[1:1000,],gamma=0.001,cost=100)> s.preds <- predict(sv,Tdata.train[1001:2000,])> sigs.svm <- trading.signals(s.preds,0.1,-0.1)> true.sigs <- trading.signals(Tdata.train[1001:2000,'T.ind.GSPC'],0.1,-0.1)> sigs.PR(sigs.svm,true.sigs)    precision      recalls         0.3 0.016666667b         NaN 0.000000000s+b       0.3 0.008902077

注意上面的e1071R添加包需要进行安装。

　　在这个例子中，我们使用了函数svm()，除了参数gamma和cost外，大部分的参数我们都采用默认值。在本节例子中，该函数使用了一个径向基核函数。

K(x,y)=exp(−γ×||x−y||2)

　　其中，γ是一个用户参数，在上面的函数调用中，它的值设置为0.001。
　　参数cost给出违反边际所引入的算式 。根据上面的代码得出的运行结果，我们可以看到svm的决策精度比nnet好许多。 尽管回溯精确度比较低。

　　下一步，我们考虑分类任务，这次使用的是R的kernlab添加包，代码如下：

>install.packages("kernlab")>library(kernlab)>data <- cbind(signals=signals,Tdata.train[,-1])>ksv <- ksvm(signals ~ .,data[1:1000,],C=10)>ks.preds <- predict(ksv,data[1001:2000,])>sigs.PR(ks.preds,data[1001:2000,1])    precision    recalls   0.1592920 0.2000000b   0.2795699 0.1656051s+b 0.1943574 0.1839763

　　使用kernlab添加包中的函数ksvm()，其中的参数C用来指定违反约束的不同损失。参数的默认值为1。除此之外，其他参数使用默认值。

　　支持向量机分类的结果并不如支持向量机回归的结果好。这并不意味这我们声明这是我们用支持向量机技术所能获得的最好结果。这些都只是说明如何在R中使用这些模型技术的简单例子。

3.4.2.3 多元自适应回归样条

　　多元自适应回归样条(MARS)是自适应回归模型（Hastie and Tibshirani）的一个例子。一个多元自适应回归样条模型具有以下一般形式：

mars(x)=c0+∑i=1kciBi(x)

其中ci是常数，Bi(x)是基函数。
　　基函数可以有多种不同的表现形式，从简单常亮到描述两个或多个变量相互关系的函数。但是，最常见的基函数就是所谓的铰链函数，有如下形式：

H[−(xi−t)]=max(0,t−xi)H[+(xi−t)]=max(0,xi−t)

　　其中xi是一个预测变量，t是一个预测变量的界限值。
　　在R中已经至少有两个添加包实现了多元自适应回归样条模型。添加包mda包含了mars()，该函数实现了多元自适应回归样条方法。添加包earth中的函数earth()同样也实现了多元自适应回归样本方法。从建模函数提供的基于公式的接口来看，函数earth()有着能够遵循更标准的R架构的优势。同时earth还实现了其他的添加包所具备的功能。下面是应用函数earth()进行回归分析的代码。

> library(earth)载入需要的程辑包：plotmo载入需要的程辑包：plotrix载入需要的程辑包：TeachingDemos> e <- earth(Tform,Tdata.train[1:1000,])> e.preds <- predict(e,Tdata.train[1001:2000,])> sigs.e <- trading.signals(e.preds,0.1,-0.1)> true.sigs <- trading.signals(Tdata.train[1001:2000,'T.ind.GSPC'],0.1,-0.1)> sigs.PR(sigs.e,true.sigs)    precision    recalls   0.2991453 0.1944444b   0.4375000 0.2675159s+b 0.3615023 0.2284866

　　得到的结果和支持向量机的归回结果相差不大，决策精确度约为30%，回溯精确度要更低一些。
多元回归样条只适用于回归问题

3.5 从预测到实践

这一届描述如何应用上节的模型所得到的预测信号。给定模型输出的一组信号，可以很许多方式将它们运用到市场的交易决策。

3.5.1 如何应用预测模型

在本案例中，我们假设在期货市场进行交易。期货市场是在合约基础上进行交易的，合约规定在未来某个确定的时间、以未来市场决定的价格买入或者卖出商品。这些合约的技术细节不在此进行讨论。但是，从客观的角度看，这意味着我们的交易系统将可以采取两种交易寸头：多头头寸和空头头寸。多头头寸
是指在时刻t、以价格p买入商品，随后在时刻t+x卖出。当交易者预期未来价格上涨时，这样的头寸对交易者来说是有意义的，进行这种交易使他获得利润。做空头头寸时，交易者在时刻t、以价格p卖出证券，同事他们有义务在未来买回同样的证券。由于可以借入证券的交易模式，这种空头头寸也是可以在时刻t后的某个时刻买入并偿还所借证券，从而获利。简略地说，当认为价格会下降时开空头仓位，认为价格会上涨时多开头仓位。