1
$\begingroup$

Assuming $0<q<1$. I built these two functions

A[k_,q_]:=Sum[PDF[BinomialDistribution[i, q], k]*PDF[ZipfDistribution[n], i], {i, k, Infinity}]

and

B[k_,q_]:=Piecewise[{
{(1/Zeta[n + 1])*PolyLog[n + 1, 1 - q], k == 0},
{(1/(Zeta[n + 1]*k!))*(q/(1 - q))^k*Sum[CoefficientList[Product[x - j, {j, 1, k - 1}],x][[-i]]* PolyLog[n + i - k, 1 - q], {i, 1, k}], k > 0}}]

Functions A and B are mathematically equivalent for $k\geq1$. (Function B is faster than function A.)

When I compute $A[25,0.6]$ and $B[25,0.6]$ with $n=3$, I start to get a little difference between the results $(5.61124*10^{-7},5.61132*10^{-7})$. So, the question is, why is this difference? Does anyone know how to avoid this round-off problem?

Improve this question
$\endgroup$
5
  • 2
    $\begingroup$ Mor me: B[25, 0.6] give me: 11.09 not: 5.61124*10^-7 ? $\endgroup$
    – Mariusz Iwaniuk
    Commented Aug 14, 2020 at 9:50
  • $\begingroup$ @MariuszIwaniuk It looks like there are some characters missing in the definition, note the big empty spaces in e.g. n 1. $\endgroup$
    – Lukas Lang
    Commented Aug 14, 2020 at 11:32
  • $\begingroup$ It's normal for different forms of the same mathematical expression to give different numerical approximations, especially when including things like infinite sums. For your particular case, you can use exact numbers: A[25, 6/10] == B[25, 6/10] // FullSimplify and N[A[25, 6/10], 20]. In general, see the many questions on this site about increasing precision of the various Mathematica functions, such as Sum $\endgroup$
    – Lukas Lang
    Commented Aug 14, 2020 at 11:37
  • $\begingroup$ @MariuszIwaniuk I'm sorry, you are right. I just edited the post with the correction. $\endgroup$
    – Slash020
    Commented Aug 14, 2020 at 13:54
  • 1
    $\begingroup$ Exact solvers (like Sum) sometimes have difficulty with inexact input (e.g. 0.6 instead of 6/10). The number 5.61124*10^-7 is the result for A[] returned by NSum. Possibly Sum failed and called NSum as a backup. $\endgroup$
    – Michael E2
    Commented Aug 14, 2020 at 14:11

1 Answer 1

Reset to default
1
$\begingroup$

This calculates the symbolic infinite sum A[] for a given integer n and caches the result. The result is used to compute subsequent values of A[k, q, n].

ClearAll[A, B];
A[k_, q_, n_Integer?Positive] := (
  DownValues[A] = Prepend[
    DownValues[A],
    HoldPattern[A[kk_, qq_, n]] :> Evaluate@
      Sum[
       PDF[BinomialDistribution[i, qq], kk]*
        PDF[ZipfDistribution[n], i],
       {i, kk, Infinity}]];
  A[k, q, n])

B[k_, q_, n_Integer?Positive] := 
  Piecewise[{{(1/Zeta[n + 1])*PolyLog[n + 1, 1 - q], 
     k == 0}, {(1/(Zeta[n + 1]*k!))*(q/(1 - q))^k*
      Sum[CoefficientList[Product[x - j, {j, 1, k - 1}], x][[-i]]*
        PolyLog[n + i - k, 1 - q], {i, 1, k}], k > 0}}];

Difference on floating-point input:

A[25, 0.6, 3] - B[25, 0.6, 3]

(*  3.27923*10^-12  *)

The exact diffference is zero (it takes a few seconds to evaluate the HypergeometricPFQ that appear in A[k, q, 3] on exact input):

Simplify[A[25, 6/10, 3] - B[25, 6/10, 3]]  

(* 0  *)

N[A[25, 6/10, 3] - B[25, 6/10, 3], 20]

N::meprec: Internal precision limit $MaxExtraPrecision = 50.` reached while evaluating....

(*  0.*10^-66  *)
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.