An exact finite mixture model follows from the binomial theorem:
Define positive mixture weights
Here is the Beta function, so the component probability distributions are . All coefficients are computable from the observed integer counts.
Use inverse transform sampling with the supplied to choose : writing and , take . Then recycle its position inside that interval,
Recycling a uniform random variable after discrete sampling makes uniform independently of , since . It is also independent of all supplied random variables with gamma distributions. Transform these by the probability integral transform:
They are independent random variables with uniform distributions because has a unit-rate exponential distribution.
Conditional on , put and . Form the values and return their th order statistic. This uses at most the supplied random variables with gamma distributions, since . The uniform order statistic formula gives
Averaging these component probability density functions with weights recovers exactly the normalized posterior density. Unused inputs can be discarded. When , , , and the output is , as required. Interval endpoints and ties have probability zero; any consistent convention there is harmless. Exactness here concerns the mathematical algorithm; finite-precision arithmetic can be stabilized with logarithmic mixture weights.
Take observation periods to have positive lengths , and let be the Markov kernel of a sweep which first draws given the current and then draws the new vector . Write . Conditional on , independence and the gamma distribution mean and variance give
Averaging over the full conditional distribution of preserves this bound, uniformly in the current state. Therefore . This is the mechanism in bounded conditional moments imply a geometric drift.
Choose , , and . Outside , ; inside , . Hence
For completeness, is a small set, even though it approaches the boundary of the positive orthant. Choose small enough that . On , . The conditional gamma distribution of has fixed shape and rate in a compact positive interval. Its probability density function therefore has a positive common lower bound on . For those , the product probability density function of similarly has a positive common lower bound on , since every shape is positive and every rate lies in a compact positive interval. Integrating over gives on , where is the uniform distribution on and .
The transition probability density function is positive throughout the positive orthant, giving an irreducible Markov chain. Since , the box minorization includes starts in the same box, which gives an aperiodic Markov chain. The proper prior distributions and bounded Poisson distribution likelihood have positive finite evidence, so the invariant posterior distribution is proper. Together with the geometric drift condition, these hypotheses imply geometric ergodicity. In particular, the argument proves the required drift rather than assuming that every Gibbs sampler is geometrically ergodic.
The positive-period convention matters. If zero lengths are allowed, the printed assertion can fail: take , , . Then and , so
This example has a proper posterior distribution but cannot satisfy a finite quadratic geometric drift condition. Thus the proof establishes the intended claim for positive observation periods, not an unrestricted extension to zero exposure.
Write . Mean-field variational inference minimizes the reverse Kullback-Leibler divergence
Each is a probability density function; without further parametric restrictions, the optimization is over all such product probability distributions. Equivalently it maximizes the evidence lower bound for any unnormalized posterior density . Its difference from is exactly the Kullback-Leibler divergence.
Fix and write
Assume and that the displayed expected values and objective decomposition are well defined. The optimal factor is
Indeed, the part of the objective depending on is
The remaining term is fixed. Gibbs inequality gives nonnegativity of the Kullback-Leibler divergence, with equality precisely at almost everywhere. This proves global optimality of that coordinate update; it does not assert global optimality of a sequence of coordinate updates for the full nonconvex product-family problem.
The printed upper index in the list of remaining coordinates is inconsistent with the dimension ; the natural interpretation is . Also, if the exponential expression has zero or infinite normalizing constant, or the expected values are undefined, the usual coordinate formula needs additional support or integrability hypotheses. A restricted parametric factor family need not contain this unrestricted optimal factor.
Parameterize the initial product normal distribution by , using
Run the given recursion with independent random variables with uniform distributions to obtain . Because the noise laws do not depend on , the chain rule and automatic differentiation give reparameterization gradients through the entire simulated trajectory. The variational objective is , where and .
One must account for the output differential entropy, not just optimize . Simulation and differentiability of do not themselves supply an evaluable . If retains the reversible Markov chain assumption from part (b), a concrete algorithm avoids a transition-density oracle using reversible density-ratio propagation. Write , with an evaluable unnormalized posterior density, and define
Detailed balance implies
For example, integrate and replace by to obtain the first identity. This also proves absolute continuity of measures of with respect to when . An everywhere finite differentiable log-posterior has positive density, so an ordinary product normal distribution has this domination.
For each outer sample , independently run trajectories of length starting at , with endpoints . These are reverse trajectories because of detailed balance; the same recursion simulates them. Estimate by
Average over an outer batch. Compute it stably by a logarithmic sum of exponentials. Differentiate through both outer and reverse trajectories and through the explicit product-normal probability density function . Only the given state gradients of and , together with elementary normal-density derivatives, are needed. Apply stochastic gradient descent to and , refreshing the independent noise batches. The logarithmic scale parameters keep all variances positive.
This is nested Monte Carlo optimization of the exact objective in the limit of adequately increasing inner and outer sample sizes. For fixed , is generally biased, even though is an unbiased estimator. Finite logarithmic moments, suitable uniform integrability, and justified differentiation under the integral sign are needed for consistency and gradient interchange. Under those conditions increasing and controlling the optimization error yields a justified numerical approximation; nonconvexity prevents a generic guarantee of a global optimum.
There is a real assumption issue: part (c) redefines without explicitly restating reversibility. The algorithm above uses the reversible interpretation inherited from part (b). For an arbitrary differentiable recursion, one additionally needs access to output densities or a simulatable reverse kernel and regularity of the objective. Differentiability alone is insufficient: is differentiable and, with a standard-normal target, makes a point mass for , so every candidate has infinite Kullback-Leibler divergence. Thus a universal finite-objective algorithm cannot be justified solely from the literal differentiability assumptions. An invertible recursion with computable density transformations offers a different sufficient route, but invertibility is not stated here.
The most plausible cause is strong dependence between the latent variables and covariance hyperparameters in their joint posterior distribution. In the usual conditionally independent Gaussian process classification model, put , , and . The probit regression likelihood is , while has a multivariate normal distribution with covariance matrix .
For distinct inputs and positive length scales, is positive definite. Its Gaussian probability density function contributes
The unit-rate exponential prior distribution on therefore gives
Its scale is tightly tied to the current magnitude and shape of . Changing changes the eigenvectors and eigenvalues of its covariance matrix; a latent vector typical under the old covariance may be very atypical under a substantially different one. Holding fixed can thus restrict the range of a hyperparameter update, and holding the hyperparameters fixed restricts the next latent-function update. Even exact sampling of each full conditional distribution can move slowly along a narrow ridge of the joint posterior distribution.
Binary observations provide limited information about the absolute size of large correctly signed latent values. This can leave substantial scale uncertainty and strengthen the dependence. Long correlation lengths can also make the covariance matrix nearly singular and the latent coordinates highly dependent. These are plausible explanations for a large mixing time, not proofs that every data set causes slow convergence. A blocked Gibbs sampler removes dependence within its chosen blocks, but not dependence between them. Repeated identical inputs should be represented by one shared latent value; otherwise the nominal Gaussian matrix is singular and the displayed inverse formula must be replaced by a representation on its support.
Let and . With conditional independence of the binary observations given the latent function, define
The prior distribution is specified on the precision , not on . Its Jacobian determinant gives
and the marginal posterior density is proportional to . Working instead with precision or logarithmic coordinates is possible, provided all corresponding Jacobian determinants are included.
Choose an importance sampling density covering the support of , draw independently from it, and use the nonnegative unbiased likelihood estimator
Each summand integrates to , proving unbiasedness. A simple always-valid choice is the latent Gaussian prior distribution, , when the covariance matrix is nonsingular. Then the estimator is just the average of the bounded values . For duplicate inputs, sample only the distinct latent values and multiply all likelihood contributions at each common input. A Gaussian approximation to the latent posterior distribution, or its mixture with the Gaussian prior distribution, can reduce relative variance while keeping complete support. Its approximation does not replace the exact weights.
The pseudo-marginal Metropolis–Hastings algorithm stores both the current and its estimator randomness ; write their sampling law as . Propose , independently generate , compute , and accept the pair with
Here denotes the proposal probability density function from to . On rejection retain both the old parameter and the old likelihood estimate. Initialize with a positive estimate. One must not independently refresh the denominator after rejection.
Correctness follows from the extended posterior density
The pair proposal is . Its Metropolis–Hastings acceptance probability simplifies to the displayed ratio, since the auxiliary densities cancel. Hence detailed balance holds for . Integrating over gives by unbiasedness, so the parameter marginal distribution is exactly the intended posterior distribution for any fixed positive . This is exactness of the stationary target, not independent exact samples or a guarantee of rapid convergence. Large relative variance of the likelihood estimates can produce long holding times, so marginalizing does not by itself guarantee an improvement in mixing time.
Use the extended target probability density function
The desired stationary distribution for positions is ; the full position-momentum joint probability distribution is , not alone. Gaussian momentum refreshment preserves , since it redraws its momentum marginal distribution independently while keeping the position fixed.
Let and be the surrogate leapfrog integration map. The proposal in the paper is : its notation includes the final momentum flip. Reversibility gives , so . Both and preserve volume. Thus is an involutive Metropolis proposal, and the required acceptance probability is
Equivalently, in terms of the surrogate Hamiltonian ,
Only evaluations of the true target probability density function are needed at endpoints; the trajectory uses the surrogate gradient. Target and surrogate normalizing constants cancel.
For , the accepted flux satisfies
Changing variables has unit absolute Jacobian determinant, so this identity proves detailed balance for accepted moves. The rejection mass stays at the same point and is also reversible. The accept/reject step therefore preserves . Since momentum refreshment also preserves , their composition and either phase of the alternating process preserve . Projecting onto positions proves the position stationary distribution is exactly the original target. Stationarity alone does not establish uniqueness or convergence from every start; those require additional irreducible Markov chain and aperiodic Markov chain hypotheses.
There is a genuine defect in the printed smoothing claim. For the positive, nondifferentiable Laplace distribution on the real line, direct minimization gives, for every ,
For the minimizer is ; for it is ; at zero both signs minimize. Hence normalization leaves , still nondifferentiable at zero. The displayed construction does not in general produce the asserted smooth surrogate. For example, the also nondifferentiable target makes that infimum for every , since the negative quartic term dominates the quadratic penalty. The invariance proof above is valid conditional on actually having a usable smooth surrogate, as the subsequent algorithm assumes.
A corrected sufficient construction is Moreau smoothing of a negative log-density. For a proper convex function with sequential lower semicontinuity , set
The Moreau envelope is differentiable, with ; require to be integrable to normalize it. The signs differ from the printed formula. For the Laplace distribution, this gives
a genuinely differentiable potential with integrable exponential tails. Using that surrogate with the boxed acceptance probability still targets the original Laplace distribution.
For this stellar population, let denote initial mass in units of the solar mass, and let denote present stellar age. Constant formation of equal numbers of stars per unit time makes the stellar age a uniform distribution on Gyr, so is uniform on . The normalized initial mass function has probability density function for , since . Consequently
These formulae have the stated age and mass domains; outside them the relevant cumulative fractions saturate at zero or one. The probability integral transform makes uniform on : gives . The time-independent initial mass function and constant number formation rate give independence of and . All fractions here count objects, including white dwarfs, using the stipulated stellar evolution law.
A red giant has
Thus the mass boundaries are and for positive . Equivalently, for the red giant region runs from to the age cap , and for it runs from to . There are no red giants with . Boundaries have zero probability and their endpoint convention does not affect the fractions. In the uniform square the red giant region is , with triangle vertices , and . The white dwarf region is the triangle .
Figure 1.
Mass-age regions and their uniform-coordinate images
. The right panel magnifies the evolved part of the unit square. All of the remaining region at larger X is main sequence. Details.
At fixed , define the conditional probabilities of a red giant, white dwarf and main sequence star by , and . Their interval widths are
Integration over the uniform distribution of age gives the individual-star fractions
The systems form a coeval binary population: the two binary star components have the same age. Their masses are independent, so has uniform probability density function on and their states have conditional independence given . Unconditional independence of their states would be incorrect: older systems make both evolved states more likely. The law of total probability now gives
The subtraction removes the double counting of systems containing two red giants.
For , and . Iterating this Gauss continued-fraction map on an irrational number never reaches zero, and its continued fraction is
The branches have inverse , for . The change of variables formula for a probability density function gives, almost everywhere,
The endpoints of these branches are a countable set of probability zero. Thus and have the same distribution; this Gauss measure is an invariant measure for the Gauss continued-fraction map.
For , the Markov property and the definition of transition intensity give
Summing over the possible pre-jump states gives the density of the expected Markov-chain entry count:
For a finite-state continuous-time Markov chain, the chance of two or more jumps in this short interval is . There is no contribution: diagonal entries of the infinitesimal generator encode exit rates, not jumps into the same state. With repeated entries possible, this density is an expected event rate rather than the probability density function of a single normalized random time.
For , the nearest neighbour distance is the th order statistic of . Each has continuous distribution function on . The probability integral transform makes independent and identically distributed random variables with the uniform distribution on . Monotonicity of gives , a uniform order statistic. Strict monotonicity is not needed at this stage.
For , counting observations below gives
Differentiating and cancelling consecutive terms gives the Beta distribution:
The Gamma function identity identifies this with the stated normalization. Ratios of Beta function integrals yield the moments
This is the probability content of a nearest neighbour ball identity, specialized to intervals on the real line.
Now write . The Lipschitz bound gives
Thus
Strict positivity of the probability density function makes continuous and strictly increasing from zero to one, so its inverse function is defined for every .
For the requested conditional calculation, suppose . At the preceding bound gives
Therefore , and substituting in the approximation proves the inverse probability content bound for a Lipschitz density:
The same argument works, without , on the meaningful local range with .
Since the nearest neighbour density estimator obeys , its conditional bound follows from the Beta distribution moments:
Hence, under the printed location condition,
The conditional inverse function bound also bounds , so the expected value in this calculation is finite under those formal hypotheses.
There is, however, a genuine flaw in the original PDF: the printed condition has no admissible points when is strictly positive on all of . A probability density function on all of cannot have Lipschitz constant zero, so . Its Lipschitz bound forces
Integrating this triangular lower envelope gives the Lipschitz density height bound
The inequality is strict because the probability density function has positive mass outside the finite interval. Thus everywhere. The last two conditional conclusions above are valid implications but vacuous as printed. Dropping strict positivity gives only , with equality forcing the entire probability density function to be the normalized triangular envelope. The local inverse bound does not justify replacing the printed condition by its reverse in the final expected value bound, since ranges over all of .
For independent and identically distributed random variables with distribution function , let be the sample maximum. In extreme value theory, belongs to the maximum domain of attraction of a nondegenerate distribution function when there exist and such that
at every continuity point of . Thus the normalized sample maxima have convergence in distribution to .
A positive measurable function has regular variation at infinity with index , written , if
The case defines a slowly varying function; equivalently a regularly varying function is with a slowly varying function.
Write for the survival function and for the right endpoint of a distribution. For , use the following standard extreme value theory normalization conventions:
They are respectively the Fréchet distribution, the negative Weibull distribution, and the Gumbel distribution for maxima. The negative Weibull distribution is supported to the left of its finite endpoint; it is not the usual positive Weibull distribution.
The necessary and sufficient maximum domain of attraction criteria are:
  • For the Fréchet distribution,
    Equivalently for every .
  • For the negative Weibull distribution,
    Equivalently as , for every .
  • For the Gumbel distribution,
    Here is a Gumbel auxiliary function, and can be finite or infinite. For each fixed , the shifted argument lies below eventually. This exponential tail-ratio condition, rather than regular variation with a fixed finite index, characterizes this case.
For the corresponding sufficient hazard function conditions, assume is absolutely continuous near , its survival function tends to zero there, and its probability density function is positive there. Write . The Von Mises conditions for extreme values give
The last condition additionally assumes that the reciprocal hazard function is continuously differentiable near the endpoint. It concerns the derivative of the reciprocal hazard function, not the derivative of the hazard function itself. Equivalently, when is differentiable, it is .
To see the connection, taking an integral of the hazard function gives
The first two limits make the logarithmic tail ratios over multiplicative or endpoint-distance intervals tend to or , giving the required regular variation. For the Gumbel distribution, set . The derivative condition makes uniformly on bounded -intervals, so
The endpoint assumptions ensure these shifts are admissible: at an infinite endpoint , and at a finite endpoint . These are sufficient smooth-tail conditions; the tail-ratio characterizations above are necessary and sufficient without differentiability assumptions. For a source using the same reciprocal-hazard convention, see Smith and Weissman's university lecture notes, equations (2.18)–(2.24).
For the Gamma distribution with integer shape and rate one, let . The base case is , and integration by parts gives
Induction proves
This survival function has the asymptotic form ; for the relation is exact. Put . For every fixed real , and . The proposed centering satisfies
so
Only the highest-degree term contributes to the limit, with the same reasoning covering . The Taylor expansion now gives
because . Therefore the Gumbel limit for gamma maxima is
Equivalently, has convergence in distribution to the standard Gumbel distribution, with scale .
Suppose the inclined population samples its arguments of periapsis and nodal orientations, or that sufficiently slow precession supplies this sampling. Near the planet the vertical amplitude is . For uniformly distributed vertical phase, has probability density function , so the midplane volume probability is .
In the genuinely three-dimensional regime , the geometric collision cross-section is . The relative velocity differs from the coplanar value only at order for the specified high-speed encounters. Therefore the inclined comet collision time is
At fixed the leading scaling is , with no leading dependence when . The displayed factor retains the first dependence on . Also in this averaging convention.
Inclination alone does not specify a collision time for one fixed orbit. If its orbital nodes miss the planet's orbit, a tilted Kepler orbit can have no collisions at all. The finite result is an orientation-averaged population result, or requires precession or scattering to sample those orientations. When , the three-dimensional estimate crosses over to the coplanar result instead of tending to zero.
Integrating the continuity equation over state space and applying the divergence theorem gives
Assume the probability density function and its derivatives decay sufficiently rapidly at infinity, or impose a zero-flux boundary condition. The boundary integral then vanishes, so is constant. Since it is initially one, the distribution remains normalized for all time.
Write , the absolute Jacobian determinant. Under the usual meaning of a map on domains, the injective area formula gives, for any nonnegative measurable ,
This version of the change of variables formula does not require to be nonsingular everywhere. A map is a locally Lipschitz function, and the formula applies by exhaustion if the domains are unbounded. For non-open , one needs the corresponding extension or regularity interpretation of the stated assumption.
First suppose almost everywhere with respect to Lebesgue measure. For each bounded nonnegative measurable ,
Taking indicator functions proves .
Conversely, the area formula shows that the measure with probability density function has total mass one and satisfies . The assumed equality then implies , since the measurable bijection has a measurable inverse on these domains. Equivalently, apply both pushforward measure identities to for each Borel set . Uniqueness of the Radon-Nikodym derivative yields
Define ; the printed quotient is undefined at the origin, but this continuous extension changes no transport cost because . In polar coordinates, the two probability density functions give
Both angular distributions are uniform, with independence of angle and radius. The proposed transport map preserves the angle and sends to . Therefore
and preservation of the angle proves . As a local check using the Jacobian determinant, its radial and tangential derivatives for are and , respectively, so and .
Now use the convex function
It is convex because the Euclidean norm is convex and is increasing and convex on . It is differentiable, including at zero, and
Its graph transport plan lies in the graph of , so the Knott–Smith optimality criterion proves quadratic optimality. Since that plan is induced by a map, part 1(b) proves optimality for the Monge optimal transport problem as well. The minimum provides a useful independent check:
Let be the probability density function of the unknown and that of the noise, both relative to Lebesgue measure. By independence of random variables, the conditional distribution of the measurement given has density . Thus its likelihood function and the joint density are
The Bayesian model evidence is the marginal data density
Bayes' theorem gives the posterior density solving this Bayesian inverse problem:
This formula applies when . Since by Tonelli theorem, these conditions hold almost everywhere under the marginal data law. Values of a conditional distribution at exceptional data are not determined by the joint law. The posterior mean or a maximum a posteriori estimate can then give a point estimate, while the full posterior density describes the remaining uncertainty.
Use the Hellinger distance convention in the original PDF:
The TeX transcription incorrectly places the integral outside the outer square root in the definition; the displayed formula for the squared distance and the original PDF both use the convention above. Although Lebesgue measure is not a probability measure, we may use it as the common dominating measure: the integral is unchanged by replacing it with any equivalent dominating probability measure.
For unit-variance normal distributions, put and . The product of the square roots of their probability density functions is
The final factor integrates to one, so
For the bound at arbitrary , use the given formula for the Hellinger distance between normal distributions and write
Since for ,
Therefore one valid choice is . In particular this bound is uniform when the two standard deviations stay bounded away from zero.
Use the shape-rate convention for the Gamma distribution. Define
The zero null space assumption makes a positive-definite matrix, and is the unique ordinary least squares solution. The Gaussian likelihood contributes , which must be retained because is unknown. Multiplying it by the improper prior and the Gamma distribution gives the normal-gamma posterior with a flat prior:
Here the constant density of the improper prior cancels on normalization. Put
The normal equation makes the residual orthogonal to the range of . The orthogonal decomposition by a closed subspace into that range and its orthogonal complement gives
Integration over by the Gaussian integral, followed by integration over using the Gamma distribution, yields
Thus the posterior distribution is proper, despite the improper prior; here follows from the zero null space assumption.
Completing the square in and collecting the powers of give the two conditional distributions:
Explicitly their probability density functions are
For a joint maximum a posteriori estimate relative to , first fix . The unique maximizing is . The remaining log density is plus a constant. Consequently, if ,
The stated assumptions do not always imply . If , the supremum occurs as but is not attained on ; if , the density is unbounded there. In either case there is no joint maximum a posteriori estimate in the stated parameter space. Calling zero a boundary mode does not make it an admissible positive precision.
If “MAP estimators” instead means separate modes of the marginal distributions, integration gives
Thus the marginal maximum a posteriori estimate of is always , and that of is when , with the same nonattainment at zero when . Joint and marginal maximization are different operations; the usual joint interpretation gives the preceding boxed pair.
The displayed density is that of a gamma distribution with shape and rate one. Put
By the additivity of independent gamma distributions with a common rate,
Equivalently, the independence of random variables makes the Laplace transform of a nonnegative random variable of the sum equal
Therefore the required probability density function is