Algorithm / formula / method to determine optimal weightings given expected return, % of volume and slippage

Question

Please bear with me - I know I'm supposed to do this with a bunch of Greek letters but I don't know how so I'll just describe the data I have and what I am trying to do.

I have an expected return for the next month for a universe of say 200 stocks. I also have the volume over the past month for each stock. Finally I have an estimate for the slippage I will experience, given the % of the daily volume that I try to trade.

For example, stock XYZ has an expected return of 2% over the next month, the highest in the universe. So I would want to put lots of weight on this stock. But it only trades about $100k per day and my estimated is slippage is say 5 bps for 0.1% of daily volume, 25bps for 1% of daily volume, 50bps for 2% of daily volume, etc. As more weight is added to this stock, the slippage adjusted expected return drops.

Then I have stock ABC with an expected return of -0.5% over the next month, at the low end of the universe and negative return so I would short it. And it trades $100m per day, so slippage would be much less. I'd want to short this stock and the expected return would not drop much due to slippage given the big liquidity.

If my goal is to maximize the absolute dollar return on my portfolio given an equity amount of say 250k, is there an algorithm / methodology I can use to optimize my portfolio that incorporate the expected return and the slippage and spits out dollar neutral weightings for each stock in my universe?

It seems to me like it has to be an iterative calculation, because each weighting affects all the other weightings?

Answer 1

Let $k$ be the funding amount (in \$), let $N$ be the number of stocks, let $V$ be the vector of expected daily volumes $v_i$, let $R$ be the vector of expected returns $r_i$, let $W$ be a vector of relative position weights $w_i$, summing to $1$ in absolute values, let $w < 0$ denote short positions and let $f$ be a function mapping from a position of size $wk$ to a slippage factor, given corresponding $v$. For simplicity and demonstration, $f$ is assumed to be symmetric around the origin (here $f(w, v) = \beta\frac{wk}{v}$, aligning with the provided estimated slippages for $\beta=2.5$) and the slippage to be additive to the expected return (however not changing the sign of the return which slippage was accounted for).

Then the objective is to maximize the total return

$$\sum_{i=1}^N\left(w_i(r_i - f(w_i, v_i))\right)$$

subject to

$$\sum_{i=1}^N|w_i| \leq 1$$

(assuming no leverage).

The following is an implementation in Python that optimizes the weights $W$ (staring at zero weights) with respect to the objective using SLSQP with SciPy, here simulating $N = 200$ random expected returns $R$ in range $[-2\%, 2\%]$ and expected daily volumes in range $[100k, 100M]$, assuming a funding of $\\\$250k$.

The cases of assigning random weights (yielding an approximately zero return) and of using uniform weights with signs equal to the expected return (yielding an approximately $1\%$ return) are intended as benchmarks. The optimization converges and the optimal weights yield a return fairly close to the $2\%$.

from scipy.optimize import minimize, Bounds

np.random.seed(42)

k = 2.5e5
N = 200
beta = 2.5
V = np.random.uniform(1e5, 1e8, N)
R = np.random.uniform(-0.02, 0.02, N)

def slippage(x, V, beta=beta):
    return beta * x * k / V

def objective(x, R=R, V=V, k=k):
    slips = slippage(x=x, V=V)
    rets = R - slips
    rets[np.sign(rets) != np.sign(R)] = 0
    return -sum(x * rets)

def constraint(x):
    return 1 - sum(np.abs(x))

random_W = np.random.uniform(-1, 1, N)
random_W /= sum(np.abs(random_W))
uniform_W = np.sign(R) / N

sol = minimize(
    fun=objective,
    x0=np.full(N, 0),
    method='SLSQP',
    bounds=Bounds(-1, 1),
    constraints=[{'type': 'ineq', 'fun': constraint}],
    options={'maxiter': 690}
)
assert sol.success
assert sum(np.abs(sol.x)) < 1 + 1e-6
assert -1 < min(sol.x) and max(sol.x) < 1
best_W = sol.x

for name, W in [('random', random_W), ('uniform', uniform_W), ('optimal', best_W)]:
    ret = -objective(W)
    print(f'{name} weights return: {int(ret * k)}$ ({ret * 100:.3f}%)')```

#The code prints the following: 

random  weights return: -482$ (-0.193%)
uniform weights return: 2485$ ( 0.994%)
optimal weights return: 4549$ ( 1.820%)