Liquidity risk: optimal execution

Interactive exploration of optimal liquidation strategies under nonlinear market impact

A trader who needs to liquidate a large position faces a fundamental dilemma. Selling everything at once minimizes exposure to adverse price movements but incurs severe transaction costs, because trading a large volume widens the bid-ask spread through temporary market impact. Selling gradually reduces transaction costs but exposes the trader to the risk that prices move unfavorably while the position is still held (see Hull 2023, chap. 21).

This page explores the optimal execution model presented in Hull (2023, chap. 21), which uses a nonlinear spread function to capture how transaction costs grow with trade size. The general framework for balancing execution costs against market risk was introduced by Almgren and Chriss (2001), who derive closed-form solutions for the linear impact case.

The model

A trader holds a position of \(V\) units to be liquidated over \(n\) days. Each day \(i\), the trader sells \(q_i\) units, and the remaining position is:

\[ x_i = V - \sum_{j=1}^{i} q_j, \quad i = 1, \ldots, n, \qquad x_0 = V \]

The bid-ask spread

The cash bid-ask spread depends on how many units are traded. Following Hull (2023, chap. 21), we model the spread as a nonlinear function of trade size:

\[ p(q) = a + b\, e^{c\,q} \]

the spread widens exponentially as the daily trading volume \(q\) increases. The parameters are:

• \(a\): the base spread when very few units are traded
• \(b\) and \(c\): control how rapidly the spread grows with volume

Each trade of \(q_i\) units costs \(q_i \times p(q_i)/2\) in transaction costs (half the spread per unit, times the number of units).

The objective

The trader chooses \(q_1, \ldots, q_n\) to minimize the cost at a given confidence level:

\[ \min_{\{q_i\}} \quad \underbrace{\lambda \, \sigma \sqrt{\sum_{i=1}^{n} x_i^2}}_{\text{Market risk}} + \sum_{i=1}^{n} \underbrace{q_i \, \frac{p(q_i)}{2}}_{\text{Transaction cost}} \]

subject to \(\sum q_i = V\) and \(q_i \geq 0\), where:

• \(\sigma\) is the standard deviation of the daily price change per unit
• \(\lambda = \Phi^{-1}(\text{confidence level})\) converts the desired confidence into a number of standard deviations

The first term is the market risk component — the Value-at-Risk of the price changes that occur while the position is held. The second term is the expected transaction cost from bid-ask spreads.

Note

Two extreme strategies

• Sell uniformly (\(q_i = V/n\)): minimizes transaction costs (by exploiting the convexity of \(q \times p(q)\)) but carries the most market risk.
• Sell immediately (\(q_1 = V\), \(q_i = 0\) for \(i > 1\)): eliminates market risk but incurs very large transaction costs due to the high spread on a massive trade.

The optimal strategy lies between these extremes, and its shape depends on the confidence level. At the 50th percentile (\(\lambda = 0\)), the trader ignores risk and sells uniformly. As the confidence level rises, the strategy becomes more front-loaded — selling more in early days to reduce exposure.

Note

The linear impact case

When the spread function is linear in \(q\), Almgren and Chriss (2001) show that the optimal trajectory has an elegant closed-form solution involving hyperbolic functions: \(x_j \propto \sinh(\kappa(T - t_j))\), where the parameter \(\kappa\) captures the trade-off between risk and cost. With the nonlinear exponential spread from Hull (2023), numerical optimization is required.

Interactive exploration

Adjust the parameters below to see how the optimal execution strategy responds to changes in position size, market conditions, and risk tolerance.

Tip

How to experiment

1. Raise the confidence level from 50% (uniform selling) toward 99% to see the strategy become more front-loaded.
2. Increase volatility \(\sigma\) to make market risk more expensive, pushing the trader to sell faster.
3. Increase the spread growth rate \(c\) to make large trades costlier, pushing the trader to spread sales more evenly.
4. Check the efficient frontier tab to see the full risk–cost trade-off.

Position

viewof V_val = Inputs.range([10, 500], {
  label: "Position V (million units)",
  step: 10,
  value: 100
})

viewof n_days = Inputs.range([2, 20], {
  label: "Liquidation horizon n (days)",
  step: 1,
  value: 5
})

Market parameters

viewof sigma = Inputs.range([0.01, 1], {
  label: "Daily price volatility σ ($/unit)",
  step: 0.01,
  value: 0.10
})

viewof a_val = Inputs.range([0, 1], {
  label: "Base spread a ($)",
  step: 0.01,
  value: 0.10
})

viewof b_val = Inputs.range([0.01, 0.5], {
  label: "Spread coefficient b ($)",
  step: 0.01,
  value: 0.05
})

viewof c_val = Inputs.range([0.005, 0.2], {
  label: "Spread growth rate c",
  step: 0.005,
  value: 0.03
})

Risk tolerance

viewof conf_pct = Inputs.range([50, 99.9], {
  label: "Confidence level (%)",
  step: 0.5,
  value: 95
})

// --- Inverse normal CDF (Abramowitz & Stegun 26.2.23) ---
qnorm = p => {
  if (p <= 0) return -Infinity
  if (p >= 1) return Infinity
  if (p === 0.5) return 0
  const u = p < 0.5 ? p : 1 - p
  const t = Math.sqrt(-2 * Math.log(u))
  const c0 = 2.515517, c1 = 0.802853, c2 = 0.010328
  const d1 = 1.432788, d2 = 0.189269, d3 = 0.001308
  const z = t - (c0 + c1 * t + c2 * t * t) / (1 + d1 * t + d2 * t * t + d3 * t * t * t)
  return p < 0.5 ? -z : z
}

lambdaV = qnorm(conf_pct / 100)

// --- Optimizer: BFGS with analytical gradients and softmax parameterization ---
optimizeTrade = (V, n, a, b, c, sig, lam) => {
  // λ ≈ 0 → uniform is optimal
  if (lam < 1e-12) return new Array(n).fill(V / n)
  if (n === 1) return [V]

  const pf = q => a + b * Math.exp(c * q)
  const dpf = q => b * c * Math.exp(c * q)

  function softmax(theta) {
    const mx = Math.max(...theta)
    const ex = theta.map(t => Math.exp(t - mx))
    const sm = ex.reduce((s, e) => s + e, 0)
    return ex.map(e => V * e / sm)
  }

  function evalObjGrad(theta) {
    const q = softmax(theta)
    const x = new Array(n)
    let cum = 0
    for (let i = 0; i < n; i++) { cum += q[i]; x[i] = V - cum }

    const sumX2 = x.reduce((s, xi) => s + xi * xi, 0)
    const denom = Math.sqrt(sumX2 + 1e-12)

    // suffix sums of x
    const sx = new Array(n)
    sx[n - 1] = x[n - 1]
    for (let i = n - 2; i >= 0; i--) sx[i] = sx[i + 1] + x[i]

    // gradient w.r.t. q
    const gq = new Array(n)
    for (let j = 0; j < n; j++) {
      gq[j] = -lam * sig * sx[j] / denom + (pf(q[j]) + q[j] * dpf(q[j])) / 2
    }

    // weighted mean for softmax chain rule
    let gBar = 0
    for (let j = 0; j < n; j++) gBar += q[j] * gq[j]
    gBar /= V

    // gradient w.r.t. theta
    const grad = new Array(n)
    for (let j = 0; j < n; j++) grad[j] = q[j] * (gq[j] - gBar)

    // objective value
    let tx = 0
    for (let i = 0; i < n; i++) tx += q[i] * pf(q[i]) / 2
    const obj = lam * sig * denom + tx

    return { obj, grad, q }
  }

  function bfgs(theta0) {
    let theta = [...theta0]
    let { obj, grad } = evalObjGrad(theta)

    // Initialize inverse Hessian as identity
    let H = Array.from({length: n}, (_, i) => {
      const row = new Array(n).fill(0)
      row[i] = 1
      return row
    })

    for (let iter = 0; iter < 500; iter++) {
      // Convergence check
      const gNorm = Math.sqrt(grad.reduce((s, g) => s + g * g, 0))
      if (gNorm < 1e-12) break

      // Search direction d = -H * grad
      const d = new Array(n).fill(0)
      for (let i = 0; i < n; i++)
        for (let j = 0; j < n; j++)
          d[i] -= H[i][j] * grad[j]

      // Directional derivative
      const dg = d.reduce((s, di, i) => s + di * grad[i], 0)
      if (dg >= 0) {
        // Reset inverse Hessian to identity
        H = Array.from({length: n}, (_, i) => {
          const row = new Array(n).fill(0)
          row[i] = 1
          return row
        })
        continue
      }

      // Backtracking line search (Armijo condition)
      let alpha = 1.0
      let newTheta, newObj, newGrad
      for (let ls = 0; ls < 40; ls++) {
        newTheta = theta.map((t, i) => t + alpha * d[i])
        const res = evalObjGrad(newTheta)
        newObj = res.obj
        newGrad = res.grad
        if (newObj <= obj + 1e-4 * alpha * dg) break
        alpha *= 0.5
      }

      // s = x_{k+1} - x_k, y = g_{k+1} - g_k
      const s = newTheta.map((t, i) => t - theta[i])
      const y = newGrad.map((g, i) => g - grad[i])
      const ys = y.reduce((sum, yi, i) => sum + yi * s[i], 0)

      // BFGS update of inverse Hessian (skip if curvature condition not met)
      if (ys > 1e-12) {
        const Hy = new Array(n).fill(0)
        for (let i = 0; i < n; i++)
          for (let j = 0; j < n; j++)
            Hy[i] += H[i][j] * y[j]
        const yHy = y.reduce((sum, yi, i) => sum + yi * Hy[i], 0)
        const rho = 1 / ys
        for (let i = 0; i < n; i++)
          for (let j = 0; j < n; j++)
            H[i][j] += -rho * (Hy[i] * s[j] + s[i] * Hy[j]) + rho * (1 + yHy * rho) * s[i] * s[j]
      }

      theta = newTheta
      obj = newObj
      grad = newGrad
    }

    return { obj, q: evalObjGrad(theta).q }
  }

  // Multiple starting points
  const starts = [
    new Array(n).fill(0),
    Array.from({length: n}, (_, i) => 2 * (1 - 2 * i / Math.max(n - 1, 1))),
    Array.from({length: n}, (_, i) => 5 * (1 - 2 * i / Math.max(n - 1, 1))),
    Array.from({length: n}, (_, i) => 10 * (1 - 2 * i / Math.max(n - 1, 1)))
  ]

  let bestQ = new Array(n).fill(V / n)
  let bestVal = Infinity

  for (const th0 of starts) {
    const res = bfgs(th0)
    if (res.obj < bestVal) { bestVal = res.obj; bestQ = [...res.q] }
  }
  return bestQ
}

// --- Cost metrics for a given strategy ---
evalStrategy = (q, V, a, b, c, sig, lam) => {
  const n = q.length
  const pf = qq => a + b * Math.exp(c * qq)
  const x = new Array(n)
  let cum = 0
  for (let i = 0; i < n; i++) { cum += q[i]; x[i] = V - cum }
  const sumX2 = x.reduce((s, xi) => s + xi * xi, 0)
  const mktSD = sig * Math.sqrt(sumX2)
  const mktVaR = lam * mktSD
  let txCost = 0
  for (let i = 0; i < n; i++) txCost += q[i] * pf(q[i]) / 2
  const totalVaR = mktVaR + txCost
  let avgTime = 0
  for (let i = 0; i < n; i++) avgTime += (i + 1) * q[i] / V
  return { x, mktSD, mktVaR, txCost, totalVaR, avgTime }
}

optQ = optimizeTrade(V_val, n_days, a_val, b_val, c_val, sigma, lambdaV)
uniformQ = new Array(n_days).fill(V_val / n_days)
optStats = evalStrategy(optQ, V_val, a_val, b_val, c_val, sigma, lambdaV)
uniStats = evalStrategy(uniformQ, V_val, a_val, b_val, c_val, sigma, lambdaV)

confLevels = [50, 75, 90, 95, 99]
multiConf = {
  return confLevels.map(cl => {
    const lam = qnorm(cl / 100)
    const q = optimizeTrade(V_val, n_days, a_val, b_val, c_val, sigma, lam)
    const stats = evalStrategy(q, V_val, a_val, b_val, c_val, sigma, lam)
    return { cl, lam, q, stats }
  })
}

// --- Efficient frontier: sweep lambda from 0 to ~4 ---
frontier = {
  const pts = []
  const nPts = 40
  const lamMax = 4
  for (let i = 0; i <= nPts; i++) {
    const lam = (i / nPts) * lamMax
    const q = optimizeTrade(V_val, n_days, a_val, b_val, c_val, sigma, lam)
    const s = evalStrategy(q, V_val, a_val, b_val, c_val, sigma, lam)
    pts.push({ lam, txCost: s.txCost, mktSD: s.mktSD })
  }
  return pts
}

fmtM = (v, d) => v.toFixed(d === undefined ? 2 : d)

Trading strategy

tradeData = {
  const data = []
  for (let i = 0; i < n_days; i++) {
    data.push({ day: i + 1, fraction: optQ[i] / V_val, strategy: "Optimal" })
    data.push({ day: i + 1, fraction: uniformQ[i] / V_val, strategy: "Uniform" })
  }
  return data
}

Plot.plot({
  height: 320,
  marginLeft: 50,
  marginRight: 20,
  fx: { label: "Day", padding: 0.2, sort: (a, b) => a - b },
  x: { axis: null, padding: 0.1 },
  y: {
    label: "Fraction of position sold",
    grid: true,
    domain: [0, Math.max(...tradeData.map(d => d.fraction)) * 1.1]
  },
  color: {
    legend: true,
    domain: ["Optimal", "Uniform"],
    range: ["#2f71d5", "#bbb"]
  },
  marks: [
    Plot.barY(tradeData, {
      fx: "day",
      x: "strategy",
      y: "fraction",
      fill: "strategy",
      fillOpacity: 0.8,
      tip: true,
      title: d => `Day ${d.day} (${d.strategy}): ${(d.fraction * 100).toFixed(2)}% of position`
    })
  ]
})

The uniform strategy (gray) sells the same fraction each day. The optimal strategy (blue) is front-loaded — it sells more in early days to reduce exposure to market risk, accepting higher transaction costs on those larger trades. The degree of front-loading increases with the confidence level.

Position trajectory

trajData = {
  const data = []
  // time 0: full position
  data.push({ day: 0, fraction: 1, strategy: "Optimal" })
  data.push({ day: 0, fraction: 1, strategy: "Uniform" })
  let cumOpt = 0, cumUni = 0
  for (let i = 0; i < n_days; i++) {
    cumOpt += optQ[i]; cumUni += uniformQ[i]
    data.push({ day: i + 1, fraction: (V_val - cumOpt) / V_val, strategy: "Optimal" })
    data.push({ day: i + 1, fraction: (V_val - cumUni) / V_val, strategy: "Uniform" })
  }
  return data
}

Plot.plot({
  height: 320,
  marginLeft: 50,
  marginRight: 20,
  x: { label: "Day", grid: false },
  y: { label: "Fraction of position remaining", domain: [0, 1], grid: true },
  color: {
    legend: true,
    domain: ["Optimal", "Uniform"],
    range: ["#2f71d5", "#bbb"]
  },
  marks: [
    Plot.line(trajData, {
      x: "day",
      y: "fraction",
      stroke: "strategy",
      strokeWidth: d => d.strategy === "Optimal" ? 2.5 : 1.5,
      strokeDasharray: d => d.strategy === "Uniform" ? "6,4" : undefined
    }),
    Plot.dot(trajData.filter(d => d.strategy === "Optimal"), {
      x: "day",
      y: "fraction",
      fill: "#2f71d5",
      r: 4,
      tip: true,
      title: d => `Day ${d.day}: ${(d.fraction * 100).toFixed(1)}% remaining`
    })
  ]
})

The position decays linearly under the uniform strategy (dashed gray) but follows a concave path under the optimal strategy (solid blue) — the trader reduces the position rapidly at first, then more slowly as the remaining exposure shrinks.

Cost analysis

html`<table class="table" style="width: 100%;">
<thead><tr>
  <th>Metric</th>
  <th style="text-align:right;">Optimal</th>
  <th style="text-align:right;">Uniform</th>
</tr></thead>
<tbody>
<tr>
  <td>Transaction cost ($M)</td>
  <td style="text-align:right;">${fmtM(optStats.txCost)}</td>
  <td style="text-align:right;">${fmtM(uniStats.txCost)}</td>
</tr>
<tr>
  <td>Market risk std. dev. ($M)</td>
  <td style="text-align:right;">${fmtM(optStats.mktSD)}</td>
  <td style="text-align:right;">${fmtM(uniStats.mktSD)}</td>
</tr>
<tr>
  <td>Market risk VaR at ${fmtM(conf_pct, 1)}% ($M)</td>
  <td style="text-align:right;">${fmtM(optStats.mktVaR)}</td>
  <td style="text-align:right;">${fmtM(uniStats.mktVaR)}</td>
</tr>
<tr style="font-weight: bold; border-top: 2px solid #666;">
  <td>Total cost at ${fmtM(conf_pct, 1)}% confidence ($M)</td>
  <td style="text-align:right;">${fmtM(optStats.totalVaR)}</td>
  <td style="text-align:right;">${fmtM(uniStats.totalVaR)}</td>
</tr>
<tr>
  <td>Average selling time (days)</td>
  <td style="text-align:right;">${fmtM(optStats.avgTime)}</td>
  <td style="text-align:right;">${fmtM(uniStats.avgTime)}</td>
</tr>
</tbody></table>

<p style="margin-top: 0.5rem; font-size: 0.9rem; color: #666;">
  λ = Φ⁻¹(${fmtM(conf_pct, 1)}%) = ${fmtM(lambdaV, 4)}.
  The optimal strategy reduces total cost by <strong>$${fmtM(uniStats.totalVaR - optStats.totalVaR)}M</strong>
  compared to uniform selling (a ${fmtM((1 - optStats.totalVaR / uniStats.totalVaR) * 100, 1)}% improvement).
</p>`

The optimal strategy has higher transaction costs but lower market risk than the uniform strategy. The net effect is a lower total cost at the chosen confidence level, because the risk reduction more than offsets the extra transaction costs.

Note

How market risk VaR is computed

Each day \(i\), the remaining position \(x_i\) is exposed to a random price change with standard deviation \(\sigma\). Assuming daily price changes are independent and normally distributed, the variance of the total P&L from price movements is \(\sigma^2 \sum_{i=1}^{n} x_i^2\), so the standard deviation is \(\sigma \sqrt{\sum x_i^2}\). The market risk VaR at a given confidence level is then \(\lambda \, \sigma \sqrt{\sum x_i^2}\), where \(\lambda = \Phi^{-1}(\text{confidence level})\).

Across confidence levels

multiConfData = {
  const data = []
  for (const mc of multiConf) {
    for (let i = 0; i < n_days; i++) {
      data.push({
        day: i + 1,
        fraction: mc.q[i] / V_val,
        conf: `${mc.cl}%`
      })
    }
  }
  return data
}

Plot.plot({
  height: 340,
  marginLeft: 50,
  marginRight: 20,
  x: { label: "Day" },
  y: {
    label: "Fraction of position sold",
    grid: true,
    domain: [0, Math.max(...multiConfData.map(d => d.fraction)) * 1.1]
  },
  color: { legend: true },
  marks: [
    Plot.line(multiConfData, {
      x: "day",
      y: "fraction",
      stroke: "conf",
      strokeWidth: 2,
      marker: "dot"
    })
  ]
})

html`<table class="table" style="width: 100%;">
<thead><tr>
  <th>Confidence</th>
  <th style="text-align:right;">λ</th>
  <th style="text-align:right;">Tx cost ($M)</th>
  <th style="text-align:right;">Market risk σ ($M)</th>
  <th style="text-align:right;">Total cost ($M)</th>
  <th style="text-align:right;">Avg. selling time</th>
</tr></thead>
<tbody>
${multiConf.map(mc => `<tr>
  <td>${mc.cl}%</td>
  <td style="text-align:right;">${fmtM(mc.lam, 3)}</td>
  <td style="text-align:right;">${fmtM(mc.stats.txCost)}</td>
  <td style="text-align:right;">${fmtM(mc.stats.mktSD)}</td>
  <td style="text-align:right;">${fmtM(mc.stats.totalVaR)}</td>
  <td style="text-align:right;">${fmtM(mc.stats.avgTime)} days</td>
</tr>`).join("")}
</tbody></table>`

At 50% confidence (\(\lambda = 0\)), the trader only cares about expected costs and sells uniformly. As the confidence level rises, the strategy becomes progressively more front-loaded, selling a larger fraction of the position in the first few days. The average selling time decreases from \(\frac{n+1}{2}\) (uniform) toward 1 (sell everything immediately).

Efficient frontier

currentFrontierPt = ({ mktSD: optStats.mktSD, txCost: optStats.txCost, label: `${fmtM(conf_pct, 0)}%` })
uniformFrontierPt = ({ mktSD: uniStats.mktSD, txCost: uniStats.txCost, label: "Uniform" })

// VaR tangent line through optimal point
varTangent = {
  const xMax = frontier[0].mktSD * 1.15
  return [
    { mktSD: 0, txCost: optStats.totalVaR },
    { mktSD: xMax, txCost: optStats.totalVaR - lambdaV * xMax }
  ].filter(d => d.txCost >= 0)
}

Plot.plot({
  height: 360,
  marginLeft: 60,
  marginRight: 20,
  x: {
    label: "Market risk std. dev. ($M)",
    grid: true,
    domain: [0, frontier[0].mktSD * 1.15]
  },
  y: {
    label: "Expected transaction cost ($M)",
    grid: true
  },
  marks: [
    // Efficient frontier curve
    Plot.line(frontier, {
      x: "mktSD",
      y: "txCost",
      stroke: "#4682b4",
      strokeWidth: 2.5
    }),
    // VaR iso-cost line (tangent)
    Plot.line(varTangent, {
      x: "mktSD",
      y: "txCost",
      stroke: "#ff7f0e",
      strokeWidth: 1.5,
      strokeDasharray: "6,4"
    }),
    // Current strategy point
    Plot.dot([currentFrontierPt], {
      x: "mktSD",
      y: "txCost",
      fill: "#d62728",
      r: 6,
      tip: true,
      title: d => `Confidence: ${d.label}\nTx cost: $${d.txCost.toFixed(2)}M\nMkt risk σ: $${d.mktSD.toFixed(2)}M`
    }),
    Plot.text([currentFrontierPt], {
      x: "mktSD", y: "txCost", text: "label",
      fill: "#d62728", fontWeight: "bold", fontSize: 11, dy: -12
    }),
    // Uniform strategy point
    Plot.dot([uniformFrontierPt], {
      x: "mktSD",
      y: "txCost",
      fill: "#999",
      r: 5,
      symbol: "square"
    }),
    Plot.text([uniformFrontierPt], {
      x: "mktSD", y: "txCost", text: "label",
      fill: "#999", fontSize: 11, dy: 14
    })
  ]
})

Each point on the efficient frontier (blue curve) represents the optimal strategy for a different level of risk aversion. The bottom-right is the uniform strategy (minimum transaction cost, maximum market risk). Moving left along the curve, the trader sells more aggressively, reducing market risk but increasing transaction costs. The red dot marks the optimal strategy for the current confidence level; the dashed orange line shows its VaR iso-cost contour.

Solving in R

The optimization can also be solved in R using the optim function. The code below runs directly in the browser, you can edit the parameters and re-run each block.

Parameters

Set up the problem inputs. Edit these to explore different scenarios.

Spread function

The nonlinear cash bid-ask spread from Hull (2023):

Softmax parameterization

To enforce the constraints \(q_i > 0\) and \(\sum q_i = V\), we optimize over unconstrained parameters \(\theta_i\) and recover the trade quantities via softmax:

\[ q_j = V \frac{e^{\theta_j}}{\sum_k e^{\theta_k}} \]

This guarantees feasibility for any values of \(\theta\).

Objective function

The objective combines market risk and transaction cost. It takes the unconstrained \(\theta\) vector as input and returns the total cost at the chosen confidence level.

Optimization

We run BFGS from multiple starting points to avoid local minima. Each starting point uses a different degree of front-loading (from uniform to aggressive).

Results

Compute cost metrics and display the optimal strategy.

References

Almgren, Robert, and Neil Chriss. 2001. “Optimal Execution of Portfolio Transactions.” Journal of Risk 3 (2): 5–39.

Hull, John. 2023. Risk Management and Financial Institutions. 6th ed. New Jersey: John Wiley & Sons.