The kernels submodule

This submodule contains the Random-Walk Metropolis, Metropolis-Adjusted Langevin Algorithm (MALA) and simplified manifold MALA kernels and a base kernel base class. Those are responsible for the actual steps of MCMC.

class mhn.mcmc.kernels.Kernel(grad_and_log_likelihood: tuple[Callable[[ndarray], tuple[ndarray, float]]], log_prior: Callable[[ndarray], float], shape: tuple[int, int], rng: Generator | None = None)

Base class for kernels used in MCMC sampling.

grad_and_log_likelihood

function that takes a flattened (!) log_theta matrix and returns a tuple of the gradient of the log-likelihood and the log-likelihood itself.

Important: Here, the likelihood should not be normalized by the number of samples like in Schill et al. (2019).

Type:

tuple[Callable[[np.ndarray], tuple[np.ndarray, float]]]

log_prior

Callable[[np.ndarray], float]: function that takes a flattened (!) log_theta matrix and returns the log-prior. This corresponds to the negative penalty term times lambda, multiplied by the number of samples, if applicable.

shape

tuple[int, int]: the shape of the log_theta matrix, i.e. (n_events + 1, n_events) for oMHN and (n_events, n_events) for cMHN.

rng

np.random.Generator | None = None: random number generator to be used for sampling. If None, a new RNG will be created.

class mhn.mcmc.kernels.MALAKernel(step_size: float, grad_and_log_likelihood: tuple[Callable[[ndarray], tuple[ndarray, float]]], log_prior: Callable[[ndarray], float], log_prior_grad: Callable[[ndarray], ndarray], shape: tuple[int, int], rng: Generator | None = None)

Class for kernel used in the Metropolis-Adjusted Langevin Algorithm. Extends Kernel.

step_size

step size for the MALA sampling.

Type:

float

grad_and_log_likelihood

function that takes a flattened (!) log_theta matrix and returns a tuple of the gradient of the log-likelihood and the log-likelihood itself.

Important: Here, the likelihood should not be normalized by the number of samples like in Schill et al. (2019)

Type:

tuple[Callable[[np.ndarray], tuple[np.ndarray, float]]]

log_prior

Callable[[np.ndarray], float]: function that takes a flattened (!) log_theta matrix and returns the log-prior. This corresponds to the negative penalty term times lambda, multiplied by the number of samples, if applicable.

log_prior_grad

Callable[[np.ndarray], np.ndarray]: function that takes a flattened (!) log_theta matrix and returns the gradient of the log-prior.

shape

tuple[int, int]: the shape of the log_theta matrix, i.e. (n_events + 1, n_events) for oMHN and (n_events, n_events) for cMHN.

rng

np.random.Generator | None = None: random number generator to be used for sampling. If None, a new RNG will be created.

get_params(initial_step: ndarray) → MALAResult

Get the results for a given step.

Parameters:

initial_step (np.ndarray) – The step for which to get the results.

Returns:

The results for the given step.

Return type:

MALAResult

log_accept(prev_step: ndarray, prev_step_res: MALAResult, new_step: ndarray, new_step_res: MALAResult) → float

Give acceptance ratio of accepting the new step.

Parameters:

• prev_step (np.ndarray) – the previous step

• prev_step_res (MALAResult) – Results of the previous step.

• new_step (np.ndarray) – the new step

• new_step_res (MALAResult) – Results of the new step.

Returns:

the acceptance ratio

Return type:

float

log_q(theta_proposed: ndarray, mu: ndarray) → float

Compute the logarithm of the proposal distribution density q(theta_new | theta) for the MMALA algorithm. This is according to https://en.wikipedia.org/wiki/Multivariate_normal_distribution#Density_function.

Parameters:

• theta_new (np.ndarray) – New proposed theta

• G (np.ndarray) – Metric tensor w.r.t. the old theta

• det_sqrt_G (float) – Square root of the determinant of the

• theta (the old)

• mu (np.ndarray) – Mean of the proposal distribution w.r.t.

• theta

Returns:

The logarithm of the proposal distribution q(theta_new | theta) density

Return type:

float

one_step(prev_step: ndarray, prev_step_res: MALAResult, return_info: bool = False) → tuple[ndarray, MALAResult] | tuple[ndarray, MALAResult, float, int]

Perform one MALA step.

Parameters:

• prev_step (np.ndarray) – the previous step

• prev_step_res (MALAResult) – Results of the previous step.

• return_info (bool, optional) – Whether to return the acceptance ratio and acceptance indicator. Defaults to False.

Returns:

tuple[np.ndarray, MALAResult] | tuple[np.ndarray, MALAResult, float, int]: The new step and

its results, and optionally acceptance ratio and acceptance indicator.

propose(prev_step: ndarray, prev_step_res: MALAResult) → tuple[ndarray, MALAResult]

Propose a new step based on the previous step and its results.

Parameters:

• prev_step (np.array) – the previous step

• prev_step_res (MALAResult) – Results of the previous step.

Returns:

the new step and its results.

Return type:

tuple[np.ndarray, MALAResult]

class mhn.mcmc.kernels.MALAResult(log_likelihood: float, log_prior: float, mu: np.ndarray)

Results for one smMALAKernel step.

Parameters:

• log_likelihood (float) – log-likelihood of the current step

• log_prior (float) – log-prior of the current step

• mu (np.ndarray) – mean of the proposal distribution at the current step

log_likelihood: float

Alias for field number 0

log_prior: float

Alias for field number 1

mu: ndarray

Alias for field number 2

class mhn.mcmc.kernels.RWMKernel(grad_and_log_likelihood: tuple[Callable[[ndarray], ndarray], Callable[[ndarray], float]], log_prior: Callable[[ndarray], float], step_size: float, shape: tuple[int, int], rng: Generator | None = None)

Class for kernel used in simplified Random-Walk Metropolis. Extends Kernel.

grad_and_log_likelihood

function that takes a flattened (!) log_theta matrix and returns a tuple of the gradient of the log-likelihood and the log-likelihood itself.

Important: Here, the likelihood should not be normalized by the number of samples like in Schill et al. (2019)

Type:

tuple[Callable[[np.ndarray], tuple[np.ndarray, float]]]

log_prior

Callable[[np.ndarray], float]: function that takes a flattened (!) log_theta matrix and returns the log-prior. This corresponds to the negative penalty term times lambda, multiplied by the number of samples, if applicable.

scale

scale of the normal proposal distribution.

Type:

float

shape

tuple[int, int]: the shape of the log_theta matrix, i.e. (n_events + 1, n_events) for oMHN and (n_events, n_events) for cMHN.

rng

np.random.Generator | None = None: random number generator to be used for sampling. If None, a new RNG will be created.

get_params(initial_step)

Get the results for a given step.

Parameters:

initial_step (np.ndarray) – The step for which to get the results.

Returns:

The results for the given step.

Return type:

RWMResult

log_accept(prev_step: ndarray, prev_step_res: RWMResult, new_step: ndarray, new_step_res: RWMResult) → float

Give acceptance ratio of accepting the new step. This is given by p(new_step | D) pr(new_step) / p(prev_step|D) pr(prev_step) where p(.| D) is the likelihood and pr(.) is the prior.

Parameters:

• prev_step (np.ndarray) – the previous step

• prev_step_res (RWMResult) – Results of the previous step.

• new_step (np.ndarray) – the new step

• new_step_res (RWMResult) – Results of the new step.

Returns:

the acceptance ratio

Return type:

float

one_step(prev_step: ndarray, prev_step_res: RWMResult, return_info: bool = False)

Perform one RWM step.

Parameters:

• prev_step (np.ndarray) – the previous step

• prev_step_res (RWMResult) – Results of the previous step.

• return_info (bool, optional) – Whether to return the acceptance ratio and acceptance indicator. Defaults to False.

Returns:

tuple[np.ndarray, RWMResult] | tuple[np.ndarray, RWMResult, float, int]: The new step and its results, and optionally acceptance ratio and acceptance indicator.

propose(prev_step: ndarray, prev_step_res: RWMResult)

Propose a new step based on the previous step and its results. The proposal is made according to a normal distribution centered at the previous step with standard deviation step_size.

Parameters:

• prev_step (np.array) – the previous step

• prev_step_res (RWMResult) – Results of the previous step.

Returns:

the new step and its results.

Return type:

tuple[np.ndarray, RWMResult]

class mhn.mcmc.kernels.RWMResult(log_likelihood: float, log_prior: float)

Results for one RWMKernel step.

Parameters:

• log_likelihood (float) – log-likelihood of the current step

• log_prior (float) – log-prior of the current step

log_likelihood: float

Alias for field number 0

log_prior: float

Alias for field number 1

class mhn.mcmc.kernels.smMALAKernel(step_size: float, grad_and_log_likelihood: tuple[Callable[[ndarray], tuple[ndarray, float]]], log_prior: Callable[[ndarray], float], log_prior_grad: Callable[[ndarray], ndarray], log_prior_hessian: Callable[[ndarray], ndarray], shape: tuple[int, int], use_cuda: bool = False, rng: Generator | None = None)

Class for kernel used in simplified manifold MALA. Extends Kernel.

step_size

step size for the smMALA sampling. This is the \(\epsilon^2/4\) in the proposal distribution.

Type:

float

grad_and_log_likelihood

function that takes a flattened (!) log_theta matrix and returns a tuple of the gradient of the log-likelihood and the log-likelihood itself.

Important: Here, the likelihood should not be normalized by the number of samples like in Schill et al. (2019)

Type:

tuple[Callable[[np.ndarray], tuple[np.ndarray, float]]]

log_prior

Callable[[np.ndarray], float]: function that takes a flattened (!) log_theta matrix and returns the log-prior. This corresponds to the negative penalty term times lambda, multiplied by the number of samples, if applicable.

log_prior_grad

Callable[[np.ndarray], np.ndarray]: function that takes a flattened (!) log_theta matrix and returns the gradient of the log-prior.

log_prior_hessian

Callable[[np.ndarray], np.ndarray]: function that takes a flattened (!) log_theta matrix and returns the Hessian of the log-prior.

shape

tuple[int, int]: the shape of the log_theta matrix, i.e. (n_events + 1, n_events) for oMHN and (n_events, n_events) for cMHN.

use_cuda

bool = False: whether to use CUDA for the computation of the Fisher information matrix. Only applicable if mhn.cuda_available() says that CUDA is available.

rng

np.random.Generator | None = None: random number generator to be used for sampling. If None, a new RNG will be created.

get_params(initial_step: ndarray) → smMALAResult

Get the results for a given step.

Parameters:

initial_step (np.ndarray) – The step for which to get the results.

Returns:

The results for the given step.

Return type:

smMALAResult

log_accept(prev_step: ndarray, prev_step_res: smMALAResult, new_step: ndarray, new_step_res: smMALAResult) → float

Give acceptance ratio of accepting the new step.

This is given by

\(\frac{p(\theta' | D) q(\theta | \theta')\pi(\theta')}{p(\theta|D) q(\theta' | \theta) \pi(\theta)}\)

where \(p(.| D)\) is the likelihood, \(\pi(.)\) is the prior and \(q(.|.)\) is the proposal distribution.

Parameters:

• prev_step (np.ndarray) – the previous step

• prev_step_res (smMALAResult) – Results of the previous step.

• new_step (np.ndarray) – the new step

• new_step_res (smMALAResult) – Results of the new step.

Returns:

the acceptance ratio

Return type:

float

log_q(theta_proposed: ndarray, G: ndarray, det_sqrt_G: float, mu: ndarray) → float

Compute the logarithm of the proposal distribution density

q(theta_new | theta) for the smMALA algorithm. This is according to https://en.wikipedia.org/wiki/Multivariate_normal_distribution #Density_function.

Parameters:

• theta_proposed (np.ndarray) – New proposed theta

• G (np.ndarray) – Metric tensor w.r.t. the old theta

• det_sqrt_G (float) – Square root of the determinant of the metric tensor w.r.t. the old theta

• mu (np.ndarray) – Mean of the proposal distribution w.r.t. the old theta

Returns:

The logarithm of the proposal distribution

q(theta_new | theta) density

Return type:

float

one_step(prev_step: ndarray, prev_step_res: smMALAResult, return_info: bool = False) → tuple[ndarray, smMALAResult] | tuple[ndarray, smMALAResult, float, int]

Perform one smMALA step.

Parameters:

• prev_step (np.ndarray) – the previous step

• prev_step_res (smMALAResult) – Results of the previous step.

• return_info (bool, optional) – Whether to return the acceptance ratio and acceptance indicator. Defaults to False.

Returns:

tuple[np.ndarray, smMALAResult] | tuple[np.ndarray, smMALAResult, float, int]: The new step

and its results, and optionally acceptance ratio and acceptance indicator.

propose(prev_step: ndarray, prev_step_res: smMALAResult) → tuple[ndarray, smMALAResult]

Propose a new step based on the previous step and its results. This is done according to

\(\mathcal N\bigg(\theta + \frac{\epsilon}{2} G(\theta)^{-1} \frac{\partial \log p(\theta)}{\partial \theta}, \epsilon G(\theta)^{-1}\bigg)\)

where

\(G(\theta) = I(\theta) - H(\log \pi(\theta))\)

with \(I\) the Fisher information matrix and \(H\) the Hessian of the log-prior.

Parameters:

• prev_step (np.array) – the previous step

• prev_step_res (smMALAResult) – Results of the previous step.

Returns:

the new step and its results.

Return type:

tuple[np.ndarray, smMALAResult]

class mhn.mcmc.kernels.smMALAResult(log_likelihood: float, log_prior: float, gradient: np.ndarray, G: np.ndarray, cholesky: np.ndarray, mu: np.ndarray, det_sqrt: float)

Results for one smMALAKernel step.

Parameters:

• log_likelihood (float) – log-likelihood of the current step

• log_prior (float) – log-prior of the current step

• gradient (np.ndarray) – gradient of the log-posterior at the current step

• G (np.ndarray) – metric tensor at the current step

• cholesky (np.ndarray) – Cholesky decomposition of the metric tensor at the current step

• mu (np.ndarray) – mean of the proposal distribution at the current step

• det_sqrt (float) – square root of the determinant of the metric tensor at the current step

G: ndarray

Alias for field number 3

cholesky: ndarray

Alias for field number 4

det_sqrt: float

Alias for field number 6

gradient: ndarray

Alias for field number 2

log_likelihood: float

Alias for field number 0

log_prior: float

Alias for field number 1

mu: ndarray

Alias for field number 5