RadonLab

Monte Carlo Greeks for discontinuous payoffs, done correctly and benchmarked against closed-form references.

When an option's payoff is discontinuous (a digital that pays a fixed amount if the underlying finishes above a strike, a barrier that knocks out on touch), the textbook pathwise (infinitesimal perturbation) estimator for its sensitivities quietly breaks: the derivative of the payoff is zero almost everywhere, so the estimate collapses toward zero. radonlab implements the established remedies, tests each one against a known analytic answer, and benchmarks them side by side so you can see which method to use and what it costs.

Pure NumPy/SciPy. Apache-2.0. No copyleft dependencies.

Companion C++ implementation

A header-only C++20 port lives at github.com/quantsingularity/radonlab-cpp. The two libraries share the same methods and the same closed-form references, and the C++ analytics are cross-checked against this library's values (which are verified to machine precision by complex-step differentiation), so the two agree to 1e-9 across languages.

	Python (radonlab)	C++ (radonlab-cpp)
Role	Reference implementation; carries the research extensions	Header-only production port
Validated by	Complex-step differentiation, to machine precision	Cross-checked against this library, to 1e-9
Extras	Merton jump-diffusion, path-dependent Malliavin delta	Adjoint algorithmic differentiation, exact second-order Greeks, threaded estimator

What's inside

All estimators return a value together with its per-path standard error.

Estimator	Differentiates	Notes
Likelihood-ratio (LRM)	the density, never the payoff	Unbiased for digitals and other discontinuous payoffs. Covers delta, vega, gamma.
Smoothing	a $C^1$ approximation $f_h$ of bandwidth $h$	A tunable bias/variance trade-off.
Conditional Monte Carlo	integrates out the barrier-crossing event with the Brownian-bridge non-crossing probability	Unbiased for the continuously-monitored price at any number of steps, and Lipschitz in spot, so its delta comes out by the pathwise method.
Pathwise / finite-difference	baselines	Kept so their failure modes on discontinuous payoffs are visible and measurable.

Closed-form Black-Scholes analytics (the ground truth) are provided for European calls/puts, cash-or-nothing and asset-or-nothing digitals, and continuously-monitored down-and-out / down-and-in calls. Every explicit Greek formula is checked against complex-step differentiation of the price to machine precision.

Install

pip install -e ".[benchmark,test]"

Package	Needed for
numpy	core (required)
scipy	core (required)
pandas	benchmark table
matplotlib	benchmark plots

Quick start

from radonlab import GeometricBrownianMotion, CashOrNothingCall, Greek
from radonlab import likelihood_ratio, pathwise, analytics

model = GeometricBrownianMotion(S0=100, r=0.03, sigma=0.20, T=1.0)
payoff = CashOrNothingCall(strike=100)

# Ground truth
ref = analytics.cash_or_nothing_call(100, 100, 0.03, 0.20, 1.0, greek="delta")

# Likelihood ratio: unbiased for the digital
lrm = likelihood_ratio(model, payoff, Greek.DELTA, n_paths=1_000_000)
print(lrm)                    # delta=... +/- ...  [likelihood-ratio]
print(lrm.confidence_interval(0.95))

# Naive pathwise: collapses to ~0, badly biased
pw = pathwise(model, payoff, Greek.DELTA, n_paths=1_000_000)
print(pw.value, "vs analytic", ref)

The math, briefly

Under geometric Brownian motion, the terminal spot is

$$ S_T = S_0 \exp!\left(\left(r - q - \tfrac{1}{2}\sigma^2\right)T + \sigma\sqrt{T},Z\right), \qquad Z \sim \mathcal{N}(0,1), $$

and a price is $P = e^{-rT},\mathbb{E}[f(S_T)]$.

Method	Identity	Behaviour on a jump in $f$
Likelihood ratio	$\dfrac{\partial P}{\partial \theta} = e^{-rT},\mathbb{E}!\left[f(S_T)\cdot \dfrac{\partial}{\partial \theta}\log p(S_T;\theta)\right]$	Unbiased: $f$ is never differentiated.
Pathwise	$\dfrac{\partial P}{\partial \theta} = e^{-rT},\mathbb{E}!\left[f'(S_T)\cdot \dfrac{\partial S_T}{\partial \theta}\right]$	Degenerate: $f'(S_T)=0$ almost everywhere for a digital.
Conditional MC (barrier)	multiplies the terminal payoff by $\displaystyle\prod_i\left(1-\exp!\left(\dfrac{-2(X_i-b)(X_{i+1}-b)}{\sigma^2\Delta t}\right)\right)$	Exact bridge probability of not breaching $b=\ln B$ between simulated points.

Writing $Z = \dfrac{\ln S_T - \mu}{\sigma\sqrt{T}}$, the likelihood-ratio score functions are:

Greek	Score function
Delta	$\dfrac{Z}{S_0,\sigma\sqrt{T}}$
Vega	$\dfrac{Z^2-1}{\sigma} - \sqrt{T},Z$
Gamma	$\dfrac{Z^2 - \sigma\sqrt{T},Z - 1}{S_0^2,\sigma^2,T}$

Full references are in each module's docstring: Broadie and Glasserman (1996), Glasserman (2003), Reiner and Rubinstein (1991), and the smoothed-perturbation and Malliavin literature.

Benchmark

Reproduce with python examples/benchmark_digital_delta.py (1,000,000 paths, seed 0). Delta of an at-the-money cash-or-nothing call, $S_0=K=100$, $r=3%$, $\sigma=20%$, $T=1$. Analytic reference delta = 0.019333. Ranked by standard error, the metric that actually reflects an estimator's accuracy:

Method	Value	Std. error	Abs. error vs analytic	Biased?
Likelihood-ratio	0.01939	2.84e-05	5.9e-05	No
Smoothing (h = default)	0.01921	4.56e-05	1.2e-04	No
Finite-difference (1% bump)	0.01934	9.49e-05	8.5e-06	No
Pathwise	0.00000	0.00e+00	1.9e-02	Yes

Two things to read off it. First, pathwise is 100% wrong: it returns exactly zero, because the digital's derivative is zero almost everywhere. Second, among the unbiased methods the likelihood-ratio estimator has the lowest variance (about 3.3 times smaller standard error than the 1% finite difference) and needs no bump to tune.

Finite difference forces a bad choice: sweeping the bump size in the same run shows the standard error exploding as the bump shrinks, while a large bump introduces its own bias:

Relative bump	Value	Std. error	Abs. error
1e-01	0.018620	2.36e-05	7.1e-04 (bias)
1e-02	0.019311	9.49e-05	2.2e-05
1e-03	0.019161	3.04e-04	1.7e-04
1e-04	0.017274	9.15e-04	2.1e-03

The likelihood-ratio method sidesteps that trade-off entirely. This is the whole reason the library exists.

Testing

pytest

The suite validates the analytics by complex-step differentiation, confirms each Monte Carlo estimator agrees with the analytic Greek within a few standard errors, and documents the pathwise failure on digitals as an explicit test.

Scope and honesty

This release covers the Black-Scholes / GBM world and, in radonlab.merton, the Merton jump-diffusion model with a closed-form digital reference: it shows the likelihood-ratio delta weight is unchanged by jumps and validates the estimator against the closed form. It is a rigorous, tested core rather than a kitchen sink: the estimators and their validation are real, and where a method is expected to be biased (pathwise on a digital), the library shows it rather than hiding it. It also includes a genuinely path-dependent frontier estimator in radonlab.asian: the likelihood-ratio / Malliavin delta of a geometric Asian digital, validated against the closed form, with a weight built from the whole path. Stochastic volatility, arithmetic-Asian and lookback payoffs, and a JIT/vectorized backend are natural next steps.

License

Apache-2.0. See LICENSE.