CoCalc -- script.workshop2010.R

GitHub Repository: braverock/portfolioanalytics
Path: blob/master/sandbox/script.workshop2010.R
¹⁴³³ views
1
### For R/Finance workshop on Portfolio Analytics
2
# Chicago, 16 April 2010
3
# Peter Carl and Brian Peterson
4

5
### Load the necessary packages
6
# Include optimizer and multi-core packages
7
library(PortfolioAnalytics)
8
require(xts)
9
require(DEoptim)
10
require(doMC)
11
registerDoMC()
12
require(TTR)
13

14
### Load the data
15
# Monthly total returns of four asset-class indexes
16
data(indexes)
17
#only look at 2000 onward
18
indexes<-indexes["2000::"]
19

20
#'## Review the data
21
postscript(file="WeightsVsRiskReview.eps", height=6, width=5, paper="special", horizontal=FALSE, onefile=FALSE)
22
# Generate charts to show 12m sma returns and CVaR by asset
23
charts.BarVaR(indexes[,1:4], p=(1-1/12), clean='boudt', show.cleaned=TRUE, methods=c("ModifiedVaR","ModifiedES"), colorset=rainbow6equal, cex.axis=1)
24
# Weights v %Contrib to CVaR
25
table = as.matrix(ES(indexes[,1:4], weights=rep(1/4,4), portfolio_method="component", p=(1-1/12))$pct_contrib_MES)
26
table = cbind(rep(1/4,4),table)
27
colnames(table) = c("Weights", "%Contrib to CVaR")
28
plot(table, ylim=c(-0.1,1), xlim=c(0,1), col=1:4, main="Weight and Contribution to Risk")
29
text(table[,1],table[,2],rownames(table), pos=4, cex = 0.8, col=1:4)
30
abline(a=0,b=1, col="darkgray", lty="dotted")
31
dev.off()
32

33
### Create a benchmark using equal-weighted portfolio returns
34
# Rebalance an equal-weight portfolio quarterly
35
dates=c(as.Date("1999-12-31"),time(indexes[endpoints(indexes, on="quarters")]))
36
weights = xts(matrix(rep(1/4,length(dates)*4), ncol=4), order.by=dates)
37
colnames(weights)= colnames(indexes[,1:4])
38
EqWgt = Return.rebalancing(indexes[,1:4],weights)
39

40
# Chart EqWgt Results
41
postscript(file="EqWgtPlot1.eps", height=6, width=5, paper="special", horizontal=FALSE, onefile=FALSE)
42
charts.PerformanceSummary(EqWgt, main="Eq Wgt Portfolio", methods=c("ModifiedVaR", "ModifiedES"), p=(1-1/12), clean='boudt', show.cleaned=TRUE, gap=36, colorset=bluefocus, lwd=3)
43
dev.off()
44

45
### Chart the VaR Sensitivity
46
postscript(file="VaRSense.eps", height=6, width=6, paper="special", horizontal=FALSE, onefile=FALSE)
47
layout(matrix(c(1,2,3,4), nrow=2))
48
for(i in 1:4)
49
  chart.VaRSensitivity(indexes[,i], methods=c("ModifiedES", "HistoricalES", "GaussianES"), colorset=bluemono, lwd=c(3,2,2), clean="boudt", main=paste("VaR Sensitivity for", colnames(indexes)[i], sep=" "), cex=.8)
50
dev.off()
51

52
# chart the assets
53
postscript(file="assetReturns.eps", height=6, width=5, paper="special", horizontal=FALSE, onefile=FALSE)
54
charts.BarVaR(indexes[,1:4], methods=c("ModifiedVaR", "ModifiedES"), colorset=rep("black",4), clean="boudt", show.clean=TRUE)
55
dev.off()
56

57
#'# EXAMPLE 1: Constrained Mean-CVaR Portfolio
58
### Show an example of a constraint set
59
aConstraintObj <- constraint(assets = colnames(indexes[,1:4]),
60
min = .05, # minimum position weight
61
max = c(.85,.5,.5,.3), #1, # maximum position weight
62
min_sum=0.99, # minimum sum of weights must be equal to 1-ish
63
max_sum=1.01, # maximum sum must also be about 1
64
weight_seq = generatesequence()) # possible weights for random or brute force portfolios
65

66
### Add a return objective to a constraint
67
# Create a small weighted annualized trailing-period mean wrapper function
68
pamean <- function(n=12, R, weights, geometric=TRUE)
69
{ sum(Return.annualized(last(R,n), geometric=geometric)*weights) }
70

71
# Portfolio annualized exponential moving average monthly return
72
aConstraintObj <- add.objective(constraints=aConstraintObj,
73
type="return",
74
name="pamean",
75
enabled=TRUE,
76
multiplier=-1, # to maximize this objective using a minimizer
77
arguments = list(n=12))
78

79
### Add a risk objective to a constraint
80
aConstraintObj <- add.objective(aConstraintObj,
81
type="risk", # the kind of objective this is
82
name="CVaR", # the function to minimize
83
enabled=TRUE, # enable or disable the objective
84
arguments=list(p=(1-1/12), clean="boudt")
85
)
86

87
### Use the Random Portfolios engine
88
# Evaluate the constraint object with Random Portfolios
89
rndResult<-optimize.portfolio(R=indexes[,1:4],
90
constraints=aConstraintObj,
91
optimize_method='random',
92
search_size=1000, trace=TRUE, verbose=TRUE)
93

94
# Chart the results
95
postscript(file="rpPerformancePlot1.eps", height=6, width=5, paper="special", horizontal=FALSE, onefile=FALSE)
96
charts.RP(rndResult, risk.col="CVaR", return.col="pamean", main="Constrained Mean-CVaR", neighbors=25)
97
dev.off()
98

99
# Chart the weights versus contribution to risk
100
postscript(file="WeightsVsRisk1.eps", height=6, width=5, paper="special", horizontal=FALSE, onefile=FALSE)
101
table = ES(indexes[,1:4],weights=rndResult$weights, portfolio_method="component", p=(1-1/12))$pct_contrib_MES
102
table = cbind(rndResult$weights,table)
103
colnames(table) = c("Weights", "%Contrib to CVaR")
104
plot(table, ylim=c(0,1), xlim=c(0,1), col=rainbow6equal, main="Weight and Contribution to Risk")
105
text(table[,1],table[,2],rownames(table), pos=4, cex = 1, col=rainbow6equal)
106
abline(a=0,b=1, col="darkgray", lty="dotted")
107
dev.off()
108

109
## Evaluate Constrained Mean-CVaR through time
110
#on one line for easy cut/paste/editing in interactive mode
111
registerDoMC()
112
rndResults<-optimize.portfolio.rebalancing(R=indexes[,1:4], constraints=aConstraintObj, optimize_method="random", trace=TRUE, rebalance_on='quarters', trailing_periods=NULL, training_period=36, search_size=1000)
113

114
# Chart the cumulative returns
115
weights1=extractWeights.rebal(rndResults)
116
Ex1=Return.rebalancing(indexes[,1:4], weights1)
117
index(weights1)<-as.Date(index(weights1))+1
118
colnames(Ex1)="Mean CVaR"
119
results = cbind(EqWgt,Ex1)
120
postscript(file="ReturnsEx1.eps", height=6, width=5, paper="special", horizontal=FALSE, onefile=FALSE)
121
charts.PerformanceSummary(results[,2:1], main="Constrained Mean-CVaR", colorset=bluefocus[-3], lwd=c(3,2), method=c("ModifiedVaR","ModifiedES"), p=(1-1/12), gap=36)
122
dev.off()
123

124
# Chart the weights and contribution to risk through time
125
# First, extract weights and calculate ES
126
numColumns = length(rndResults[[1]]$weights)
127
numRows = length(rndResults)
128
contrib1 <- matrix(nrow=numRows, ncol=numColumns)
129
for(i in 1:numRows){
130
    todate=paste("::",as.Date(names(rndResults)[i]), sep="")
131
    contrib1[i,]=ES(indexes[todate,1:4], portfolio_method="component", weights=rndResults[[i]]$weights, clean="boudt", p=(1-1/12))$pct_contrib_MES
132
}
133
colnames(contrib1) = names(unlist(rndResults[[1]]$weights))
134
contrib1<-xts(contrib1,order.by=index(weights1))
135

136
# test
137
postscript(file="InSampleStats.eps", height=6, width=5, paper="special", horizontal=FALSE, onefile=FALSE)
138

139
  retrisk1 <- matrix(nrow=numRows, ncol=2)
140
  for(i in 1:numRows)
141
    retrisk1[i,] = unlist(rndResults[[i]]$objective_measures)
142
  rownames(retrisk1) = names(rndResults)
143
  colnames(retrisk1) = c("pamean", "CVaR")
144
 chart.TimeSeries(retrisk1, legend="topright", main="In Sample Estimates")
145
dev.off()
146
# end test
147

148
postscript(file="WeightsContribEx1.eps", height=6, width=5, paper="special", horizontal=FALSE, onefile=FALSE)
149
layout(rbind(1, 2, 3), height = c(3, 3, 1.2), width = 1)
150
par(mar = c(2, 4, 4, 2) + 0.1)
151
PerformanceAnalytics:::chart.StackedBar.xts(weights1, main="Mean-CVaR Weights", legend.loc=NULL, cex.axis=1, colorset=bluemono[c(-2,-4,-6)], space=0, border="darkgray")
152
par(mar = c(2, 4, 4, 2) + 0.1)
153
PerformanceAnalytics:::chart.StackedBar.xts(contrib1, main="Mean-CVaR Risk Contribution", legend.loc=NULL, cex.axis=1, colorset=bluemono[c(-2,-4,-6)], space=0, border="darkgray")
154
plot.new()
155
par(oma = c(0, 0, 0, 0), mar=c(2,4,4,2)) #c(bottom, left, top, right)
156
legend("top", legend = colnames(weights1), fill = bluemono[c(-2,-4,-6)], ncol = 4, box.col="darkgray", border.col="white", cex=1)
157
dev.off()
158

159
#'# EXAMPLE 2: Mean-CVaR Risk Limit Portfolio
160
### Add a risk contribution constraint
161
# No more than 40% of the risk may be contributed by any one asset
162
# Reset the position constraints
163
aConstraintObj$max <-rep(1,4)
164
names(aConstraintObj$max) <- names(aConstraintObj$min)
165

166
aConstraintObj <- add.objective(aConstraintObj, type="risk_budget", name="CVaR",  enabled=TRUE, min_prisk=-Inf, max_prisk=.4, arguments = list(clean='boudt', method="modified",p=(1-1/12)))
167

168
rndResult2<-optimize.portfolio(R=indexes[,1:4],
169
constraints=aConstraintObj,
170
optimize_method='random',
171
search_size=1000, trace=TRUE, verbose=TRUE)
172

173
### Chart the results
174
postscript(file="rpPerformancePlot2.eps", height=6, width=5, paper="special", horizontal=FALSE, onefile=FALSE)
175
charts.RP(rndResult2, risk.col="CVaR", return.col="pamean", main="Mean-CVaR With Risk Limits", neighbors=25)
176
dev.off()
177

178
# Chart the weights versus contribution to risk
179
postscript(file="WeightsVsRisk2.eps", height=6, width=5, paper="special", horizontal=FALSE, onefile=FALSE)
180
table = cbind(rndResult2$weights, rndResult2$objective_measures$CVaR$pct_contrib_MES)
181
colnames(table) = c("Weights", "%Contrib to CVaR")
182
plot(table, ylim=c(0,1), xlim=c(0,1), col=rainbow6equal, main="Weight and Contribution to Risk")
183
text(table[,1],table[,2],rownames(table), pos=4, cex = 1, col=rainbow6equal)
184
abline(a=0,b=1, col="darkgray", lty="dotted")
185
dev.off()
186

187
### Evaluate Mean-CVaR Risk Limit through time
188
#on one line for easy cut/paste/editing in interactive mode
189
registerDoMC()
190
rndResults2<-optimize.portfolio.rebalancing(R=indexes[,1:4], constraints=aConstraintObj, optimize_method="random", trace=TRUE, rebalance_on='quarters', trailing_periods=NULL, training_period=36, search_size=3000)
191

192
### Chart the results
193
# Chart the cumulative returns
194
weights2=extractWeights.rebal(rndResults2)
195
Ex2=Return.rebalancing(indexes[,1:4], weights2)
196
index(weights2)<-as.Date(index(weights2)) + 1
197
colnames(Ex2) = "Risk Limit"
198
results = cbind(results, Ex2) 
199
postscript(file="ReturnsEx2.eps", height=6, width=5, paper="special", horizontal=FALSE, onefile=FALSE)
200
charts.PerformanceSummary(results[,3:1], main="Mean-CVaR Risk Limit", colorset=bluefocus[-3], lwd=c(3,2,2), method=c("ModifiedVaR","ModifiedES"), p=(1-1/12), gap=36)
201
dev.off()
202

203
# Chart the weights and contribution to risk through time
204
numColumns = length(rndResults2[[1]]$weights)
205
numRows = length(rndResults2)
206
contrib2 <- matrix(nrow=numRows, ncol=numColumns)
207
for(i in 1:numRows)
208
  contrib2[i,] = unlist(rndResults2[[i]]$objective_measures$CVaR$pct_contrib_MES)
209
colnames(contrib2) = names(unlist(rndResults2[[1]]$objective_measures$CVaR$pct_contrib_MES))
210
contrib2<-xts(contrib2,order.by=index(weights1))
211

212
# op <- par(no.readonly = TRUE)
213
postscript(file="WeightsContribEx2.eps", height=6, width=5, paper="special", horizontal=FALSE, onefile=FALSE)
214
layout(rbind(1, 2, 3), height = c(3, 3, 1.2), width = 1)
215
par(mar = c(2, 4, 4, 2) + 0.1)
216
PerformanceAnalytics:::chart.StackedBar.xts(weights2, main="Risk Limit Weights", legend.loc=NULL, cex.axis=1, colorset=bluemono[c(-2,-4,-6)], space=0, border="darkgray")
217
par(mar = c(2, 4, 4, 2) + 0.1)
218
PerformanceAnalytics:::chart.StackedBar.xts(contrib2, main="Risk Limit Risk Contribution", legend.loc=NULL, cex.axis=1, colorset=bluemono[c(-2,-4,-6)], space=0, border="darkgray")
219
plot.new()
220
par(oma = c(0, 0, 0, 0), mar=c(2,4,4,2))
221
legend("top", legend = colnames(weights2), fill = bluemono[c(-2,-4,-6)], ncol = 4, box.col="darkgray", border.col="white", cex=1)
222
dev.off()
223
# par(op)
224

225
#'# EXAMPLE 3: Equal Risk Portfolio
226
### Constraints for an Equal risk contribution portfolio
227
EqRiskConstr <- constraint(assets = colnames(indexes[,1:4]), min = 0.05, max = c(0.85,0.5,0.5,0.3), min_sum=1, max_sum=1, weight_seq = generatesequence())
228
EqRiskConstr <- add.objective(EqRiskConstr, type="risk_budget", name="CVaR",  enabled=TRUE, min_concentration=TRUE, arguments = list(clean='boudt', p=(1-1/12)))
229
EqRiskConstr <- add.objective(constraints=EqRiskConstr, type="return", name="pamean", enabled=TRUE, multiplier=0, arguments = list(n=12))
230

231
### Use DEoptim engine
232
EqRiskResultDE<-optimize.portfolio(R=indexes[,1:4], constraints=EqRiskConstr, optimize_method='DEoptim', search_size=2000, trace=TRUE, verbose=FALSE) #itermax=55, CR=0.99, F=0.5,
233

234
postscript(file="EqRiskDE.eps", height=6, width=5, paper="special", horizontal=FALSE, onefile=FALSE)
235
charts.DE(EqRiskResultDE, return.col="pamean", risk.col="CVaR")
236
dev.off()
237

238
### Evaluate through time
239
EqRiskResultDERebal<-optimize.portfolio.rebalancing(R=indexes[,1:4],
240
constraints=EqRiskConstr, optimize_method="DEoptim", trace=FALSE, rebalance_on='quarters', trailing_periods=NULL, training_period=36, itermax=75, CR=0.99, F=0.5, search_size=3000)
241

242
### Chart results
243
# Panel 1: Equal Risk Performance Summary
244
EqRiskWeights=extractWeights.rebal(EqRiskResultDERebal)
245
EqRisk=Return.rebalancing(indexes, EqRiskWeights)
246
index(EqRiskWeights)<-as.Date(index(EqRiskWeights)) + 1
247
R=cbind(EqRisk,EqWgt)
248
colnames(R)=c("Equal Risk","Equal Weight")
249
postscript(file="EqRiskPerf.eps", height=6, width=5, paper="special", horizontal=FALSE, onefile=FALSE)
250
charts.PerformanceSummary(R, methods=c("ModifiedVaR", "HistoricalVaR"), p=(1-1/12), colorset=bluefocus)
251
dev.off()
252

253
# Panel 2: Equal Risk Allocations
254
postscript(file="EqRiskBars.eps", height=6, width=5, paper="special", horizontal=FALSE, onefile=FALSE)
255
layout(rbind(1, 2, 3), height = c(3, 3, 1.2), width = 1)
256
par(mar = c(2, 4, 4, 2) + 0.1)
257
PerformanceAnalytics:::chart.StackedBar.xts(EqRiskWeights, main="Equal Risk Weights", legend.loc=NULL, cex.axis=1, colorset=bluemono[c(-2,-4,-6)], space=0, border="darkgray")
258
par(mar = c(2, 4, 4, 2) + 0.1)
259
### @TODO: Make this an extract function or calculate in the optim results
260
EqRiskPercContrCVaR=matrix(nrow=nrow(EqRiskWeights), ncol=ncol(EqRiskWeights))
261
for(i in 1:nrow(EqRiskWeights)){
262
    dates = paste(index(indexes)[1], index(EqRiskWeights)[i], sep="::")
263
    EqRiskPercContrCVaR[i,] = ES(indexes[dates,1:4], weights=EqRiskWeights[i,], p=(1-1/12), portfolio_method="component")$pct_contrib_MES
264
}
265
colnames(EqRiskPercContrCVaR) = names(EqRiskResultDERebal[[1]]$weights)
266
EqRiskPercContrCVaR<-xts(EqRiskPercContrCVaR,order.by=index(weights1))
267
par(mar = c(2, 4, 4, 2) + 0.1)
268
PerformanceAnalytics:::chart.StackedBar.xts(EqRiskPercContrCVaR, main="Equal Risk Risk Contribution", legend.loc=NULL, cex.axis=1, colorset=bluemono[c(-2,-4,-6)], space=0, border="darkgray")
269
plot.new()
270
par(oma = c(0, 0, 0, 0), mar=c(2,4,4,2))
271
legend("top", legend = colnames(weights1), fill = bluemono[c(-2,-4,-6)], ncol = 4, box.col="darkgray", border.col="white", cex=1)
272
dev.off()
273

274
#'## extended Equal Risk example
275
#' shart with the prior example
276
CDDConstr<-EqRiskConstr
277
#' turn back on the return objective
278
CDDConstr$objectives[[2]]$multiplier = -1
279
#' turn off risk_budget objective
280
CDDConstr$objectives[[1]]$multiplier = 0
281
#' add CDD objective
282
CDDConstr <- add.objective(CDDConstr, type="risk", name="CDD",  enabled=TRUE, arguments = list(p=(1-1/12)))
283

284
### Use DEoptim engine
285
CDDResultDE2<-optimize.portfolio(R=indexes[,1:4], constraints=CDDConstr, optimize_method='DEoptim', search_size=3000, trace=TRUE, verbose=FALSE, itermax=75, CR=0.99, F=0.5)
286
postscript(file="CDDDE2.eps", height=6, width=5, paper="special", horizontal=FALSE, onefile=FALSE)
287
charts.DE(CDDResultDE2, return.col="pamean", risk.col="CDD")
288
dev.off()
289

290
### Evaluate through time
291
CDDResultDERebal2<-optimize.portfolio.rebalancing(R=indexes[,1:4],
292
        constraints=CDDConstr, optimize_method="DEoptim", trace=FALSE, rebalance_on='quarters', trailing_periods=NULL, training_period=36, itermax=75, CR=0.99, F=0.5, search_size=3000)
293

294
### Chart results
295
# Panel 1: CDD Performance Summary
296
CDDweights=extractWeights.rebal(CDDResultDERebal2)
297
CDDRet=Return.rebalancing(indexes, CDDweights)
298
index(CDDweights)<-as.Date(index(CDDweights)) + 1
299
R=cbind(CDDRet,EqWgt)
300
colnames(R)=c("CDD","Equal Weight")
301
postscript(file="CDDPerf.eps", height=6, width=5, paper="special", horizontal=FALSE, onefile=FALSE)
302
charts.PerformanceSummary(R, methods=c("ModifiedVaR", "HistoricalVaR"), p=(1-1/12), colorset=bluefocus)
303
dev.off()
304

305
# Panel 2: CDD Allocations
306
postscript(file="CDDBars2.eps", height=6, width=5, paper="special", horizontal=FALSE, onefile=FALSE)
307
layout(rbind(1, 2, 3), height = c(3, 3, 1.2), width = 1)
308
par(mar = c(2, 4, 4, 2) + 0.1)
309
PerformanceAnalytics:::chart.StackedBar.xts(CDDweights, main="CDD w/ Return Obj. Weights", legend.loc=NULL, cex.axis=1, colorset=bluemono[c(-2,-4,-6)], space=0, border="darkgray")
310
par(mar = c(2, 4, 4, 2) + 0.1)
311
### @TODO: Make this an extract function or calculate in the optim results
312
CDDPercContrCVaR2=matrix(nrow=nrow(CDDweights), ncol=ncol(CDDweights))
313
for(i in 1:nrow(CDDweights)){
314
    dates = paste(index(indexes)[1], index(CDDweights)[i], sep="::")
315
    CDDPercContrCVaR2[i,] = ES(indexes[dates,1:4], weights=CDDweights[i,], p=(1-1/12), portfolio_method="component")$pct_contrib_MES
316
}
317
colnames(CDDPercContrCVaR2) = names(unlist(CDDResultDERebal2[[1]]$weights))
318
CDDPercContrCVaR2<-xts(CDDPercContrCVaR2,order.by=index(weights1))
319
par(mar = c(2, 4, 4, 2) + 0.1)
320
PerformanceAnalytics:::chart.StackedBar.xts(CDDPercContrCVaR2, main="CDD w/ Ret. Obj. Risk Contribution", legend.loc=NULL, cex.axis=1, colorset=bluemono[c(-2,-4,-6)], space=0, border="darkgray")
321
plot.new()
322
par(oma = c(0, 0, 0, 0), mar=c(2,4,4,2))
323
legend("top", legend = colnames(weights1), fill = bluemono[c(-2,-4,-6)], ncol = 4, box.col="darkgray", border.col="white", cex=1)
324
dev.off()
325

326

327

328
#'## APPENDIX EXAMPLES? ###
329
stopifnot(isTRUE(extended_ex),"halting")
330

331
## Markowitz-like constrained mean-variance
332
MeanVarConstr <- constraint(assets = colnames(indexes[,1:4]),
333
min = 0, max = 1, min_sum=0.99, max_sum=1.01, 
334
weight_seq = generatesequence())
335

336
MeanVarConstr <- add.objective(constraints=MeanVarConstr, type="return", name="mean", enabled=TRUE, multiplier=-1, arguments = list())
337

338
MeanVarConstr <- add.objective(MeanVarConstr,
339
type="risk", # the kind of objective this is
340
name="sd",
341
enabled=TRUE, # enable or disable the objective
342
arguments= list())
343

344
MeanVarResultRP<-optimize.portfolio(R=indexes[,1:4],
345
constraints=MeanVarConstr,
346
optimize_method='random',
347
search_size=1000,
348
trace=TRUE, verbose=TRUE)
349
# Chart the results
350
postscript(file="mvRP.eps", height=6, width=5, paper="special", horizontal=FALSE, onefile=FALSE)
351
charts.RP(MeanVarResultRP, risk.col="sd", return.col="mean", main="Mean Variance", neighbors=25)
352
dev.off()
353

354
MeanVarResultRPRebal<-optimize.portfolio.rebalancing(R=indexes[,1:4],
355
constraints=MeanVarConstr, 
356
optimize_method="random",
357
trace=FALSE, 
358
rebalance_on='quarters', 
359
trailing_periods=NULL, 
360
training_period=36, 
361
search_size=1000)
362

363
## Risk budget
364
RiskBudget.constr <- constraint(assets = colnames(indexes[,1:4]), min = 0, max = 1, min_sum=1, max_sum=1, weight_seq = generatesequence())
365

366
RiskBudget.constr <- add.objective(RiskBudget.constr, type="risk_budget", name="CVaR",  enabled=TRUE, min_prisk=-Inf, max_prisk=.4, arguments = list(clean='boudt', method="modified",p=.95))
367

368
RiskBudget.constr <- add.objective(constraints=RiskBudget.constr, type="risk", name="CVaR", enabled=TRUE, arguments = list(method="modified", portfolio_method="single", enabled=TRUE, p=.95, clean="boudt"))
369

370
RiskBudget.optimResult<-optimize.portfolio(R=indexes[,1:4], constraints=RiskBudget.constr,
371
optimize_method='DEoptim', search_size=1000, verbose=TRUE)
372

373
## Equal risk
374
# First method: minimize risk concentration across assets
375
EqRiskConstr1 <- constraint(assets = colnames(indexes[,1:4]), min = 0.05, max = c(0.85,0.5,0.5,0.3), min_sum=0.99, max_sum=1.01, weight_seq = generatesequence())
376
EqRiskConstr1 <- add.objective(EqRiskConstr1, type="risk_budget", name="CVaR",  enabled=TRUE, min_concentration=TRUE, arguments = list(clean='boudt', p=(1-1/12)))
377
EqRiskResultDE1<-optimize.portfolio(R=indexes[,1:4], constraints=EqRiskConstr1, optimize_method='DEoptim', search_size=2000, trace=TRUE, verbose=FALSE) #itermax=55, CR=0.99, F=0.5,
378

379
# Second method: force risk levels between limits
380
EqRiskConstr <- constraint(assets = colnames(indexes[,1:4]), min = 0.05, max = c(0.85,0.5,0.5,0.3), min_sum=1, max_sum=1, weight_seq = generatesequence())
381
EqRiskConstr <- add.objective(EqRiskConstr, type="risk_budget", name="CVaR", enabled=TRUE, arguments = list(p=(1-1/12), clean="boudt"), min_prisk=0.24, max_prisk=0.26)
382
EqRiskResultDE<-optimize.portfolio(R=indexes[,1:4], constraints=EqRiskConstr, optimize_method='DEoptim', search_size=2000, trace=TRUE, verbose=FALSE)
383

384
#EqRiskResultRP<-optimize.portfolio(R=indexes[,1:4], constraints=EqRiskConstr, optimize_method='random', search_size=1000, trace=TRUE, verbose=TRUE)
385
# Chart the results
386
# postscript(file="EqRiskRP.eps", height=6, width=5, paper="special", horizontal=FALSE, onefile=FALSE)
387
# charts.RP(EqRiskResultRP, risk.col="CVaR", return.col="mean", main="Mean Variance", neighbors=25)
388
# dev.off()
389
EqRiskResultDERebal<-optimize.portfolio.rebalancing(R=indexes[,1:4],
390
constraints=EqRiskConstr, optimize_method="DEoptim", trace=FALSE, rebalance_on='quarters', trailing_periods=NULL, training_period=36, itermax=45, CR=0.99, F=0.5, search_size=1000)
391

392
# Panel 1: Equal Risk Performance Summary
393
EqRiskWeights=extractWeights.rebal(EqRiskResultDERebal)
394
EqRisk=Return.rebalancing(indexes, EqRiskWeights)
395
index(EqRiskWeights)<-as.Date(index(EqRiskWeights)) + 1
396
R=cbind(EqRisk,EqWgt)
397
colnames(R)=c("Equal Risk","Equal Weight")
398
charts.PerformanceSummary(R, methods=c("ModifiedVaR", "HistoricalVaR"), p=(1-1/12), colorset=redfocus)
399

400
# Panel 2: Equal Risk Allocations
401
chart.StackedBar(EqRiskWeights)
402
## @TODO: Make this an extract function or calculate in the optim results
403
EqRiskPercContrCVaR=matrix(nrow=nrow(EqRiskWeights), ncol=ncol(EqRiskWeights))
404
for(i in 1:nrow(EqRiskWeights)){
405
    dates = paste(index(indexes)[1], index(EqRiskWeights)[i], sep="::")
406
    EqRiskPercContrCVaR[i,] = ES(indexes[dates,1:4], weights=EqRiskWeights[i,], p=(1-1/12), portfolio_method="component")$pct_contrib_MES
407
}
408
colnames(EqRiskPercContrCVaR) = names(unlist(EqRiskResultDERebal[[1]]$weights))
409
rownames(EqRiskPercContrCVaR) = colnames(weights1)
410
chart.StackedBar(EqRiskPercContrCVaR)
411

412
## Return target with risk budget
413

414
## You want to do what?
415

416

417
### ---------------------- Scratch area ---------------------- ###
418
# Single period DEOptim results
419
assets.pema=matrix(nrow=1,ncol=4)
420
# for(i in 1:4){ assets.pema[,i] = paEMA(n=12,indexes[,i],1) }
421
for(i in 1:4){ assets.pema[,i] = Return.annualized(last(indexes[,i],12),geometric=TRUE) }
422
colnames(assets.pema)=colnames(indexes[,1:4])
423
rownames(assets.pema)="Tr 12m Ann Return"
424
assets.CVaR=ES(indexes[,1:4], invert=FALSE)
425
assets = rbind(assets.CVaR,assets.pema)
426
pdf()
427
plot(t(assets))
428
points(optimResult$objective_measures[1], optimResult$objective_measures[2])
429
text(t(assets), colnames(assets))
430
text(optimResult$objective_measures[1], optimResult$objective_measures[2], "Optimal")
431
dev.off()
432
paEMA <- function(n=10, R, weights, ...)
433
{# call Exponential Moving Average from TTR, return the last observation
434
    sum((12*last(apply(R,2,FUN=TTR::EMA,n=n)))*weights)
435
}
436
#     xtract = extractStats.rp(rndResult)
437
#     plot(xtract[,6],xtract[,7], xlab="CVaR", ylab="pEMA", col="lightgray", main="Random Portfolios")
438
#     points(xtract[1,6],xtract[1,7], col="orange", pch=16) # equal weighted (seed)
439
#     points(rndResult$constrained_objective[[1]], rndResult$constrained_objective[[2]], col="red", pch=16)
440

441
## Use DEoptim as an engine
442
optimResult<-optimize.portfolio(R=indexes[,1:4],
443
constraints=aConstraintObj,
444
optimize_method='DEoptim',
445
itermax=45, CR=0.99, F=0.5,
446
search_size=2000 #, verbose=TRUE
447
)
448

449
### Run it several times
450
optimResultList <- foreach(ii=iter(1:20),# find 20 sol'ns
451
.errorhandling='pass') %dopar%
452
optimize.portfolio(R=indexes[,1:4], aConstraintObj,
453
optimize_method='DEoptim', trace=TRUE,
454
itermax=25, CR=0.99, F=0.5, search_size=1000)
455

456
rndResultList <- foreach(ii=iter(1:20),# 120 sol'ns
457
.errorhandling='pass') %dopar%
458
optimize.portfolio(R=indexes[,1:4], aConstraintObj,
459
optimize_method='random', trace=TRUE,
460
search_size=1000, verbose=TRUE)
461

462
### Evaluate through time
463
#on one line for easy cut/paste/editing in interactive mode
464
registerDoMC()
465
rndResults<-optimize.portfolio.rebalancing(R=indexes[,1:4], constraints=aConstraintObj, optimize_method="random", trace=TRUE, rebalance_on='quarters', trailing_periods=NULL, training_period=36, search_size=1000)
466

467
# on multiple, for demo
468
rndResults<-optimize.portfolio.rebalancing(R=indexes[,1:4],
469
constraints=aConstraintObj, # our constraints object
470
optimize_method="random", # allows kitchen sink sol'ns
471
trace=FALSE, # set verbosity for tracking
472
rebalance_on='quarters', # any xts 'endpoints'
473
trailing_periods=NULL, # calculation from inception
474
training_period=36, # starting period for calculation
475
search_size=1000) # how many portfolios to test
476

477
Product

Resources

Company
