2. API Documentation

2.1. Overview

This package provides tools to state and numerically integrate Ito Stochastic Differential Equations (SDEs), including equations with time-dependent parameters, time-dependent correlations, and stochastic jumps, and to compute with, and extract statistics from, their realized paths.

Package contents:

  1. A set of tools to ease computations with stochastic processes, as obtained from numerical integration of the corresponding SDE, is provided via the process and montecarlo classes (see Infrastructure):
    • The process class, a subclass of numpy.ndarray representing a sequence of values in time, realized in one or several paths. Algebraic manipulations and ufunc computations are supported for instances that share the same timeline, or are constant, and comply with numpy broadcasting rules. Interpolation along the timeline is supported via callability of process instances. Process-specific functionalities, such as averaging and indexing along time or across paths, are delegated to process-specific methods, attributes and properties (no overriding of numpy.ndarray operations).
    • The montecarlo class, as an aid to cumulate the results of several Monte Carlo simulations of a given stochastic variable, and to extract summary estimates for its probability distribution function and statistics.
  2. Numerical realizations of the differentials commonly found as stochasticity sources in SDEs, are provided via the source class and its subclasses, with or without memory of formerly invoked realizations (see Stochasticity Sources).
  3. A general framework for stochastic step by step simulations, and for numerical SDE integration, is provided via the paths_generator class, and its cooperating subclasses integrator, SDE and SDEs (see SDE Integration Framework). The full API allows for extensive customization of preprocessing, post-processing, stochasticity sources instantiation and handling, integration algorithms etc. The integrate decorator provides a simple and concise interface to handle standard use cases, via Euler-Maruyama integration.
  4. Several preset stochastic processes are provided, including lognormal, Ornstein-Uhlenbeck, Hull-White n-factor, Heston, and jump-diffusion processes (see Stochastic Processes). Each process consists of a process generator class, a subclass of integrator and SDE, named with a _process suffix, and a definition of the underlying SDE, a subclass of SDE or SDEs, named with a _SDE suffix.
  5. Several analytical results relating to the preset stochastic processes are made available, as a general reference and for testing purposes (see Analytical Results). They are limited to the case of constant process parameters, and with some further limitations on the parameters’ domains. Function arguments are consistent with those of the corresponding processes. Suffixes _pdf, _cdf and _chf stand respectively for probability distribution function, cumulative probability distribution function, and characteristic function. Black-Scholes formulae for the valuation of call and put options have been included (with prefix bs).
  6. As an aid to interactive and notebook sessions, shortcuts are provided for stochasticity sources and preset processes (see Shortcuts). Shortcuts have been wrapped as “kfuncs”, objects with managed keyword arguments that simplify interactive workflow when frequent parameters tuning operations are needed (see kfunc decorator documentation). Analytical results are wrapped as kfuncs as well.

For all sources and processes, values can take any shape, scalar or multidimensional. Correlated multivariate stochasticity sources are supported. Poisson jumps are supported, and may be compounded with any random variable supported by scipy.stats. Time-varying process parameters (correlations, intensity of Poisson processes, volatilities etc.) are allowed whenever applicable. process instances act as valid stochasticity source realizations (as does any callable object complying with a source protocol), and may be passed as a source specification when computing the realization of a given process.

Computations are fully vectorized across paths, providing an efficient infrastructure for simulating a large number of process realizations. Less so, for large number of time steps: integrating 100 time steps across one million paths takes seconds, one million time steps across 100 paths takes minutes.

2.2. Infrastructure

process([t, x, v, c, dtype]) Array representation of a process (a subclass of numpy.ndarray).
piecewise([t, x, v, dtype, mode]) Return a process that interpolates to a piecewise constant function.
montecarlo([sample, axis, bins, range, use, …]) Summary statistics of Monte Carlo simulations.

2.3. Stochasticity Sources

source(*[, paths, vshape, dtype]) Base class for stochasticity sources.
wiener_source(*[, paths, vshape, dtype, …]) dw, a source of standard Wiener process (Brownian motion) increments.
poisson_source(*[, paths, vshape, dtype, lam]) dn, a source of Poisson process increments.
cpoisson_source(*[, paths, vshape, dtype, …]) dj, a source of compound Poisson process increments (jumps).
odd_wiener_source(*[, paths, vshape]) dw, a source of standard Wiener process (Brownian motion) increments with antithetic paths exposing opposite increments (averages exactly to 0 across paths).
even_poisson_source(*[, paths, vshape]) dn, a source of Poisson process increments with antithetic paths exposing identical increments.
even_cpoisson_source(*[, paths, vshape]) dj, a source of compound Poisson process increments (jumps) with antithetic paths exposing identical increments.
true_source(*[, paths, vshape, dtype, rtol, …]) Base class for stochasticity sources with memory.
true_wiener_source(*[, paths, vshape, …]) dw, source of standard Wiener process (brownian motion) increments with memory.
true_poisson_source(*[, paths, vshape, …]) dn, a source of Poisson process increments with memory.
true_cpoisson_source(*[, paths, vshape, …]) dj, a source of compound Poisson process increments (jumps) with memory.
norm_rv([a, b]) Normal distribution with mean a and standard deviation b, possibly time-dependent.
uniform_rv([a, b]) Uniform distribution between a and b, possibly time-dependent.
exp_rv([a]) Exponential distribution with scale a, possibly time-dependent.
double_exp_rv([a, b, pa]) Double exponential distribution, with scale a with probability pa, and -b with probability (1 - pa), possibly time-dependent.
rvmap(f, y) Map f to random variates of distribution y, possibly time-dependent.

2.4. SDE Integration Framework

paths_generator(*[, paths, xshape, wshape, …]) Step by step generation of stochastic simulations across multiple paths, intended for subclassing.
integrator(*[, paths, xshape, wshape, …]) Step by step numerical integration of Ito Stochastic Differential Equations (SDEs), intended for subclassing.
SDE(*[, paths, vshape, dtype, steps, i0, …]) Class representation of a user defined Stochastic Differential Equation (SDE), intended for subclassing.
SDEs(*[, paths, vshape, dtype, steps, i0, …]) Class representation of a user defined system of Stochastic Differential Equations (SDEs), intended for subclassing.
integrate([sde, q, sources, log, addaxis]) Decorator for Ito Stochastic Differential Equation (SDE) integration.

2.5. Stochastic Processes

wiener_process([paths, vshape, dtype, …]) Wiener process (Brownian motion) with drift.
lognorm_process([paths, vshape, dtype, …]) Lognormal process.
ornstein_uhlenbeck_process([paths, vshape, …]) Ornstein-Uhlenbeck process (mean-reverting Brownian motion).
hull_white_process([paths, vshape, dtype, …]) F-factors Hull-White process (sum of F correlated mean-reverting Brownian motions).
hull_white_1factor_process([paths, vshape, …]) 1-factor Hull-White process (F=1 Hull-White process with F-index collapsed to a scalar).
cox_ingersoll_ross_process([paths, vshape, …]) Cox-Ingersoll-Ross mean reverting process.
full_heston_process([paths, vshape, dtype, …]) Heston stochastic volatility process (returns both process and volatility).
heston_process([paths, vshape, dtype, …]) Heston stochastic volatility process (stores and returns process only).
jumpdiff_process([paths, vshape, dtype, …]) Jump-diffusion process (lognormal process with compound Poisson logarithmic jumps).
merton_jumpdiff_process([paths, vshape, …]) Merton jump-diffusion process (jump-diffusion process with normal jump size distribution).
kou_jumpdiff_process([paths, vshape, dtype, …]) Double exponential (Kou) jump-diffusion process (jump-diffusion process with double exponential jump size distribution).
wiener_SDE(*[, paths, vshape, dtype, steps, …]) SDE for a Wiener process (Brownian motion) with drift.
lognorm_SDE(*[, paths, vshape, dtype, …]) SDE for a lognormal process with drift.
ornstein_uhlenbeck_SDE(*[, paths, vshape, …]) SDE for an Ornstein-Uhlenbeck process.
hull_white_SDE(*[, paths, vshape, dtype, …]) SDE for an F-factors Hull White process.
cox_ingersoll_ross_SDE(*[, paths, vshape, …]) SDE for a Cox-Ingersoll-Ross mean reverting process.
full_heston_SDE(*[, paths, vshape, dtype, …]) SDE for a Heston stochastic volatility process.
heston_SDE(*[, paths, vshape, dtype, steps, …]) SDE for a Heston stochastic volatility process.
jumpdiff_SDE(*[, paths, vshape, dtype, …]) SDE for a jump-diffusion process (lognormal process with compound Poisson logarithmic jumps).
merton_jumpdiff_SDE(*[, paths, vshape, …]) SDE for a Merton jump-diffusion process.
kou_jumpdiff_SDE(*[, paths, vshape, dtype, …]) SDE for a double exponential (Kou) jump-diffusion process.

2.6. Analytical Results

wiener_mean(t, *[, x0, mu, sigma]) Mean of values at time t of a Wiener process (as per the wiener_process class) with time-independent parameters.
wiener_var(t, *[, x0, mu, sigma]) Variance of values at time t of a Wiener process (as per the wiener_process class) with time-independent parameters.
wiener_std(t, *[, x0, mu, sigma]) Standard deviation of values at time t of a Wiener process (as per the wiener_process class) with time-independent parameters.
wiener_pdf(t, x, *[, x0, mu, sigma]) Probability distribution function of values at time t of a Wiener process (as per the wiener_process class) with time-independent parameters, evaluated at x.
wiener_cdf(t, x, *[, x0, mu, sigma]) Cumulative probability distribution function of values at time t of a Wiener process (as per the wiener_process class) with time-independent parameters, evaluated at x.
wiener_chf(t, u, *[, x0, mu, sigma]) Characteristic function of the probability distribution of values at time t of a Wiener process (as per the wiener_process class) with time-independent parameters, evaluated at u.
lognorm_mean(t, *[, x0, mu, sigma]) Mean of values at time t of a lognormal process (as per the lognorm_process class) with time-independent parameters.
lognorm_var(t, *[, x0, mu, sigma]) Variance of values at time t of a lognormal process (as per the lognorm_process class) with time-independent parameters.
lognorm_std(t, *[, x0, mu, sigma]) Standard deviation of values at time t of a lognormal process (as per the lognorm_process class) with time-independent parameters.
lognorm_pdf(t, x, *[, x0, mu, sigma]) Probability distribution function of values at time t of a lognormal process (as per the lognorm_process class) with time-independent parameters, evaluated at x.
lognorm_cdf(t, x, *[, x0, mu, sigma]) Cumulative probability distribution function of values at time t of a lognormal process (as per the lognorm_process class) with time-independent parameters, evaluated at x.
lognorm_log_chf(t, u, *[, x0, mu, sigma]) Characteristic function of the probability distribution of values at time t of the logarithm of a lognormal process (as per the lognorm_process class) with time-independent parameters, evaluated at u.
oruh_mean(t, *[, x0, theta, k, sigma]) Mean of values at time t of an Ornstein-Uhlenbeck process (as per the ornstein_uhlenbeck_process class) with time-independent parameters.
oruh_var(t, *[, x0, theta, k, sigma]) Variance of values at time t of an Ornstein-Uhlenbeck process (as per the ornstein_uhlenbeck_process class) with time-independent parameters.
oruh_std(t, *[, x0, theta, k, sigma]) Standard deviation of values at time t of an Ornstein-Uhlenbeck process (as per the ornstein_uhlenbeck_process class) with time-independent parameters.
oruh_pdf(t, x, *[, x0, theta, k, sigma]) Probability distribution function of values at time t of an Ornstein-Uhlenbeck process (as per the ornstein_uhlenbeck_process class) with time-independent parameters, evaluated at x.
oruh_cdf(t, x, *[, x0, theta, k, sigma]) Cumulative probability distribution function of values at time t of an Ornstein-Uhlenbeck process (as per the ornstein_uhlenbeck_process class) with time-independent parameters, evaluated at x.
hw2f_mean(t, *[, x0, theta, k, sigma, rho]) Mean of values at time t of a Hull-White 2-factors process (as per the hull_white_process class) with time-independent parameters.
hw2f_var(t, *[, x0, theta, k, sigma, rho]) Variance of values at time t of a Hull-White 2-factors process (as per the hull_white_process class) with time-independent parameters.
hw2f_std(t, *[, x0, theta, k, sigma, rho]) Standard deviation of values at time t of a Hull-White 2-factors process (as per the hull_white_process class) with time-independent parameters.
hw2f_pdf(t, x, *[, x0, theta, k, sigma, rho]) Probability distribution function of values at time t of a Hull-White 2-factors process (as per the hull_white_process class) with time-independent parameters, evaluated at x.
hw2f_cdf(**args) Cumulative probability distribution function of values at time t of a Hull-White 2-factors process (as per the hull_white_process class) with time-independent parameters, evaluated at x.
cir_mean(t, *[, x0, theta, k, xi]) Mean of values at time t of a Cox-Ingersoll-Ross process (as per the cox_ingersoll_ross_process class) with time-independent parameters.
cir_var(t, *[, x0, theta, k, xi]) Variance of values at time t of a Cox-Ingersoll-Ross process (as per the cox_ingersoll_ross_process class) with time-independent parameters.
cir_std(t, *[, x0, theta, k, xi]) Standard deviation of values at time t of a Cox-Ingersoll-Ross process (as per the cox_ingersoll_ross_process class) with time-independent parameters.
cir_pdf(t, x, *[, x0, theta, k, xi]) Probability distribution function of values at time t of a Cox-Ingersoll-Ross process (as per the cox_ingersoll_ross_process class) with time-independent parameters, evaluated at x.
heston_log_mean(t, *[, x0, mu, sigma, y0, …]) Mean of the logarithm of values at time t of a Heston process (as per the full_heston_process class) with time-independent parameters.
heston_log_var(**args) Variance of the logarithm of values at time t of a Heston process (as per the full_heston_process class) with time-independent parameters.
heston_log_std(t, *[, x0, mu, sigma, y0, …]) Standard deviation of the logarithm of values at time t of a Heston process (as per the full_heston_process class) with time-independent parameters.
heston_log_pdf(t, logx, *[, x0, mu, sigma, …]) Probability distribution function of values at time t of the logarithm of a Heston process, (as per the full_heston_process class) with time-independent parameters, evaluated at logx.
heston_log_chf(t, u, *[, x0, mu, sigma, y0, …]) Characteristic function of the probability distribution of values at time t of the logarithm of a Heston process (as per the full_heston_process class) , with time-independent parameters, evaluated at u.
mjd_log_pdf(t, logx, *[, x0, mu, sigma, …]) Probability distribution function of values at time t of the logarithm of a Merton jump-diffusion process (as per the merton_jumpdiff_process class), with time-independent parameters, evaluated at logx.
mjd_log_chf(t, u, *[, x0, mu, sigma, lam, a, b]) Characteristic function of the probability distribution of values at time t of the logarithm of a Merton jump-diffusion process (as per the merton_jumpdiff_process class), with time-independent parameters, evaluated at u.
kou_mean(t, *[, x0, mu, sigma, lam, a, b, pa]) Mean of values at time t of a double exponential (Kou) jump-diffusion process (as per the kou_jumpdiff_process class) with time-independent parameters.
kou_log_pdf(t, logx, *[, x0, mu, sigma, …]) Probability distribution function of values at time t of the logarithm of a double-exponential (Kou) jump-diffusion process (as per the kou_jumpdiff_process class), with time-independent parameters, evaluated at logx.
kou_log_chf(t, u, *[, x0, mu, sigma, lam, …]) Characteristic function of the probability distribution of values at time t of the logarithm of a Kou jump-diffusion process, (as per the kou_jumpdiff_process class) with time-independent parameters, evaluated at u.
bsd1d2(k, t, *[, x0, r, q, sigma]) Black-Scholes d1 and d2 coefficients.
bscall(k, t, *[, x0, r, q, sigma]) Black-Scholes call option value.
bscall_delta(k, t, *[, x0, r, q, sigma]) Black-Scholes call option delta.
bsput(k, t, *[, x0, r, q, sigma]) Black-Scholes put option value.
bsput_delta(k, t, *[, x0, r, q, sigma]) Black-Scholes put option delta.

2.7. Shortcuts

Stochasticity sources and preset processes may be addressed using the following shortcuts:

Full name Shortcut
wiener_source dw
poisson_source dn
cpoisson_source dj
odd_wiener_source odd_dw
even_poisson_source even_dn
even_cpoisson_source even_dj
true_wiener_source true_dw
true_poisson_source true_dn
true_cpoisson_source true_dj
wiener_process wiener
lognorm_process lognorm
ornstein_uhlenbeck_process oruh
hull_white_process hwff
hull_white_1factor_process hw1f
cox_ingersoll_ross_process cir
full_heston_process heston_xy
heston_process heston
jumpdiff_process jumpdiff
merton_jumpdiff_process mjd
kou_jumpdiff_process kou

Shortcuts have been wrapped as “kfuncs”, objects with managed keyword arguments (see kfunc decorator documentation below).

Analytical results are named according to the shortcut of the corresponding process (e.g. lognorm_mean, lognorm_cdf etc. from the lognorm shortcut) and are wrapped as kfuncs as well.

kfunc([f, nvar]) Decorator to wrap classes or functions as objects with managed keyword arguments.
iskfunc(cls_or_object) Tests if the given class or instance has been wrapped as a kfunc.

2.8. Testing

Tests have been set up within the numpy.testing framework. To launch tests, invoke sdepy.test() or sdepy.test('full').

The testing subpackage sdepy.tests was written in pursuit of the following goals:

  • Maximize case coverage, by exposing the package functions and methods to a plurality of different input shapes, values, data types etc., and of different combinations thereof, as may be encountered in practice.
  • Provide a quantitative validation of the algorithms, functions and processes covered in sdepy.
  • Keep dependencies of the test code on the adopted testing framework to a bare minimum.

Most often, a number of testing cases is declared as a list or lists of classes and inputs, a general testing procedure is set up, and the latter is iteratively applied to the former. Unfortunately, all this resulted in a thinly documented (if at all), hard to read, and hard to maintain testing code base - sorry about that.

The quantitative validation of the package, via tests marked as 'slow' and 'quant', is done in two steps:

  • To validate a sdepy release, tests are run with both 100 and 100000 paths against a fixed random seed. Numerical integration results for the mean, standard deviation, probability distribution, and/or characteristic function are compared against their exact values computed analytically from the process parameters. Comparisons are then plotted and visually inspected, and the occasional larger than usual deviation is manually checked to be statistically acceptable, i.e. only so few standard deviations off the mark. The plots and the average and maximum errors are recorded in png and text files located in the ./tests/cfr directory, relative to the package home directory where sdepy.__file__ is located.
  • Each time sdepy.test('full') is invoked, to keep testing times manageable and the testing procedure uninvasive, tests are run with 100 paths against the same fixed random seed, without plotting or storing results. The realized errors are then compared and checked against the expected errors, as distributed with the package and stored in the ./tests/cfr directory.

Note that the tests rely on the reproducibility of expected errors, once random numbers have been seeded with np.random.seed(), across platforms and versions of Python, NumPy and SciPy.

In order to reproduce the full tests and inspect the graphs, change the following configuration settings in the file of the sdepy._config subpackage (private, not part of the API, may change in the future):

PLOT = True
SAVE_ERRORS = True
QUANT_TEST_MODE = 'HD'

With these settings, tests are run with 100000 paths, and realized errors and plots are stored in the ./tests/cfr directory. In case some tests fail, to carry out the whole procedure and get the failing errors and plots, set in the same configuration file:

QUANT_TEST_FAIL = False