[Math] Recovering the probability mass function from the characteristic function of a discrete probability distribution using Mathematica

fourier analysismathematicaprobability

I would like to recover the probability mass function (pmf) from the characteristic function (CF) of a discrete probability distribution using Mathematica.

Ideally, I'd like to do calculations like this example (not necessarily this simple). To compute the pmf fX+X[x] of the sum of two iid discrete uniform distributions X with support S={1,…,6}, it seems reasonable to try something like this:

Subscript[\[CurlyPhi], X + X][t_] := 
 CharacteristicFunction[DiscreteUniformDistribution[{1, 6}], t]^2

Subscript[f, X + X][x_] := 
 InverseFourierTransform[Subscript[\[CurlyPhi], X + X][t], t, x]

but InverseFourierTransform doesn't recover the pmf from the CF of discrete distributions (it does recover pdfs from the CFs of continuous distributions). Am I using the wrong function? (There are other candidates, but InverseFourier only works on lists of numbers, and InverseZTransform doesn't seem to work here either.) Am I forgetting to set some necessary options to the inverse function? Or is there just no built-in to recover the pmf from the CF of a discrete probability distribution?

Best Answer

Here's one way to solve the problem. The trick is to substitute -I*Log[x] for t in the characteristic function $\varphi[t]$. This turns it into a generating function. Then extract the coefficient of $x^n$ from the resulting series using the SeriesCoefficient built-in. The coefficient of $x^n$ is the value of the pmf f at n. Here's a simple example:

        (* Characteristic function (CF) of the discrete probability
           distribution for summing two 6-sided dice *)
In[1]:= \[CurlyPhi][t_] :=
          CharacteristicFunction[DiscreteUniformDistribution[{1, 6}], t]^2

        (* Recover the probability mass function f from this CF *)
In[2]:= f[n_] := SeriesCoefficient[\[CurlyPhi][-I*Log[x]], {x, 0, n}]

        (* Show the range of the pmf f over its domain, {2,..,12} *)
In[3]:= f[#] & /@ Range[2, 12]
Out[3]= {1/36, 1/18, 1/12, 1/9, 5/36, 1/6, 5/36, 1/9, 1/12, 1/18, 1/36}

        (* Note that applying f to a number outside its domain (like 1)
           produces 0 *)
In[4]:= f[1]
Out[4]= 0

        (* Now use this technique to show the range of the pmf of summing three
           6-sided dice *)
In[5]:= Map[
          Function[n, 
            SeriesCoefficient[
              CharacteristicFunction[
                DiscreteUniformDistribution[{1, 6}], -I*Log[x]]^3, {x, 0, n}]], 
                  Range[3, 18]]
Out[5]= {1/216, 1/72, 1/36, 5/108, 5/72, 7/72, 25/216, 1/8, 1/8, 25/216, 7/72,
          5/72, 5/108, 1/36, 1/72, 1/216}

        (* And verify that the values in the range of this pmf sum to 1 *)
In[6]:= Fold[Plus, 0, %]
Out[6]= 1
Related Question