{'fraction_non_alphanumeric': 0.03686507568019315, 'fraction_numerical': 0.011003807224440524, 'mean_word_length': 4.674130663856691, 'pattern_counts': {'":': 0, '<': 1, '<?xml version=': 0, '>': 1, 'https://': 0, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 2, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
Local Control Regression: Improving the Least Squares Monte Carlo Method for Portfolio Optimization

The least squares Monte Carlo algorithm has become popular for solving portfolio optimization problems. A simple approach is to approximate the value functions on a discrete grid of portfolio weights, then use control regression to generalize the discrete estimates. However, the classical global control regression can be expensive and inaccurate. To overcome this difficulty, we introduce a local control regression technique, combined with adaptive grids. We show that choosing a coarse grid for local regression can produce sufficiently accurate results.


Introduction
The least squares Monte Carlo (LSMC) algorithm, originally introduced by [3], [15] and [18] for the pricing of American options, has become a popular tool for solving stochastic control problems in many fields, such as exotic option pricing, inventory management, project valuation, portfolio optimization and many others. The main strength of the LSMC algorithm is its ability to handle multivariate stochastic state variables, by combining the Monte Carlo simulation of the state variables with the regression estimates of the conditional expectations in the dynamic programming equations. The purpose of this paper is to improve the control discretization approach of the LSMC algorithm for solving general portfolio optimization problems, including general investment objective functions, general intermediate costs and general asset dynamics. Before describing the details of the method, we briefly review the literature on the LSMC algorithm applied to portfolio optimization, with an emphasis on how portfolio weights are handled and how the maximization over portfolio weights is performed. The Taylor expansion method in [1] was the first attempt to solve the portfolio optimization problems using the LSMC algorithm. The authors determine a semi-closed form solution by deriving the first-order condition of the Taylor series expansion of the investor's future value function, see for other examples, [19], [10], and [6]. Another approach to derive and solve the first-order condition is to use [11]'s regress-later technique, see for example, [4] and [5]. However, these two first-order condition approaches, relying on analytical derivation, are restricted in the range of applications. The simplest, most straightforward and most general approach to deal with the portfolio weights in the LSMC algorithm is to discretize them, perform one regression per discretized portfolio weight, and compute the optimal allocation for each Monte Carlo path by naive exhaustive search, see, for example, [14], [6] and [20]. Although this approach is stable and accurate, it is computationally very expensive if a fine grid is used, especially for multi-dimensional problems. To improve the computational efficiency of the control discretization approach, control regression is a possible solution. The simplest approach is to regress state variables and portfolio weights at once. For example, in [7] and [8], the value functions are regressed on the simulated exogenous state variables, discretized wealth levels and discretized portfolio weights, and then the first-order condition is derived and solved. In [20], a control regression approach is implemented via the control randomization approach of [12]: the portfolio weights are randomized and the wealth variable is correspondingly computed in the forward simulation, afterwards the value functions are recursively regressed on the exogenous state variables, the randomized portfolio weights and wealth. However, as shown in [20], this approach may not be as accurate as the grid search approach, especially when the time horizon is long or the payoff function is highly nonlinear. Moreover, as discussed in [13], regressing on state variables and portfolio weights at once requires a very large regression, which can make it computationally difficult. Another strand of control regression approach is introduced in [13]. The authors first regress the conditional expectations given in the first-order condition for each portfolio weight on a coarse grid, then regress the resulting regression coefficients on these portfolio weights, obtaining a hybrid generalized conditional expectation estimate expressed in both portfolio weight and state variables while keeping the size of the coarse grid of portfolio weights manageable. The regression basis of this approach is usually chosen to be polynomial (linear in [13], quadratic in [16] and [9]), in order to keep the maximization over the portfolio weights simple (constrained optimization solver in [13], analytical solution in [16] and [9]). One such attempt is in [2]. Their approach is to sample the portfolio weights and the asset prices using low-discrepancy Sobol sequences, then use a global regression to obtain a conditional expectation estimate for every possible portfolio weight, and finally solve for the optimal portfolio weight by gradient optimization. Building upon the aforementioned approaches, we propose a hybrid framework that adapts [20]'s combined method of control discretization and control randomization to incorporate the control regression approach of [2]. Our control regression approach differs from [2] in two main aspects, namely the use of local regression instead of global regression, and the use of adaptive refinement grids to target the region of optimal weights instead of gradient optimization. Our numerical experiments show that local control regression is more accurate and more efficient than global control regression, and that using a very coarse grid for the local control regression can produce sufficiently accurate results.

Problem Description
We consider a finite horizon portfolio allocation problem with one risk-free asset and d risky assets available. Let T be the investment horizon. We assume the portfolio can be rebalanced at any discrete time t n from the equally-spaced time grid Π = {t 0 = 0, t 1 , . . . , t N = T }. Denote by r f the return of the risk-free asset over a single period, by {r tn } 0≤n≤N = {r i tn } 1≤i≤d 0≤n≤N the asset returns , by the compounding factors of the asset returns, and by {S tn } 0≤n≤N = {S i tn } 1≤i≤d 0≤n≤N the asset prices. Let q tn = (q i tn ) 1≤i≤d 0≤n≤N −1 describe the number of units held in each risky asset and let {q f tn } 0≤n≤N −1 denote the amount allocated in the risk-free cash. Finally, let {W tn } 0≤n≤N denote the portfolio value (wealth) process. Let × and ÷ denote the element-wise multiplication and division between two vectors. The asset prices evolve as and the wealth process evolves as Let F = {F t } 0≤t≤T be the filtration generated by all the state variables. At any rebalancing time t n ∈ Π, the optimization problem is Here, {Z tn } 0≤n≤N is the vector of return predictors which drive the asset price dynamics {S tn } 0≤n≤N . Thus, the variables {Z tn } 0≤n≤N form the exogenous state variables in our problem. The investment objective is to maximize the expected utility E[U (W t N )] of the investor's final wealth. The objective is maximized over the portfolio allocation weights in the risky assets α tn = (α i tn ) 1≤i≤d 0≤n≤N −1 , with the allocation in the risk-free asset being equal to α f tn = 1 − 1≤i≤d α i tn . The relationship between portfolio weights α tn and portfolio positions q tn is given by α tn × W tn = S tn × q tn . The set A ⊆ R d of admissible strategies is assumed time-invariant for each discrete rebalancing time.

Least Squares Monte Carlo
In this section, we describe how to use the LSMC method to solve the dynamic portfolio optimization problem. The objective function in Eq. (2.3) can be formulated as a discrete-time dynamic programming principle, The conventional LSMC scheme simulates the state variables forward, then approximates the value functions recursively backwards in time by least squares regressions. The difficulty here stems from the conditioning wealth variable W tn , because it requires the information of all the previous portfolio decisions {α t n } 0≤n ≤n which are not known yet in the backward dynamic programming loop. To overcome this difficulty, a popular approach is to approximate the conditional expectations in Eq. (2.3) on a discrete grid of wealth levels, and then interpolate them on these discrete wealth levels, see for example, [1], [7] and [8].
The discretization approach can be computationally demanding if a fine grid is required, especially when there are multiple endogenous state variables. For example, when market impact is accounted for, the asset prices become endogenous state variables, which can make the discretization approach computationally infeasible. Instead, we use a control randomization approach that simulates random portfolio weights {α m tn } 1≤m≤M 0≤n≤N −1 in the forward loop in addition to the return predictors {Z m tn } 1≤m≤M 0≤n≤N , the asset returns {r m tn } 1≤m≤M 0≤n≤N , and the asset prices {S m tn } 1≤m≤M 0≤n≤N . Simulations for the wealth variable {W m tn } 1≤m≤M 0≤n≤N can be subsequently computed using the randomized portfolio weights. These portfolio weights are uniformly drawn from the admissible set A, so as to ensure that the simulated sample covers the whole space of possible portfolio weights. After all the state variables have been simulated, the next task is to approximate the conditional expectations in the dynamic programming formula in Eq. (3.1) by regression. A possible approach, called control regression, is to regress on the randomized portfolio weights as well as the return predictors, and then compute the optimal decisions by solving the first-order conditions. Such a control regression scheme, combined with control randomization, is analyzed in [12], and an application of this approach for solving dynamic portfolio optimization problem is given in [20]. However, regressing on the portfolio weights and the state variables at once requires a very large regression, which can make it computationally difficult, especially for high dimensional portfolios. Instead, we resort to control discretization, which has been shown to be more accurate than control regression in [20]. The control space is discretized as A ≈ A disc = {a 1 , ..., a J }. At time t n , assuming the mappinĝ v tn+1 : (z, w) →v tn+1 (z, w) has been estimated, the dynamic programming principle in Eq.
where CV j tn denotes the continuation value function of choosing the portfolio weight a j at time t n . To evaluate CV j tn for each a j ∈ A disc , we first update the portfolio weights α m tn = a j for each Monte Carlo path m = 1, ..., M . Then, the wealth variables at time t n+1 can be recomputed based on the portfolio weight a j , is the recomputed wealth at time t n+1 conditioning on the information at time t n at the j th node of the control grid. Let {ψ k (z, w)} 1≤k≤K be the basis functions of the state variables for the regression. For each portfolio decision a j ∈ A disc , the corresponding CV j tn is approximated by least squares minimization: Finally, CV j tn is parametrized asĈ and the optimal decision policy iŝ where the maximization is performed by an exhaustive grid search. As tested in [20], the control discretization approach is numerically very stable, even when transaction costs and price impact are taken into account. However, the exhaustive grid search is computationally very expensive, especially for high dimensional portfolios. In the next section, we will introduce a local control regression method to reduce the computational burden of control discretization.

Control Regression and Maximization
After the value functions have been estimated on the discrete grid A disc , the next step is to use control regression to project these discrete estimates onto a continuous space of control. This section describes in detail how to use control regression to improve the estimates. For notional convenience, we denote bŷ CV j =ĈV j (z, w) the estimated continuation value of choosing the portfolio weight a j at an arbitrary time and for an arbitrary input of state variables (z, w). Suppose we have already implemented the grid search approach to estimate the conditional valuesĈV j on a coarse discrete grid A disc , and have estimated the optimal control by the following equationα = arg max (4.1) The usual approach to generalize these conditional value estimates on the discrete grid A disc is to use global regression, see for example [13] and [2]. However, it can be hard to find a regression basis that fits the global grid well, even when using advanced regression techniques. As an example, Figure 4.1 illustrates the global regression approach for a portfolio with risk-free cash and two risky assets (d = 2). Moreover, even if adequate regression basis can be found, the global regression can be computationally very expensive. In fact, it is unnecessary to project the discrete value functions onto the whole global control space, because what matters for maximizing the value functions is the accuracy on the small local region containing the optimum.

Local control regression
In order to improve the numerical accuracy and efficiency of control regression, we propose to perform a local regression around the optimal estimateα. We construct a local grid L disc (α) ⊂ A disc that contains the neighbors ofα : where δ is the mesh size of the discrete grid A disc . Then, we regress the estimated conditional valuesĈV j with respect to the points in the local grid L disc (α), obtaining the parametric conditional value estimateŝ Φ (a) over the whole continuous local control hypercube 3) The main advantage of local regression is a better goodness-of-fit, as the local grid L disc (α) contains only very few points and are therefore easier to fit even with simple bases such as local second-order polynomial bases. In addition, a simple second-order polynomial is a sensible candidate for such a local parametric regression, as the local grid only contains a maximum of three points in each dimension and the maximum lies at either the middle point or the boundary of the local coarse grid. Furthermore, to guarantee the robustness of local regression, we perform local Ridge regressions with second-order polynomial bases.

Maximization by adaptive grids
The next step is to improve the discrete optimal control estimate in Eq. (4.1) by solving the continuous equation in Eq. (4.4). Classical approaches to do so are to use gradient optimization which requires either analytically or numerically compute the gradients. Instead, we perform the maximization using a gradientfree adaptive refinement strategy.
Here is a simple illustration of this adaptive grid method. Suppose we are given a portfolio of a risk-free asset and a risky asset. We denote by {A p } 0≤p≤P the sequence of adaptive grids and denote by {δ p } 0≤p≤P the sequence as adaptive mesh sizes such that δ p+1 = δ p /2. We define the initial adaptive grid as A 0 = A disc and define the initial control mesh size as δ 0 = δ. Suppose the optimal controlα 0 has been estimated on the space A 0 = A disc . Then, the updated adaptive grid is defined as A 1 = {α 0 − δ 0 /2,α 0 ,α 0 + δ 0 /2} and the updated estimate of optimal control is given byα 1 = arg max a∈A1Φ (a). Iteratively, the adaptive grid will be A p = {α p−1 − δ/2 p ,α p−/1 ,α p−1 + δ/2 p }. For every refined adaptive grid, there is a maximum of two new points in each dimension of the control. This hierarchical structure guarantees a limited number of operations will be performed in the algorithm. Below is an example of our adaptive grid refinement strategy starting with a δ = 1/4 control mesh for A disc : One can see that it takes only five iterations to reach a control mesh precision higher than 1/100, and only a maximum of two extra evaluation points are required at each iteration. This adaptive grid method can be described as a modified multivariate bisection search for the local maximum of a function. The combination of local regression and adaptive grids is described in Algorithm 1. In our experience, both gradient optimization and adaptive grids can reach the optimum but the latter is more efficient.

Numerical Experiments
In this section, we validate the numerical method proposed in the previous sections by applying it to a multiperiod portfolio optimization problem as formulated in Section 2. We obtain monthly close prices from October 2007 to January 2016 of 13 financial assets (ETF tickers: BND, SPY, EFA, EEM, GLD, BWX, SLV, USO, UUP, FXE, FXY, and FXA) from Yahoo finance. We calibrate a first-order vector autoregressive model to log-returns (log S t −log S t−1 ) to obtain the return dynamics. We fix the annual interest rate on the cash account to 4.5%. We consider a three-dimensional portfolio optimization with the cash component, BND and SPY. The rest of the 11 financial assets are used as exogenous return predictors for BND and SPY. In the numerical tests, we use a sample of M = 10 4 Monte Carlo paths. We assume transaction costs of 0.3% proportional to the portfolio turnover. For liquidity costs and market impact which depend on the transaction volume, we describe the price movement during the transaction by the marginal supply-demand curve (MSDC) calibrated by [17] to the European medium-and large-cap equities. The MSDC reads when q > 0 (position increase) when q < 0 (position reduction) , where S B denotes the best bid price, S A denotes the best ask price, ∆q denotes the position change of the investor at time t and k ∈ R + is called the liquidity risk factor. By integrating with respect to ∆q , we obtain an explicit parametric form for the liquidity costs, . Following the calibration results from [17], we set k = 8 × 10 −6 . The calibrated bid-ask spread S A − S B is about 10 −5 , which, according to our tests, is small enough to have virtually no impact on the optimal portfolio allocation. For convenience, we simply set S B = S A . The permanent market impact can be deemed proportional to the temporary peak of the MSDC. We assume 2/3 for the proportional rate for the market impact. For the details of how to incorporate these switching costs into the LSMC method, we refer to [20]. We use the constant relative risk aversion (CRRA) utility, i.e., U (w) = w 1−γ /(1 − γ). for the optimization objective. We report our numerical results in terms of monthly adjusted certainty equivalent returns (CER) given by The magnitude of monthly returns is usually less than one percent, thus we report CER in basis points to make comparisons easier.

Local regression v.s. global regression
As explained in Section 4.1, the optimal portfolio allocation computed on a discrete grid can be improved by control regression. Here, we compare local control regression and global control regression in terms of accuracy and efficiency, given that both use adaptive grids to find the optimum.  Table 5.1 for 4 th -order global polynomial regression and δ = 1/8 suggests that high-order polynomial bases may produce unreliable results.

Mesh size
We then test the sensitivity of the accuracy (Table 5.3) and efficiency (Table 5.4) with respect to the mesh size of the control grid. Table 5.3 shows that a coarse grid (δ = 1/8) combined with local control regression is able to produce accurate results: the errors in CER and initial portfolio allocation are negligible compared to a fine grid (δ = 1/32). Table 5.4 reports the runtime of one time step iteration of the backward dynamic programming for different control meshes. The main message coming from these results is that although we still do not completely get rid of the curse of dimensionality, the presented technique allows us to solve portfolio optimization problems of much greater size (more than ten risky assets) than what is usually considered in the literature (less than three risky assets).
It is important to note that, in order to properly compare the runtime, the code is written based on naive "for loop" for each Monte Carlo path and the computational speed can be greatly improved by vectorization or by using lower level programming languages. In particular, the code is implemented in Python 3.4.3 on a single processor Intel Core i7 2.2 GHz. This table compares the monthly adjusted CER (in basis points) and the initial optimal portfolio allocation using different sizes of control mesh (δ = 1/2, 1/4, 1/8, 1/16, 1/32) under different investment horizons (N = 3, 6, 12 months) and different risk-aversion parameters of the CRRA utility (γ = 5, 10, 15).

Conclusion
This paper designs an improved LSMC method for solving stochastic control problems. This LSMC modification can be summarized as follows: we first estimate the value functions on a discrete grid of controls, then generalize the discrete estimates by local control regression, and finally determine the optimal (continuous) control by an adaptive refinement algorithm. We apply the method to a dynamic multivariate portfolio optimization problem in the presence of transaction costs, liquidity costs and market impact. We numerically show that the local regression of continuation estimates with respect to portfolio weights is more accurate and more efficient than the classical global control regression. We show that the combination of local control regression and adaptive grids significantly improves the computational efficiency without sacrificing accuracy, and that accurate results can be obtained from very coarse grids.