Barndorff-Nielson & Shephard process¶

quantflow.sp.bns.BNS `pydantic-model` ¶

Bases: StochasticProcess1D

Barndorff-Nielson & Shephard (BNS) stochastic volatility model.

This is a stochastic volatility model where the variance process is given by a non-Gaussian Ornstein-Uhlenbeck process driven by a pure-jump Lévy process. The BNS model is defined by the following system of SDEs:

\[\begin{equation} \begin{aligned} dx_t &= \sqrt{v_t} dw_t + \rho dz_{\kappa t} \\ dv_t &= -\kappa v_t dt + dz_{\kappa t} \end{aligned} \end{equation}\]

The model is flexible and can capture various stylized facts of financial markets, such as volatility clustering and leverage effects.

This implementation uses a GammaOU process for the variance, which is a common choice in the BNS model.

Fields:

variance_process (GammaOU)
rho (float)

variance_process `pydantic-field` ¶

variance_process

Variance process

rho `pydantic-field` ¶

rho = 0

Correlation between variance and price processes, in (-1, 1)

create `classmethod` ¶

create(vol, kappa, decay, rho)

Convenience constructor for BNS process with parameters of the variance process

Source code in quantflow/sp/bns.py

@classmethod
def create(cls, vol: float, kappa: float, decay: float, rho: float) -> Self:
    """Convenience constructor for BNS process with parameters
    of the variance process
    """
    return cls(
        variance_process=GammaOU.create(rate=vol * vol, kappa=kappa, decay=decay),
        rho=rho,
    )

characteristic_exponent ¶

characteristic_exponent(t, u)

Characteristic exponent of the BNS process with Gamma-OU variance.

\[\begin{equation} \phi(t, u) = B v_0 - \lambda \left[ \frac{\kappa t \, A}{\beta - A} + \frac{\beta}{\beta - A} \log\frac{\beta - i u \rho + B}{\beta - i u \rho} \right] \end{equation}\]

with

\[\begin{equation} \begin{aligned} A &= i u \rho - \frac{u^2}{2 \kappa}, \\ B &= \frac{u^2 (1 - e^{-\kappa t})}{2 \kappa}. \end{aligned} \end{equation}\]

where \(v_0\) is the initial variance, \(\kappa\) is the mean-reversion speed, \(\rho\) is the leverage parameter, and \((\lambda, \beta)\) are the intensity and exponential-jump rate of the background driving Lévy process.

Source code in quantflow/sp/bns.py

def characteristic_exponent(self, t: FloatArrayLike, u: Vector) -> Vector:
    r"""Characteristic exponent of the BNS process with Gamma-OU variance.

    \begin{equation}
        \phi(t, u) = B v_0
            - \lambda \left[
                \frac{\kappa t \, A}{\beta - A}
                + \frac{\beta}{\beta - A}
                \log\frac{\beta - i u \rho + B}{\beta - i u \rho}
            \right]
    \end{equation}

    with

    \begin{equation}
        \begin{aligned}
            A &= i u \rho - \frac{u^2}{2 \kappa}, \\
            B &= \frac{u^2 (1 - e^{-\kappa t})}{2 \kappa}.
        \end{aligned}
    \end{equation}

    where $v_0$ is the initial variance, $\kappa$ is the mean-reversion speed,
    $\rho$ is the leverage parameter, and $(\lambda, \beta)$ are the intensity
    and exponential-jump rate of the background driving Lévy process.
    """
    return self._characteristic_exponent(t, u, 1.0)

sample ¶

sample(n, time_horizon=1, time_steps=100)

PARAMETER	DESCRIPTION
`n`	Number of paths TYPE: `int`
`time_horizon`	Time horizon TYPE: `float` DEFAULT: `1`
`time_steps`	Number of time steps to arrive at horizon TYPE: `int` DEFAULT: `100`

Source code in quantflow/sp/bns.py

def sample(
    self,
    n: Annotated[int, Doc("Number of paths")],
    time_horizon: Annotated[float, Doc("Time horizon")] = 1,
    time_steps: Annotated[
        int, Doc("Number of time steps to arrive at horizon")
    ] = 100,
) -> Paths:
    return self.sample_from_draws(Paths.normal_draws(n, time_horizon, time_steps))

sample_from_draws ¶

sample_from_draws(draws, *args)

PARAMETER	DESCRIPTION
`draws`	Pre-drawn standard normal increments for the Brownian motion TYPE: `Paths`
`*args`	Optional pre-drawn subordinator paths; new draws are generated if omitted TYPE: `Paths` DEFAULT: `()`

Source code in quantflow/sp/bns.py

def sample_from_draws(
    self,
    draws: Annotated[
        Paths,
        Doc("Pre-drawn standard normal increments for the Brownian motion"),
    ],
    *args: Annotated[
        Paths,
        Doc(
            "Optional pre-drawn subordinator paths; "
            "new draws are generated if omitted"
        ),
    ],
) -> Paths:
    path_dz = (
        args[0]
        if args
        else self.variance_process.bdlp.sample(
            draws.samples,
            self.variance_process.kappa * draws.t,
            draws.time_steps,
        )
    )
    v = self._variance_path(path_dz)
    # Price: x_t = integral sqrt(v) dW + rho * z_{kappa t}
    diffusion = np.sqrt(v[:-1] * draws.dt) * draws.data[:-1]
    paths = np.zeros_like(draws.data)
    paths[1:] = np.cumsum(diffusion, axis=0) + self.rho * path_dz.data[1:]
    return Paths(t=draws.t, data=paths)

characteristic ¶

characteristic(t, u)

Characteristic function at time t for a given input parameter u

The characteristic function represents the Fourier transform of the probability density function

\[\begin{equation} \Phi = {\mathbb E} \left[e^{i u x_t}\right] = e^{-\phi(t, u)} \end{equation}\]

where \(\phi\) is the characteristic exponent, which can be more easily computed for many processes.

PARAMETER	DESCRIPTION
`t`	Time horizon TYPE: `FloatArrayLike`
`u`	Characteristic function input parameter TYPE: `Vector`

Source code in quantflow/sp/base.py

def characteristic(
    self,
    t: Annotated[FloatArrayLike, Doc("Time horizon")],
    u: Annotated[Vector, Doc("Characteristic function input parameter")],
) -> Vector:
    r"""Characteristic function at time `t` for a given input parameter `u`

    The characteristic function represents the Fourier transform of the
    probability density function

    \begin{equation}
        \Phi = {\mathbb E} \left[e^{i u x_t}\right] = e^{-\phi(t, u)}
    \end{equation}

    where $\phi$ is the characteristic exponent, which can be more easily
    computed for many processes.
    """
    return np.exp(-self.characteristic_exponent(t, u))

convexity_correction ¶

convexity_correction(t)

Convexity correction for the process

Source code in quantflow/sp/base.py

def convexity_correction(self, t: FloatArrayLike) -> Vector:
    """Convexity correction for the process"""
    return -self.characteristic_exponent(t, complex(0, -1)).real

analytical_std ¶

analytical_std(t)

Analytical standard deviation of the process at time t

This has a closed form solution if the process has an analytical variance

Source code in quantflow/sp/base.py

def analytical_std(self, t: FloatArrayLike) -> FloatArrayLike:
    """Analytical standard deviation of the process at time `t`

    This has a closed form solution if the process has an analytical variance
    """
    return np.sqrt(self.analytical_variance(t))

analytical_mean ¶

analytical_mean(t)

Analytical mean of the process at time t

Implement if available

Source code in quantflow/sp/base.py

def analytical_mean(self, t: FloatArrayLike) -> FloatArrayLike:
    """Analytical mean of the process at time `t`

    Implement if available
    """
    raise NotImplementedError

analytical_variance ¶

analytical_variance(t)

Analytical variance of the process at time t

Implement if available

Source code in quantflow/sp/base.py

def analytical_variance(self, t: FloatArrayLike) -> FloatArrayLike:
    """Analytical variance of the process at time `t`

    Implement if available
    """
    raise NotImplementedError

analytical_pdf ¶

analytical_pdf(t, x)

Analytical pdf of the process at time t

Implement if available

Source code in quantflow/sp/base.py

def analytical_pdf(self, t: FloatArrayLike, x: FloatArrayLike) -> FloatArrayLike:
    """Analytical pdf of the process at time `t`

    Implement if available
    """
    raise NotImplementedError

analytical_cdf ¶

analytical_cdf(t, x)

Analytical cdf of the process at time t

Implement if available

Source code in quantflow/sp/base.py

def analytical_cdf(self, t: FloatArrayLike, x: FloatArrayLike) -> FloatArrayLike:
    """Analytical cdf of the process at time `t`

    Implement if available
    """
    raise NotImplementedError

marginal ¶

marginal(t)

Marginal distribution of the process at time t

Source code in quantflow/sp/base.py

def marginal(self, t: FloatArrayLike) -> StochasticProcess1DMarginal:
    """Marginal distribution of the process at time `t`"""
    return StochasticProcess1DMarginal(process=self, t=t)

domain_range ¶

domain_range()

Source code in quantflow/sp/base.py

def domain_range(self) -> Bounds:
    return default_bounds()

frequency_range ¶

frequency_range(std, max_frequency=None)

Maximum frequency when calculating characteristic functions

Source code in quantflow/sp/base.py

def frequency_range(self, std: float, max_frequency: float | None = None) -> Bounds:
    """Maximum frequency when calculating characteristic functions"""
    if max_frequency is None:
        max_frequency = np.sqrt(40 / std / std)
    return Bounds(0, max_frequency)

support ¶

support(mean, std, points)

Support of the process at time t

Source code in quantflow/sp/base.py

def support(self, mean: float, std: float, points: int) -> FloatArray:
    """Support of the process at time `t`"""
    bounds = self.domain_range()
    start = float(sigfig(bound_from_any(bounds.lb, mean - std)))
    end = float(sigfig(bound_from_any(bounds.ub, mean + std)))
    return np.linspace(start, end, points + 1)

quantflow.sp.bns.BNS2 `pydantic-model` ¶

Bases: StochasticProcess1D

Two-factor Barndorff-Nielson & Shephard stochastic volatility model.

The original multi-factor BNS extension drives a single log-price with a single Brownian motion against a convex combination of independent Gamma-OU variances. With weight \(w \in [0, 1]\) for the first factor:

\[\begin{equation} \begin{aligned} \sigma^2_t &= w\, v^1_t + (1 - w)\, v^2_t \\ dx_t &= \sigma_t\, dw_t + \rho_1\, dz^1_{\kappa_1 t} + \rho_2\, dz^2_{\kappa_2 t} \\ dv^i_t &= -\kappa_i v^i_t\, dt + dz^i_{\kappa_i t}, \quad i = 1, 2 \end{aligned} \end{equation}\]

A fast and a slow factor combined this way add flexibility to the term structure of volatility while retaining the analytic tractability of BNS.

Fields:

bns1 (BNS)
bns2 (BNS)
weight (float)

bns1 `pydantic-field` ¶

bns1

First BNS variance factor

bns2 `pydantic-field` ¶

bns2

Second BNS variance factor

weight `pydantic-field` ¶

weight = 0.5

Weight \(w\) of the first variance factor in the convex combination; the second factor receives weight \(1 - w\)

characteristic_exponent ¶

characteristic_exponent(t, u)

Characteristic exponent as the sum of two weighted BNS exponents.

Conditional on the variance paths the diffusion is Gaussian with variance \(w \int v^1_s\,ds + (1-w) \int v^2_s\,ds\). Independence of the two BDLPs then factorises the unconditional expectation into a product, giving

\[\begin{equation} \phi_{x_t,u} = \phi^{(1)}_{x_t,u}\big|_{u^2 \mapsto w u^2} + \phi^{(2)}_{x_t,u}\big|_{u^2 \mapsto (1 - w) u^2} \end{equation}\]

where the substitution applies only to the diffusion term (\(u^2\)) and leaves the leverage term (\(i u \rho_i\)) unchanged.

PARAMETER	DESCRIPTION
`t`	Time horizon or array of evaluation times TYPE: `FloatArrayLike`
`u`	Characteristic exponent argument TYPE: `Vector`

Source code in quantflow/sp/bns.py

def characteristic_exponent(
    self,
    t: Annotated[FloatArrayLike, Doc("Time horizon or array of evaluation times")],
    u: Annotated[Vector, Doc("Characteristic exponent argument")],
) -> Vector:
    r"""Characteristic exponent as the sum of two weighted BNS exponents.

    Conditional on the variance paths the diffusion is Gaussian with variance
    $w \int v^1_s\,ds + (1-w) \int v^2_s\,ds$. Independence of the two BDLPs
    then factorises the unconditional expectation into a product, giving

    \begin{equation}
        \phi_{x_t,u} = \phi^{(1)}_{x_t,u}\big|_{u^2 \mapsto w u^2}
            + \phi^{(2)}_{x_t,u}\big|_{u^2 \mapsto (1 - w) u^2}
    \end{equation}

    where the substitution applies only to the diffusion term ($u^2$) and
    leaves the leverage term ($i u \rho_i$) unchanged.
    """
    w = self.weight
    return self.bns1._characteristic_exponent(
        t, u, w
    ) + self.bns2._characteristic_exponent(t, u, 1.0 - w)

sample ¶

sample(n, time_horizon=1, time_steps=100)

PARAMETER	DESCRIPTION
`n`	Number of sample paths TYPE: `int`
`time_horizon`	Time horizon TYPE: `float` DEFAULT: `1`
`time_steps`	Number of discrete time steps TYPE: `int` DEFAULT: `100`

Source code in quantflow/sp/bns.py

def sample(
    self,
    n: Annotated[int, Doc("Number of sample paths")],
    time_horizon: Annotated[float, Doc("Time horizon")] = 1,
    time_steps: Annotated[int, Doc("Number of discrete time steps")] = 100,
) -> Paths:
    return self.sample_from_draws(Paths.normal_draws(n, time_horizon, time_steps))

sample_from_draws ¶

sample_from_draws(draws, *args)

PARAMETER	DESCRIPTION
`draws`	Single Brownian motion driving both factors TYPE: `Paths`
`*args`	Optional pre-drawn BDLP paths for bns1 and bns2 (in that order) TYPE: `Paths` DEFAULT: `()`

Source code in quantflow/sp/bns.py

def sample_from_draws(
    self,
    draws: Annotated[Paths, Doc("Single Brownian motion driving both factors")],
    *args: Annotated[
        Paths,
        Doc("Optional pre-drawn BDLP paths for bns1 and bns2 (in that order)"),
    ],
) -> Paths:
    time_horizon = draws.t
    time_steps = draws.time_steps
    n = draws.samples
    path_dz1 = (
        args[0]
        if len(args) > 0
        else self.bns1.variance_process.bdlp.sample(
            n, self.bns1.variance_process.kappa * time_horizon, time_steps
        )
    )
    path_dz2 = (
        args[1]
        if len(args) > 1
        else self.bns2.variance_process.bdlp.sample(
            n, self.bns2.variance_process.kappa * time_horizon, time_steps
        )
    )
    v1 = self.bns1._variance_path(path_dz1)
    v2 = self.bns2._variance_path(path_dz2)
    w = self.weight
    sigma2 = w * v1 + (1.0 - w) * v2
    diffusion = np.sqrt(sigma2[:-1] * draws.dt) * draws.data[:-1]
    paths = np.zeros_like(draws.data)
    paths[1:] = (
        np.cumsum(diffusion, axis=0)
        + self.bns1.rho * path_dz1.data[1:]
        + self.bns2.rho * path_dz2.data[1:]
    )
    return Paths(t=time_horizon, data=paths)

characteristic ¶

characteristic(t, u)

Characteristic function at time t for a given input parameter u

The characteristic function represents the Fourier transform of the probability density function

\[\begin{equation} \Phi = {\mathbb E} \left[e^{i u x_t}\right] = e^{-\phi(t, u)} \end{equation}\]

where \(\phi\) is the characteristic exponent, which can be more easily computed for many processes.

PARAMETER	DESCRIPTION
`t`	Time horizon TYPE: `FloatArrayLike`
`u`	Characteristic function input parameter TYPE: `Vector`

Source code in quantflow/sp/base.py

def characteristic(
    self,
    t: Annotated[FloatArrayLike, Doc("Time horizon")],
    u: Annotated[Vector, Doc("Characteristic function input parameter")],
) -> Vector:
    r"""Characteristic function at time `t` for a given input parameter `u`

    The characteristic function represents the Fourier transform of the
    probability density function

    \begin{equation}
        \Phi = {\mathbb E} \left[e^{i u x_t}\right] = e^{-\phi(t, u)}
    \end{equation}

    where $\phi$ is the characteristic exponent, which can be more easily
    computed for many processes.
    """
    return np.exp(-self.characteristic_exponent(t, u))

convexity_correction ¶

convexity_correction(t)

Convexity correction for the process

Source code in quantflow/sp/base.py

def convexity_correction(self, t: FloatArrayLike) -> Vector:
    """Convexity correction for the process"""
    return -self.characteristic_exponent(t, complex(0, -1)).real

analytical_std ¶

analytical_std(t)

Analytical standard deviation of the process at time t

This has a closed form solution if the process has an analytical variance

Source code in quantflow/sp/base.py

def analytical_std(self, t: FloatArrayLike) -> FloatArrayLike:
    """Analytical standard deviation of the process at time `t`

    This has a closed form solution if the process has an analytical variance
    """
    return np.sqrt(self.analytical_variance(t))

analytical_mean ¶

analytical_mean(t)

Analytical mean of the process at time t

Implement if available

Source code in quantflow/sp/base.py

def analytical_mean(self, t: FloatArrayLike) -> FloatArrayLike:
    """Analytical mean of the process at time `t`

    Implement if available
    """
    raise NotImplementedError

analytical_variance ¶

analytical_variance(t)

Analytical variance of the process at time t

Implement if available

Source code in quantflow/sp/base.py

def analytical_variance(self, t: FloatArrayLike) -> FloatArrayLike:
    """Analytical variance of the process at time `t`

    Implement if available
    """
    raise NotImplementedError

analytical_pdf ¶

analytical_pdf(t, x)

Analytical pdf of the process at time t

Implement if available

Source code in quantflow/sp/base.py

def analytical_pdf(self, t: FloatArrayLike, x: FloatArrayLike) -> FloatArrayLike:
    """Analytical pdf of the process at time `t`

    Implement if available
    """
    raise NotImplementedError

analytical_cdf ¶

analytical_cdf(t, x)

Analytical cdf of the process at time t

Implement if available

Source code in quantflow/sp/base.py

def analytical_cdf(self, t: FloatArrayLike, x: FloatArrayLike) -> FloatArrayLike:
    """Analytical cdf of the process at time `t`

    Implement if available
    """
    raise NotImplementedError

marginal ¶

marginal(t)

Marginal distribution of the process at time t

Source code in quantflow/sp/base.py

def marginal(self, t: FloatArrayLike) -> StochasticProcess1DMarginal:
    """Marginal distribution of the process at time `t`"""
    return StochasticProcess1DMarginal(process=self, t=t)

domain_range ¶

domain_range()

Source code in quantflow/sp/base.py

def domain_range(self) -> Bounds:
    return default_bounds()

frequency_range ¶

frequency_range(std, max_frequency=None)

Maximum frequency when calculating characteristic functions

Source code in quantflow/sp/base.py

def frequency_range(self, std: float, max_frequency: float | None = None) -> Bounds:
    """Maximum frequency when calculating characteristic functions"""
    if max_frequency is None:
        max_frequency = np.sqrt(40 / std / std)
    return Bounds(0, max_frequency)

support ¶

support(mean, std, points)

Support of the process at time t

Source code in quantflow/sp/base.py

def support(self, mean: float, std: float, points: int) -> FloatArray:
    """Support of the process at time `t`"""
    bounds = self.domain_range()
    start = float(sigfig(bound_from_any(bounds.lb, mean - std)))
    end = float(sigfig(bound_from_any(bounds.ub, mean + std)))
    return np.linspace(start, end, points + 1)

Barndorff-Nielson & Shephard process¶

quantflow.sp.bns.BNS pydantic-model ¶

variance_process pydantic-field ¶

rho pydantic-field ¶

create classmethod ¶

characteristic_exponent ¶

sample ¶

sample_from_draws ¶

characteristic ¶

convexity_correction ¶

analytical_std ¶

analytical_mean ¶

analytical_variance ¶

analytical_pdf ¶

analytical_cdf ¶

marginal ¶

domain_range ¶

frequency_range ¶

support ¶

quantflow.sp.bns.BNS2 pydantic-model ¶

bns1 pydantic-field ¶

bns2 pydantic-field ¶

weight pydantic-field ¶

characteristic_exponent ¶

sample ¶

sample_from_draws ¶

characteristic ¶

convexity_correction ¶

analytical_std ¶

analytical_mean ¶

analytical_variance ¶

analytical_pdf ¶

analytical_cdf ¶

marginal ¶

domain_range ¶

frequency_range ¶

support ¶

quantflow.sp.bns.BNS `pydantic-model` ¶

variance_process `pydantic-field` ¶

rho `pydantic-field` ¶

create `classmethod` ¶

quantflow.sp.bns.BNS2 `pydantic-model` ¶

bns1 `pydantic-field` ¶

bns2 `pydantic-field` ¶

weight `pydantic-field` ¶