# Gamma processes

$\renewcommand{\var}{\operatorname{Var}} \renewcommand{\dd}{\mathrm{d}} \renewcommand{\bb}{\mathbb{#1}} \renewcommand{\vv}{\boldsymbol{#1}} \renewcommand{\rv}{\mathsf{#1}} \renewcommand{\gvn}{\mid} \renewcommand{\Ex}{\mathbb{E}} \renewcommand{\Pr}{\mathbb{P}}$

Gamma processes provide the classic subordinator models, i.e. non-decreasing Lévy processes. By “gamma process” in fact I mean specifically a Lévy process with gamma increments. Gamma processes are a natural model for spiky things

Other processes that happen to have gamma marginals, e.g. are, rather, extensions and are handled separately..

Ground zero for these processes appears to be , and then the weaponisation of these processes to construct the Dirichlet prior in .

Tutorial introductions to gamma processes can be found in . Existence proofs etc are deferred to those sources. You could also see Wikipedia, although that article is not particularly helpful.

The marginal distribution of the an increment of duration $$t$$ is given by the gamma distribution, which we had better cover first.

## Gamma distribution

Let us take a brief divergence into the gamma distribution which is the increment distribution the gamma process (Indeed, all divisible distributions induce Lévy processes.)

The density $$g(x;t,\alpha, \lambda )$$ of the univariate gamma is

$g(x; \alpha, \lambda)= \frac{ \lambda^{\alpha} }{ \Gamma (\alpha) } x^{\alpha\,-\,1}e^{-\lambda x}, x\geq 0.$ This is the shape-rate parameterisation, with rate $$\lambda$$ and shape $$\alpha,$$. We can think of the gamma distribution as the distribution at time 1 of a gamma process.

If $$\rv{g}\sim \operatorname{Gamma}(\alpha, \lambda)$$ then $$\bb E(\rv{g})=\alpha/\lambda$$ and $$\var(\rv{g})=\alpha/\lambda^2.$$

We use various facts about the gamma distribution which quantify its divisibility properties.

1. If $$\rv{g}_1\sim \operatorname{Gamma}(\alpha_1, \lambda),\,\rv{g}_2\sim \operatorname{Gamma}(\alpha_2, \lambda),$$ and $$\rv{g}_1\perp \rv{g}_2,$$ then $$\rv{g}_1+\rv{g}_2\sim \operatorname{Gamma}(\alpha_1+\alpha_2, \lambda)$$ (additive rule)
2. If $$\rv{g}\sim \operatorname{Gamma}(\alpha, \lambda)$$ then $$c \rv{g}\sim \operatorname{Gamma}(\alpha, \lambda/c)$$ (multiplicative rule)
3. If $$\rv{g}_1\sim \operatorname{Gamma}(\alpha_1, \lambda)\perp \rv{g}_1\sim \operatorname{Gamma}(\alpha_2, \lambda)$$ then $$\frac{\rv{g}_1}{\rv{g}_1+\rv{g}_2}\sim \operatorname{Beta}(\alpha_1, \alpha_2)$$ independent of $$\rv{g}_1+\rv{g}_2$$ (stick-breaking rule)

### Moments

Also note that the moment generating function of the gamma distribution is

$\bb{E}[\exp(\rv{g} s)]=\left(1-{\frac {s}{\lambda }}\right)^{-\alpha }{\text{ for }}s<\lambda$ which gives us expressions for all moments fairly easily;

\begin{aligned} \bb{E}[\rv{g}]&=\left.\frac{\dd}{\dd s}\bb{E}[\exp(\rv{g} t)]\right|_{s=0}\\ &=\left.\frac{\dd}{\dd s}\left(1-{\frac {s}{\lambda }}\right)^{-\alpha }\right|_{s=0}\\ &=\left.\frac{\alpha}{\lambda}\left(1-{\frac {s}{\lambda }}\right)^{-\alpha -1}\right|_{s=0}\\ &=\alpha/\lambda\\ \bb{E}[\rv{g}^2] &=\left.\frac{\dd}{\dd s}\frac{\alpha}{\lambda}\left(1-{\frac {s}{\lambda }}\right)^{-\alpha -1}\right|_{s=0}\\ &=\left.\frac{\alpha^2+\alpha}{\lambda^2}\left(1-{\frac {s}{\lambda }}\right)^{-\alpha -2}\right|_{s=0}\\ &=\frac{\alpha(\alpha+1)}{\lambda^2}\\ \bb{E}[\rv{g}^3] &=\left.\frac{\dd}{\dd s}\frac{\alpha^2+\alpha}{\lambda^2}\left(1-{\frac {s}{\lambda }}\right)^{-\alpha -2}\right|_{s=0}\\ &=\frac{\alpha(\alpha+1)(\alpha+2)}{\lambda^3}\left(1-{\frac {s}{\lambda }}\right)^{-\alpha -2}\\ &=\frac{\alpha(\alpha+1)(\alpha+2)}{\lambda^3}\\ \bb{E}[\rv{g}^4] &=\frac{\alpha(\alpha+1)(\alpha+2)(\alpha+3)}{\lambda^4}\\ &\dots\\ \bb{E}[\rv{g}^n] &=\frac{\langle \alpha \rangle_{n}}{\lambda^n}\\ \end{aligned}

Here $$\langle \alpha \rangle_{n}:=\frac{\Gamma(\alpha+n)}{\Gamma(\alpha)}$$ is the rising factorial.

## Multivariate gamma distribution with dependence

I am not sure what general correlations are possible here, but one obvious possibility is to choose a transform matrix $$M$$ with non-negative entries. Then the RV $$\{M\rv{g}(t)\}$$ is still marginally a gamma variate, but the components of the vector are no longer independent. Is this the most general possible gamma distribution? What is the covariance structure of that process? 🏗

## Gamma superpositions

If we wish to know the distribution of the sum of a set of scaled gamma random variables, we can use moment-generating-function approaches . It does not come out so simply as in the Gaussian case, with a recursive coefficient definition.

## The Gamma process

The univariate gamma process $$\{\rv{g}(t;\alpha,\lambda)\}_t$$ is an independent-increment process, with time index $$t$$ and parameters by $$\alpha, \lambda.$$ We assume it is started at $$\rv{g}(0)=0$$.

The marginal density $$g(x;t,\alpha, \lambda )$$ of the process at time $$t$$ is a gamma RV, specifically, $g(x;t, \alpha, \lambda) =\frac{ \lambda^{\alpha t} } { \Gamma (\alpha t) } x^{\alpha t\,-\,1}e^{-\lambda x}, x\geq 0.$ That is, $$\rv{g}(t) \sim \operatorname{Gamma}(\alpha(t_{i+1}-t_{i}), \lambda)$$.

which corresponds to increments per unit time in terms of $$\bb E(\rv{g}(1))=\alpha/\lambda$$ and $$\var(\rv{g}(1))=\alpha/\lambda^2.$$

Note that if $$\alpha t=1,$$ then $$\rv{g}(t;\alpha ,\lambda )\sim \operatorname{Exp}(\lambda).$$

This leads to a method for simulating a path of a gamma process at a sequence of increasing times, $$\{t_1, t_2, t_3, t_L\}.$$ Given $$\rv{g}(t_1;\alpha, \lambda),$$ we know that the increments are distributed as independent variates $$\rv{g}_i:=\rv{g}(t_{i+1})-\rv{g}(t_{i})\sim \operatorname{Gamma}(\alpha(t_{i+1}-t_{i}), \lambda)$$. Presuming we may simulate from the Gamma distribution, it follows that

$\rv{g}(t_i)=\sum_{j < i}\left( \rv{g}(t_{i+1})-\rv{g}(t_{i})\right)=\sum_{j < i} \rv{g}_j.$

A standard $$d$$-dimensional gamma process is the concatenation of $$d$$ independent univariate gamma processes.

## Gamma bridge

Consider a univariate gamma process, $$\rv{g}(t)$$ with $$\rv{g}(0)=0.$$ The gamma bridge, analogous to the Brownian bridge, is the gamma process conditional upon attaining a fixed the value $$S=\rv{g}(1)$$ at terminal time $$1.$$ We write $$\rv{g}_{S}:=\{\rv{g}(t)\mid \rv{g}(1)=S\}_{0< t < 1}$$ for the paths of this process.

We can simulate from the gamma bridge easily. Given the increments of the process are independent, if we have a gamma process $$\rv{g}$$ on the index set $$[0,1]$$ such that $$\rv{g}(1)=S$$, then we can simulate from the bridge paths which connect these points at intermediate time $$t,\, 0<t<1$$ by recalling that we have known distributions for the increments; in particular $$\rv{g}(t)\sim\operatorname{Gamma}(\alpha, \lambda)$$ and $$\rv{g}(1)-\rv{g}(t)\sim\operatorname{Gamma}(\alpha (1-t), \lambda)$$ and these increments, as increments over disjoints sets, are themselves independent. Then, by the stick breaking rule,

$\frac{\rv{g}(t)}{\rv{g}(1)}\sim\operatorname{Beta}(\alpha t, \alpha(1-t))$ independent of $$\rv{g}(1).$$ We can therefore sample from a path of the bridge $$\rv{g}_{S}(t)$$ for some $$t< 1$$ by simulating $$\rv{g}_{S}(t)=B S,$$ where $$B\sim \operatorname{Beta}(\alpha (t),\alpha (1-t)).$$

## Time-warped gamma process

walks us through the mechanics of (deterministically) time-warping Gamma processes, which ends up being not too unpleasant. Predictable stochastic time-warps look like they should be OK. See for an application. Why bother? Linear superpositions of gamma processes can be hard work, and sometime the generalisation from time-warping can come out nicer. 🏗

## Matrix gamma processes

I am tempted to look at matrix-valued extensions, much as the Wishart distribution is a matrix valued extension of the gamma distribution. “Wishart processes” are indeed a thing but the common definition, unlike the common definition of the gamma process, is not necessarily (element-wise, or in any other sense) monotonic. That is, it generalises the square Bessel process, which is indeed marginally gamma distributed (specifically $$\chi^2$$ distributed) but also non-monotonic, so it is not a natural extension of what we have here.

## Centred gamma process

Define $$H_t:=\rv{g}_t-\alpha t /\lambda.$$ Then we have a mean-zero process, in fact, a martingale, since we have subtracted the compensator from it.

We call the marginal distribution at time $$t=1$$ a centred gamma distribution, and will recycle the letter $$H=G-\alpha /\lambda$$. As a linear transform of a random variable, the MGF of the centred gamma distribution is (ignoring questions of region of convergence)…

\begin{aligned} \bb{E}[\exp(H s)]&=\exp (\alpha s/\lambda)\bb{E}[\exp(G s)]\\ &=\exp (\alpha s/\lambda)\left(1-{\frac {s}{\lambda }}\right)^{-\alpha }\\ &=\exp (\alpha s/\lambda - \alpha \log(1-s/\lambda))\end{aligned}

## As a Lévy process

For arguments $$x, t>0$$ and parameters $$\alpha, \lambda>0,$$ we have the increment density as simply a gamma density:

$p_{X}(t, x)=\frac{\lambda^{\alpha t} x^{\alpha t-1} \mathrm{e}^{-x \lambda}}{ \Gamma(\alpha t)},$

This gives us a spectrally positive Lévy measure

$\pi_{X}(x)=\frac{\alpha}{x} \mathrm{e}^{-\lambda x}$

and Laplace exponent

$\Phi_{X}(z)=\alpha \ln (1+ z/\lambda), z \geq 0.$

That is, the Poisson rate, with respect to time $$t$$ of jumps whose size is in the range $$[x, x+dx)$$, is $$\pi(x)dx.$$ We think of this as an infinite superposition of Poisson processes driving different sized jumps, where the jumps are mostly tiny. This is how I think about Lévy process theory, at least.