Entropy: The Markov Ordering Approach

The focus of this article is on entropy and Markov processes. We study the properties of functionals which are invariant with respect to monotonic transformations and analyze two invariant"additivity"properties: (i) existence of a monotonic transformation which makes the functional additive with respect to the joining of independent systems and (ii) existence of a monotonic transformation which makes the functional additive with respect to the partitioning of the space of states. All Lyapunov functionals for Markov chains which have properties (i) and (ii) are derived. We describe the most general ordering of the distribution space, with respect to which all continuous-time Markov processes are monotonic (the {\em Markov order}). The solution differs significantly from the ordering given by the inequality of entropy growth. For inference, this approach results in a convex compact set of conditionally"most random"distributions.


A Bit of History: Classical Entropy
Two functions, energy and entropy, rule the Universe. In 1865 R. Clausius formulated two main laws [1]: 1. The energy of the Universe is constant.
2. The entropy of the Universe tends to a maximum.
The universe is isolated. For non-isolated systems energy and entropy can enter and leave, the change in energy is equal to its income minus its outcome, and the change in entropy is equal to entropy production inside the system plus its income minus outcome. The entropy production is always positive.
Entropy was born as a daughter of energy. If a body gets heat ∆Q at the temperature T then for this body dS = ∆Q/T . The total entropy is the sum of entropies of all bodies. Heat goes from hot to cold bodies, and the total change of entropy is always positive.
Ten years later J.W. Gibbs [2] developed a general theory of equilibrium of complex media based on the entropy maximum: the equilibrium is the point of the conditional entropy maximum under given values of conserved quantities. The entropy maximum principle was applied to many physical and chemical problems. At the same time J.W. Gibbs mentioned that entropy maximizers under a given energy are energy minimizers under a given entropy.
The classical expression p ln p became famous in 1872 when L. Boltzmann proved his H-theorem [3]: the function decreases in time for isolated gas which satisfies the Boltzmann equation (here f (x, v) is the distribution density of particles in phase space, x is the position of a particle, v is velocity). The statistical entropy was born: S = −kH. This was the one-particle entropy of a many-particle system (gas).
In 1902, J.W. Gibbs published a book "Elementary principles in statistical dynamics" [4]. He considered ensembles in the many-particle phase space with probability density ρ(p 1 , q 1 , . . . p n , q n ), where p i , q i are the momentum and coordinate of the ith particle. For this distribution, S = −k ρ(p 1 , q 1 , . . . p n , q n ) ln(ρ(p 1 , q 1 , . . . p n , q n ))dq 1 . . . dq n dp 1 . . . dp n Gibbs introduced the canonical distribution that provides the entropy maximum for a given expectation of energy and gave rise to the entropy maximum principle (MaxEnt). The Boltzmann period of history was carefully studied [5]. The difference between the Boltzmann entropy which is defined for coarse-grained distribution and increases in time due to gas dynamics, and the Gibbs entropy, which is constant due to dynamics, was analyzed by many authors [6,7]. Recently, the idea of two functions, energy and entropy which rule the Universe was implemented as a basis of two-generator formalism of nonequilibrium thermodynamics [8,9].
In information theory, R.V.L. Hartley (1928) [10] introduced a logarithmic measure of information in electronic communication in order "to eliminate the psychological factors involved and to establish a measure of information in terms of purely physical quantities". He defined information in a text of length n in alphabet of s symbols as H = n log s.
In 1948, C.E. Shannon [11] generalized the Hartley approach and developed "a mathematical theory of communication", that is information theory. He measured information, choice and uncertainty by the entropy: Here, p i are the probabilities of a full set of n events ( n i=1 p i = 1). The quantity S is used to measure of how much "choice" is involved in the selection of the event or of how uncertain we are of the outcome. Shannon mentioned that this quantity form will be recognized as that of entropy, as defined in certain formulations of statistical mechanics. The classical entropy (1), (2) was called the Boltzmann-Gibbs-Shannon entropy (BGS entropy). (In 1948, Shannon used the concave function (2), but under the same notation H as for the Boltzmann convex function. Here we use H for the convex H-function, and S for the concave entropy.) In 1951, S. Kullback and R.A. Leibler [12] supplemented the BGS entropy by the relative BGS entropy, or the Kullback-Leibler divergence between the current distribution P and some "base" (or "reference") distribution Q: The Kullback-Leibler divergence is always non-negative D KL (P Q) ≥ 0 (the Gibbs inequality). It is not widely known that this "distance" has a very clear physical interpretation. This function has been well known in physical thermodynamics since 19th century under different name. If Q is an equilibrium distribution at the same temperature as P has, then where F is free energy and T is thermodynamic temperature. In physics, F = U − T S, where physical entropy S includes an additional multiplier k, the Boltzmann constant. The thermodynamic potential −F/T has its own name, Massieu function. Let us demonstrate this interpretation of the Kullback-Leibler divergence. The equilibrium distribution Q provides the conditional entropy (2) maximum under a given expectation of energy i u i q i = U and the normalization condition i q i = 1.
After the classical work of Zeldovich (1938Zeldovich ( , reprinted in 1996), the expression for free energy in the "Kullback-Leibler form" where c i is concentration and c * i (T ) is the equilibrium concentration of the ith component, is recognized as a useful instrument for the analysis of kinetic equations (especially in chemical kinetics [14,15]).
Each given positive distribution Q could be represented as an equilibrium Boltzmann distribution for given T > 0 if we take u i = −kT log q i + u 0 for an arbitrary constant level u 0 . If we change the order of arguments in the Kullback-Leibler divergence then we get the relative Burg entropy [16,17]. It has a much more exotic physical interpretation: for a current distribution P we can define the "auxiliary energy" functional U P for which P is the equilibrium distribution under a given temperature T . We can calculate the auxiliary free energy of any distribution Q and this auxiliary energy functional: F P (Q). (Up to an additive constant, for P = P * this F P (Q) turns into the classical free energy, F * P (Q) = F (Q).) In particular, we can calculate the auxiliary free energy of the physical equilibrium, F P (P * ). The relative Burg entropy is D KL (P * P ) = F P (P * ) − F P (P ) kT This functional should also decrease in any Markov process with given equilibrium P * .
Information theory developed by Shannon and his successors focused on entropy as a measure of uncertainty of subjective choice. This understanding of entropy was returned from information theory to statistical mechanics by E.T. Jaynes as a basis of "subjective" statistical mechanics [18,19]. He followed Wigner's idea "entropy is an antropocentric concept". The entropy maximum approach was declared as a minimization of the subjective uncertainty. This approach gave rise to a MaxEnt "anarchism". It is based on a methodological hypothesis that everything unknown could be estimated by the principle of the entropy maximum under the condition of fixed known quantities. At this point the classicism in entropy development changed to a sort of scientific modernism. The art of model fitting based on entropy maximization was developed [20]. The principle of the entropy maximum was applied to plenty of problems: from many physical problems [21], chemical kinetics and process engineering [15] to econometrics [22,23] and psychology [24]. Many new entropies were invented and now one has rich choice of entropies for fitting needs [25]. The most celebrated of them are the Rényi entropy [26], the Burg entropy [16,17], the Tsallis entropy [27,28] and the Cressie-Read family [29,30]. The nonlinear generalized averaging operations and generalized entropy maximization procedures were suggested [31].
Following this impressive stream of works we understand the MaxEnt approach as conditional maximization of entropy for the evaluation of the probability distribution when our information is partial and incomplete. The entropy function may be the classical BGS entropy or any function from the rich family of non-classical entropies. This rich choice causes a new problem: which entropy is better for a given class of applications?
The MaxEnt "anarchism" was criticized many times as a "senseless fitting". Arguments pro and contra the MaxEnt approach with non-classical entropies (mostly the Tsallis entropy [27]) were collected by Cho [32]. This sometimes "messy and confusing situation regarding entropy-related studies has provided opportunities for us: clearly there are still many very interesting studies to pursue" [33].

Key Points
In this paper we do not pretend to invent new entropies. (There appear new functions as limiting cases of the known entropy families, but this is not our main goal). Entropy is understood in this paper as a measure of uncertainty which increases in Markov processes. In our paper we consider a Markov process as a semigroup on the space of positive probability distributions. The state space is finite. Generalizations to compact state spaces are simple. We analyze existent relative entropies (divergences) using several simple ideas: 1. In Markov processes probability distributions P (t) monotonically approach equilibrium P * : divergence D(P (t) P * ) monotonically decrease in time.
2. In most applications, conditional minimizers and maximizers of entropies and divergences are used, but the values are not. This means that the system of level sets is more important than the functions' values. Hence, most of the important properties are invariant with respect to monotonic transformations of entropy scale.

3.
The system of level sets should be the same as for additive functions: after some rescaling the divergences of interest should be additive with respect to the joining of statistically independent systems. 4. The system of level sets should after some rescaling the divergences of interest should have the form of a sum (or integral) over where the function f is the same for all states. In information theory, divergences of such form are called separable, in physics the term trace-form functions is used The first requirement means that if a distribution becomes more random then it becomes closer to equilibrium (Markov process decreases the information excess over equilibrium). For example, classical entropy increases in all Markov processes with uniform equilibrium distributions. This is why we may say that the distribution with higher entropy is more random, and why we use entropy conditional maximization for the evaluation of the probability distribution when our information is partial and incomplete.
It is worth to mention that some of the popular Bregman divergences, for example, the squared Euclidean distance or the Itakura-Saito divergence, do not satisfy the first requirement (see Section 4.3).
The second idea is just a very general methodological thesis: to evaluate an instrument (a divergence) we have to look how it works (produces conditional minimizers and maximizers). The properties of the instrument which are not related to its work are not important. The requirements number three and four are, in their essence, conditions for separation of variables in conditional minimization (maximization) problems.
The number three allows to separate variables if the system consists of independent subsystems, the number four relates to separation of variables for partitions of the space of probability distributions.
Amongst a rich world of relative entropies and divergences, only two families meet these requirements. Both were proposed in 1984. The Cressie-Read (CR) family [29,30]: and the convex combination of the Burg and Shannon relative entropies proposed in [34] and further analyzed in [35,36]: When λ → 0 the CR divergence tends to the KL divergence (the relative Shannon entropy) and when λ → −1 it turns into the Burg relative entropy. The Tsallis entropy was introduced four years later [27] and became very popular in thermostatistics (there are thousands of works that use or analyze this entropy [37]). The Tsallis entropy coincides (up to a constant multiplier λ + 1) with the CR entropy for λ ∈] − 1, ∞[ and there is no need to study it separately (see discussion in Section 2.2).
A new problem arose: which entropy is better for a specific problem? Many authors compare performance of different entropies and metrics for various problems (see, for example, [39,40]). In any case study it may be possible to choose "the best" entropy but in general we have no sufficient reasons for such a choice. We propose a possible way to avoid the choice of the best entropy.
Let us return to the idea: the distribution Q is more random than P if there exists a continuous-time Markov process (with given equilibrium distribution P * ) that transforms P into Q. We say in this case that P and Q are connected by the Markov preorder with equilibrium P * and use notation P ≻ 0 P * Q. The Markov order ≻ P * is the transitive closure of the Markov preorder.
If a priori information gives us a set of possible distributions W then the conditionally "maximally random distributions" (the "distributions without additional information", the "most indefinite distributions" in W ) should be the minimal elements in W with respect to Markov order. If a Markov process (with equilibrium P * ) starts at such a minimal element P then it cannot produce any other distribution from W because all distributions which are more random that P are situated outside W . In this approach, the maximally random distributions under given a priori information may be not unique. Such distributions form a set which plays the same role as the standard MaxEnt distribution. For the moment based a priori information the set W is an intersection of a linear manifold with the simplex of probability distributions, the set of minimal elements in W is also polyhedron and its description is available in explicit form. In low-dimensional case it is much simpler to construct this polyhedron than to find the MaxEnt distributions for most of specific entropies.

Structure of the Paper
The paper is organized as follows. In Section 2 we describe the known non-classical divergences (relative entropies) which are the Lyapunov functions for the Markov processes. We discuss the general construction and the most popular families of such functions. We pay special attention to the situations, when different divergences define the same order on distributions and provide the same solutions of the MaxEnt problems (Section 2.2). In two short technical Sections 2.3 and 2.4 we present normalization and symmetrization of divergences (similar discussion of these operations was published very recently [38]. The divergence between the current distribution and equilibrium should decrease due to Markov processes. Moreover, divergence between any two distributions should also decrease (the generalized data processing Lemma, Section 3).
Definition of entropy by its properties is discussed in Section 4. Various approaches to this definition were developed for the BGS entropy by Shannon [11], [41] and by other authors for the Rényi entropy [43,44], the Tsallis entropy [42], the CR entropy and the convex combination of the BGS and Burg entropies [46]. Csiszár [45] axiomatically characterized the class of Csiszár-Morimoto divergences (formula (6) below). Another characterization of this class was proved in [46] (see Lemma 1, Section 4.3 below).
From the celebrated properties of entropy [47] we selected the following three: 1. Entropy should be a Lyapunov function for continuous-time Markov processes; 2. Entropy is additive with respect to the joining of independent systems; 3. Entropy is additive with respect to the partitioning of the space of states (i.e., has the trace-form).
To solve the MaxEnt problem we have to find the maximizers of entropy (minimizers of the relative entropy) under given conditions. For this purpose, we have to know the sublevel sets of entropy, but not its values. We consider entropies with the same system of sublevel sets as equivalent ones (Section 2.2). From this point of view, all important properties have to be invariant with respect to monotonic transformations of the entropy scale. Two last properties from the list have to be substituted by the following: 2'. There exists a monotonic transformation which makes entropy additive with respect to the joining of independent systems (Section 4.1); 3'. There exists a monotonic transformation which makes entropy additive with respect to the partitioning of the space of states (Section 4.2).
These properties imply specific separation of variables for the entropy maximization problems (Sections 4.1 and 4.2). Several "No More Entropies" Theorems are proven in Section 4.3: if an entropy has properties 1, 2' and 3' then it belongs to one of the following one-parametric families: to the Cressie-Read family, or to a convex combination of the classical BGS entropy and the Burg entropy (may be, after a monotonic transformation of scale). It seems very natural to consider divergences as orders on distribution spaces (Section 5.1), the sublevel sets are the lower cones of this orders. For several functions, H 1 (P ), . . . , H k (P ) the sets {Q | H i (P ) > H i (Q) for all i} give the simple generalization of the sublevel sets. In Section 5 we discuss the more general orders in which continuous time Markov processes are monotone, define the Markov order and fully characterize the local Markov order. The Markov chains with detailed balance define the Markov order for general Markov chains (Section 5.2). It is surprising that there is no necessity to consider other Markov chains for the order characterization (Section 5.2) because no reversibility is assumed in this analysis.
In Section 6.1 we demonstrate how is it possible to use the Markov order to reduce the uncertainty in the standard settings when a priori information is given about values of some moments. Approaches to construction of the most random distributions are presented in Section 6.2.
Various approaches for the definition of the reference distribution (or the generalized canonical distribution) are compared in Section 7.
In Conclusion we briefly discuss the main results.

Csiszár-Morimoto Functions H h
During the time of modernism plenty of new entropies were proposed. Esteban and Morales [25] attempted to systemize many of them in an impressive table. Nevertheless, there are relatively few entropies in use now. Most of the relative entropies have the form proposed independently in 1963 by I. Csiszar [49] and T. Morimoto [48]: where h(x) is a convex function defined on the open (x > 0) or closed x ≥ 0 semi-axis. We use here notation H h (P P * ) to stress the dependence of H h both on p i and p * i . These relative entropies are the Lyapunov functions for all Markov chains with equilibrium P * = (p * i ). Moreover, they have the relative entropy contraction property given by the generalized data processing lemma (Section 3.2 below).
For h(x) = x log x this function coincides with the Kullback-Leibler divergence from the current distribution p i to the equilibrium p * i . Some practically important functions h have singularities at 0. For example, if we take h(x) = − log x, then the correspondent H h is the relative Burg entropy

Required Properties of the Function h(x)
Formally, h(x) is an extended real-valued proper convex function on the closed positive real half-line [0, ∞[. An extended real-valued function can take real values and infinite values ±∞. A proper function has at least one finite value. An extended real valued function on a convex set U is called convex if its epigraph is a convex set [50]. For a proper function this definition is equivalent to the Jensen inequality Not everywhere differentiable functions h(x) are often used, for example, h(x) = |x − 1|. Nevertheless, it is convenient to consider the twice differentiable on ]0, ∞[ functions h(x) and to produce a non-smooth h(x) (if necessary) as a limit of smooth convex functions. We use widely this possibility.
The Most Popular Divergences H h (P P * ) The quantity −H h is the number of non-zero probabilities p i and does not depend on P * . Sometimes it is called the Hartley entropy.
this is the l 1 -distance between P and P * .
this is the usual Kullback-Leibler divergence or the relative BGS entropy; this is the relative Burg entropy. It is obvious that this is again the Kullback-Leibler divergence, but for another order of arguments.
5. Convex combinations of h = x ln x and h = − ln x also produces a remarkable family of divergences: The convex combination of divergences becomes a symmetric functional of (P, P * ) for β = 1/2. There exists a special name for this case, "Jeffreys' entropy".
This is the quadratic term in the Taylor expansion of the relative Botzmann-Gibbs-Shannon entropy, D KL (P P * ), near equilibrium. Sometimes, this quadratic form is called the Fisher entropy.
This is the CR family of power divergences [29,30]. For this family we use notation H CR λ . If λ → 0 then H CR λ → D KL (P P * ), this is the classical BGS relative entropy; if λ → −1 then H CR λ → D KL (P * P ), this is the relative Burg entropy.
8. For the CR family in the limits λ → ±∞ only the maximal terms "survive". Exactly as we get the limit l ∞ of l p norms for p → ∞, we can use the root (λ(λ + 1)H CR λ ) 1/|λ| for λ → ±∞ and write in these limits the divergences: The existence of two limiting divergences H CR ±∞ seems very natural: there may be two types of extremely non-equilibrium states: with a high excess of current probability p i above p * i and, inversely, with an extremely small current probability p i with respect to p * i .
9. The Tsallis relative entropy [27] corresponds to the choice h For this family we use notation H Ts α .

Rényi Entropy
The Rényi entropy of order α > 0 is defined [26] as It is a concave function, and H R α (P ) → S(P ) when α → 1, where S(P ) is the Shannon entropy. When α → ∞, the Rényi entropy has a limit H ∞ (X) = − log max i=1,...n p i , which has a special name "Min-entropy".
It is easy to get the expression for a relative Rényi entropy H R α (P P * ) from the requirement that it should be a Lyapunov function for any Markov chain with equilibrium P * : For the Min-entropy, the correspondent divergence (the relative Min-entropy) is is a Lyapunov function for any Markov chain with equilibrium P * , hence, the relative Min-entropy is also the Lyapunov functional.

Entropy Level Sets
A level set of a real-valued function f of is a set of the form : where c is a constant (the "level"). It is the set where the function takes on a given constant value. A sublevel set of f is a set of the form A superlevel set of f is given by the inequality with reverse sign: The intersection of the sublevel and the superlevel sets for a given value c is the level set for this level. In many applications, we do not need the entropy values, but rather the order of these values on the line. For any two distributions P, Q we have to compare which one is closer to equilibrium P * , i.e., to answer the question: which of the following relations is true: H(P P * ) > H(Q P * ), H(P P * ) = H(Q P * ) or H(P P * ) < H(Q P * )? To solve the MaxEnt problem we have to find the maximizers of entropy (or, in more general settings, the minimizers of the relative entropy) under given conditions. For this purpose, we have to know the sublevel sets, but not the values. All the MaxEnt approach does not need the values of the entropy but the sublevel sets are necessary.
Let us consider two functions, φ and ψ on a set U. For any V ⊂ U we can study conditional minimization problems φ(x) → min and ψ(x) → min, x ∈ V . The sets of minimizers for these two problems coincide for any V ⊂ U if and only if the functions φ and ψ have the same sets of sublevel sets. It should be stressed that here just the sets of sublevel sets have to coincide without any relation to values of level.
Let us compare the level sets for the Rényi, the Cressie-Read and the Tsallis families of divergences (for α − 1 = λ and for all values of α). The values of these functions are different, but the level sets are the same (outside the Burg limit, where α → 0): for α = 0, 1 For α → 1 all these divergences turn into the Shannon relative entropy. Hence, if α = 0 then for any P , P * , Q, Q * the following equalities A, B, C are equivalent, A⇔B⇔C: This equivalence means that we can select any of these three divergences as a basic function and consider the others as functions of this basic one.
For any α ≥ 0 and λ = α + 1 the Rényi, the Cressie-Read and the Tsallis divergences have the same family of sublevel sets. Hence, they give the same maximizers to the conditional relative entropy minimization problems and there is no difference which entropy to use.
The CR family has a more convenient normalization factor 1/λ(λ+1) which gives a proper convexity for all powers, both positive and negative, and provides a sensible Burg limit for λ → −1 (in contrary, when α → 0 both the Rényi and Tsallis entropies tend to 0).
When α < 0 then for the Tsallis entropy function h = (x α −x) α−1 loses convexity, whereas for the Cressie-Read family convexity persists for all values of λ. The Rényi entropy also loses convexity for α < 0. Neither the Tsallis, nor the Rényi entropy were invented for use with negative α.
There may be a reason: for α < 0 the function x α is defined for x > 0 only and has a singularity at x = 0. If we assume that the divergence should exist for all non-negative distributions, then the cases α ≤ 0 should be excluded. Nevertheless, the Burg entropy which is singular at zeros is practically important and has various attractive properties. The Jeffreys' entropy (the symmetrized Kullback-Leibler divergence) is also singular at zero, but has many important properties. We can conclude at this point that it is not obvious that we have to exclude any singularity at zero probability. It may be useful to consider positive probabilities instead ("nature abhors a vacuum").
Finally, for the MaxEnt approach (conditional minimization of the relative entropy), the Rényi and the Tsallis families of divergences (α > 0) are particular cases of the Cressie-Read family because they give the same minimizers. For α ≤ 0 the Rényi and the Tsallis relative entropies lose their convexity, while the Cressie-Read family remains convex for λ ≤ −1 too.

Minima and normalization
For a given P * , the function H h achieves its minimum on the hyperplane i p i = i p * i =const at equilibrium p * i , because at this point . Therefore, we can always assume that the function h(x) achieves its minimal value at point x = 1, and this value is zero. For this purpose, one should just transform h: This normalization transforms The normalized version of any divergence H h (P P * ) could be produced by the normalization and does not need separate discussion.

Symmetrization
Another technical issue is symmetry of a divergence. If h(x) = x ln x then both H h (P P * ) (the KL divergence) and H h (P * P ) (the relative Burg entropy) are the Lyapunov functions for the Markov chains, and H h ( The transformation (20) is an involution: The fixed points of this involution are such functions h(x) that H h (P P * ) is symmetric with respect to transpositions of P and P * . There are many such h(x). An example of symmetric H h (P P * ) gives the choice h(x) = − √ x: Essentially (up to a constant addition and multiplier) this function coincides with a member of the CR family, H CR − 1 2 (12), and with one of the Tsallis relative entropies H Ts 1 2 (15). The involution (20) is a linear operator, hence, for any convex h(x) we can produce its symmetrization:

Lyapunov Functionals for Markov Chains
Let us consider continuous time Markov chains with positive equilibrium probabilities p * j . The dynamics of the probability distribution p i satisfy the Master equation (the Kolmogorov equation): where coefficients q ij (i = j) are non-negative. For chains with a positive equilibrium distribution p * j another equivalent form is convenient: where p * i and q ij are connected by identity The time derivative of the Csiszár-Morimoto function H h (p) (6) due to the Master equation is To prove this formula, it is worth to mention that for any n numbers for all positive x, y then h(x) is convex on R + . Therefore, if for some function h(x) H h (p) is the Lyapunov function for all the Markov chains with equilibrium P * then h(x) is convex on R + . The Lyapunov functionals H h do not depend on the kinetic coefficients q ij directly. They depend on the equilibrium distribution p * which satisfies the identity (23). This independence of the kinetic coefficients is the universality property.

"Lyapunov Divergences" for Discrete Time Markov Chains
The Csiszár-Morimoto functions (6) are also Lyapunov functions for discrete time Markov chains. Moreover, they can serve as a "Lyapunov distances" [51] between distributions which decreases due to time evolution (and not only the divergence between the current distribution and equilibrium). In more detail, let A = (a ij ) be a stochastic matrix in columns: The ergodicity contraction coefficient for A is a number α(A) [52,53]: Let us consider in this subsection the normalized Csiszár- Theorem about relative entropy contraction. (The generalized data processing Lemma.) For each two probability positive distributions P, Q the divergence H h (P Q) decreases under action of stochastic matrix A [54,55]: The generalizations of this theorem for general Markov kernels seen as operators on spaces of probability measures was given by [56]. The shift in time for continuous-time Markov chain is a column-stochastic matrix, hence, this contraction theorem is also valid for continuous-time Markov chains.
The question about a converse theorem arises immediately. Let the contraction inequality hold for two pairs of positive distributions (P, Q) and (U, V ) and for all H h : Could we expect that there exists such a stochastic matrix A that U = AP and V = AP ? The answer is positive: The converse generalized data processing lemma. Let the contraction inequality (27) hold for two pairs of positive distributions (P, Q) and (U, V ) and for all normalized H h . Then there exists such a column-stochastic matrix A that U = AP and V = AQ [54].
This means that for the system of inequalities (27) (for all normalized convex functions h on ]0, ∞[) is necessary and sufficient for existence of a (discrete time) Markov process which transform the pair of positive distributions (P, Q) in (U, V ). It is easy to show that for continuous-time Markov chains this theorem is not valid: the attainable regions for them are strictly smaller than the set given by inequalities (27) and could be even non-convex (see [60] and Section 8.1 below).

Additivity Property
The additivity property with respect to joining of subsystems is crucial both for the classical thermodynamics and for the information theory.
Let us consider a system which is result of joining of two subsystems. A state of the system is an ordered pair of the states of the subsystems and the space of states of the system is the Cartesian product of the subsystems spaces of state. For systems with finite number of states this means that if the states of subsystems are enumerated by indexes j and k then the states of the system are enumerated by pairs jk. The probability distribution for the whole system is p jk , and for the subsystems the probability distributions are the marginal distributions q j = k p jk , r k = j p jk .
The additive functions of state are defined for each state of the subsystems and for a state of the whole system they are sums of these subsystem values: where v j and w k are functions of the subsystems state.
In classical thermodynamics such functions are called the extensive quantities. For expected values of additive quantities the similar additivity condition holds: Let us consider these expected values as functionals of the probability distributions: and w({r k }). Then the additivity property for the expected values reads: where q j and the r k are the marginal distributions. Such a linear additivity property is impossible for non-linear entropy functionals, but under some independence conditions the entropy can behave as an extensive variable.
Let P be a product of marginal distributions. This means that the subsystems are statistically independent: p jk = q j r k . Assume also that the distribution P * is also a product of marginal distributions p * jk = q * j r * k . Then some entropies reveal the additivity property with respect to joining of independent systems.

The BGS relative entropy
. It is obvious that a convex combination of the Shannon and Burg entropies has the same additivity property.

The Rényi entropy H
For α → ∞ the Min-entropy also inherites this property.
This property implies the separation of variables for the entropy maximization problems. Let functionals u 1 ({p jk }), . . . u m ({p jk }) be additive (28) (29) and let the relative entropy H(P P * ) be additive with respect to joining of independent systems. Then the solution to the problem is p min jk = q min j r min k , where q min j , r min k are solutions of partial problems: and for some redistribution of the additive functionals values Let us call this property the separation of variables for independent subsystems. Neither the CR, nor the Tsallis divergences families have the additivity property. It is proven [46] that a function H h has the additivity property if and only if it is a convex combination of the Shannon and Burg entropies. See also Theorem 3 in Appendix.
Nevertheless, they have the property of separation of variables for independent subsystems because of the coincidence of the level sets with the additive function, the Rényi entropy (for all α > 0).
The Tsallis entropy family has absolutely the same property of separation of variables as the Rényi entropy. To extend this property of the Rényi Tsallis entropies for negative α, we have to change there min to max.
For the CR family the result sounds even better: because of better normalization, the separation of variables is valid for H CR λ → min problem for all values λ ∈] − ∞, ∞[.

Separation of Variables for Partition of the State Space
Another important property of separation of variables is valid for all divergences which have the form of a sum of convex functions f (p i , p * i ). Let the set of states be divided into two subsets, I 1 and I 2 , and let the functionals u 1 , . . . u m be linear. We represent each probability distribution as a direct sum P = P 1 ⊕ P 2 , where P 1,2 are restrictions of P on I 1,2 .
Let us consider the problem H(P P * ) → min subject to conditions u i (P ) = U i for a set of linear functionals u i (P ). The solution P min to this problem has a form P min = P min 1 ⊕ P min 2 , where P 1,2 are solutions to the problems H(P 1,2 P * 1,2 ) → min subject to conditions u i (P 1,2 ) = U 1,2 i and i∈I 1,2 p 1,2 i = π 1,2 for some redistribution of the linear functionals values, U i = U 1 i + U 2 i , and of the total probability, 1 = π 1 + π 2 (π 1,2 ≥ 0) . Again, the solution to the divergence minimization problem is composed from solutions of the partial maximization problems. Let us call this property the separation of variables for incompatible events (because I 1 ∩ I 2 = ∅).
This property is trivially valid for the Tsallis family (for α > 0, and for α < 0 with the change of minimization to maximization) and for the CR family. For the Rényi family it also holds (for α > 0, and for α < 0 with the change from minimization to maximization), because the Rényi entropy is a function of those trace-form entropies, their level sets coincide.
Again, a simple check shows that this separation of variables property holds also for the convex combination of Shannon's and Burg's entropies, βD KL (P P * ) + (1 − β)D KL (P * P ).
The question arises: is there any new divergence that has the following three properties: (i) the divergence H(P P * ) should decrease in Markov processes with equilibrium P * , (ii) for minimization problems the separation of variables for independent subsystems holds and (iii) the separation of variables for incompatible events holds. A new divergence means here that it is not a function of a divergence from the CR family or from the convex combination of the Shannon and the Burg entropies.
The answer is: no, any divergence which has these three properties and is defined and differentiable for positive distributions is a function of H h for h(x) = p α or h(x) = βx ln x − (1 − β) ln x. If we relax the differentiability property, then we have to add to the CR family the CR analogue of min-entropy The limiting case for the CR family for λ → −∞ is less known but is also a continuous and piecewise differential Lyapunov function for the Master equation: Both properties of separation of variables are based on the specific additivity properties: additivity with respect to the composition of independent systems and additivity with respect to the partitioning of the space of states. Separation of variables can be considered as a weakened form of additivity: not the minimized function should be additive but there exists such a monotonic transformation of scale after which the function becomes additive (and different transformations may be needed for different additivity properties).

"No More Entropies" Theorems
The classical Shannon work included the characterization of entropy by its properties. This meant that the classical notion of entropy is natural, and no more entropies are expected. In the seminal work of Rényi, again the characterization of entropy by its properties was proved, and for this, extended family the no more entropies theorem was proved too. In this section, we prove the next no more entropies theorem, where two one-parametric families are selected as sensible: the CR family and the convex combination of Shannon's and Burg's entropies. They are two branches of solutions of the correspondent functional equation and intersect at two points: Shannon's entropy (λ = 1 in the CR family) and Burg's entropy (λ = 0). We consider entropies as equivalent if their level sets coincide. In that sense, the Rényi entropy and the Tsallis entropy (with α > 0) are equivalent to the CR entropy with α − 1 = λ, λ > −1.
Following Rényi, we consider entropies of incomplete distributions: p i ≥ 0, i p i ≤ 1. The divergence H(P P * ) is a C 1 smooth function of a pair of positive generalized probability distributions P = (p i ), p i > 0 and P * = (p * i ), p * i > 0, i = 1, . . . n. The following 3 properties are required for characterization of the "natural" entropies.
1. To provide the separation of variables for incompatible events together with the symmetry property we assume that the divergence is separable, possibly, after a scaling transformation: there exists such a function of two variables f (p, p * ) and a monotonic function of one variable φ(x) that H(P P * ) = φ( i f (p i , p * i )). This formula allows us to define H(P P * ) for all n.
2. H(P P * ) is a Lyapunov function for the Kolmogorov equation (22) for any Markov chain with equilibrium P * . (One can call these functions the universal Lyapunov functions because they do not depend on the kinetic coefficients directly, but only on the equilibrium distribution P * .) 3. To provide separation of variables for independent subsystems we assume that H(P P * ) is additive (possibly after a scaling transformation): there exists such a function of one variable ψ(x) that the function ψ(H(P P * )) is additive for the union of independent subsystems: if P = (p ij ), p ij = q j r j , p * ij = q * j r * j , then ψ(H(P P * )) = ψ(H(Q Q * )) + ψ(H(R R * )). In a paper [46] this family was identified as the Tsallis relative entropy with some abuse of language, because in the Tsallis entropy the case with α < 0 is usually excluded.
First of all, let us prove that any function which satisfies the conditions 1 and 2 is a monotone function of a Csiszár-Morimoto function (6) for some convex smooth function h(x). This was mentioned in 2003 by P. Gorban [46]. Recently, a similar statement was published by S. Amari (Theorem 1 in [57]). Proof. Let us consider a Markov chain with two states. For such a chain If H is a Lyapunov function thenḢ ≤ 0 and the following inequality holds: We can consider p 1 , p 2 as independent variables from an open triangle D = {(p 1 , p 2 ) | p 1,2 > 0, p 1 + p 2 < 1}. For this purpose, we can include the Markov with two states into a chain with three states and q 3i = q i3 = 0.
This lemma has important corollaries about many popular divergences H(P (t) P * ) which are not Lyapunov functions of Markov chains. This means that there exist such distributions P 0 and P * and a Markov chain with equilibrium distribution P * that due to the Kolmogorov equations if P (0) = P 0 . This Markov process increases divergence between the distributions P, P * (in a vicinity of P 0 ) instead of making them closer. For example, Corollary 1. The following Bregman divergences [58] are not universal Lyapunov functions for Markov chains: These divergences violate the requirement: due to the Markov process distributions always monotonically approach equilibrium. (Nevertheless, among the Bregman divergences there exists a universal Lyapunov function for Markov chains, the Kulback-Leibler divergence.) We place the proof of Theorem 1 in Appendix.
Remark. If we relax the requirement of smoothness and consider in conditions of Theorem 1 just continuous functions, then we have to add to the answer the limit divergences,

Entropy: a Function or an Order?
Theorem 1 gives us all of the divergences for which (i) the Markov chains monotonically approach their equilibrium, (ii) the level sets are the same as for a separable (sum over states) divergence and (iii) the level sets are the same as for a divergence which is additive with respect to union of independent subsystems.
We operate with the level sets and their orders, compare where the divergence is larger (for monotonicity of the Markov chains evolution), but the values of entropy are not important by themselves. We are interested in the following order: P precedes Q with respect to the divergence H ... (P P * ) if there exists such a continuous curve P (t) (t ∈ [0, 1]) that P (0) = P , P (1) = Q and the function H(t) = H ... (P (t) P * ) monotonically decreases on the interval t ∈ [0, 1]. This property is invariant with respect to a monotonic (increasing) transformation of the divergence. Such a transformation does not change the conditional minimizers or maximizers of the divergence.
There exists one important property that is not invariant with respect to monotonic transformations. The increasing function F (H) of a convex function H(P ) is not obligatorily a convex function.
Nevertheless, the sublevel sets given by inequalities H(P ) ≤ a coincide with the sublevel sets F (H(P )) ≤ F (a). Hence, sublevel sets for F (H(P )) remain convex.

The Jensen inequality
is not invariant with respect to monotonic transformations. Instead of them, there appears the max form analogue of the Jensen inequality: This inequality is invariant with respect to monotonically increasing transformations and it is equivalent to convexity of sublevel sets.

Proposition 1. All sublevel sets of a function H on a convex set V are convex if and only if for any two points P, Q ∈ V and every θ ∈ [0, 1] the inequality (31) holds.
It seems very natural to consider divergences as orders on distribution spaces, and discuss only properties which are invariant with respect to monotonic transformations. From this point of view, the CR family appears absolutely naturally from the additivity (ii) and the "sum over states" (iii) axioms, as well as the convex combination βD KL (P P * ) + (1 − β)D KL (P * P ) (α ∈ [0, 1]), and in the above property context there are no other smooth divergences.

Description of Markov Order
The CR family and the convex combinations of Shannon's and Burg relative entropies are distinguished families of divergences. Apart from them there are many various "divergences", and even the Csiszár-Morimoto functions (6) do not include all used possibilities. Of course, most users prefer to have an unambiguous choice of entropy: it would be nice to have "the best entropy" for any class of problems. But from some point of view, ambiguity of the entropy choice is unavoidable. In this section we will explain why the choice of entropy is necessarily non unique and demonstrate that for many MaxEnt problems the natural solution is not a fixed distribution, but a well defined set of distributions.
The most standard use of divergence in many application is as follows: 1. On a given space of states an "equilibrium distribution" P * is given. If we deal with the probability distribution in real kinetic processes then it means that without any additional restriction the current distribution will relax to P * . In that sense, P * is the most disordered distribution. On the other hand, P * may be considered as the "most disordered" distribution with respect to some a priori information.
2. We do not know the current distribution P , but we do know some linear functionals, the moments u(P ).
3. We do not want to introduce any subjective arbitrariness in the estimation of P and define it as the "most disordered" distribution for given value u(P ) = U and equilibrium P * . That is, we define P as solution to the problem: Without the condition u(P ) = U the solution should be simply P * . Now we have too many entropies and do not know what is the optimal choice of H ... and what should be the optimal estimate of P . In this case the proper question may be: which P could not be such an optimal estimate? We can answer the exclusion question. Let for a given P 0 the condition hold, u(P 0 ) = U. If there exists a Markov process with equilibrium P * such that at point P 0 due to the Kolmogorov equation (22) dP dt = 0 and d(u(P )) dt = 0 then P 0 cannot be the optimal estimate of the distribution P under condition u(P ) = U. The motivation of this approach is simple: any Markov process with equilibrium P * increases disorder and brings the system "nearer" to the equilibrium P * . If at P 0 it is possible to move along the condition plane towards the more disordered distribution then P 0 cannot be considered as an extremely disordered distribution on this plane. On the other hand, we can consider P 0 as a possible extremely disordered distribution on the condition plane, if for any Markov process with equilibrium P * the solution of the Kolmogorov equation (22) P (t) with initial condition P (0) = P 0 has no points on the plane u(P ) = U for t > 0.
Markov process here is considered as a "randomization". Any set C of distributions can be divided in two parts: the distributions which retain in C after some non-trivial randomization and the distributions which leave C after any non-trivial randomization. The last are the maximally random elements of C: they cannot become more random and retain in C. Conditional minimizers of relative entropies H h (P P * ) in C are maximally random in that sense.
There are too many functions H h (P P * ) for effective description of all their conditional minimizers. Nevertheless, we can describe the maximally random distributions directly, by analysis of Markov processes.
To analyze these properties more precisely, we need some formal definitions.

Definition 1. (Markov preorder).
If for distributions P 0 and P 1 there exists such a Markov process with equilibrium P * that for the solution of the Kolmogorov equation with P (0) = P 0 we have P (1) = P 1 then we say that P 0 and P 1 are connected by the Markov preorder with equilibrium P * and use notation P 0 ≻ 0 P * P 1 .

Definition 2.
Markov order is the closed transitive closure of the Markov preorder. For the Markov order with equilibrium P * we use notation P 0 ≻ P * P 1 .
For a given P * = (p * i ) and a distribution P = (p i ) the set of all vectors v with coordinates where p * i and q ij ≥ 0 are connected by identity (23) is a closed convex cone. This is a cone of all possible time derivatives of the probability distribution at point P for Markov processes with equilibrium P * = (p * i ). For this cone, we use notation Q (P,P * ) Definition 3. For each distribution P and a n-dimensional vector ∆ we say that ∆ < (P,P * ) 0 if ∆ ∈ Q (P,P * ) . This is the local Markov order.

Proposition 2. Q (P,P * ) is a proper cone, i.e., it does not include any straight line.
Proof. To prove this proposition its is sufficient to analyze the formula for entropy production (for example, in form (24)) and mention that for strictly convex h (for example, for traditional x ln x or (x − 1) 2 /2) dH h /dt = 0 if and only if dP/dt = 0. If the cone Q (P,P * ) includes both vectors x and −x (x = 0 it means that there exist Markov chains with equilibrium P * and with opposite time derivatives at point P . Due to the positivity of entropy production (24) this is impossible.
The connection between the local Markov order and the Markov order gives the following proposition, which immediately follows from definitions. Proposition 3. P 0 ≻ P * P 1 if and only if there exists such a continuous almost everywhere differentiable curve P (t) in the simplex of probability distribution that P (0) = P 0 , P (1) = P 1 and for all t For our purposes, the following estimate of the Markov order through the local Markov order is important.
This proposition follows from the characterization of the local order and detailed description of the cone Q (P (t),P * ) (Theorem 2 below).
Let us recall that a convex pointed cone is a convex envelope of its extreme rays. A ray with directing vector x is a set of points λx (λ ≥ 0). We say that l is an extreme ray of Q if for any u ∈ l and any x, y ∈ Q, whenever u = (x + y)/2, we must have x, y ∈ l. To characterize the extreme rays of the cones of the local Markov order Q (P,P * ) we need a graph representation of the Markov chains. We use the notation A i for states (vertices), and designate transition from state A i to state A j by an arrow (edge) A i → A j . This transition has its transition intensity q ji (the coefficient in the Kolmogorov equation (21)).

Lemma 2. Any extreme ray of the cone Q (P,P * ) corresponds to a Markov process which transition graph is a simple cycle
where k ≤ n, all the indices i 1 , . . . i k are different, and transition intensities for a directing vector of such an extreme ray q i j+1 i j may be selected as 1/p * i j : (here we use the standard convention that for a cycle q i k+1 i k = q i 1 i k ).
Proof. First of all, let us mention that if for three vectors x, y, u ∈ Q (P,P * ) we have u = (x + y)/2 then the set of transitions with non-zero intensities for corresponding Markov processes for x and y are included in this set for u (because negative intensities are impossible). Secondly, just by calculation of the free variables in the equations (23) (with additional condition) we find that the the amount of non-zero intensities for a transition scheme which represents an extreme ray should be equal to the amount of states included in the transition scheme. Finally, there is only one scheme with k vertices, k edges and a positive equilibrium, a simple oriented cycle.

Theorem 2.
Any extreme ray of the cone Q (P,P * ) corresponds to a Markov process whose transition graph is a simple cycle of the length 2: A i ⇄ A j . A transition intensities q ij , q ji for a directing vector of such an extreme ray may be selected as Proof. Due to Lemma 2, it is sufficient to prove that for any distribution P the right hand side of the Kolmogorov equation (22) (under the standard convention regarding cyclic order). Other coordinates of v are zeros. Let us find the minimal value of p i j /p * i j and rearrange the indices by a cyclic permutation to put this minimum in the first place: The vector v is a sum of two vectors: a directing vector for the cycle A i 2 → . . . A i k → A i 2 of the length n − 1 with transition intensities given by formula (34) (under the standard convention about the cyclic order for this cycle) and a vector where v 2 is the directing vector for a cycle of length 2, A i 1 ⇄ A i 2 which can have only two non-zero coordinates: The coefficient in front of v 2 is positive because

be omitted). A conic combination of conic combinations is a conic combination again.
It is quite surprising that the local Markov order and, hence, the Markov order also are generated by the reversible Markov chains which satisfy the detailed balance principle. We did not include any reversibility assumptions, and studied the general Markov chains. Nevertheless, for the study of orders, the system of cycles of length 2 all of which have the same equilibrium is sufficient.

Combinatorics of Local Markov Order
Let us describe the local Markov order in more detail. First of all, we represent kinetics of the reversible Markov chains. For each pair A i , A j (i = j) we select an arbitrary order in the pair and write the correspondent cycle of the length 2 in the form A i ⇆ A j . For this cycle we introduce the directing vector γ ij with coordinates γ ij k = −δ ik + δ jk (36) where δ ik is the Kronecker delta. This vector has the ith coordinate −1, the jth coordinate 1 and other coordinates are zero. Vectors γ ij are parallel to the edges of the standard simplex in R n . They are antisymmetric in their indexes: γ ij = −γ ji . We can rewrite the Kolmogorov equation in the form where i = j, each pair is included in the sum only once (in the preselected order of i, j) and The coefficient r ji ≥ 0 satisfies the detailed balance principle: We use the three-value sign function: With this function we can rewrite Equation (37) again as follows: The non-zero coefficients r ji may be arbitrary positive numbers. Therefore, using Theorem 2, we immediately find that the cone of the local Markov order at point P is where cone{} stands for the conic hull. The number sign p i   In Figure 1 we represent compartments and cones of the local Markov order for the Markov chains with three states, A 1,2,3 . The reversible Markov chain consists of three reversible transitions The topology of the partitioning of the standard simplex into compartments and the possible values of the cone Q (P,P * ) do not depend on the position of the equilibrium distribution P * .
Let us describe all possible compartments and the correspondent local Markov order cones. For every natural number k ≤ n − 1 the k-dimensional compartments are numerated by surjective functions σ : {1, 2, . . . , n} → {1, 2, . . . , k + 1}. Such a function defines the partial ordering of quantities p j p * j inside the compartment: Let us use for the correspondent compartment notation C σ and for the Local Markov order cone Q σ . Let k i be a number of elements in preimage of i (i = 1, . . . , k): k i = |{j | σ(j) = i}|. It is convenient to represent surjection σ as a tableau with k rows and k i cells in the ith row filled by numbers from {1, 2, . . . , n}. First of all, let us draw diagram, that is a finite collection of cells arranged in left-justified rows. The ith row has k i cells. A tableau is obtained by filling cells with numbers {1, 2, . . . , n}. Preimages of i are located in the ith row. The entries in each row are increasing. (This is convenient to avoid ambiguity of the representation of the surjection σ by the diagram.) Let us use for tableaus the same notation as for the corresponding surjections.
Let a tableau A have k rows. We say that a tableau B follows A (and use notation A → B) if B has k − 1 rows and B can be produced from A by joining of two neighboring rows in A (with ordering the numbers in the joined row). For the transitive closure of the relation → we use notation ⇛.
Here r∂U stands for the "relative boundary" of a set U in the minimal linear manifold which includes U.
The following Proposition characterizes the local order cone through the surjection σ. It is sufficient to use in definition of Q σ (40) vectors γ ij (36) with i and j from the neighbor rows of the diagram (see Figure 1).

Proposition 6.
For a given surjection σ compartment C σ and cone Q σ have the following description: Compartment C σ is defined by equalities p i where i, j belong to one row of the tableau σ and where j is situated in a row one step down from i in the tableau (σ(j) = σ(i) + 1).
Cone Q σ is a conic hull of k−1 i=1 k i k i+1 vectors γ ij . For these vectors, j is situated in a row one step down from i in the tableau. Extreme rays of Q σ are products of the positive real half-line on vectors γ ij (43).
Each compartment has the lateral faces and the base. We call the face a lateral face, if its closure includes the equilibrium P * . The base of the compartment belongs to a border of the standard simplex of probability distributions.
To enumerate all the lateral faces of a k-dimensional compartment C σ of codimension s (in C σ ) we have to take all subsets with s elements in {1, 2, . . . , k}. For any such a subset J the correspondent k − s-dimensional lateral face is given by additional equalities p i The 1-dimensional lateral faces (extreme rays) of compartment C σ are given by selection of one number from {1, 2, . . . , k} (this number is the complement of J). For this number r, the correspondent 1-dimensional face is a set parameterized by a positive number a ∈]1, a r ], a r = 1/ σ(i)≤r p * i : The compartment C σ is the interior of the k-dimensional simplex with vertices P * and v r (r = 1, 2, . . . k). The vertex v r is the intersection of the correspondent extreme ray (45) with the border of the standard simplex of probability distributions: P = v r if The base of the compartment C σ is a k − 1-dimensional simplex with vertices v r (r = 1, 2, . . . k).
It is necessary to stress that we use the reversible Markov chains for construction of the general Markov order due to Theorem 2.

Conditionally Extreme Distributions in Markov Order
The Markov order can be used to reduce the uncertainty in the standard settings. Let the plane L of the known values of some moments be given: u i (P ) = U i on L. Assume also that the "maximally disordered" distribution (equilibrium) P * is known and we assume that the probability distribution is P * if there is no restrictions. Then, the standard way to evaluate P for given moment conditions u i (P ) = U i is known: just to minimize H ... (P P * ) under these conditions. For the Markov order we also can define the conditionally extreme points on L. First of all, it is obvious that in the case when all the moments u i (P ) are just some of the values p i , then there exists only one extreme point of the Markov order on L, and this point is, at the same time, the conditional minimum on L of all Csiszár-Morimoto functions H h (P ) (6) (see, for example, Figure 2). This situation is unstable, and for a small perturbation of L the set of extreme points of the Markov order on L includes the intersection of L with one of compartments (Figure 3a). For the Markov chains with three states, each point of this intersection is a conditional minimizer of one of the CR divergences (see Fig. 3a). Such a situation persists for all L in general positions (Figure 3b). The extreme points of the family βD KL (P P * ) + (1 − β)D KL (P * P ) form an interval which is strictly inside the interval of the extreme points of the Markov order on L. For higher dimensions of L ∩ Σ n the Markov order on L also

How to Find the Most Random Distributions?
Let the plane L of the known values of some moments be given: u i (P ) = j u i j p j = U i (i = 1, . . . m) on L. For a given divergence H(P P * ) we are looking for a conditional minimizer P : We can assume that H(P P * ) is convex. Moreover, usually it is one of the Csiszár-Morimoto functions (6). This is very convenient for numerical minimization because the matrix of second derivatives is diagonal. Let us introduce the Lagrange multipliers µ i (i = 1, . . . m) and write the system of equations (µ 0 is the Lagrange multiplier for the total probability identity j p j = 1 : Here we have n + m + 1 equations for n + m + 1 unknown variables (p j , µ i , µ 0 ). Usually H is a convex function with a diagonal matrix of second variables and the method of choice for solution of this equation (48) is the Newton method. On the l + 1st iteration to find P l+1 = P l + ∆P we have to solve the following system of linear equations For a diagonal matrix of the second derivatives the first n equations can be explicitly resolved. If for the solution of this system (49) the positivity condition p l j + ∆p j > 0 does not hold (for some of j) then we should decrease the step, for example by multiplication ∆P := θ∆P , where For initial approximation we can take any positive normalized distribution which satisfies the conditions u i (P ) = U i (i = 1, . . . m).
For the Markov orders the set of conditionally extreme distributions consists of intersections of L with compartments.
Here we find this set for one moment condition of the form u(P ) = j u j p j = U. First of all, assume that U = U * , where U * = u(P * ) = j u j p * j (if U = U * then equilibrium is the single conditionally extreme distribution). In this case, the set of conditionally extreme distributions is the intersection of the condition hyperplane with the closure of one compartment and can be described by the following system of equations and inequalities (under standard requirements p i ≥ 0, i p i = 1 ): To find this solution it is sufficient to study dynamics of u(P ) due to equations (37) and to compare it with dynamics of u(P ) due to a model systemṖ = P * − P . This model system is also a Markov chain and, therefore, P * − P ∈ Q (P,P * ) . Equations and inequalities (50) mean that the set of conditionally extreme distributions is the intersection of the condition hyperplane with the closure of compartment C. In C, numbers p i p * i have the same order on the real line as numbers u i (U − U * ) have, these two tuples of numbers correspond to the same tableau σ and C = C σ .
For several linearly independent conditions there exists a condition plane L: Let us introduce the m-dimensional space T with coordinates u i . Operator u(P ) = (u i (P )) maps the distribution space into T and the affine manifold L (51) maps into a point with coordinates u i = U i . If P * ∈ L then the problem is trivial and the only extreme distribution of the Markov order on L is P * . Let us assume that P * / ∈ L.
For each distribution P ∈ L we can study the possible direction of motions of projection distributions onto T due to the Markov processes.
First of all, let us mention that if u(γ ij ) = 0 then the transitions A i ⇆ A j move the distribution along L. Hence, for any conditionally extreme distribution P ∈ L this transition A i ⇆ A j should be in equilibrium and the partial equilibrium condition holds: p i Let us consider processes with u(γ ij ) = 0. If there exists a convex combination (39) of vectors u(γ ij )sign p i p * i − p j p * j (u(γ ij ) = 0) that is equal to zero then P cannot be an extreme distribution of the Markov order on L.
These two conditions for vectors γ ij with u(γ ij ) = 0 and for the set of vectors with non-zero projection on the condition space define the extreme distributions of the Markov order on the condition plane L for several conditions.

Reference Distributions for Main Divergences
A system with equilibrium P * is given and expected values of some variables u j (P ) = U j are known. We need to find a distribution P with these values u j (P ) = U j and is "the closest" to the equilibrium distribution under this condition.
This distribution parameterized through expectation values is often called the reference distribution or generalized canonical distribution. After Gibbs and Jaynes, the standard statement of this problem is an optimization problem: H(P P * ) → min, u j (P ) = U j for appropriate divergence H(P P * ). If the number of conditions is m then this optimization problem can be often transformed into m + 1 equations with m + 1 unknown Lagrange multipliers.
In this section, we study the problem of the generalized canonical distributions for single condition For the Csiszár-Morimoto functions H h (P P * ) We assume that the function h ′ (x) has an inverse function g: g(h ′ (x)) = x for any x ∈]0, ∞[. The method of Lagrange multipliers gives for the generalized canonical distribution: As a result, we get the final expression for the distribution and equations for Lagrange multipliers µ 0 and µ: If the image of h ′ (x) is the whole real line (h ′ (]0, ∞[) = R) then for any real number y the value g(y) ≥ 0 is defined and there exist no problems about positivity of p i due to (54).
). Therefore, g(x) = exp x and for the generalized canonical distribution we get As a result, we get one equation for µ and an explicit expression for µ 0 through µ. These µ 0 and µ have the opposite sign comparing to (5) just because the formal difference between the entropy maximization and the relative entropy minimization. Equation (55) is essentially the same as (5).
For the Burg entropy h ′ (x) = − 1 x , g(x) = − 1 x too and For the Lagrange multipliers µ 0 , µ we have a system of two algebraic equations For the convex combination of the BGS and Burg entropies h ′ (x) = β ln x − 1−β x (0 < β < 1), and the function x = g(y) is a solution of a transcendent equation Such a solution exists for all real y because this h ′ (x) is a (monotonic) bijection of ]0, ∞[ on the real line.
Solution to Equation (58) can be represented through a special function, the Lambert function [62]. This function is a solution to the transcendent equation we w = z and is also known as W function, Ω function or modified logarithm lmz [36]. Below we use the main branch w = lmz for which lmz > 0 if z > 0 and lm0 = 0. Let us write (58) in the form where δ = (1 − β)/β, Λ = −y/β. Then Another equivalent representation of the solution gives Indeed, let us take z = δ/x and calculate exponent of both sides of (59). After simple transformations, we obtain ze z = δe Λ .
The identity lma = ln a − ln lma is convenient for algebraic operations with this function. Many other important properties are collected in [62].
The generalized canonical distribution for the convex combination of the BGS and Burg divergence is [36] where Λ i = − 1 β (µ 0 + u i µ), δ = (1 − β)/β and equations (54) hold for the Lagrange multipliers. For small 1 − β (small addition of the Burg entropy to the BGS entropy) we have For λ = 1 (a quadratic divergence) we easily get linear equations and explicit solutions for µ 0 and µ. If λ = 1 2 then equations for the Lagrange multipliers (54) become quadratic and also allow explicit solution. The same is true for λ = 1 3 and 1 4 but explicit solutions to the correspondent cubic or quartic equations are too cumbersome.
We studied the generalized canonical distributions for one condition u(P ) = U and main families of entropies. For the BGS entropy, the method of Lagrange multipliers gives one transcendent equation for the multiplier µ 1 and explicit expression for µ 0 as a function of µ 1 (55). In general, for functions H h , the method gives a system of two equations (54). For the Burg entropy this is a system of algebraic equation (57). For a convex combination of the BGS and the Burg entropies the expression for generalized canonical distribution function includes the special Lambert function (60). For the CR family the generalized canonical distribution is presented by formula (61). for several values of λ it can be represented in explicit form. The Tsallis entropy family is a subset of the CR family (up to constant multipliers).

Polyhedron of Generalized Canonical Distributions for the Markov Order
The set of the most random distributions with respect to the Markov order under given condition consists of those distributions which may be achieved by randomization which has the given equilibrium distribution and does not violate the condition.
In the previous section, this set was characterized for a single condition i p i u i = U, U = U * by a system of inequalities and equations (50). It is a polyhedron that is an intersection of the closure of one compartment with the hyperplane of condition. Here we construct the dual description of this polyhedron as a convex envelope of the set of extreme points (vertices).
The Krein-Milman theorem gives general backgrounds of such a representation of convex compact sets in locally convex topological vector spaces [63]: a compact convex set is the closed convex hull of its extreme points. (An extreme point of a convex set K is a point x ∈ K which cannot be represented as an average x = 1 2 (y + z) for y, z ∈ K, y, z = x.) Let us assume that there are k+1 ≤ n different numbers in the set of numbers u i (U −U * ). There exists the unique surjection σ : {1, 2, . . . n} → {1, 2, . . . k + 1} with the following properties: The polyhedron of generalized canonical distributions is the intersection of the condition plane i p i u i = U with the closure of C σ .
This closure is a simplex with vertices P * and v r (r = 1, 2, . . . k) (46). The vertices of the intersection of this simplex with the condition hyperplane belong to edges of the simplex, hence we can easily find all of them: the edge [x, y] has nonempty intersection with the condition hyperplane if either u(x) ≥ U&u(y) ≤ U or u(x) ≤ U&u(y) ≥ U. This intersection is a single point P if u(x) = u(y): If u(x) = u(y) then the intersection is the whole edge, and the vertices are x and y. For example, if U is sufficiently close to U * then the intersection is a simplex with k vertices w r (r = 1, 2, . . . k). Each w r is the intersection of the edge [P * , v r ] with the condition hyperplane.
Let us find these vertices explicitly. We have a system of two equations a i, σ(i)≤r Position of the vertex w r on the edge [P * , v r ] is given by the following expressions If b ≥ 0 for all r then the polyhedron of generalized canonical distributions is a simplex with vertices w r . If the solution becomes negative for some r then the set of vertices changes qualitatively and some of them belong to the base of C σ . For example, in Figure 3a the interval of the generalized canonical distribution (1D polyhedron) has vertices of two types: one belongs to the lateral face, another is situated on the basement of the compartment. In Figure 3b both vertices belong to the lateral faces. Vertices w r on the edges [P * , v r ] have very special structure: the ratio p i /p * i can take for them only two values, it is either a or b.
Another form for representation of vertices w r (64) can be found as follows. w r belongs to the edge [P * , v r ], hence, w r = λP * + (1 − λ)v r for some λ ∈ [0, 1]. Equation for the value of λ follows from the condition u(w r ) = U: λU * + (1 − λ)u(v r ) = U. Hence, we can use (62) with x = P * , y = v r .
For sufficiently large value of U − U * for some of these vertices b loses positivity, and instead of them the vertices on edges [v r , v q ] (46) appear.
There exists a vertex on the edge then his vertex has the form P = λv r + (1 − λ)v q and for λ the condition u(P ) = U gives (62) belongs to the condition plane and the extreme distributions are u(v r ) u(v q ).
For each of v r the ratio p i /p * i can take only two values: a r or 0. Without loss of generality we can assume that q > r. For a convex combination λv r + (1 − λ)v q (1 > λ > 0) the ratio p i /p * i can take three values: λa r + (1 − λ)a q (for σ(i) ≤ r), (1 − λ)a q (for r < σ(i) ≤ q) and 0 (for σ(i) > q).
The case when a vertex is one of the v r is also possible. In this case, there are two possible values of p i /p * i , it is either a r or 0. All the generalized canonical distributions from the polyhedron are convex combinations of its extreme points (vertices). If the set of vertices is {w r }, then for any generalized canonical distributions P = λ i w i (λ i ≥ 0, i λ i = 1). The vertices can be found explicitly. Explicit formulas for the extreme generalized canonical distributions are given in this section: (64) and various applications of (62). These formulas are based on the description of compartment C σ given in Proposition 7 and Equation (46).

Continuous Time Kinetics
We have to discuss the history of the Markov order in the wider context of orders, with respect to which the solutions of kinetic equations change monotonically in time. The Markov order is a nice and constructive example of such an order and at the same time the prototype of all of them (similarly the Master Equation is a simple example of kinetic equations and, at the same time, the prototype of all kinetic equations).
The idea of orders and attainable domains (the lower cones of these orders) in phase space was developed in many applications: from biological kinetics to chemical kinetics and engineering. A kinetic model includes information of various levels of detail and of variable reliability. Several types of building block are used to construct a kinetic model. The system of these building blocks can be described, for example, as follows: 1. The list of components (in chemical kinetics) or populations (in mathematical ecology) or states (for general Markov chains); 2. The list of elementary processes (the reaction mechanism, the graph of trophic interactions or the transition graph), which is often supplemented by the lines or surfaces of partial equilibria of elementary processes; 3. The reaction rates and kinetic constants.
We believe that the lower level information is more accurate and reliable: we know the list of component better than the mechanism of transitions, and our knowledge of equilibrium surfaces is better than the information about exact values of kinetic constants.
It is attractive to use the more reliable lower level information for qualitative and quantitative study of kinetics. Perhaps, the first example of such a analysis was performed in biological kinetics.
In 1936, A.N. Kolmogorov [64] studied the dynamics of a pair of interacting populations of prey (x) and predator (y) in general form:ẋ = xS(x, y),ẏ = yW (x, y) under monotonicity conditions: ∂S(x, y)/∂y < 0, ∂W (x, y)/∂y < 0. The zero isoclines, the lines at which the rate of change for one population is zero (given by equations S(x, y) = 0 or W (x, y) = 0), are graphs of two functions y(x). These isoclines divide the phase space into compartments (generically with curvilinear borders). In every compartment the angle of possible directions of motion is given (compare to Figure 1).
Analysis of motion in these angles gives information about dynamics without an exact knowledge of the kinetic constants. The geometry of the zero isoclines intersection together with some monotonicity conditions give important information about the system dynamics [64] without exact knowledge of the right hand sides of the kinetic equations.
This approach to population dynamics was further developed by many authors and applied to various problems [65,66]. The impact of this work on population dynamics was analyzed by K. Sigmund in review [67].
It seems very attractive to use an attainable region instead of the single trajectory in situations with incomplete information or with information with different levels of reliability. Such situations are typical in many areas of engineering. In 1964, F. Horn proposed to analyze the attainable regions for chemical reactors [68]. This approach was applied both to linear and nonlinear kinetic equations and became popular in chemical engineering. It was applied to the optimization of steady flow reactors [69], to batch reactor optimization by use of tendency models without knowledge of detailed kinetics [70] and for optimization of the reactor structure [71]. Analysis of attainable regions is recognized as a special geometric approach to reactor optimization [72] and as a crucially important part of the new paradigm of chemical engineering [73]. Plenty of particular applications was developed: from polymerization [74] to particle breakage in a ball mill [75]. Mathematical methods for study of attainable regions vary from the Pontryagin's maximum principle [76] to linear programming [77], the Shrink-Wrap algorithm [78] and convex analysis. Attainable regions from an initial distribution a 0 for a linear system with three components A 1 , A 2 , A 3 in coordinates c 1 , c 2 (concentrations of A 1 , A 2 ) (c 3 = const − c 1 − c 2 ) [60]: for a full mechanism A 1 ⇄ A 2 ⇄ A 3 ⇄ A 1 (outlined region), for a two-step mechanism A 1 ⇄ A 2 , A 1 ⇄ A 3 (horizontally shaded region) and for a two-step mechanism A 1 ⇄ A 2 , A 2 ⇄ A 3 (vertically shaded region). Equilibrium is a * . The dashed lines are partial equilibria. The connection between attainable regions, thermodynamics and stoichiometric reaction mechanisms was studied by A.N. Gorban in the 1970s. In 1979, he demonstrated how to utilize the knowledge about partial equilibria of elementary processes to construct the attainable regions [60].
He noticed that the set (a cone) of possible direction for kinetics is defined by thermodynamics and the reaction mechanism (the system of the stoichiometric equation of elementary reactions).
Thermodynamic data are more robust than the reaction mechanism and the reaction rates are known with lower accuracy than the stoichiometry of elementary reactions. Hence, there are two types of attainable regions. The first is the thermodynamic one, which use the linear restrictions and the thermodynamic functions [79]. The second is generated by thermodynamics and stoichiometric equations of elementary steps (but without reaction rates) [60,80].
It was demonstrated that the attainable regions significantly depend on the transition mechanism ( Figure 4) and it is possible to use them for the mechanisms discrimination [81].
Already simple examples demonstrate that the sets of distributions which are accessible from a given initial distribution by Markov processes with equilibrium are, in general, non-convex polytopes [60,82] (see, for example, the outlined region in Figure 4, or, for particular graphs of transitions, any of the shaded regions there). This non-convexity makes the analysis of attainability for continuous time Markov processes more difficult (and also more intriguing).
Beyond these two distinguished one-parametric families there is the whole world of the Csiszár-Morimoto Lyapunov functionals for the Master equation (6). These functions monotonically decrease along any solution of the Master equation. The set of all these functions can be used to reduce the uncertainty by conditional minimization: for each h we could find a conditional minimizer of H h (p).
Most users prefer to have an unambiguous choice of entropy: it would be nice to have "the best entropy" for any class of problems. But from a certain point of view, ambiguity of the entropy choice is unavoidable, and the choice of all conditional optimizers instead of a particular one is a possible way to avoid an arbitrary choice. The set of these minimizers evaluates the possible position of a "maximally random" probability distribution. For many MaxEnt problems the natural solution is not a fixed distribution, but a well defined set of distributions.
The task to minimize functions H h (p) which depend on a functional parameter h seems too complicated. The Markov order gives us another way for the evaluation of the set of possible "maximally random" probability distribution, and this evaluation is, in some sense, the best one. We defined the Markov order, studied its properties and demonstrated how it can be used to reduce uncertainty.
It is quite surprising that the Markov order is generated by the reversible Markov chains which satisfy the detailed balance principle. We did not include any reversibility assumptions and studied the general Markov chains. There remain some questions about the structure and full description of the global Markov order. Nevertheless, to find the set of conditionally extreme ("most random") probability distributions, we need the local Markov order only. This local order is fully described in Section 5.2 and has a very clear geometric structure. For a given equilibrium distribution P * , the simplex of probability distributions is divided by n(n − 1)/2 hyperplanes of "partial equilibria" (this terminology comes from chemical kinetics [60,61]): p i p * i = p j p * j (there is one hyperplane for each pair of states (i, j)). In each compartment a cone of all possible time derivatives of the probability distribution is defined as a conic envelope of n(n − 1)/2 vectors (40). The extreme rays of this cone are explicitly described in Proposition 6 (43). This cone defines the local Markov order. When we look for conditionally extreme distributions, this cone plays the same role as a hyperplane given by entropy growth condition ( dS/dt > 0) in the standard approach.
For the problem of the generalized canonical (or reference) distribution the Markov order gives a polyhedron of the extremely disordered distributions. The vertices of that polyhedron can be computed explicitly.
The construction of efficient algorithms for numerical calculation of conditionally extreme compacts in high dimensions is a challenging task for our future work as well as the application of this methodology to real life problems.