Maximum-entropy module

In maxent.f90, all subroutines and functions are contained within the maxent module. The basis functions are computed in function phimaxent(), where Eq. (4) is used. To this end, the Lagrange multipliers are determined using the convex optimization algorithms that are presented in Section 2. Most of the functions and subroutines in this module are self-explanatory. We elaborate on the computation of the gradient of $F({\lambda}) := \ln Z({\lambda})$, and also on the Hessian ${H}$, which is required in the Newton method. On taking the gradient (with respect to ${\lambda}$) of $F({\lambda})$, we obtain the negative of the left-hand side of the linear reproducing conditions in Eq. (3c):

$${g} = - \sum_{i=1}^n \phi_i \tilde x^i,$$ (6)

where $\tilde x^i = x^i - x$ ( $i = 1, 2, \ldots, n$). This is coded in function dfunc(). The expressions for the Hessian of $F({\lambda})$ are given in Reference [11] (matrix form) as well as in the Appendix of Reference [12]. By definition, the components of the Hessian are:

$$H_{rs}({\lambda}) = \dfrac{\partial^2 F({\lambda})}{\partial \lambda_r \partial \lambda_s}, \quad (r,s = 1,\ldots,d).$$ (7)

The Hessian matrix (see Eq. (45) in Reference [11]) in three dimensions ($d = 3$) is:

$$\hspace*{-0.13in} {H} = \begin{bmatrix}\! \! {<\!\tilde x_1^2\!> ... ...! \! {<\!\tilde x_3^2\!> } \! - \! {<\!\tilde x_3\!> }^2 \! \! \! \end{bmatrix}$$ (8)

where ${<\!\cdot\!> }$ is the expectation operator, which for a scalar-valued function $u({x})$ is:

$${<\!u({x})\!> } = \sum_{i=1}^n \phi_i({x}) u_i , \quad u_i = u({x}^i).$$ (9)

The Hessian matrix and its inverse are computed in function hessian(flag) (flag is an optional Boolean argument) and function invhessian(flag), respectively.

Let ${\lambda} = {\lambda}^*$ denote the converged solution for the Lagrange multipliers and $\phi_i^*$ the corresponding basis function solution for the $i$th node. Since ${<\!\tilde x_r\!> }^* = 0$ ($r = 1$-$3$), the Hessian (hessian(.true.) is the call) is

$${H}^* = \begin{bmatrix}{<\!\tilde x_1^2\!> }^* & {<\!\tilde x_1 \... ...trix} , \quad {H}^* = \sum_{k=1}^n \phi_k^* \tilde {x}^k \otimes \tilde {x}^k .$$ (10)

The expressions for the derivatives of the basis functions are given below. We adopt the notations and approach presented in Arroyo and Ortiz [12]; in the interest of space, just the final results are indicated. We can write Eq. (4) as

$$\phi_i^*({x};{\lambda}^*) = \dfrac{\exp [{f_i^*({x};{\lambda}^*)}... ..., \quad f_i^*({x};{\lambda}^*) = \ln w_i({x}) - {\lambda^* \cdot \tilde {x}^i},$$ (11)

where ${\lambda}^*$ is implicitly dependent on ${x}$. On using Eq. (11), we have

and therefore the gradient of $\phi_i^*$ is

$${\nabla} \phi_i^* = \phi_i^* \left\{ \tilde {x}^i \cdot \left[ ({... ...c{{\nabla}w_i}{w_i} - \sum_{j=1}^n \phi_j^* \dfrac{{\nabla}w_j}{w_j} \right\} .$$ (13)

If the prior $w_i({x})$ is a Gaussian radial basis function (see Reference [13]), then $w_i({x})=\exp(-\beta \vert{x}^i - {x}\vert^2)$ and Eq. (13) reduces to ${\nabla}\phi_i^* = \phi_i^* ({H}^*)^{-1} \cdot \tilde {x}^i$! This result appears in the Appendix of Reference [12]. In general, ${\nabla} \phi_i^*$ depends on ${\nabla}w_i$ through the expression given in Eq. (13). The gradient of the basis functions is computed in function dphimaxent().

On taking the gradient of Eq. (13), we obtain the following expression for the second derivatives of the max-ent basis functions:

$$\begin{split}{\nabla} {\nabla} \phi_i^* = & \dfrac{{\nabla}... ... \left( \dfrac{{\nabla}w_j}{w_j} \right) \right\} , \end{split}$$ (14)

where on using the identity ${H}^* \cdot ({H}^*)^{-1} = {I}$, the gradient of the inverse of the Hessian can be written as

$${\nabla} \bigl( ({H}^*)^{-1} \bigr) = - ({H}^*)^{-1} \cdot {\nabl... ...}^* = \sum_{k=1}^n \tilde {x}^k \otimes \tilde {x}^k \otimes {\nabla} \phi_k^*,$$ (15)

and therefore

$$\tilde {x}^i \cdot {\nabla} \bigl( ({H}^*)^{-1} \bigr) = - \tilde... ...{\nabla} \phi_k^* \right), \quad \tilde v^k = \tilde {x}^k \cdot ({H}^*)^{-1} .$$ (16)

The term $- \tilde {x}^i \cdot {\nabla} \bigl( ({H}^*)^{-1} \bigr) \cdot {A}^*$ in Eq. (14) is given by

$$- \tilde {x}^i \cdot {\nabla} \bigl( ({H}^*)^{-1} \bigr) \cdot {A... ... \otimes {\nabla} \phi_k^* \right), \quad \tilde w^k = \tilde v^k \cdot {A}^* ,$$ (17)

and the term $- \tilde {x}^i \cdot ({H}^*)^{-1} \cdot {\nabla} {A}^*$ is:

The Hessian of the basis functions is computed in function ddphimaxent(). It satisfies the conditions $\sum_i {\nabla} {\nabla} \phi_i^* = {0}$ and $\sum_i {\nabla} {\nabla} \phi_i^* \otimes \tilde {x}^i = {0}$.