The parma package

Description

Portfolio Allocation and Risk Management. Models and Methods for scenario and moment based optimization of portfolios.

Details

The portfolio allocation and risk managament applications (parma) package contains a unique set of methods and models for the optimal allocation of capital in financial portfolios. It uniquely represents certain discontinuous problems using their smooth approximation counterparts and implements fractional based programming for the direct optimization of risk-to-reward ratios. In combination with the rmgarch package, it enables the confident solution to scenario based optimization problems using such risk and deviation measures as Mean Absolute Deviation (MAD), Variance (EV), Minimax, Conditional Value at Risk (CVaR), Conditional Drawdown at Risk (CDaR) and Lower Partial Moments (LPM). In addition, it implements moment based optimization for use with the quadratic EV problem, and a higher moment CARA utility expansion using the coskewness and cokurtosis matrices generated from the GO-GARCH with affine GH or NIG distributions. Benchmark relative optimization (tracking error) is also implemented as are basic mixed integer cardinality constraints. Finally, for non-convex problem formulations such as the upper to lower partial moments function, global optimization methods using a penalty based method are available. The key functions in the package are parmaspec which defines the optimization setup, and parmasolve which solves the problem given a chosen representation and solver. A portfolio frontier function is implemented in parmafrontier, utility optimization in parmautility and a custom translation of the cmaes global optimization solver of Hansen (2006) with full features is implemented in cmaes.

How to cite this package

Whenever using this package, please cite as

@Manual{Galanos_2014,
 author       = {Alexios Galanos and Bernhard Pfaff},
 title        = {{parma}: Portfolio Allocation and Risk Management Applications.},
 year         = {2014},
 note 	      = {R package version 1.5-1.},}

License

The releases of this package is licensed under GPL version 3.

Author(s)

Alexios Galanos and Bernhard Pfaff

References

Charnes, A. and Cooper, W. 1962, Programming with linear fractional functionals, Naval Research Logistics Quarterly, 9, 181–186.
Dinkelbach, W. 1967, On nonlinear fractional programming, Management Science, 13(7), 492–498.
Fishburn, P.C. 1977, Mean-risk analysis with risk associated with below-target returns, The American Economic Review, 67(2), 116-126.
Galanos, A. 2012, Higher Moment Models for Risk and Portfolio Management, Thesis (submitted) Cass Business School.
Hansen, N. 2006, The CMA Evolution Strategy: A Comparing Review, Towards a New Evolutionary Computation (Studies in Fuzziness and Soft Computing), 192, 75–102.
Holthausen, D. 1981, A risk-return model with risk and return measured as deviations from a target return, The American Economic Review, 71, 182–188.
Konno, H. and Yamazaki, H. 1991, Mean-absolute deviation portfolio optimization model and its applications to Tokyo stock market, Management Science, 37(5), 519–531.
Markowitz, H. 1952, Portfolio selection, The Journal of Finance, 7(1), 77–91.
Rockafellar, R.T. and Uryasev, S. and Zabarankin, M., 2006, Generalized deviations in risk analysis, Finance and Stochastics, 10(1), 51–74.
Stoyanov, S.V. and Rachev, S.T. and Fabozzi, F.J. 2007, Optimal financial portfolios, Applied Mathematical Finance, 14(5), 401–436.

Second-order Cone Programming

Description

The function solves second-order cone problem by primal-dual interior point method. It is a wrapper function to the C-routines written by Lobo, Vandenberghe and Boyd (see reference below).

Usage

Socp(f, A, b, C, d, N,
     x = NULL, z = NULL, w = NULL, control = list())

Arguments

Details

The primal formulation of an SOCP is given as:

minimise f' * x

subject to

||A_i*x + b_i|| <= c_i' * x + d_i

for i = 1,\ldots, L. Here, x is the (n \times 1) vector to be optimised. The dual form of an SOCP is expressed as:

maximise \sum_{i = 1}^L -(b' * z_i + d_i * w_i)

subject to

\sum_{i = 1}^L (A_i' * z_i + c_i * w_i) = f

and

||z_i || = w_i

for i = 1,\ldots, L, given strictly feasible primal and dual initial points.

The algorithm stops, if one of the following criteria is met:

1. abs.tol – maximum absolute error in objective function; guarantees that for any x: abs(f'*x - f'*x\_opt) <= abs\_tol.

2. rel.tol – maximum relative error in objective function; guarantees that for any x: abs(f'*x - f'*x\_opt)/(f'*x\_opt) <= rel\_tol (if f'*x\_opt > 0). Negative value has special meaning, see target next.

3. target – if rel\_tol<0, stops when f'*x < target or -b'*z >= target.

4. max.iter – limit on number of algorithm outer iterations. Most problems can be solved in less than 50 iterations. Called with max_iter = 0 only checks feasibility of x and z, (and returns gap and deviation from centrality).

5. The target value is reached. rel\_tol is negative and the primal objective p is less than the target.

Value

A list-object with the following elements:

Note

This function has been ported from the Rsocp package contained in the Rmetrics-Project on R-Forge. In contrast to the former implementation, allowance is made for specifying more than one cone constraint.

Author(s)

Bernhard Pfaff

References

Lobo, M. and Vandenberghe, L. and Boyd, S., SOCP: Software for Second-order Cone Programming, User's Guide, Beta Version, April 1997, Stanford University.

See Also

SocpPhase1, SocpPhase2, SocpControl

Control Variables for Socp

Description

This function returns a list object of control parameters that are passed down to the C-function SOCP. It's default values are used in Socp.

Usage

SocpControl(abs.tol = 1e-18, rel.tol = 1e-16, target = 0,
            max.iter = 500, Nu = 10, out.mode = 0, BigM.K = 2,
            BigM.iter = 5)

Arguments

Details

For details about these control parameters, the reader is referred to the reference below, in particular sections 2.7, 2.8 and 4.3 to 4.5. A pdf-version of the user's guide is shipped in the packages doc subdirectory.

Value

A list object with the control parameters.

Note

This function has been ported from the Rsocp package contained in the Rmetrics-Project on R-Forge.

Author(s)

Bernhard Pfaff

References

Lobo, M. and Vandenberghe, L. and Boyd, S., SOCP: Software for Second-order Cone Programming, User's Guide, Beta Version, April 1997, Stanford University.

See Also

Socp

SOCP: Initialising objective variable x in primal form

Description

This function determines values for x, whence they have not been specified by the user. Here, a feasibility problem is solved first and its solution is then used as an initial point for the original problem.

Usage

SocpPhase1(f, A, b, N, control)

Arguments

Details

The finding of an initial point x is described in the user's guide, sectionb 2.8.

Value

A vector with the initial point for x.

Note

This function has been ported from the Rsocp package contained in the Rmetrics-Project on R-Forge.

Author(s)

Bernhard Pfaff

References

Lobo, M. and Vandenberghe, L. and Boyd, S., SOCP: Software for Second-order Cone Programming, User's Guide, Beta Version, April 1997, Stanford University.

See Also

Socp, SocpPhase2, SocpControl

SOCP: Initialising objective variable z in dual form

Description

This function determines values for z, whence they have not been specified by the user.

Usage

SocpPhase2(f, A, b, N, x, control)

Arguments

Value

A vector with the initial point for z (dual form of SOCP).

Note

This function has been ported from the Rsocp package contained in the Rmetrics-Project on R-Forge.

Author(s)

Bernhard Pfaff

References

Lobo, M. and Vandenberghe, L. and Boyd, S., SOCP: Software for Second-order Cone Programming, User's Guide, Beta Version, April 1997, Stanford University.

See Also

Socp, SocpPhase1, SocpControl

The Covariance Matrix Adaptation Evolution Strategy (cmaes) Solver

Description

The direct translation of the Hansen's cmaes matlab code v3.60.

Usage

cmaes(pars, fun, lower = rep(0, length(pars)), upper = rep(1, length(pars)), 
insigma = 1, ctrl = cmaes.control(), ...)

cmaes.control(
options = list(StopFitness = -Inf, MaxFunEvals = Inf,
MaxIter = '1e3*(N+5)^2/sqrt(popsize)', StopFunEvals = Inf,
StopIter = Inf, TolX = '1e-11*max(insigma)', TolUpX = '1e3*max(insigma)',
TolFun = 1e-12, TolHistFun = 1e-13, StopOnStagnation = TRUE,
StopOnWarnings = TRUE, StopOnEqualFunctionValues = '2 + N/3',
DiffMaxChange = Inf, DiffMinChange = 0, WarnOnEqualFunctionValues = FALSE,
EvalParallel = FALSE, EvalInitialX = TRUE, Restarts = 0,
IncPopSize = 2, PopSize = '4 + floor(3*log(N))', ParentNumber = 'floor(popsize/2)',
RecombinationWeights = c("superlinear", "linear", "constant"), 
DiagonalOnly = '0*(1+100*N/sqrt(popsize))+(N>=1000)', 
CMA = TRUE, Seed = 'as.integer(Sys.time())', DispFinal = TRUE, 
DispModulo = 100, Warnings = FALSE), 
CMA = list(cs = '(mueff+2)/(N+mueff+3)', 
damps = '1 + 2*max(0,sqrt((mueff-1)/(N+1))-1) + cs',
ccum = '(4 + mueff/N) / (N+4 + 2*mueff/N)', ccov1 = '2 / ((N+1.3)^2+mueff)',
ccovmu = '2 * (mueff-2+1/mueff) / ((N+2)^2+mueff)', active = 0))

Arguments

Details

This solver has been translated from the matlab version created by Nikolaus Hansen and available on his website http://www.cmap.polytechnique.fr/~nikolaus.hansen/cmaes_inmatlab.html. There is also a cmaes on CRAN but this does not offer the same level of options and flexibility that the matlab version offers. For more details on what the options mean and generally how the cmaes solver works, consult the relevant website and literature.

Author(s)

Alexios Galanos

References

Hansen, N. 2006, The CMA Evolution Strategy: A Comparing Review, Towards a New Evolutionary Computation (Studies in Fuzziness and Soft Computing), 192, 75–102.

Examples

## Not run: 
ctrl = cmaes.control()
ctrl$options$StopOnWarnings = FALSE
ctrl$cma$active = 1
ctrl$options$TolFun = 1e-12
ctrl$options$DispModulo=100
ctrl$options$Restarts = 0
ctrl$options$MaxIter = 3000
ctrl$options$TolUpX = 5
ctrl$options$PopSize = 300
test1 = cmaes(rnorm(10), fun = parma:::fsphere, 
lower = -Inf, upper = Inf, insigma = 1, ctrl = ctrl)
test2 = cmaes(rnorm(10), fun = parma:::frosenbrock, 
lower = -Inf, upper = Inf, insigma = 1, ctrl = ctrl)

ctrl = cmaes.control()
ctrl$options$StopOnWarnings = FALSE
ctrl$cma$active = 1
ctrl$options$TolFun = 1e-12
ctrl$options$DispModulo=100
ctrl$options$Restarts = 0
ctrl$options$MaxIter = 3000
ctrl$options$PopSize = 400
test3 = cmaes(rep(1, 10), fun = parma:::frastrigin10, 
lower = -50, upper = 50, insigma = 1, ctrl = ctrl)

## End(Not run)

NLP custom constraint functions

Description

Provides a number of custom constraints and their jacobians for use with the NLP representation (both minimum risk and the fractional problem).

Usage

ineqfun.turnover.min(w, optvars, uservars)
ineqjac.turnover.min(w, optvars, uservars)
ineqfun.bsturnover.min(w, optvars, uservars)
ineqjac.bsturnover.min(w, optvars, uservars)
ineqfun.turnover.opt(w, optvars, uservars)
ineqjac.turnover.opt(w, optvars, uservars)
ineqfun.bsturnover.opt(w, optvars, uservars)
ineqjac.bsturnover.opt(w, optvars, uservars)
ineqfun.variance.opt(w, optvars, uservars)
ineqjac.variance.opt(w, optvars, uservars)
ineqfun.variance.min(w, optvars, uservars)
ineqjac.variance.min(w, optvars, uservars)

Arguments

Details

Functions preceded be “ineqfun” denote the functions which may be passed to the ineqfun option in the parmaspec function, while “eqfun” to the eqfun option. Functions which include “jac” in the first part of the name denote the jacobians of the equivalent constraints and should be passed to the ineqgrad or eqgrad options in parmaspec. Functions ending with “min” denote formulations for use with the minrisk type problems while functions ending with “opt” denote the fractional risk formulation.
For the simple turnover constraint, the uservars list must contain an entry called ‘wold’ denoting the previous vector of weights with which the comparison will be made. Additionally, an entry called ‘turnover’ is required which denotes the (positive) value for the maximum turnover.
For the buy and sell turnover constraint (bsturnover), there should instead be (positive) ‘buyturnover’ and ‘sellturnover’ entries in the uservars list. Finally note that when using this type of constraints in a fractional programming setup, care should be taken that the combination of bounds, turnover limits and the forecast return vector do not result in a negative expected return in which case the problem is not solvable.
The variance constraint allows the targeting of a maximum acceptable variance. The extra arguments which must be passed to the uservars list are ‘Cov’ for the asset covariance matrix and ‘varbound’ representing the maximum acceptable upper variance.

Value

Used internally by the NLP solver. The fun return a scalar, while the jac return a matrix with n.cols equal to the length of the decision vector (which may be greater than the length of the weights as in the fractional problem which contains the fractional multiplier as well as other formulations which have additional decision variables).

Author(s)

Alexios Galanos

15 Exchange Traded Funds (ETFs)

Description

The xts dataset consists of the adjusted daily closing prices of 15 iShare ETFs for the period 2003-05-28 to 2012-06-01 (2272 periods) representing a selection of US style and international equity benchmark. The iShare series are IWF (Russell 1000 Growth Index), IWD (Russell 1000 Value Index), IWO (Russell 2000 Growth Index), IWN (Russell 2000 Value Index), EEM (MSCI Emerging Markets Index), TLT (Barclays 20+ Year T-Bond), EWC (MSCI Canada Index), EWA (MSCI Australia Index), EWJ (MSCI Japan Index), EWG (MSCI Germany Index), EWL (MSCI Switzerland Index), EWQ (MSCI France Index), EWU (MSCI UK Index), EPP (MSCI Pacific ex-Japan), EZA (MSCI South Africa Index).

Usage

data(etfdata)

Format

An xts matrix containing 2272x15 observations.

Source

Yahoo Finance

Class "parmaPort"

Description

The parma optimal portfolio class.

Objects from the Class

Objects can be created by calls to parmasolve..

Slots

solution:

Object of class "vector" The list with the optimal values.

model:

Object of class "vector" A list with details of the risk model.

Methods

show

signature(object = "parmaPort"): Summary.

weights

signature(object = "parmaPort"): Extracts the optimal weights.

tictoc

signature(object = "parmaPort"): Extracts the time elapsed to solve the problem.

checkarbitrage

signature(object = "parmaPort"): Extracts the arbitrage check on the scenario.

parmarisk

signature(object = "parmaPort"): Extracts the expected risk of the optimized portfolio.

parmareward

signature(object = "parmaPort"): Extracts the expected reward of the optimized portfolio.

parmastatus

signature(object = "parmaPort"): Solver termination code for the LP and NLP solvers.

Author(s)

Alexios Galanos

Examples

showClass("parmaPort")

Class "parmaSpec"

Description

Object returned from calling parmaspec.

Objects from the Class

Objects can be created by calls of the form new("parmaSpec").

Slots

model:

A list with details of the risk and optimization model.

modeldata:

A list with the data.

constraints:

A list with details on the optimization constraints.

Methods

parmasolve

signature(spec = "parmaSpec"): Solves for the optimal weights.

show

signature(object = "parmaSpec"): Summary method.

parmaset<-

signature(object = "parmaSpec"): Set a specified ‘arg’ from the list of arguments defined in the parmaspec to some ‘value’ and return a modified specification.

parmaget

signature(object = "parmaSpec"): Get a specified ‘arg’ from the list of arguments defined in the parmaspec.

Author(s)

Alexios Galanos

Examples

showClass("parmaSpec")

Efficient Frontier Generator

Description

Solves for the portfolios on the efficient frontier given a specification object.

Usage

parmafrontier(spec, n.points = 100, miny = NULL, maxy = NULL, type = NULL, 
solver = NULL, solver.control = list(), parma.control = list(ubounds = 10000, 
mbounds = 1e+05, penalty = 10000), cluster = NULL)

Arguments

Details

If using a cluster object, care should be taken since some of the LP based problems are quite memory intensive. Though some care some gone into using slam matrices and forcing garbage collection in intermediate setup steps of these problems, there is no guarantee that memory will not be quickly depleted for large problems using many cores/sockets.

Value

A matrix object with columns for the weights, risk, reward and termination status of the solver (not available for the QP solver), and n.points rows. The user should investigate the non-converged solutions as well as remove any NA based rows (non-converged with error).

Author(s)

Alexios Galanos

Portfolio Allocation Model Solver

Description

Solves for the optimal weights given parmaSpec defined model.

Usage

parmasolve(spec, type = NULL, solver = NULL, solver.control = list(), x0 = NULL, 
w0 = NULL, parma.control = list(ubounds = 1e4, mbounds = 1e5, 
penalty = 1e4, eqSlack = 1e-05), ...)

Arguments

Details

In most of the cases, the intersection of objective and constraints will define whether a problem is LP, MILP, QP, MIQP, QCQP, NLP, MINLP or GNLP. However, there are cases when the problem can be solved by more than one type of solver so that the ‘type’ option allows the user to choose between the options. For the GNLP type, the ‘solver’ allows for a choice of “cmaes” or “crs”, while the LP type there is a choice of “glpk” or “symphony”. Future development will likely expand on these choices.

Value

A parmaPort object containing details of the PARMA optimized portfolio.

Author(s)

Alexios Galanos

Portfolio Allocation Model Specification

Description

Defines the type model and method for optimization using either a scenario or covariance matrix.

Usage

parmaspec(scenario = NULL, probability = NULL, S = NULL, Q = NULL, qB = NULL, 
		benchmark = NULL, benchmarkS = NULL, forecast = NULL, target = NULL, 
		targetType =  c("inequality", "equality"), 
		risk = c("MAD", "MiniMax", "CVaR", "CDaR", "EV", "LPM", "LPMUPM"), 
		riskType = c("minrisk", "optimal", "maxreward"), riskB = NULL,
		options = list(alpha = 0.05, threshold = 999, moment = 1, lmoment=1,
		umoment=1, lthreshold = -0.01, uthreshold = 0.01), 
		LB = NULL, UB = NULL, budget = 1, leverage = NULL, 
		ineqfun = NULL, ineqgrad = NULL, eqfun = NULL, eqgrad = NULL, 
		uservars = list(), ineq.mat = NULL, ineq.LB = NULL, 
		ineq.UB = NULL, eq.mat = NULL, eqB = NULL, max.pos = NULL, 
		asset.names = NULL, ...)

Arguments

Details

The parmaspec method is the entry point for specifying and solving portfolio problems in the parma package. Currently 7 measures of risk are supported, 3 based on tail measures: Conditional Value at Risk (CVaR), Conditional Drawdown at Risk (CDaR) and Lower Partial Moments (LPM), and 3 based on the Lp-Norm: Mean Absolute Deviation (L_1, MAD), Mean Variance (L_2, EV) and MiniMax (L_inf, Minimax). The LPMUPM measure is the ratio of lower to upper partial moments, a non convex measure discussed in Holthausen (1981). Additionally, the problems may be solved based on minimization of risk subject to a target return, else on the optimal risk-reward ratio using fractional programming (see references), thus avoiding the estimation of the entire frontier. Problems are classified and solved according to whether they can be formulated as Linear (LP), Mixed Integer LP (MILP), Quadratic (QP), Mixed Integer Quadratic (MIQP), Second Order Cone Programming (SOCP), Non-Linear (NLP), Mixed Integer NLP (MINLP) and Global NLP (GNLP). This in turn depends on the intersection of objectives and constraints. It is possible that a problem may be solved both as an LP and NLP (or QP and NLP), and this can be defined during the solver stage (parmasolve). Because all NLP models, make use of analytical derivatives, the results should be the same for any formulation chosen, and considerations such as memory usage should guide the choice of formulation (with some LP models being particularly expensive). Not all problem types are supported, but this might change subject to the availability of solvers in R which can deal with these specific types e.g. MINLP and MIQP. The parmaspec also allows the input of a benchmark so that benchmark relative optimization is carried out. User defined equality and inequality functions for NLP problems need to be properly defined to be accepted by the model, and their analytic gradients also provided, unless the problem is solved as a GNLP in which case a derivative free penalty function is used. These custom constraint functions should be provided in a list, and should take as arguments the vector of decision variables ‘w’, an argument called ‘optvars’ which is used by the program internally, and an argument called ‘uservars’ which is a list with optional user defined values for the constraints. The examples in the inst folder provide some guidance, and the user is left to his own devices to study the underlying workings of the program to understand how to supply these. Finally, the NLP functions which are known to be discontinuous because of the presence of functions such as the min and max, have been re-written to take advantage of smooth approximations to such functions, details of which may be founds in the vignette. The package support for GNLP is based on a choice of the cmaes solver of the cmaes package (which is not production level) or the crs solver of the nloptr package which may be defined in solver.control option of the parmasolve method with named argument ‘solver’. High quality GNLP solvers are not available in R and as such support for these types of problems is experimental at best and your mileage will vary. The problems which must be solved as GNLP include the ‘LPMUPM’ measure, all problems with risk type ‘optimal’ AND cardinality constraints (‘max.pos’), and all problems with custom NLP constraints without derivatives, non-convex inequalities or non-affine equalities.

Value

A parmaSpec object containing details of the PARMA specification.

Author(s)

Alexios Galanos

Utility Based Optimization

Description

Utility based portfolio optimization using either Taylor series expansion of utility function with moments or scenario based.

Usage

parmautility(U = c("CARA", "Power"), method = c("moment", "scenario"),
		scenario = NULL, M1 = NULL, M2 =  NULL, M3 = NULL, M4 = NULL, RA = 1, 
		budget = 1, LB = rep(0, length(M1)), UB = rep(1, length(M1)))

Arguments

Details

The function currently only implements the CARA moment based approach, but will be expanded in the future. The moment approach can take as inputs either M1 and M2 (2-moment approximation), or M1, M2, M3 and M4 (4-moment approximation). Not many models generate M3 and M4, but the “gogarch” model with manig or magh distribution will.

Value

A parmaPort object containing details of the PARMA optimized portfolio.

Author(s)

Alexios Galanos

References

Galanos, A. and Rossi, E. and Urga, G. 2012, Independent Factor Autoregressive Conditional Density Model submitted-TBA

Portfolio Risk Measures

Description

Calculates a given portfolio risk/deviation measure given a set of weights and matrix of returns, possible representing a forecast scenario.

Usage

riskfun(weights, Data, risk = c("mad", "ev", "minimax", "cvar", "cdar", "lpm"), 
benchmark = NULL, alpha = 0.05, moment = 1, threshold = 0, VaR = NULL, DaR = NULL)

Arguments

Details

A simple utility function for the calculation and understanding of some of the risk and deviation measures implemented in the package.

Value

A numeric value representing the risk/deviation measure.

Author(s)

Alexios Galanos