I tried

Expectation[-q*Log[q], q \[Distributed] BetaDistribution[a, b]]

and got

(a (HarmonicNumber[a] - HarmonicNumber[a + b]))/(a + b)

(I get the same result using integration.) This seems pretty different from the expression in Wikipedia, even after substituting the identity (mentioned at MathWorld that H[n] = gamma + PolyGamma[0,n+1], where gamma is Euler's constant and PolyGamma[0, n+1] = Digamma[n+1].

Any ideas which identity I'm missing?


P.S. I trust the Wikipedia answer, I have seen this result before, but I want to re-derive it with Mathematica as a sanity check, before I go on to derive other expressions of interest.

  • 2
    $\begingroup$ You could plot the Wikipedia result and yours for a number of parameter values to see whether they might be the same after all. $\endgroup$ – Sjoerd C. de Vries Nov 6 '13 at 7:02

Your definition of entropy is incorrect. It's $E(-\ln(P(x)))$, with $E$ the expectation operator and $P$ the probability mass function of the random variable $x$. I believe you may have been mixing up a few things.

The formal definition of the expectation is $E(x)=\int{x P(x)dx}$. I assume that you had this in mind and you further confused your random variable $q$ with its PDF to write

Expectation[-q*Log[q], q \[Distributed] BetaDistribution[a, b]] 

for the entropy.

Instead this should have been:

Expectation[-Log[PDF[BetaDistribution[a, b], q]], q \[Distributed] BetaDistribution[a, b]]

Unfortunately, I couldn't prod Mathematica to show the equality between Wikipedia's result and this formal definition symbolically. However, numerically it works out quite well.

This is Wikipedia's version:

Mathematica graphics

In Mathematica terms:

Ψ[n_] := HarmonicNumber[n - 1] - EulerGamma

entropyBeta[a_, b_] = 
  FullSimplify[Log[Beta[a, b]] - (a - 1) Ψ[a] - (b - 1) Ψ[b] + (a + b - 2) Ψ[a + b],     
    Assumptions -> {a > 0, b > 0}];


  Plot[entropyBeta[a, b], {b, 0.01, 5}],
    NExpectation[-Log[PDF[BetaDistribution[a, b], q]], q \[Distributed] BetaDistribution[a, b]], 
   {b, 0.1, 5, .1}]
 {a, 0.01, 5}]

Mathematica graphics

As you can see, they match very well.

  • $\begingroup$ One spurious bracket removed. $\endgroup$ – Sjoerd C. de Vries Oct 23 '19 at 19:34

dist = BetaDistribution[a, b];

The constraints on the parameters are

assume = DistributionParameterAssumptions[dist];

The integral needs a little help:

(entropy[BetaDistribution[a, b]] = Assuming[assume && 0 < x < 1,
       -Log[PDF[dist, x]]*PDF[dist, x] // Simplify //
        FunctionExpand, {x, 0, 1}] // FullSimplify] /.
    Gamma[x_ + y_] :> Gamma[x]*Gamma[y]/Beta[x, y]) // TraditionalForm

enter image description here


Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

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