6
$\begingroup$

I have a probability mass function that should sum to unity algebraically, but neither FullSimplify nor FunctionExpand can deal with the Hypergeometric2F1 therein.

myPMF[x_, p_, k_] := (   p^k *(1 - p)^(x - k) 
                      + (1 - p)^k * p^(x - k)  ) Binomial[ x - 1, k - 1] // Simplify    
(*  coin toss that ends when either k heads or k  tails appears ;
    x = total number of flips;
    k <= x <= 2 k  
*)

Assuming[ 0 < p < 1, Element[x | n, Integers]
   Sum[ myPMF[x, p, n], {x, n, 2 n - 1}] // FullSimplify 
 ]    (* expecting exact 1 like the numerical verification below ;
        FunctionExpand basically does nothing to HyperGeometric *)

start = 16; stop = 40; (* this takes maybe 6 seconds *)
Table[
       Sum[ myPMF[x, p, n], {x, n, 2 n - 1} ]
 , {n, start, stop} ] // Simplify

The exact sum can be done with partial numerical inputs as the code shows.

How do I carry out the algebraic sum and arrive at unity? Perhaps a different way of using FunctionExpand?

Thank you.

Improve this question
$\endgroup$
3
  • $\begingroup$ Your function can be simplified to 2 - (BetaRegularized[1 - p, n, n] + BetaRegularized[p, n, n]). $\endgroup$
    – J. M.'s missing motivation ♦
    Commented Oct 9, 2016 at 5:17
  • $\begingroup$ @J.M. This is interesting. Thanks. Are you suggesting that I define my probability mass function with BetaRegularized to begin with, or that I can manipulate the current form with built-in commands into the BetaRegularized version? If the latter, how would I go about doing that? $\endgroup$
    – Lee David Chung Lin
    Commented Oct 9, 2016 at 5:58
  • $\begingroup$ Nono, your PMF is as it should be, but its form was the tip-off that the summation would involve the regularized incomplete beta function at some point. Admittedly, Mathematica had to be helped a bit here; I'll write a detailed answer later, but if you want to try things out yourself in the meantime, the Pfaff transformation is particularly useful. $\endgroup$
    – J. M.'s missing motivation ♦
    Commented Oct 9, 2016 at 6:06

1 Answer 1

Reset to default
6
$\begingroup$

Here is the "no thinking needed" way to demonstrate this:

FunctionExpand[DifferenceRootReduce[
               Sum[(p^n (1 - p)^(x - n) + (1 - p)^n p^(x - n)) Binomial[x - 1, n - 1],
                   {x, n, 2 n - 1}] // FunctionExpand, n]]
   1

If, like me, you felt that that was unsatisfactory, here is a way that requires one to assist Mathematica in a few places, with the guidance of knowing a few hypergeometric identities.

For reference, this is the expression we get from directly evaluating the sum:

expr = Sum[(p^n (1 - p)^(x - n) + (1 - p)^n p^(x - n)) Binomial[x - 1, n - 1],
           {x, n, 2 n - 1}] // Simplify
   2 - ((1 - p) p)^n Binomial[2 n - 1, n - 1]
       (Hypergeometric2F1[1, 2 n, 1 + n, 1 - p] + Hypergeometric2F1[1, 2 n, 1 + n, p])

As I had noted in the comments, the Pfaff transformation is one of the first things to try; making this replacement gives

expr /. Hypergeometric2F1[a_, b_, c_, z_] :>
        (1 - z)^(c - a - b) Hypergeometric2F1[c - a, c - b, c, z] // FullSimplify
   2 - n (Beta[1 - p, n, n] + Beta[p, n, n]) Binomial[2 n - 1, n - 1]

where we now encounter the incomplete beta function. (Alternatively, one could use this identity.) A further simplification can be done using the identity

n Binomial[2 n - 1, n - 1] == 1/Beta[n, n] // FullSimplify

which yields

2 - (Beta[1 - p, n, n] + Beta[p, n, n])/Beta[n, n]

(or alternatively, 2 - (BetaRegularized[1 - p, n, n] + BetaRegularized[p, n, n])).

Finally, one can use the reflection formula for the incomplete beta function to get the desired result.

Improve this answer
$\endgroup$

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge you have read our privacy policy.

Not the answer you're looking for? Browse other questions tagged or ask your own question.