Monte Carlo Option Pricing Simulation

Pricing by simulation

Under the risk-neutral measure the stock follows geometric Brownian motion with drift equal to the risk-free rate r. A European option only cares about the price at maturity T, and that price can be sampled exactly in one step, with no time grid:

The price is the discounted mean payoff over n simulated values, max(S_T − K, 0) for a call and max(K − S_T, 0) for a put:

where s is the sample standard deviation of the payoffs. The drift is r, not the stock's expected return: pricing by replication makes the real-world drift irrelevant, one of the central results of Black and Scholes. The histogram shows the lognormal spread of S_T; only the orange bars contribute to the price, weighted by the payoff line. Out-of-the-money options have few paths that pay, which is why their relative error is large — raise the strike of a call to see it.

The Black–Scholes check

For this model the expectation has a closed form. With Φ the standard normal distribution function,

The lab evaluates Φ through the complementary error function (a power series below 2.5 and a continued fraction above, accurate to about 10⁻¹⁵), so the reference price is exact to far more digits than the simulation can reach. The Monte Carlo error then behaves like any sample mean: it shrinks as 1/√n, and the 95% interval should contain the Black–Scholes price for about 19 seeds in 20. Φ(d₂) is also the risk-neutral probability that a call finishes in the money (Φ(−d₂) for a put), which the last readout estimates.

Nobody needs simulation for a plain European option. The point of checking against the formula is trust: the same code, with the payoff replaced, prices options that have no formula — averages over a path (Asian options), barriers, baskets of several assets — which is why Boyle introduced the method in 1977.

Antithetic variates

If Z is standard normal, so is −Z. The antithetic estimator uses each draw twice, averaging the payoffs at Z and −Z, and treats the pair average as one sample:

Compared at equal cost — the same number of simulated prices — it wins whenever that covariance is negative, which holds for any payoff that is monotone in Z, such as a call or a put. With the default parameters and 400,000 simulated prices, it cuts the standard error by about 29% for the call and 23% for the put. For a payoff that is symmetric in Z, such as a straddle near the money, it can do worse than plain sampling. Switch the Variance reduction control and compare the standard error readout at the same seed.

Example result

For an at-the-money call (S₀ = K = 100, r = 5%, σ = 20%, T = 1 year) with 10,000 paths and seed 42, the Monte Carlo price is 10.6060 with a standard error of 0.1468, a 95% interval from 10.318 to 10.894. The Black–Scholes price is 10.4506, so the error is 0.1554, and 56.8% of the simulated prices finish in the money.

Implementation

The page runs simulation.js from the package below; the Python file is a line-by-line port with the same random-number generator, so a given seed produces the same result in both languages. The excerpts show the core of each.

// Risk-neutral GBM, sampled exactly at maturity: S_T = S₀·exp((r − σ²/2)T + σ√T·Z).
// The price is the discounted mean payoff, e^(−rT)·E[payoff(S_T)].
export const payoff = (type, s, k) => (type === 'put' ? Math.max(k - s, 0) : Math.max(s - k, 0));

export function simulate(st, count, rng, onPrice) {
  const { p, drift, vol } = st;
  const end = Math.min(st.m, st.i + count);
  while (st.i < end) {
    const z = rng.normal();
    const s1 = p.s0 * Math.exp(drift + vol * z);
    let y = payoff(p.option_type, s1, p.strike);
    let prices = [s1];
    if (p.variance_reduction === 'antithetic') {
      // Reuse −Z: the pair average has lower variance when the payoff is monotone in Z.
      const s2 = p.s0 * Math.exp(drift - vol * z);
      y = (y + payoff(p.option_type, s2, p.strike)) / 2;
      prices = [s1, s2];
    }
    for (const s of prices) { st.prices++; if (payoff(p.option_type, s, p.strike) > 0) st.itm++; if (onPrice) onPrice(s); }
    st.i++;
    const d = y - st.mean; // Welford running mean / variance of the (pair-averaged) payoff
    st.mean += d / st.i;
    st.m2 += d * (y - st.mean);
  }
}

export function summarize(st) {
  const price = st.df * st.mean;
  const se = st.i > 1 ? st.df * Math.sqrt(st.m2 / (st.i - 1) / st.i) : null;
  const { p } = st;
  const bs = blackScholes(p.option_type, p.s0, p.strike, p.rate, p.sigma, p.maturity);
  return {
    price,
    std_error: se,
    ci_low: se === null ? null : price - 1.96 * se,
    ci_high: se === null ? null : price + 1.96 * se,
    bs_price: bs,
    error: Math.abs(price - bs),
    itm_fraction: st.prices ? st.itm / st.prices : 0,
    simulated_prices: st.prices,
    covers: se === null ? null : Math.abs(price - bs) <= 1.96 * se,
  };
}

def payoff(kind, s, k):
    return max(k - s, 0) if kind == "put" else max(s - k, 0)


def price_option(p, rng, every=None, rows=None):
    """Discounted mean payoff under risk-neutral GBM, S_T = S0*exp((r - s^2/2)T + s*sqrt(T)*Z)."""
    kind, s0, k = p["option_type"], p["s0"], p["strike"]
    anti = p["variance_reduction"] == "antithetic"
    m = p["paths"] // 2 if anti else p["paths"]
    drift = (p["rate"] - p["sigma"] * p["sigma"] / 2) * p["maturity"]
    vol = p["sigma"] * math.sqrt(p["maturity"])
    df = math.exp(-p["rate"] * p["maturity"])
    mean = m2 = 0.0
    itm = prices = 0
    for i in range(1, m + 1):
        z = rng.normal()
        s1 = s0 * math.exp(drift + vol * z)
        y = payoff(kind, s1, k)
        sims = [s1]
        if anti:  # reuse -Z: lower variance when the payoff is monotone in Z
            s2 = s0 * math.exp(drift - vol * z)
            y = (y + payoff(kind, s2, k)) / 2
            sims = [s1, s2]
        for s in sims:
            prices += 1
            if payoff(kind, s, k) > 0:
                itm += 1
        d = y - mean  # Welford running mean / variance
        mean += d / i
        m2 += d * (y - mean)
        if rows is not None and (i % every == 0 or i == m):
            se = df * math.sqrt(m2 / (i - 1) / i) if i > 1 else ""
            rows.append([prices, df * mean, se])
    se = df * math.sqrt(m2 / (m - 1) / m) if m > 1 else None
    return df * mean, se, itm / prices if prices else 0, prices

Run the full model from a terminal:

node simulation.js --option-type call --s0 100 --seed 42
python simulation.py --option-type call --s0 100 --seed 42

Download the model

Assumptions & limitations

• Educational model, not investment advice. The page demonstrates a numerical method. It does not value any real contract or recommend any trade.
• Black–Scholes world. Constant interest rate and volatility, lognormal prices, no dividends, no transaction costs, continuous trading. Real volatility changes over time and differs by strike (the volatility smile), and real returns have fatter tails than the lognormal model.
• European exercise only. The payoff depends on S_T alone. American options, which can be exercised early, need regression-based methods such as Longstaff–Schwartz on top of the simulation.
• Normal-theory interval. The 95% interval uses the central limit theorem; for deep out-of-the-money options with few paying paths the payoff distribution is very skewed and the interval is less reliable.
• Pseudo-random normals. Box–Muller draws from a seeded generator; see pseudo-random numbers and seeds.

Related simulations

• OS-017Geometric Brownian motionThe standard stock-price model: log-normal paths whose mean and typical outcome drift apart.
• OS-002Monte Carlo integrationEstimate an integral as the average of f at random points, with an error bar you can compute.
• OS-013Brownian motionParticles jiggle with Gaussian steps; their mean squared displacement grows as 4Dt, as Einstein predicted.
• OS-011Central limit theoremAverage n draws from almost any distribution and the means turn normal — at a speed you can measure.

Background: How Accurate Is a Monte Carlo Estimate? · Pseudo-Random Numbers and Seeds in Simulation

References

1. Black, F. & Scholes, M. (1973). The pricing of options and corporate liabilities. Journal of Political Economy, 81(3), 637–654. doi.org/10.1086/260062
2. Boyle, P. P. (1977). Options: A Monte Carlo approach. Journal of Financial Economics, 4(3), 323–338. doi.org/10.1016/0304-405X(77)90005-8
3. Glasserman, P. (2003). Monte Carlo Methods in Financial Engineering. Springer.