Sampling Random Portfolios Under Linear Constraints
Summary
The document frames constrained random portfolio generation as sampling weights that satisfy linear equality and inequality conditions. Equality constraints can be maintained by taking steps within the null space of the equality matrix, for example through a Monte Carlo random walk. Inequality constraints restrict the feasible region, and the question’s proposed approach reflects steps that cross its boundaries.
The response notes that reflection can be costly because a simplex has many faces to check, and walks near corners may require many rejected or corrected moves. It mentions adaptive step sizes as a possible aid, while cautioning that this can bias samples. For an unconstrained unit simplex, it gives a direct construction: sort uniform random values between zero and one, append the endpoints, and use successive gaps as weights; those weights sum to one and can be scaled for a general simplex. This does not solve arbitrary intersections of constraints, and the discussion distinguishes continuous weights from the harder integer-weight case without specifying a complete general-purpose sampler.
Key ideas
- Equality constraints can be preserved by moving within the equality matrix’s null space.
- Inequality constraints define boundaries that a random walk must handle to stay feasible.
- Reflection can become computationally expensive as the number of boundaries and corner interactions grows.
- Sorted uniform draws and their successive gaps provide a direct way to sample weights on the unit simplex.
- Adaptive step sizing may help avoid boundary crossings but can bias the resulting samples.
Tags
Full text
# Is creating constrained random portfolios a hard problem?
# Is creating constrained random portfolios a hard problem?
Creating random portfolios with weights $x_i$ can be thought of as sampling from the surface of a simplex given by $$Ex = f$$ and $$Ax \le b$$ Where $E$ and $A$ are constraint matrices for equality and inequality constraints, and $f$ and $b$ are solutions for some portfolio you want to match on, $x_0$. Adjusting the previous equations to reflect this gives us: $$Ex = f = Ex_0$$ and $$Ax \le b = Ax_0.$$
The equality constraints are relatively easy to satisfy, the solution is simply $x = x_0 + Z r$, where $Zr$ is in the null-space of $E$, and this can be the basis of a Monte Carlo random-walk. However there is no straightforward, mathematical way to link the inequality constraints with the equality constraints. The solution seems to be to check if your walk takes you over the boundary, and then reflect back over it, however as the number of $x$'s gets very large, the number of faces to reflect over grows with $n$ and this starts to pose a problem with computation time. Is there a better way to do this algorithm, or is this problem just hard?
Edit: Here is some sample R code with the random walk and reflecting over the boundaries, only handing the inequality constraint that all $x$'s must be positive:
```
require(MASS)
getWeights <- function(Emat, x0, n, verbose = FALSE) {
Z = Null(t(Emat))
ret = matrix(0, nrow = length(x0), ncol = n + 1)
## Would it be better to use apply here?
nc = ncol(Z)
mn = mean(x0)
ret[, 1] = x0 + Z %*% rnorm(nc, 0, mn)/sqrt(nc)
k = 0
if(verbose) cat("Created Vectors: ")
if(verbose) cat(paste(k))
for (i in 2:(n + 1)) {
ret[, i] = ret[, i - 1] + Z %*% rnorm(nc, 0, mn)/sqrt(nc)
m = k + 1;
while(any(ret[, i] < 0)) {
reflection = rep(0, ncol(Emat))
reflection[which(ret[, i] < 0)] = ret[, i][which(ret[, i] < 0)]
for (j in 1:ncol(Z)) {
ret[, i] = ret[, i] - 2* Z[, j] * (reflection %*% Z[, j])/sqrt(Z[,
j] %*% Z[, j])
}
##for(i in 1:nchar(paste(k))) cat("\b")
##if(verbose) cat(paste(m))
##k = m
##m = k + 1
}
if(verbose) for(i in 1:nchar(paste(k))) cat("\b")
if(verbose) cat(paste(m))
k = m
}
ret = ret[, 2:(n + 1)]
if(verbose) cat("\n")
return(ret)
}
Emat = matrix(1, ncol = 1000, nrow = 1)
x0 = rep(1/1000, 1000)
w = getWeights(Emat, x0, 1000, TRUE)
```
As you can see this code simply goes too slow. (Do you think it would be best to try and implement this code faster, using C and multiple cores, or would it be more worthwhile changing the algorithm?)
## Answer by user1157 (score 3)
https://quant.stackexchange.com/a/10239
The approach of reflecting is expensive, since the $d$-simplex has $d$ maximal faces, all of which have to be checked for intersection at each step. Additionally, if the random walk moves into a corner, the number of moves which have to be discarded can become very high. Depending on the configuration of the constraints this could well be your best solution.
If I understood your question correctly, the question boils down to sampling from the intersection of the boundary of one simplex ($E$) and the volume of another ($A$) simultaneously.
Suggestions:
- Adaptive step-size: You can calculate the distance from the boundary of $A$ and $E$ and choose half the step size, of course this will add a bias to the generation proces, this could be a problem.
Algorithm for one simplex: There is a nice trick for the unit simplex by Donald Rubin, which is discussed on the computer science site, to sample from the unit simplex. The result can be scaled afterwards to match a general simplex:
- Create $n-1$ uniformly distributed random values between 0 and 1 and sort them. [Example: 0.1, 0.3, 0.5, 0.55]
- Add 0 and 1 to the list. [Example: 0, 0.1, 0.3, 0.5, 0.55, 1]
- Use the differences between the values, they will sum to 1. [Example: 0.1, 0.2, 0.2, 0.05, 0.45]
Edit: Strictly speaking it is not a hard problem, since it can be solved in polynomial time; it would be though if you were looking for integer weights.Shown in full with attribution under the source's licence. Licence: CC BY-SA 4.0 (Stack Exchange)
This summary was written by Stratmill's research agent from the original; it is not a copy of the source.