A Symbolic Computation Approach to a Problem involving Multivariate Poission Distributions

(Draft, do not distribute!)

Eduardo SONTAG^{1} and Doron ZEILBERGER^{1}

Introduction Eduardo, please write it. The Problem Suppose that we have n independent Poisson random variables, X_{i} (i = 1 … n) with parameters λ_{i} respectively. In other words

P r(X_{1} = k_{1}, …, X_{n} = k_{n}) = \frac{λ_{1}^{k_{1}}}{k_{1}!} \frac{λ_{2}^{k_{2}}}{k_{2}!} … \frac{λ_{n}^{k_{n}}}{k_{n}!}

Suppose that we can’t observe the X_{i}’s directly, but only a certain number, m, of linear combinations of them:

Y_{j} = ∑_{i=1}^{n} a_{ij}X_{j} , (j = 1, …, m) ,

where (a_{ij}) is a certain m × n matrix with non-negative coefficients. We are interested in the following questions:

1. Can one compute (fast!), for any given vector (b_{1}, …, b_{m}) (possibly with large coordinates), the probability

F (b_{1}, …, b_{m}) := P r(Y_{1} = b_{1}, …, Y_{m} = b_{m})

2. Can one compute (fast!), for any given vectors (b_{1}, …, b_{m}), (possibly with large coordinates) the conditional expectation

G_{i}(b_{1}, …, b_{m}) := E[X_{i} | Y_{1} = b_{1}, …, Y_{m} = b_{m}]

3. More generally, can one compute (fast!), the higher moments G_{i}^{(r)}(b_{1}, …, b_{m}) := E[X_{i}^{r} | Y_{1} = b_{1}, …, Y_{m} = b_{m}], (r ≥ 2) .

(which would immediately allow us to compute the moments about the mean.)

1 Department of Mathematics, Rutgers University (New Brunswick), Hill Center-Busch Campus, 110 Frelinghuysen

Rd., Piscataway, NJ 08854-8019, USA. [sontag, zeilberg] at math dot rutgers dot edu , http://www.math.rutgers.edu/~zeilberg . First version: xx.x, 200x. Accompanied by Maple package MVPOISSON downloadable from Zeilberger’s website. Supported in part by the NSF.

The generating function Fix a matrix A = (a_{ij}) (1 ≤ i ≤ m, 1 ≤ j ≤ n). Let f(z_{1}, …, z_{m}) be the (multivariable) generating function of F (b_{1}, …, b_{m}), in other words

f(z_{1}, …, z_{m}) = ∑_{b_{1}≥0,...,b_{m}≥0} F (b_{1}, …, b_{m})z_{1}^{b_{1}} ⋯ z_{m}^{b_{m}} .

We have:

f(z_{1}, …, z_{m}) = ∑_{b_{1}≥0,...,b_{m}≥0} F (b_{1}, …, b_{m})z_{1}^{b_{1}} ⋯ z_{m}^{b_{m}} =

= ∑_{b_{1}≥0,...,b_{m}≥0} (  ∑_{{k_{1},...,kn≥0 \atop a_{11}k_{1}+...+a_{1n}kn=b_{1},...,a_{m1}k_{1}+...+amnkn=bm}} \frac{λ_{1}^{k_{1}}}{k_{1}!} \frac{λ_{2}^{k_{2}}}{k_{2}!} … \frac{λ_{n}^{k_{n}}}{k_{n}!} ) z_{1}^{b_{1}} ⋯ z_{m}^{b_{m}}

By changing the order of summation, this equals

∑_{k_{1}≥0,...,k_{n}≥0} \frac{λ_{1}^{k_{1}}}{k_{1}!} \frac{λ_{2}^{k_{2}}}{k_{2}!} … \frac{λ_{n}^{k_{n}}}{k_{n}!} z_{1}^{a_{11}k_{1}+...+a_{1n}k_{n}} ⋯ z_{m}^{a_{m1}k_{1}+...+a_{mn}k_{n}}

= ∑_{k_{1}≥0,...,k_{n}≥0} \frac{(λ_{1}z_{1}^{a_{11}}z_{2}^{a_{21}} ⋯ z_{m}^{a_{m1}})^{k_{1}}}{k_{1}!} ⋯ \frac{(λ_{n}z_{1}^{a_{1n}}z_{2}^{a_{2n}} ⋯ z_{m}^{a_{mn}})^{k_{n}}}{k_{n}!}

= ( ∑_{k_{1}≥0} \frac{(λ_{1}z_{1}^{a_{11}}z_{2}^{a_{21}} ⋯ z_{m}^{a_{m1}})^{k_{1}}}{k_{1}!} ) ⋯ ( ∑_{k_{n}≥0} \frac{(λ_{n}z_{1}^{a_{1n}}z_{2}^{a_{2n}} ⋯ z_{m}^{a_{mn}})^{k_{n}}}{k_{n}!} )

= exp(λ_{1}z_{1}^{a_{11}}z_{2}^{a_{21}} ⋯ z_{m}^{a_{m1}}) ⋯ exp(λ_{n}z_{1}^{a_{1n}}z_{2}^{a_{2n}} ⋯ z_{m}^{a_{mn}})

= exp (λ_{1}z_{1}^{a_{11}}z_{2}^{a_{21}} ⋯ z_{m}^{a_{m1}} + ⋯ + λ_{n}z_{1}^{a_{1n}}z_{2}^{a_{21}} ⋯ z_{m}^{a_{mn}}) .

We have just derived Theorem 1: The generating function f(z_{1}, …, z_{m}) is given by

f(z) = exp (∑_{j}^{n} λ_{j} ∏_{i=1}^{m} z_{i}^{a_{ij}} ) .

So P (b_{1}, …, b_{m}) is the coefficient of z_{1}^{b_{1}} ⋯ z_{m}^{b_{m}} in the multivariable Taylor expansion about the origin of f(z_{1}, …, z_{m}). To get the factorial moments of X_{i}, E[X_{i}^{(r)}] we apply the operator λ_{i}^{r}(\frac{∂}{∂λ_{i}})^{r} to P (b_{1}, …, b_{m}) and divide by P (b_{1}, …, b_{m}). By doing the same for the generating function, we see that the numerator is λ_{i}^{r} times the coefficient of (∏_{i=1}^{m} z_{i}^{a_{ij}})^{r}, In other words we have Theorem 2: The factorial moments E[X_{j}^{(r)}] are given in terms of the P (b_{1}, …, b_{m}) by

\frac{P (b_{1} − ra_{1j}, b_{2} − ra_{2j}, …, b_{m} − ra_{mj})}{P (b_{1}, …, b_{m})}

So everything depends on a fast computation of the coefficients of f(z_{1}, …, z_{m}).