Is this implementation detail for solving maximum entropy on a computer correct? I am currently looking at a paper by Mattos and Veiga, who describe an approach to solving the maximum entropy problem subject to linear constraints:
$$\begin{aligned}
\max_{p_i} -\sum_{i=1}^N p_i \log p_i \quad \text{s.t.} & \sum_i p_i = 1\\
& \sum_i p_i g_j(x_i) = m_j,\ j = 1,\ldots M\\
& p_i \geq 0
\end{aligned}$$
(i.e. to choose the probability distribution with the maximum entropy subject to constraints on $M$ moments of the distribution). The authors make the claim (with a slight change in notation from the paper) that "Using the method of Lagrange multiplier, the MaxEnt problem (4) can be written in the
following unrestricted form" (emphasis mine)
$$\max_{p_i,\lambda_0,\lambda_j} -\sum_{i=1}^N p_i \log p_i + (\lambda_0 - 1)\left(\sum_{i=1}^N p_i - 1\right) + \sum_{j=1}^M \lambda_j\left(\sum_{i=1}^N p_i g_j(x_i) - m_j\right)$$
I am guessing that this claim is motivated by the fact that maximizing the Lagrangian has the same first order condition as the original problem, but it seems to me like the literal  equivalence they are claiming is false. Specifically, fix $p_i, \lambda_j$ such that the first and last terms are finite, but $\sum p_i > 1$. Then you can make the objective arbitrarily large by letting $\lambda_0 \to \infty$. They then go on to derive a supposed solution to the second maximization problem by first using the first order conditions to show that for fixed $\lambda_j$, we have that $p_i$ must be:
$$p_i = \frac{\exp\left(-\sum_{j=1}^M \lambda_j g_j(x_i)\right)}{\sum_{i=1}^N \exp\left(-\sum_{j=1}^M \lambda_j g_j(x_i)\right)}$$
which they plug back into the second maximization problem to obtain an unconstrained, concave problem in terms of the $\lambda_j$ only.
Given that the first step seems suspect, what am I to make of their solution? Does it actually solve the maximum entropy problem? If so, what is a rigorous proof of this fact? If not, what does their proposal end up solving? Is there some sense in which even though their algorithm does not literally solve the maximum entropy problem, they still produce something close?
 A: The solution derived in the paper turns out to give the right answer to the maximum entropy problem, but potentially for the wrong reason. Specifically, after plugging in the $p_i$ as a function of the $\lambda j$ into the objective (and imposing that both constraints are satisfied), they get that we should instead be maximizing
$$\log\left(\sum_{i=1}^N \exp\left(-\sum_{j=1}^M \lambda_j g_j(x_i)\right)\right) + \sum_{j=1}^M \lambda_j m_j$$
which they proceed to do via Newton's method. To see that this ends up giving the correct solution, notice that the above objective is strictly concave, and the first order conditions of the above maximization problem are exactly the moment conditions of the original problem:
$$\sum_{i=1}^N p_i(\lambda) g_j(x_i) = m_j$$
Therefore, the Newton's method used to solve the above optimization problem is exactly equivalent to finding the values of $\lambda$ that, when plugged into the form that the $p_i$'s must take, give the desired moments, which in turn give the solution to the original maximum entropy problem.
