# Use Mathematica to derive the probability distributions using the Principle of Maximum Entropy

I am reading the article of Deriving probability distributions using the Principle of Maximum Entropy

and I am trying to derive some of the equations in it automatically using Mathematica.

# 1. Derivation of maximum entropy probability distribution with no other constraints (uniform distribution)

First, we solve for the case where the only constraint is that the distribution is a pdf, which we will see is the uniform distribution. To maximize entropy, we want to minimize the following function: $$J(p)=\int_a^b p(x) \ln p(x) d x-\lambda_0\left(\int_a^b p(x) d x-1\right)$$ . Taking the derivative with respect ot $$p(x)$$ and setting to zero, $$\frac{\delta J}{\delta p(x)}=1+\ln p(x)-\lambda_0=0$$

Can I derive the second equation automatically using Mathematica?

if just $$J(p)=\int_a^b p(x) \ln p(x) dx$$, I can do

Needs["VariationalMethods"]
VariationalD[p[x] Log[p[x]], p[x], x]


the output is 1 + Log[p[x]] as expected.

But here $$J(p)=\int_a^b p(x) \ln p(x) d x-\lambda_0\left(\int_a^b p(x) d x-1\right)$$, how can I do that?

• Why not VariationalD[p[x] Log[p[x]] - λ0 ( p[x] - 1), p[x], x]? What's the issue?
– bmf
Jan 12 at 11:40
• @user64494 that's weird. on my laptop it returns 1 - λ0 + Log[p[x]]` which is the expected answer from the OP
– bmf
Jan 12 at 11:55
• @bmf's code works as expected also in 12.3.1 (Windows) and 13.2.0 (Wolfram Cloud). Jan 12 at 11:56
• @benjaminchanming I mean that you wrote that if $J$ has the simple form you know how to do the variational derivative in Mathematica. For the complicated form I am unsure why you got stuck, but I showed you that you can write the same command effectively and get the result that you want for $\tfrac{\delta J}{\delta \rho}$ and now you say you are unsure about the result, which is the one you wanted. Unless, I am missing something
– bmf
Jan 12 at 12:18
• @benjaminchanming, linearity is one of the basic properties of integral, and the equality holds (if both of the integrals exist). However, for ill-behaved (!!) functions, the equality might not hold, and that is why Mathematica returns False. Jan 12 at 12:42