R演示-数据提取资料

1. 获取实验数据
本次实验可以使用R软件的金融分析包quantmod中的getSymbols()函数，来获取所需的二十支股票的基本数据。

为了统一工作路径，我们首先将R的工作目录设定为D盘。

在设定好默认路径后，我们需要载入quantmod包，它可以帮助我们对金融数据进行分析和处理。

通过quantmod包中的getSymbols( )函数来下载和获取实验要求中的二十只股票从2012年03月1日到2015年03月30日的交易数据。

注意下载程序包时，需要设定CRAN网站镜像，在弹出CRAN mirror的菜单里选择距离网络代理器最近的地点，如China（Beijing4）（图1）。

运行代码如下：
#设定工作目录
setwd("D:/")
#安装quantmod程序包
install.packages("quantmod")
#读取quantmod程序包
library(quantmod)
图1下载CRAN镜像
图2安装与读取程序包
#下载股票基本交易数据，本实验所需的20支股票的数据都可以通过同样的R语句下载，需要将代码部分改成所需的股票代码并按顺序赋值（如赋值为s1、s2、…、s20）。

注意区分沪市和深市的股票，沪市股票代码后应当加上标注“.ss”，深市股票代码后应加标注“.sz”。

将R语句补充完整并在R中运行。

s1<-getSymbols("600016.ss",from="2012-03-01",to="2015-03-31",auto.assign=FALSE) s2<-getSymbols("600059.ss",from="2012-03-01",to="2015-03-31",auto.assign=FALSE) s3<-getSymbols("600060.ss",from="2012-03-01",to="2015-03-31",auto.assign=FALSE) s4<-getSymbols("600069.ss",from="2012-03-01",to="2015-03-31",auto.assign=FALSE) s5<-getSymbols("600070.ss",from="2012-03-01",to="2015-03-31",auto.assign=FALSE) s6<-getSymbols("600079.ss",from="2012-03-01",to="2015-03-31",auto.assign=FALSE) s7<-getSymbols("600080.ss",from="2012-03-01",to="2015-03-31",auto.assign=FALSE) s8<-getSymbols("600085.ss",from="2012-03-01",to="2015-03-31",auto.assign=FALSE) s9<-getSymbols("600095.ss",from="2012-03-01",to="2015-03-31",auto.assign=FALSE) s10<-getSymbols("600114.ss",from="2012-03-01",to="2015-03-31",auto.assign=FALSE) s11<-getSymbols("600123.ss",from="2012-03-01",to="2015-03-31",auto.assign=FALSE) s12<-getSymbols("600066.ss",from="2012-03-01",to="2015-03-31",auto.assign=FALSE) s13<-getSymbols("600160.ss",from="2012-03-01",to="2015-03-31",auto.assign=FALSE)
s14<-getSymbols("600436.ss",from="2012-03-01",to="2015-03-31",auto.assign=FALSE)
s15<-getSymbols("600188.ss",from="2012-03-01",to="2015-03-31",auto.assign=FALSE) s16<-getSymbols("600199.ss",from="2012-03-01",to="2015-03-31",auto.assign=FALSE) s17<-getSymbols("600218.ss",from="2012-03-01",to="2015-03-31",auto.assign=FALSE) s18<-getSymbols("600256.ss",from="2012-03-01",to="2015-03-31",auto.assign=FALSE) s19<-getSymbols("600261.ss",from="2012-03-01",to="2015-03-31",auto.assign=FALSE) s20<-getSymbols("600612.ss",from="2012-03-01",to="2015-03-31",auto.assign=FALSE)
图3下载数据
#计算月收益率。

将R语句补充完整并在R中运行。

s1.m<-monthlyReturn(s1)
s2.m<-monthlyReturn(s2)
s3.m<-monthlyReturn(s3)
s4.m<-monthlyReturn(s4)
s5.m<-monthlyReturn(s5)
s6.m<-monthlyReturn(s6)
s7.m<-monthlyReturn(s7)
s8.m<-monthlyReturn(s8)
s9.m<-monthlyReturn(s9)
s10.m<-monthlyReturn(s10)
s11.m<-monthlyReturn(s11)
s12.m<-monthlyReturn(s12)
s13.m<-monthlyReturn(s13)
s14.m<-monthlyReturn(s14)
s15.m<-monthlyReturn(s15)
s16.m<-monthlyReturn(s16)
s17.m<-monthlyReturn(s17)
s18.m<-monthlyReturn(s18)
s19.m<-monthlyReturn(s19)
s20.m<-monthlyReturn(s20)
图4计算月收益率
#将不同股票的数据合并在一起并赋值给DATA，并将数据转换成时间序列。

将R语句补充完整并在R中运行。

data<-merge(s1.m,s2.m,s3.m,s4.m,s5.m,s6.m,s7.m,s8.m,s9.m,s10.m,s11.m, s12.m,s13.m,s14.m,s15.m,s16.m,s17.m,s18.m,s19.m,s20.m)
data<-as.data.frame(data)
#将数据列用股票代码命名，将R语句补充完整并在R中运行。

colnames(data)<-c("600016","600059","600060","600069","600070","600079", "600080","600085","600095","600114","600123","600141","600160","600172","60
0188","600199","600218","600228","600256","600261")
#将数据以csv的形式存在D盘根目录下
write.csv(data,"stock.csv")
图5赋值
之后在国泰安提取债券数据，打开国泰安官网，选取数据中心-债券，在左侧编辑所需要的参数。

选定债券行情指标-行情指标-债券月交易信息表-月收盘价。

设定时间为2012-3月到2015-3月。

图6交易数据收集
输入所选取的债券的代码120605，111051，122829，112059。

图7交易数据收集
将所收集到的数据汇总到stock.csv文件中，并将股票收益率年化。

2计算回报率的均值、协方差和标准差。

计算均值时使用apply()函数。

apply()函数是通过对数组或者矩阵的一个维度使用函数生成值得到列表或者数组、向量。

第二个参数值为2表示对列向量进行操作，第三个参数为mean表示对列向量进行求平均数的操作。

计算协方差使用cov()函数。

cov()函数用以计算列与列之间的协方差。

同样，通过head()函数可以查看前四只股票的协方差矩阵，如图9。

计算标准差时使用了两个函数，分别为diag()和sqrt()。

diag()函数为求对角阵的函数，返回一个矩阵对角线的值构成的向量。

因此，diag(cov_mat)表示每只资产的方差值构成的向量。

通过函数sqrt()进行开方操作，得到各资产的标准差值。

其代码如下：
# 计算均值，并将其赋值给mean_vec
mean_vec <- apply(R, 2, mean)
#计算协方差，并将其赋值给cov_mat
cov_mat <- cov(R)
#计算标准差，并将其赋值给sd_vec
sd_vec <- sqrt(diag(cov_mat));
可以通过head()函数查看前几个数据，如输入head(mean_vec,4)，head(cov_mat,4)，head(sd_vec,4)，返回查看对应函数输出结果的前4个数据。

如图8~图10。

图8平均值图示
图9协方差矩阵图示
图10标准差图示
2.设定约束条件
本实验根据马克维茨资产组合选择模型，除了两个基本约束条件外，还需加入恒定比例的约束条件。

并且可以通过循环的方式，改变bvec值，从而得到不同期望收益下的最小方差组合，画出有效前沿。

代码如下：
#设定约束矩阵A，第一列为1，第二列为各证券均值，第三列为限制20支股票占比例70%，，第五列开始为24行以1为对角的对角矩阵
A <- cbind(rep(1,24),mean_vec,c(rep(1,20),rep(0,4)),diag(1,nrow=24))
#设定600个目标收益率，以此作为描点的基础
Mu <- seq(0.12,0.3,length=600)
#设置变量Sd用来存储资产组合收益率的标准差
Sd <- Mu
#设置变量W用来存储资产组合比例
W <- matrix(0, nrow=600, ncol=24)
3. 规划求解
本实验采用二次规划求解的方式，利用R语言中的solve.QP函数，用循环语句for()将所设定的W矩阵中每一个数当做目标收益率循环计算，得到有效前沿上的点。

并将每一次计算中的风险值和投资比例分别赋值给Sd和W。

其代码如下：
# 利用二次规划(QP)为各个目标回报率计算最优投资组合
for (i in 1:length(Mu))
{
B <- c(1, Mu[i], 0.7, rep(0,24)) # 约束向量
result <- solve.QP(Dmat=2*cov_mat,dvec=rep(0,24),Amat=A,bvec=B,meq=3) Sd[i] <- sqrt(result$value) # 最小风险
W[i,] <- result$solution # 投资比例
}
4. 通过有效前沿、最大夏普比率计算出最优投资策略
本实验设定无风险年利率为3%。

在图中标出这一点，代码如下：
mufree = 0.03
points(0,mufree,cex=3,pch=20)
图11绘制无风险年利率点# 计算夏普比率
sharpe = (Mu-mufree)/Sd
# 计算夏普比率的最大值
sharpe.max = (sharpe == max(sharpe))
#计算夏普比率最大时各资产的投资比例
W.opt <- round(W[sharpe.max,],6)
图12 找出最大夏普比率值
round()函数的作用是对数据进行四舍五入，第二个参数表示四舍五入的精度。

W.opt即为所要求的最优权重比例构成的数组。

查看W.opt的值，如图7。

图13最优投资比例
对所得结果进行加总会发现，投资比例之和不为1，并且有很多0值，这是由于round()函数设定了保留六位有效数字，很多计算出的投资比例很小，当以该精度输出时显示为0值。

5. 画出资本市场线及有效性边界
我们可以利用R中的plot()函数画出有效前沿，其形状为一条双曲线，代码如下：
plot(Sd,Mu,type="l",xlim=c(0,1.2),ylim=c(0,.6),lty=3,lwd=3,ylab="Expected Return", xlab="Risk(Standard Deviation)")
在图中画出资本市场线（CML），如图14，并且标出风险资产的最优投资组合点，即资本市场线（CML）与有效前沿曲线的切点。

代码如下：#切线投资组合
lines(c(0,2),mufree+c(0,2)*(Mu[sharpe.max]-mufree)/Sd[sharpe.max],lwd=3,lty =2,col="darkgrey")
# 风险资产的最优投资组合
points(Sd[sharpe.max],Mu[sharpe.max],cex=3,pch=20,col="red")
将风险资产的最有投资组合点标注为红色。

如图15。

图14切线投资组合
图15风险资产的最优投资组合在图中画出最小方差投资组合点。

代码如下：mvp = (Sd == min(Sd))
points(Sd[mvp],Mu[mvp],cex=3, pch=20, col="blue")
将最小方差投资组合点标注为蓝色。

如图16。

图16计算最小方程投资组合
画出有效边界。

代码如下：
ef = (Mu > Mu[mvp])
lines(Sd[ef],Mu[ef],type="l",xlim=c(0,.5),ylim=c(0,.2),lwd=3, col="black")
如图17。

图17有效边界9月数据操作同3月。

有效边界如下图18：
图17有效边界
三、实验结论或总结（可加页）
通过这次实验，基于R得出组合的投资比例，最终得出的调整情况如下表。

通过这次实验中学会了使用R进行数据提取和线性规划求解，做出相应的投资策略，对有效边界理论有了新的认识。