CoCalc -- optimize.HFindexes.R

GitHub Repository: braverock/portfolioanalytics
Path: blob/master/sandbox/symposium2013/optimize.HFindexes.R
¹⁴³³ views
1
### For Presentation at FactSet's 2013 US Investment Process Symposium
2
# November 10 - 12 , 2013
3
# Peter Carl
4

5
### Load the necessary packages
6
# Include optimizer packages
7
require(PortfolioAnalytics)
8
require(quantmod)
9
library(DEoptim)
10
library(ROI)
11
require(ROI.plugin.quadprog)
12
require(ROI.plugin.glpk)
13
# ... and multi-core packages
14
require(foreach)
15
require(doMC)
16
registerDoMC(6)
17

18
# Available on r-forge
19
# require(FactorAnalytics) # development version > build 
20

21
### Set script constants
22
runname='historical.moments'
23
rebalance_period = 'years' #'quarters' # uses endpoints identifiers from xts; how to do semi-annual?
24
clean = "boudt" # "none" 
25
permutations = 2000
26
p=1-1/12 # set confidence for VaR/mETL for monthly data
27

28
### Description
29

30
# Seven assets, in this case hedge fund indexes representing different styles and two 'types':
31
# Convertible Arbitrage   | Non-directional
32
# Equity Market Neutral   | Non-directional
33
# Fixed Income Arbitrage  | Non-directional
34
# Event Driven            | Non-directional
35
# CTA Global              | Directional
36
# Global Macro            | Directional
37
# Long/Short Equity       | Directional
38

39
# see analyze.HFindexes.R for more detail
40

41
# Set up seven objectives as buoy portfolios:
42
# - Equal contribution to...
43
#   1 Weight
44
#   2 Variance
45
#   3 Risk (mETL)
46
# - Reward to Risk ratio of...
47
#   4 Mean-Variance
48
#   5 Mean-mETL
49
# - Minimum...
50
#   6 Variance
51
#   7 Risk (mETL)
52

53
# Add constraints
54
# - Box constraints - 5% to 30%
55
# - Group constraints - Non-directional constrained to 20-70%; Directional between 10-50%
56
# - Rebalancing period - quarterly
57
# - Turnover constraints - TBD
58

59
#------------------------------------------------------------------------
60
# Set up an initial portfolio object with constraints and objectives using
61
# v2 specification
62

63
# Create initial portfolio object used to initialize ALL the bouy portfolios
64
init.portf <- portfolio.spec(assets=colnames(R), 
65
                             weight_seq=generatesequence(by=0.005)
66
)
67
# Add leverage constraint
68
init.portf <- add.constraint(portfolio=init.portf, 
69
  type="leverage"
70
)
71
# Add box constraint
72
init.portf <- add.constraint(portfolio=init.portf, 
73
                             type="box", 
74
                             min=0.05, 
75
                             max=0.3
76
)
77
# Add group constraint
78
init.portf <- add.constraint(portfolio=init.portf, type="group",
79
                  groups=list(c(1:4),
80
                              c(5:7)),
81
                  group_min=c(0.25,.05),
82
                  group_max=c(0.85,0.75)
83
)
84

85
# print(init.portf)
86
# summary(init.portf)
87

88
#------------------------------------------------------------------------
89
### Construct BUOY 1: Constrained Mean-StdDev Portfolio - using ROI
90
# Add the return and sd objectives to the constraints created above
91
MeanSD.portf <- add.objective(portfolio=init.portf,
92
  type="return", # the kind of objective this is
93
  name="mean" # name of the function
94
)
95
MeanSD.portf <- add.objective(portfolio=MeanSD.portf,
96
  type="risk", # the kind of objective this is
97
  name="var", # name of the function
98
  arguments=list(clean=clean)
99
)
100

101
### Construct BUOY 2: Constrained Mean-mETL Portfolio - using ROI
102
#@ Cannot maximize mean return per unit ETL with ROI, consider using
103
#@ random portfolios or DEoptim. - RB
104
# Add the return and mETL objectives
105
MeanmETL.portf <- add.objective(portfolio=init.portf,
106
  type="return", # the kind of objective this is
107
  name="mean" # name of the function
108
  , multiplier=-12
109
  )
110
MeanmETL.portf <- add.objective(portfolio=MeanmETL.portf,
111
  type="risk", # the kind of objective this is
112
  name="ES", # the function to minimize
113
  arguments=list(p=p, clean=clean)
114
  )
115

116
### Construct BUOY 3: Constrained Minimum Variance Portfolio - using ROI
117
# Add the variance objective
118
MinSD.portf <- add.objective(portfolio=init.portf,
119
  type="risk", # the kind of objective this is
120
  name="var", # name of the function
121
  arguments=list(p=p, clean=clean)
122
)
123

124
### Construct BUOY 4: Constrained Minimum mETL Portfolio - using ROI
125
# Add the mETL objective
126
MinmETL.portf <- add.objective(portfolio=init.portf,
127
  type="risk", # the kind of objective this is
128
  name="ES", # the function to minimize
129
  arguments=list(p=p, clean=clean)
130
)
131

132
### Construct BUOY 5: Constrained Equal Variance Contribution Portfolio - using RP
133
#@ - Add the sub-objectives first. Adding these 3 objectives means that we are
134
#@ maximizing mean per unit StdDev with equal volatility contribution portfolios. - RB
135
# Without a sub-objective, we get a somewhat undefined result, since there are (potentially) many Equal SD contribution portfolios.
136
# MRCSD.portf <- add.objective(portfolio=init.portf,
137
#                             type="risk",
138
#                             name="StdDev"
139
# ) # OR
140
# MRCSD.portf <- add.objective(portfolio=init.portf,
141
#   type="return",
142
#   name="mean"
143
# )
144
MRCSD.portf <- add.objective(portfolio=init.portf, 
145
  type="risk_budget", 
146
  name="StdDev",  
147
  min_concentration=TRUE,
148
  arguments = list(clean=clean)
149
)
150

151
# MRCSD.portf$constraints[[1]]$min_sum = 0.99 # set to speed up RP
152
# MRCSD.portf$constraints[[1]]$max_sum = 1.01
153

154
### Construct BUOY 6: Constrained Equal mETL Contribution Portfolio - using RP
155
#@ Add the sub-objectives first. These should be added to the MRCmETL portfolio.
156
#@ All objectives below mean that we are maximizing mean return per unit ES with
157
#@ equal ES contribution. - RB
158
# Without a sub-objective, we get a somewhat undefined result, since there are (potentially) many Equal SD contribution portfolios.
159
# MRCmETL.portf <- add.objective(portfolio=init.portf,
160
#                             type="risk",
161
#                             name="ES"
162
# ) # OR
163
MRCmETL.portf <- add.objective(portfolio=init.portf,
164
                            type="return",
165
                            name="mean"
166
)
167
MRCmETL.portf <- add.objective(MRCmETL.portf,
168
                              type="risk_budget",
169
                              name="ES",
170
                              min_concentration=TRUE,
171
                              arguments = list(p=p, clean=clean)
172
)
173
# Calculate portfolio variance, but don't use it in the objective; used only for plots
174
MRCmETL.portf <- add.objective(portfolio=MRCmETL.portf,
175
                                  type="risk", # the kind of objective this is
176
                                  name="StdDev", # the function to minimize
177
                                  enabled=TRUE, # enable or disable the objective
178
                                  multiplier=0, # calculate it but don't use it in the objective
179
                                  arguments=list(clean=clean)
180
)
181
# MRCmETL.portf$constraints[[1]]$min_sum = 0.99 # set to speed up RP
182
# MRCmETL.portf$constraints[[1]]$max_sum = 1.01
183

184
### Construct BUOY 7: Equal Weight Portfolio
185
# There's only one, so create a portfolio object with all the objectives we want calculated. 
186
EqWgt.portf <- portfolio.spec(assets=colnames(R))
187
EqWgt.portf <- add.constraint(portfolio=EqWgt.portf, type="leverage", min_sum=1, max_sum=1)
188
EqWgt.portf <- add.objective(portfolio=EqWgt.portf, type="return", name="mean")
189
EqWgt.portf <- add.objective(portfolio=EqWgt.portf, type="risk_budget", name="ES", arguments=list(p=p, clean=clean))
190
EqWgt.portf <- add.objective(portfolio=EqWgt.portf, type="risk_budget", name="StdDev", arguments=list(clean=clean))
191

192
### Construct BUOY 8: Inverse Volatility Portfolio
193
# There's only one, so create a portfolio object with all the objectives we want calculated. 
194
VolWgt.portf <- portfolio.spec(assets=colnames(R))
195
VolWgt.portf <- add.constraint(portfolio=VolWgt.portf, type="leverage", min_sum=0.99, max_sum=1.01)
196
VolWgt.portf <- add.objective(portfolio=VolWgt.portf, type="return", name="mean")
197
VolWgt.portf <- add.objective(portfolio=VolWgt.portf, type="risk_budget", name="ES", arguments=list(p=p, clean=clean))
198
VolWgt.portf <- add.objective(portfolio=VolWgt.portf, type="risk_budget", name="StdDev", arguments=list(clean=clean))
199

200
# REMOVED - to much to show already
201
# ### Construct RISK BUDGET Portfolio
202
# ConstrConcmETL.portf <- portfolio.spec(assets=colnames(R), 
203
#   weight_seq=generatesequence(by=0.005)
204
# )
205
# # Add leverage constraint
206
# ConstrConcmETL.portf <- add.constraint(portfolio=RiskBudget.portf, 
207
#   type="leverage", 
208
#   min_sum=0.99, # set to speed up RP, DE
209
#   max_sum=1.01
210
# )
211
# # Establish position bounds
212
# ConstrConcmETL.portf <- add.constraint(portfolio=ConstrConcmETL.portf, 
213
#   type="box", 
214
#   min=0.01, # leave relatively unconstrained
215
#   max=1.0
216
# )
217
# # Maximize mean return
218
# ConstrConcmETL.portf <- add.objective(portfolio=ConstrConcmETL.portf,
219
#   type="return", # maximize return
220
#   name="mean",
221
#   multiplier=12
222
# )
223
# # Add a risk measure
224
# # Use ES to be consistent with risk measures in other BUOY portfolios
225
# ConstrConcmETL.portf <- add.objective(portfolio=ConstrConcmETL.portf,
226
#   type="risk",
227
#   name="ETL", # using a different name to avoid clobbering slot below, workaround for bug
228
#   multiplier=1,
229
#   arguments = list(p=p, clean=clean)
230
# )
231
# 
232
# # Set contribution limits
233
# ConstrConcmETL.portf <- add.objective(portfolio=ConstrConcmETL.portf,
234
#   type="risk_budget",
235
#   name="ES",
236
#   max_prisk=0.3, # Sets the maximum percentage contribution to risk
237
#   arguments = list(p=p, clean=clean)
238
# )
239
# Calculate portfolio variance, but don't use it in the objective; used only for plots
240
# ConstrConcmETL.portf <- add.objective(portfolio=ConstrConcmETL.portf,
241
#   type="risk", # the kind of objective this is
242
#   name="StdDev", # the function to minimize
243
#   enabled=TRUE, # enable or disable the objective
244
#   multiplier=0, # calculate it but don't use it in the objective
245
#   arguments=list(clean=clean)
246
# )
247

248
#------------------------------------------------------------------------
249
### Evaluate portfolio objective objects
250
#------------------------------------------------------------------------
251
# Generate a single set of random portfolios to evaluate against all RP constraint sets
252
print(paste('constructing random portfolios at',Sys.time()))
253

254
# Modify the init.portf specification to get RP running 
255
rp.portf <- init.portf
256
# rp.portf$constraints[[1]]$min_sum = 0.99 # set to speed up RP
257
# rp.portf$constraints[[1]]$max_sum = 1.01
258
rp.portf$constraints[[1]]$min_sum = 1.00 # for more accuracy
259
rp.portf$constraints[[1]]$max_sum = 1.00
260
# rp = random_portfolios(portfolio=rp.portf, permutations=30000, max_permutations=400) # will get fewer with less accuracy
261
load(file=paste(resultsdir,'random-portfolios-2013-10-05.historical.moments.rda'))
262
rp.mean = apply(rp, 1, function(w) mean(R %*% w))
263
rp.sd = apply(rp, 1, function(x) StdDev(R=R, weights=x, p=p, clean=clean))
264
plot(rp.sd, rp.mean, col="darkgray", cex=0.5)
265

266
# REMOVED: This was fruitless
267
# rp1 = random_portfolios(portfolio=rp.portf, permutations=10000, max_permutations=400, rp_method="sample")
268
# rp1.mean = apply(rp1, 1, function(w) mean(R %*% w))
269
# rp1.sd = apply(rp1, 1, function(x) StdDev(R=R, weights=x, p=p))
270
# rp1.etl=NULL; for(i in 1:NROW(rp1)) {rp1.etl[i]=ETL(R=R, weights=as.vector(rp1[i,]), p=p, portfolio_method="component")[[1]]}
271
# plot(rp1.sd, rp1.mean, col="gray", cex=0.5)
272
# 
273
# rp2 = random_portfolios(portfolio=rp.portf, permutations=10000, max_permutations=400, rp_method="simplex", fev=2)
274
# rp2.mean = apply(rp2, 1, function(w) mean(R %*% w))
275
# rp2.sd = apply(rp2, 1, function(x) StdDev(R=R, weights=x, p=p))
276
# points(rp2.sd,rp2.mean, col="blue", cex=0.5)
277
# 
278
# rp3 = random_portfolios(portfolio=rp.portf, permutations=10000, max_permutations=400, rp_method="grid")
279
# rp3.mean = apply(rp3, 1, function(w) mean(R %*% w))
280
# rp3.sd = apply(rp3, 1, function(x) StdDev(R=R, weights=x, p=p))
281
# points(rp3.sd,rp3.mean, col="green", cex=0.5)
282

283
# print(paste('done constructing random portfolios at',Sys.time()))
284
# save(rp,file=paste(resultsdir, 'random-portfolios-', Sys.Date(), '-', runname, '.rda',sep=''))
285

286

287
start_time<-Sys.time()
288
print(paste('Starting optimization at',Sys.time()))
289

290
### Evaluate BUOY 1: Constrained Mean-StdDev Portfolio - with ROI
291
MeanSD.ROI<-optimize.portfolio(R=R,
292
  portfolio=MeanSD.portf,
293
  optimize_method='ROI',
294
  trace=TRUE
295
) 
296
plot(MeanSD.ROI, risk.col="StdDev", return.col="mean", rp=permutations, chart.assets=TRUE, main="Mean-Volatility Portfolio")
297
save(MeanSD.ROI,file=paste(resultsdir, 'MeanSD-', Sys.Date(), '-', runname, '.rda',sep='')) # Save the results
298
print(paste('Completed meanSD optimization at',Sys.time(),'moving on to meanmETL'))
299

300
### Evaluate BUOY 2: Constrained Mean-mETL Portfolio - with DE or RND
301
# Not possible with ROI - RB
302
# MeanmETL.ROI<-optimize.portfolio(R=R,
303
#   portfolio=MeanmETL.portf,
304
#   optimize_method='ROI',
305
#   trace=TRUE, verbose=TRUE
306
#   ) 
307
# 
308

309
# So use random portfolios instead
310
MeanmETL.RND<-optimize.portfolio(R=R,
311
                                 portfolio=MeanmETL.portf,
312
                                 optimize_method='random',
313
                                 rp=rp,
314
                                 trace=TRUE
315
) 
316
plot(MeanmETL.RND, risk.col="StdDev", return.col="mean", rp=permutations, chart.assets=TRUE, main="Mean-mETL Portfolio")
317
plot(MeanmETL.RND, risk.col="ES", return.col="mean", rp=permutations, chart.assets=TRUE, main="Mean-mETL Portfolio")
318
save(MeanmETL.RND,file=paste(resultsdir, 'MeanETL-RP-', Sys.Date(), '-', runname, '.rda',sep=''))
319
chart.EfficientFrontier(MeanmETL.RND)
320

321
# OR with DE optim
322
MeanmETL.DE<-optimize.portfolio(R=R,
323
  portfolio=MeanmETL.portf,
324
  optimize_method='DEoptim',
325
  search_size=20000,
326
  initialpop=rp[1:50,] # seed with a starting population that we know fits the constraint space
327
) 
328
plot(MeanmETL.DE, risk.col="StdDev", return.col="mean", chart.assets=TRUE, main="Mean-mETL Portfolio")
329
plot(MeanmETL.DE, risk.col="ES", return.col="mean", chart.assets=TRUE, main="Mean-mETL Portfolio")
330
save(MeanmETL.DE,file=paste(resultsdir, 'MeanETL-DE-', Sys.Date(), '-', runname, '.rda',sep=''))
331
print(paste('Completed meanmETL optimization at',Sys.time(),'moving on to MinSD'))
332

333

334
### Evaluate BUOY 3: Constrained Minimum Variance Portfolio - with ROI
335
MinSD.ROI<-optimize.portfolio(R=R,
336
  portfolio=MinSD.portf,
337
  optimize_method='ROI', 
338
  trace=TRUE
339
  ) # 
340
plot(MinSD.ROI, risk.col="StdDev", return.col="mean", rp=permutations, chart.assets=TRUE, main="Minimum Volatility Portfolio with ROI")
341
save(MinSD.ROI,file=paste(resultsdir, 'MinSD-', Sys.Date(), '-', runname, '.rda',sep=''))
342
print(paste('Completed MinSD optimization at',Sys.time(),'moving on to MinmETL'))
343
# OR with random portfolios
344
# MinSD.RND<-optimize.portfolio(R=R,
345
#   portfolio=MinSD.portf,
346
#   optimize_method='random',
347
#   rp=rp,
348
#   trace=TRUE, verbose=TRUE
349
#   ) # 
350
# plot(MinSD.RND, risk.col="StdDev", return.col="mean", chart.assets=TRUE, main="Minimum Volatility Portfolio with RP")
351

352
### Evaluate BUOY 4: Constrained Minimum mETL Portfolio - with ROI
353
MinmETL.ROI<-optimize.portfolio(R=R,
354
  portfolio=MinmETL.portf,
355
  optimize_method='ROI',
356
  trace=TRUE, verbose=TRUE
357
  ) 
358
plot(MinmETL.ROI, risk.col="StdDev", return.col="mean", rp=permutations, chart.assets=TRUE, main="Minimum mETL Portfolio")
359
plot(MinmETL.ROI, risk.col="ES", return.col="mean", rp=permutations, chart.assets=TRUE, main="Minimum mETL Portfolio")
360
save(MinmETL.ROI,file=paste(resultsdir, 'MinmETL-', Sys.Date(), '-', runname, '.rda',sep=''))
361
print(paste('Completed MinmETL optimization at',Sys.time(),'moving on to MRCSD'))
362

363
### Evaluate BUOY 5: Constrained Equal Variance Contribution Portfolio - with RP
364
# MRCSD.RND<-optimize.portfolio(R=R,
365
#   portfolio=MRCSD.portf,
366
#   optimize_method='random',
367
#   rp=rp,
368
#   trace=TRUE
369
#   ) 
370
# plot(MRCSD.RND, risk.col="StdDev", return.col="mean", chart.assets=TRUE, main="Equal Volatility Contribution Portfolio")
371
# chart.RiskBudget(MRCSD.RND, risk.type="percentage", neighbors=25)
372
# save(MRCSD.RND,file=paste(resultsdir, 'MRCSD.RND-', Sys.Date(), '-', runname, '.rda',sep=''))
373
# ... not a very satisfying solution
374

375
# OR DE optim - this gets very close (a nice, straight line), so use it
376
MRCSD.DE<-optimize.portfolio(R=R,
377
  portfolio=MRCSD.portf,
378
  optimize_method='DEoptim',
379
  search_size=20000,
380
  itermax=400,
381
  initialpop=rp[1:50,], # seed with a starting population that we know fits the constraint space
382
  trace=FALSE
383
  ) 
384
plot(MRCSD.DE, risk.col="StdDev", return.col="mean", chart.assets=TRUE, main="Equal Volatility Contribution Portfolio")
385
chart.RiskBudget(MRCSD.DE, risk.type="percentage", neighbors=25)
386
save(MRCSD.DE,file=paste(resultsdir, 'MRCSD.DE-', Sys.Date(), '-', runname, '.rda',sep=''))
387

388
print(paste('Completed MRCSD optimization at',Sys.time(),'moving on to MRCmETL'))
389

390
### Evaluate BUOY 6: Constrained Equal mETL Contribution Portfolio - with RP
391
MRCmETL.RND<-optimize.portfolio(R=R,
392
  portfolio=MRCmETL.portf,
393
  optimize_method='random',
394
  rp=rp,
395
  trace=TRUE
396
  ) # 
397
plot(MRCmETL.RND, risk.col="StdDev", return.col="mean", chart.assets=TRUE, main="Equal mETL Contribution Portfolio")
398
plot(MRCmETL.RND, risk.col="ES", return.col="mean", chart.assets=TRUE, main="Equal mETL Contribution Portfolio")
399
chart.RiskBudget(MRCmETL.RND, neighbors=25)
400
save(MRCmETL.RND,file=paste(resultsdir, 'MRCmETL-', Sys.Date(), '-', runname, '.rda',sep=''))
401

402
# OR DE optim - 
403
MRCmETL.DE<-optimize.portfolio(R=R,
404
  portfolio=MRCmETL.portf,
405
  optimize_method='DEoptim',
406
  search_size=20000,
407
  NP=200,
408
  initialpop=rp[1:50,], # seed with a starting population that we know fits the constraint space
409
  trace=FALSE
410
  ) 
411
plot(MRCmETL.DE, risk.col="StdDev", return.col="mean", chart.assets=TRUE, main="Equal Volatility Contribution Portfolio")
412
chart.RiskBudget(MRCmETL.DE, risk.type="percentage", neighbors=25)
413
save(MRCmETL.DE,file=paste(resultsdir, 'MRCmETL.DE-', Sys.Date(), '-', runname, '.rda',sep=''))
414

415
# # test it unconstrained:
416
# unconstr.portf <- portfolio.spec(assets=colnames(R), 
417
#                              weight_seq=generatesequence(by=0.005)
418
# )
419
# unconstr.portf <- add.constraint(portfolio=unconstr.portf, 
420
#                              type="leverage", 
421
#                              min_sum=0.99, # set to speed up RP
422
#                              max_sum=1.01
423
# )
424
# # Establish position bounds
425
# unconstr.portf <- add.constraint(portfolio=unconstr.portf, 
426
#                   type="box", 
427
#                   min=0.01, 
428
#                   max=1.0
429
# )
430
# MRCmETLun.portf <- add.objective(portfolio=unconstr.portf,
431
#                             type="return",
432
#                             name="mean"
433
# )
434
# MRCmETLun.portf <- add.objective(MRCmETL.portf,
435
#                               type="risk_budget",
436
#                               name="ES",
437
#                               min_concentration=TRUE,
438
#                               arguments = list(p=p, clean=clean)
439
# )
440
# 
441
# # ...in DE optim - 
442
# MRCmETLun.DE<-optimize.portfolio(R=R,
443
#   portfolio=MRCmETLun.portf,
444
#   optimize_method='DEoptim',
445
#   search_size=20000,
446
#   NP=200,
447
#   initialpop=rp[1:50,], # seed with a starting population that we know fits the constraint space
448
#   trace=FALSE
449
#   ) 
450

451
print(paste('Completed MRCmETL optimization at',Sys.time(),'moving on to RiskBudget'))
452

453
### Evaluate BUOY 7: Equal Weight Portfolio
454
# Calculate the objective measures for the equal weight portfolio
455
EqWgt.opt <- equal.weight(R=R, portfolio=EqWgt.portf)
456

457
### Evaluate BUOY 8: Inverse Volatility Portfolio
458
volatility.weight <- function (R, portfolio, ...) 
459
{
460
  if (!is.portfolio(portfolio)) 
461
      stop("portfolio object passed in must be of class 'portfolio'")
462
  assets <- portfolio$assets
463
  nassets <- length(assets)
464
  weights <- as.vector((1/StdDev(R))/sum(1/StdDev(R)))
465
  names(weights) <- names(assets)
466
  if (ncol(R) != nassets) {
467
      if (ncol(R) > nassets) {
468
          R <- R[, 1:nassets]
469
          warning("number of assets is less than number of columns in returns object, subsetting returns object.")
470
      }
471
      else {
472
          stop("number of assets is greater than number of columns in returns object")
473
      }
474
  }
475
  out <- constrained_objective(w = weights, R = R, portfolio = portfolio, 
476
      trace = TRUE, ...)$objective_measures
477
  return(structure(list(R = R, weights = weights, objective_measures = out, 
478
      call = match.call(), portfolio = portfolio), class = c("optimize.portfolio.invol", 
479
      "optimize.portfolio")))
480
}
481
# Calculate the objective measures for the vol weight portfolio
482
VolWgt.opt <- volatility.weight(R=R, portfolio=VolWgt.portf)
483

484

485
# REMOVED
486
# ### Evaluate Constrained Concentration to mETL Portfolio - with DE
487
# # registerDoSEQ() # turn off parallelization to keep the trace data
488
# ConstrConcmETL.DE<-optimize.portfolio(R=R,
489
#   portfolio=ConstrConcmETL.portf,
490
#   optimize_method='DEoptim',
491
#   search_size=40000,
492
#   NP=4000,
493
#   itermax=400,
494
#   trace=FALSE
495
# ) 
496
# #   list(c=0.25, # speed of crossover adaption (0,1]
497
# #   CR=0.75) # crossover probability [0,1]
498
# plot(ConstrConcmETL.DE, risk.col="StdDev", return.col="mean")
499
# plot(ConstrConcmETL.DE, risk.col="ES", return.col="mean") # several outlier portfolios
500
# chart.RiskBudget(ConstrConcmETL.DE)
501
# chart.RiskBudget(ConstrConcmETL.DE, risk.type="percentage")
502
# 
503
# save(ConstrConcmETL.DE,file=paste(resultsdir, 'ConstrConcmETL-', Sys.Date(), '-', runname, '.rda',sep=''))
504
# print(ConstrConcmETL.DE$elapsed_time)
505
print('Done with optimizations.')
506

507
#------------------------------------------------------------------------
508
### Extract data from optimizations for analysis
509
#------------------------------------------------------------------------
510
# Combine optimization objects
511
buoys <- combine.optimizations(list(MeanSD=MeanSD.ROI, MeanmETL=MeanmETL.RND, MinSD=MinSD.ROI, MinmETL=MinmETL.ROI, MRCSD=MRCSD.DE, MRCmETL=MRCmETL.DE, VolWgt=VolWgt.opt, EqWgt=EqWgt.opt))
512
# chart.Weights(buoys, plot.type="bar", ylim=c(0,1))
513
# 
514
# #@ Chart the portfolios that have mean and ES as objective measures. - RB
515
# chart.RiskReward(buoys, risk.col="ES")
516
# #@ Chart the portfolios that have mean and StdDev as objective measures. - RB
517
# chart.RiskReward(buoys, risk.col="StdDev")
518
# 
519
# #@ The MRCmETL and RB optimizations would be good to compare because they are
520
# #@ similar in that they both include component ES as an objective. - RB
521
# buoyETL <- combine.optimizations(list(MRCmETL=MRCmETL.RND, RB=RiskBudget.DE, EqWgt=EqWgt.opt))
522
# chart.RiskBudget(buoyETL, match.col="ES", risk.type="percentage", legend.loc="topright")
523
# 
524
# #@ Compare the equal weight portfolio and the equal SD contribution portfolio. - RB
525
# buoyStdDev <- combine.optimizations(list(MRCSD=MRCSD.RND, EqWgt=EqWgt.opt))
526
# chart.RiskBudget(buoyStdDev, match.col="StdDev", risk.type="absolute", legend.loc="topleft")
527

528
Wgts = extractWeights(buoys)
529

530
### Extract portfolio measures from each objective
531
# We can't just extract them, because they aren't all calculated
532
# so fill them in...
533
buoys.portfmeas =  buoys.contrib.sd = buoys.contrib.es = buoys.perc.sd = buoys.perc.es = NULL
534
for(i in 1:NROW(Wgts)){
535
  mean = sum(colMeans(R)*Wgts[i,])
536
  sd = StdDev(R, weights=Wgts[i,], portfolio_method="component", clean=clean)
537
  es = ES(R, weights=Wgts[i,], method="modified", portfolio_method="component", p=p, clean=clean)
538
  buoys.portfmeas=rbind(buoys.portfmeas, c(mean, sd[[1]][1], es[[1]][1]))
539
  buoys.contrib.sd= rbind(buoys.contrib.sd,sd[[2]])
540
  buoys.contrib.es= rbind(buoys.contrib.es,es[[2]])
541
  buoys.perc.sd = rbind(buoys.perc.sd,sd[[3]])
542
  buoys.perc.es = rbind(buoys.perc.es,es[[3]])
543
}
544
colnames(buoys.portfmeas)=c("Mean", "StdDev", "mETL")
545
rownames(buoys.portfmeas)=
546
rownames(buoys.contrib.sd)= 
547
rownames(buoys.contrib.es)= 
548
rownames(buoys.perc.sd) = 
549
rownames(buoys.perc.es) = rownames(Wgts)
550
colnames(buoys.contrib.sd)= 
551
colnames(buoys.contrib.es)= 
552
colnames(buoys.perc.sd) = 
553
colnames(buoys.perc.es) = colnames(Wgts)
554

555
# get the RP portfolios with risk and return pre-calculated
556
xtract = extractStats(MRCmETL.RND) 
557
save(xtract,file=paste(resultsdir, 'xtract-RPs-', Sys.Date(), '-', runname, '.rda',sep=''))
558
# columnnames = colnames(xtract)
559
results.names=rownames(buoys.portfmeas)
560

561
# by Asset metrics
562
assets.portfmeas=as.matrix(scatterFUN(R, FUN="mean"))
563
assets.portfmeas=cbind(assets.portfmeas, scatterFUN(R, FUN="StdDev", clean=clean))
564
assets.portfmeas=cbind(assets.portfmeas, scatterFUN(R, FUN="ES", clean=clean))
565
colnames(assets.portfmeas)=c("Mean", "StdDev", "mETL")
566
rownames(assets.portfmeas)=colnames(Wgts)
567

568

569

570
#------------------------------------------------------------------------
571
# Run select buoy optimizations through time
572
#------------------------------------------------------------------------
573
# 
574

575
# Equal Weight
576
dates=index(R[endpoints(R, on=rebalance_period)])
577
EqWgt.w = xts(matrix(rep(1/NCOL(R),length(dates)*NCOL(R)), ncol=NCOL(R)), order.by=dates)
578
colnames(EqWgt.w)= colnames(R)
579
EqWgt.R=Return.rebalancing(R, EqWgt.w)
580
chart.StackedBar(EqWgt.w, colorset=wb13color, gap=0)
581

582
VolWgt.w = NULL
583
for(i in 3:length(dates)){
584
  x = volatility.weight(R=R[paste0("::",dates[i]),], portfolio=VolWgt.portf)
585
  VolWgt.w = rbind(VolWgt.w, x$weights)
586
}
587
VolWgt.w = as.xts(VolWgt.w, order.by=dates[-1:-2])
588
VolWgt.R=Return.rebalancing(R, VolWgt.w)
589

590
# Equal SD
591
MRCSD.DE.t = optimize.portfolio.rebalancing(R=R,
592
  portfolio=MRCSD.portf, 
593
  optimize_method='DEoptim',
594
  search_size=2000,
595
  NP=200,
596
  initialpop=rp[1:50,], # seed with a starting population that we know fits the constraint space
597
  trace=FALSE,
598
  rebalance_on=rebalance_period, # uses xts 'endpoints'
599
  trailing_periods=NULL, # calculates from inception
600
  training_period=36) # starts 3 years in to the data history
601
MRCSD.w = extractWeights(MRCSD.DE.t)
602
MRCSD.gw = extractGroups(MRCSD.DE.t)
603
save(MRCSD.DE.t,file=paste(resultsdir, 'MRCSD.DE.t-', Sys.Date(), '-', runname, '.rda',sep=''))
604
chart.UnStackedBar(MRCSD.w, rotate="horizontal", colorset=wb13color, space=0, las=2)
605
MRCSD.R=Return.rebalancing(edhec.R, MRCSD.w)
606
colnames(MRCSD) = "MRCSD"
607

608
# Extract perc contrib of mES from results object
609
x=NULL
610
for(i in 1:length(names(MRCSD.DE.t)))  {
611
  x = rbind(x,MRCSD.DE.t[[i]][["objective_measures"]]$StdDev$pct_contrib_StdDev)
612
}
613
MRCSD.DE.pct_contrib_StdDev.t = as.xts(x, order.by=as.POSIXct(names(MRCSD.DE.t)))
614
chart.UnStackedBar(x.xts, rotate="horizontal", colorset=wb13color, space=0, las=2)
615

616

617
# MRC mETL
618
MRCmETL.DE.t = optimize.portfolio.rebalancing(R=R,
619
  portfolio=MRCmETL.portf, 
620
  optimize_method='DEoptim',
621
  search_size=20000,
622
  NP=200,
623
  initialpop=rp[1:50,], # seed with a starting population that we know fits the constraint space
624
  trace=FALSE,
625
  rebalance_on=rebalance_period, # uses xts 'endpoints'
626
  trailing_periods=NULL, # calculates from inception
627
  training_period=36) # starts 3 years in to the data history
628
MRCmETL.w = extractWeights(MRCmETL.DE.t)
629
MRCmETL.gw = extractGroups(MRCmETL.DE.t)
630
save(MRCmETL.DE.t,file=paste(resultsdir, 'MRCmETL.DE.t-', Sys.Date(), '-', runname, '.rda',sep=''))
631
chart.UnStackedBar(MRCmETL.w, rotate="horizontal", colorset=wb13color, space=0, las=2)
632
MRCmETL=Return.rebalancing(edhec.R, MRCmETL.w)
633
colnames(MRCmETL) = "MRCmETL"
634

635
MRCmETL.RND.t = optimize.portfolio.rebalancing(R=R,
636
  portfolio=MRCmETL.portf, 
637
  optimize_method='random',
638
  rp=rp,
639
  trace=TRUE,
640
  rebalance_on=rebalance_period, # uses xts 'endpoints'
641
  trailing_periods=NULL, # calculates from inception
642
  training_period=36) # starts 3 years in to the data history
643
MRCmETL.RND.w = extractWeights(MRCmETL.RND.t)
644
MRCmETL.gw = extractGroups(MRCmETL.RND.t)
645
save(MRCmETL.DE.t,file=paste(resultsdir, 'MRCmETL.DE.t-', Sys.Date(), '-', runname, '.rda',sep=''))
646
chart.UnStackedBar(MRCmETL.w, rotate="horizontal", colorset=wb13color, space=0, las=2)
647
MRCmETL=Return.rebalancing(edhec.R, MRCmETL.w)
648
colnames(MRCmETL) = "MRCmETL"
649
end_time<-Sys.time()
650
end_time-start_time
651

652

653
#########################################################################
654
# Optimization ends here
655
#########################################################################
656

657

658

659
Product

Resources

Company
