Functions

Documentation for all the functions in jumpdiff.

Jump-diffusion timeseries generator

jumpdiff.jd_process.jd_process(time: float, delta_t: float, a: callable, b: callable, xi: float, lamb: float, init: float = None, solver: str = 'Euler', b_prime: callable = None) → numpy.ndarray

Integrates a jump-diffusion process with drift a(x), diffusion b(x), jump amplitude xi (\(\xi\)), and jump rate lamb (\(\lambda\)).

\[\mathrm{d} X(t) = a(x,t)\;\mathrm{d} t + b(x,t)\;\mathrm{d} W(t) + \xi\;\mathrm{d} J(t),\]

with \(J\) Poisson with jump rate \(\lambda\). This integrator has both an Euler─Maruyama and a Milstein method of integration. For Milstein one has to introduce the derivative of the diffusion term b, denoted b_prime.

Parameters:

• time (float > 0) – Total integration time. Positive float or int.
• delta_t (float > 0) – Time sampling, the smaller the better.
• a (callable) –

  The drift function. Can be a function of a lambda. For an Ornstein─Uhlenbeck process with drift -2x, a takes the form

  a =  lambda x: -2x.

• b (callable) –

  The diffusion function. Can be a function of a lambda. For an Ornstein─Uhlenbeck process with diffusion 1, a takes the form

  b =  lambda x: 1.

• xi (float > 0) – Variance of the jump amplitude, which will be turned into a normal distribution like \(\mathcal{N}\)(0,√xi).
• lamb (float > 0) – Jump rate of the Poissonian jumps. This is implemented as the numpy function np.random.poisson(lam = lamb * delta_t).
• init (float (defaul None)) – Initial conditions. If None given, generates a random value from a normal distribution ~ \(\mathcal{N}\)(0,√delta_t).
• solver ('Euler' or 'Milstein' (defaul 'Euler')) – The regular Euler─Maruyama solver ‘Euler’ is the default, with an order of √delta_t. To employ a state-dependent diffusion, i.e., b(x) as a function of x, the Milstein scheme has an order of delta_t. You must introduce as well the derivative of b(x), i.e., b’(x), as the argument b_prime.

Returns:

X – Timeseries of size int(time/delta_t)

Return type:

np.array

Moments

jumpdiff.moments.moments(timeseries: numpy.ndarray, bw: float = None, bins: numpy.ndarray = None, power: int = 6, lag: list = [1], correction: bool = True, norm: bool = False, kernel: callable = None, tol: float = 1e-10, conv_method: str = 'auto', verbose: bool = False) → numpy.ndarray

Estimates the moments of the Kramers─Moyal expansion from a timeseries using a Nadaraya─Watson kernel estimator method. These later can be turned into the drift and diffusion coefficients after normalisation.

Parameters:

• timeseries (np.ndarray) – A 1-dimensional timeseries.
• bw (float) – Desired bandwidth of the kernel. A value of 1 occupies the full space of the bin space. Recommended are values 0.005 < bw < 0.4.
• bins (np.ndarray (default None)) – The number of bins for each dimension, defaults to np.array([5000]). This is the underlying space for the Kramers─Moyal conditional moments.
• power (int (default 6)) – Upper limit of the the Kramers─Moyal conditional moments to calculate. It will generate all Kramers─Moyal conditional moments up to power.
• lag (list (default 1)) – Calculates the Kramers─Moyal conditional moments at each indicated lag, i.e., for timeseries[::lag[]]. Defaults to 1, the shortest timestep in the data.
• corrections (bool (default True)) – Implements the second-order corrections of the Kramers─Moyal conditional moments directly
• norm (bool (default False)) – Sets the normalisation. False returns the Kramers─Moyal conditional moments, and True returns the Kramers─Moyal coefficients.
• kernel (callable (default None)) –

  Kernel used to convolute with the Kramers─Moyal conditional moments. To select example an Epanechnikov kernel use

  kernel = kernels.epanechnikov

  If None the Epanechnikov kernel will be used.

• tol (float (default 1e-10)) – Round to zero absolute values smaller than tol, after convolutions.
• conv_method (str (default auto)) – A string indicating which method to use to calculate the convolution. docs.scipy.org/doc/scipy/reference/generated/scipy.signal.convolve.
• verbose (bool (default False)) – If True will report on the bandwidth used.

Returns:

• edges (np.ndarray) – The bin edges with shape (D,bins.shape) of the calculated moments.
• moments (np.ndarray) – The calculated moments from the Kramers─Moyal expansion of the timeseries at each lag. To extract the selected orders of the moments, use moments[i,:,j], with i the order according to powers, j the lag (if any given).

jumpdiff.moments.corrections(m: numpy.ndarray, power: int)

The moments function will by default apply the corrections. You can turn the corrections off in that fuction by setting corrections = False.

Second-order corrections of the Kramers─Moyal coefficients (conditional moments), given by

\[\begin{split}F_1 &= M_1,\\ F_2 &= \frac{1}{2}\! \left(M_2-M_1^2\right ), \\ F_3 &= \frac{1}{6}\! \left( M_3-3M_1M_2+3M_1^3 \right ), \\ F_4 &= \frac{1}{24}\! \left(M_4-4M_1M_3+18M_1^2M_2-3M_2^2 -15M_1^4 \right) , \\ F_5 &= \frac{1}{120}\!\left(M_5 -5 M_1 M_4 +30 M_1^2 M_3 -150 M_1^3 M_2 +45 M_1 M_2^2-10 M_2 M_3 + 105 M_1^5 \right) \\ F_6 &= \frac{1} {720}\! \left (M_6 -6 M_1 M_5 + 45 M_1^2 M_4 -300 M_1^3 M_3 +1575 M_1^4 M_2-675 M_1^2 M_2^2 \right. \\ & \qquad \left . \ +180 M_1 M_2 M_3+45 M_2^3 -15 M_2 M_4-10 M_3^2 - 945 M_1^6\right ) \, ,\end{split}\]

with the prefactor the normalisation, i.e., the normalised results are the Kramers─Moyal coefficients. If norm is False, this results in the Kramers─Moyal conditional moments.

Parameters:

• (moments) (m) – The calculated conditional moments from the Kramers─Moyal expansion of the at each lag. To extract the selected orders of the moments use moments[i,:,j], with i the order according to powers, j the lag.
• power (int) – Upper limit of the Kramers─Moyal conditional moments to calculate. It will generate all Kramers─Moyal conditional moments up to power.

Returns:

F – The corrections of the calculated Kramers─Moyal conditional moments from the Kramers─Moyal expansion of the timeseries at each lag. To extract the selected orders of the moments, use F[i,:,j], with i the order according to powers, j the lag (if any introduced).

Return type:

np.ndarray

Parameters

jumpdiff.parameters.jump_amplitude(moments: numpy.ndarray, tol: float = 1e-10, full: bool = False, verbose: bool = False) → numpy.ndarray

Retrieves the jump amplitude xi (\(\xi\)) via

\[\lambda(x,t) = \frac{M_4(x,t)}{3\sigma_{\xi}^4}.\]

Take notice that the different normalisation of the moments leads to a different results.

Parameters:

• moments (np.ndarray) – Moments extracted with the function moments. Needs moments up to order 6.
• tol (float (defaul 1e-10)) – Toleration for the division of the moments.
• full (bool (defaul False)) – If True returns also the (biased) weighed standard deviation of the averaging process.
• verbose (bool (defaul True)) – Prints the result.

Returns:

xi_est – Estimator of the jump amplitude xi (\(\xi\)).

Return type:

np.ndarray

References

Anvari, M., Tabar, M. R. R., Peinke, J., Lehnertz, K., ‘Disentangling the stochastic behavior of complex time series.’ Scientific Reports, 6, 35435, 2016. doi: 10.1038/srep35435.

Lehnertz, K., Zabawa, L., and Tabar, M. R. R., ‘Characterizing abrupt transitions in stochastic dynamics.’ New Journal of Physics, 20(11):113043, 2018. doi: 10.1088/1367-2630/aaf0d7.

jumpdiff.parameters.jump_rate(moments: numpy.ndarray, xi_est: numpy.ndarray = None, tol: float = 1e-10, full: bool = False, verbose: bool = False) → numpy.ndarray

Retrieves the jump rate lamb (\(\lambda\)) via

\[\sigma_{\xi}^2 = \frac{M_6(x,t)}{5M_4(x,t)}.\]

Take notice that the different normalisation of the moments leads to a different results.

Parameters:

• moments (np.ndarray) – moments extracted with the function ‘moments’. Needs moments of order 6.
• tol (float (defaul 1e-10)) – Toleration for the division of the moments.
• full (bool (defaul False)) – If True returns also the (biased) weighed standard deviation of the averaging process.
• verbose (bool (defaul True)) – Prints the result.

Returns:

xi_est – Estimator on the jump rate lamb (\(\lambda\))

Return type:

np.ndarray

References

Anvari, M., Tabar, M. R. R., Peinke, J., Lehnertz, K., ‘Disentangling the stochastic behavior of complex time series.’ Scientific Reports, 6, 35435, 2016. doi: 10.1038/srep35435.

Lehnertz, K., Zabawa, L., and Tabar, M. R. R., ‘Characterizing abrupt transitions in stochastic dynamics.’ New Journal of Physics, 20(11):113043, 2018. doi: 10.1088/1367-2630/aaf0d7.

Q-ratio

jumpdiff.q_ratio.q_ratio(lag: numpy.ndarray, timeseries: numpy.ndarray, loc: int = None, correction: bool = False) → numpy.ndarray

q_ratio method to distinguish pure diffusion from jump-diffusion timeseries, Given by the relation of the 4th and 6th Kramers─Moyal coefficient with increasing lag

\[\begin{split}Q(x,\tau) = \frac{D_6(x,\tau)}{5 D_4(x,\tau)} = \left\{\begin{array}{ll} b(x)^2 \tau, & \text{diffusive} \\ \sigma_\xi^2(x), & \text{jumpy} \end{array}\right.\end{split}\]

Parameters:

• lag (np.ndarray of ints) – An array with the time-lag to extract the Kramers–Moyal coefficient for different lags.
• timeseries (np.ndarray) – A 1-dimensional timeseries.
• loc (float (defaul None)) – Use a particular point in space to calculate the ratio. If None given, the maximum of the probability density function is taken.
• corrections (bool (defaul False)) – Select whether to use corrective terms.

Returns:

• lag (np.ndarray of ints) – Same as input, but only lag > 0 and as ints.
• ratio (np.ndarray of len(lag)) – Ratio of the sixth-order over forth-order Kramers–Moyal coefficient.

References

Anvari, M., Tabar, M. R. R., Peinke, J., Lehnertz, K., ‘Disentangling the stochastic behavior of complex time series.’ Scientific Reports, 6, 35435, 2016. doi: 10.1038/srep35435.

Lehnertz, K., Zabawa, L., and Tabar, M. R. R., ‘Characterizing abrupt transitions in stochastic dynamics.’ New Journal of Physics, 20(11):113043, 2018. doi: 10.1088/1367-2630/aaf0d7.

Formulae

jumpdiff.formulae.m_formula(power, tau=True)

Generate the formula for the conditional moments with second-order corrections based on the relation with the ordinary Bell polynomials

\[M_n(x^{\prime},\tau) \sim (n!)\tau D_n(x^{\prime}) + \frac{(n!)\tau^2}{2} \sum_{m=1}^{n-1} D_m(x^{\prime}) D_{n-m}(x^{\prime})\]

Parameters:power (int) – Desired order of the formula. Returns:term – Expression up to given power. Return type:sympy.symbols

jumpdiff.formulae.f_formula(power)

Generate the formula for the conditional moments with second-order corrections based on the relation with the ordinary Bell polynomials

\[\begin{split}D_n(x) &= \frac{1}{\tau (n!)} \bigg[ \hat{B}_{n,1} \left(M_1(x,\tau),M_2(x,\tau),\ldots,M_{n}(x,\tau)\right) \\ &\qquad \left.-\frac{\tau}{2} \hat{B}_{n,2} \left(M_1(x,\tau),M_2(x,\tau),\ldots,M_{n-1}(x,\tau)\right)\right].\end{split}\]

Parameters:power (int) – Desired order of the formula. Returns:term – Expression up to given power. Return type:sympy.symbols

jumpdiff.formulae.f_formula_solver(power)

Generate the reciprocal relation of the moments to the Kramers─Moyal coefficients by sequential iteration.

\[\begin{split}D_n(x) &= \frac{1}{\tau (n!)} \bigg[ \hat{B}_{n,1} \left(M_1(x,\tau),M_2(x,\tau),\ldots,M_{n}(x,\tau)\right) \\ &\qquad \left.-\frac{\tau}{2} \hat{B}_{n,2} \left(M_1(x,\tau),M_2(x,\tau),\ldots,M_{n-1}(x,\tau)\right)\right].\end{split}\]

Parameters:power (int) – Desired order of the formula. Returns:term – Expression up to given power. Return type:sympy.symbols

Helping functions

Kernels function

jumpdiff.kernels.kernel(kernel_func)

Transforms a kernel function into a scaled kernel function (for a certain bandwidth bw).

Currently implemented kernels are:

Epanechnikov, Gaussian, Uniform, Triangular, Quartic.

For a good overview of various kernels see https://en.wikipedia.org/wiki/Kernel_(statistics)

jumpdiff.kernels.volume_unit_ball(dims: int) → float

Returns the volume of a unit ball in dimensions dims.

jumpdiff.kernels.epanechnikov(x: numpy.ndarray, dims: int) → numpy.ndarray

The Epanechnikov kernel in dimensions dims.

jumpdiff.kernels.gaussian(x: numpy.ndarray, dims: int) → numpy.ndarray

Gaussian kernel in dimensions dims.

jumpdiff.kernels.uniform(x: numpy.ndarray, dims: int) → numpy.ndarray

Uniform, or rectangular kernel in dimensions dims

jumpdiff.kernels.triagular(x: numpy.ndarray, dims: int) → numpy.ndarray

Triagular kernel in dimensions dims

jumpdiff.kernels.quartic(x: numpy.ndarray, dims: int) → numpy.ndarray

Quartic, or biweight kernel in dimensions dims