Abstract

Given noisy, partial observations of a time-homogeneous, finite-statespace Markov chain, conceptually simple, direct statistical inference is available, in theory, via its rate matrix, or infinitesimal generator, {mathsf {Q}}, since exp ({mathsf {Q}}t) is the transition matrix over time t. However, perhaps because of inadequate tools for matrix exponentiation in programming languages commonly used amongst statisticians or a belief that the necessary calculations are prohibitively expensive, statistical inference for continuous-time Markov chains with a large but finite state space is typically conducted via particle MCMC or other relatively complex inference schemes. When, as in many applications {mathsf {Q}} arises from a reaction network, it is usually sparse. We describe variations on known algorithms which allow fast, robust and accurate evaluation of the product of a non-negative vector with the exponential of a large, sparse rate matrix. Our implementation uses relatively recently developed, efficient, linear algebra tools that take advantage of such sparsity. We demonstrate the straightforward statistical application of the key algorithm on a model for the mixing of two alleles in a population and on the Susceptible-Infectious-Removed epidemic model.

Highlights

  • A reaction network is a stochastic model for the joint evolution of one or more populations of species

  • Sherlock is encapsulated by the number of each species that is present, and the system evolves via a set of reactions: Poisson processes whose rates depend on the current state

  • Both by way of motivation and because we shall use them later to illustrate our method, we present two examples of continuous-time Markov processes, where a finite, sparse rate matrix contains all of the information about the dynamics

Read more

Summary

Introduction

A reaction network is a stochastic model for the joint evolution of one or more populations of species. In cases where the number of states, d, is finite, direct exact likelihood-based inference is available via the exponential of the infinitesimal generator for the continuous-time Markov chain, or rate matrix, Q. As well as being directly of use for models with finite state spaces, exponentials of finite rate matrices can be used to perform inference on Markov jump processes with a countably infinite statespace; see Georgoulas et al (2017) and Sherlock and Golightly (2019) The latter uses the uniformisation and scaling and squaring algorithms as described in this article, while the former uses the less efficient but more general algorithm of Al-Mohy and Higham (2011)

Examples and motivation
Data and likelihood calculations
Likelihood for noisy and partially observed data
Statespace reduction for epidemic models
Matrix exponentiation
The uniformisation algorithm
Scaling and squaring
Improvements
Implementation
Numerical comparisons and demonstrations
Comparison with other matrix exponentiation algorithms
Bayesian inference for the Swansea measles epidemic of 2013
Discussion
Implementation details
Proof of Theorem 1
C Functions in the expQ package
D Latent-variable updates for the SIR model
E Log-likelihood R code for the Moran model
Full Text
Published version (Free)

Talk to us

Join us for a 30 min session where you can share your feedback and ask us any queries you have

Schedule a call