Binh Ho

Notes

Expository notes on probability, statistical theory, and computation, written mostly to understand things properly. Corrections are welcome.

Probability

  • Margins and dependence, separated: copulas and Sklar's theorem

    Any joint distribution with continuous margins factors uniquely into its marginal distributions and a copula, the joint law of the coordinates after each is transformed to be uniform. Sklar’s theorem is proved through the quantile transform, followed by what a copula does and does not determine.

  • Rare events that never stop: almost sure convergence versus convergence in probability

    Independent indicators with success probabilities \(1/n\) converge to zero in probability but almost surely fail to converge. The Borel–Cantelli lemmas explain exactly where the line between the two modes lies.

  • Moving mass at least cost: the Monge and Kantorovich problems of optimal transport

    Monge asked for the cheapest map carrying one distribution onto another; Kantorovich relaxed the map to a coupling and obtained a linear problem that always has a solution. The note proves existence, treats the discrete and one-dimensional cases, states Brenier’s theorem for quadratic cost, and derives the Sinkhorn algorithm for the entropically regularized problem.

  • Gaussians with a random variance: the Laplace and Student-t distributions as normal scale mixtures

    If \(Z\) is standard normal and \(V\) an independent positive variable, then \(\sqrt V Z\) is a normal scale mixture. Exponential \(V\) gives the Laplace distribution and inverse-gamma \(V\) the Student t. The representation explains their heavy tails and turns robust estimation into iteratively reweighted least squares.

  • The least committal distribution: maximum entropy under moment constraints

    Among all distributions with prescribed expectations of a few statistics, the one of largest entropy is an exponential family member, whenever such a member satisfies the constraints. We prove this with Gibbs’ inequality and give a case in which the supremum is not attained.

  • Slicing the range instead of the domain: a picture of the Lebesgue integral

    Riemann sums cut the horizontal axis into intervals; Lebesgue sums cut the vertical axis into levels and measure the set where the function exceeds each level. Made precise, the picture says that the Lebesgue integral of a nonnegative function is an ordinary integral of the decreasing function \(t \mapsto \mu(f \geq t)\).

  • Thin shells and orthogonal directions: where volume and probability live in high dimensions

    In high dimension the volume of a ball sits next to its boundary, a standard Gaussian vector has length close to \(\sqrt d\) although its density peaks at the origin, independent random directions are nearly orthogonal, and pairwise distances between random points become almost equal. One exponential inequality explains most of this.

  • Perturb and maximize: sampling a categorical distribution with Gumbel noise

    Adding independent Gumbel noise to unnormalized log-probabilities and taking the argmax gives an exact draw from the corresponding categorical distribution. The proof also shows that the maximum and the argmax are independent; the same computation, read through races of exponential clocks, gives sampling without replacement and the Gumbel-softmax relaxation.

  • Square roots of a covariance matrix: Cholesky, eigen and symmetric factors of the multivariate normal

    A multivariate normal vector can be written \(\mu + AZ\) for any matrix with \(AA^\top = \Sigma\). All such factors differ by an orthogonal matrix and give the same law, but they couple \(X\) to \(Z\) differently: the Cholesky factor encodes sequential regressions, the eigen factor the principal axes, and the SVD of any factor recovers the eigendecomposition of \(\Sigma\).

  • Tail bounds from moments and moment generating functions

    Markov, Chebyshev and higher-moment bounds, the Chernoff method, Hoeffding’s lemma and inequality, Bernstein’s inequality, and the sub-Gaussian and sub-exponential classes, with complete proofs and a comparison on binomial tails.

  • Monotone convergence, Fatou's lemma and dominated convergence

    The three limit theorems of Lebesgue integration with complete proofs, examples showing why each hypothesis is needed, and the extension of dominated convergence to sequences that converge only in measure.

  • Boole's inequality and the union bound

    The probability that at least one of countably many events occurs is at most the sum of their probabilities, whatever the dependence between them. Three proofs, the second-order Bonferroni bound, and two standard uses: multiple testing and maxima of Gaussian variables.

Stochastic processes and time series

  • Fair games and when to stop them: martingales and the optional stopping theorem

    A martingale is a process whose best forecast of the future is its present value. We prove that no betting system and no bounded stopping rule can change its mean, see exactly which conditions are needed when the stopping time is unbounded, and use the result for the gambler’s ruin and for sequential tests that remain valid under continuous monitoring.

  • The chain rule with a correction: a derivation of Itô's formula

    For a function of Brownian motion the second-order Taylor term does not vanish, because squared increments add up to elapsed time. Itô’s formula follows from this observation; a complete proof is given for twice continuously differentiable functions with bounded derivatives, followed by applications to geometric Brownian motion and the Ornstein–Uhlenbeck process.

  • From Brownian motion to Itô processes: drift, diffusion and quadratic variation

    The Wiener process, Brownian motion with drift and general Itô processes form a nested sequence of models. Quadratic variation separates their two ingredients: on a fixed time interval the diffusion coefficient is recovered exactly from a path, while the drift is not.

  • Autocovariance, autocorrelation and weak stationarity

    The second-order description of a time series: autocovariance and autocorrelation functions, weak stationarity, the examples of white noise, AR(1) and MA(1), and the sampling behaviour of the sample autocorrelation function.

  • The Kalman filter as repeated Gaussian conditioning

    For a linear Gaussian state-space model, the filtering distributions are Gaussian and their means and covariances obey the Kalman recursions. The derivation needs only two facts about Gaussian vectors.

  • Bootstrapping intermittent demand with a Markov chain

    The bootstrap of Willemain, Smart and Schwarz forecasts the distribution of intermittent demand over a lead time by combining a two-state Markov chain for demand occurrence, resampled demand sizes and a jittering step. The note gives an R implementation, compares the method with the ordinary bootstrap on simulated data, and derives the maximum likelihood estimator of a transition matrix.

Statistical theory

Bayesian statistics

  • Zero-forcing: how the reverse Kullback–Leibler divergence locks onto one mode

    Fitting a Gaussian \(q\) by minimizing \(\mathrm{KL}(q\Vert p)\) penalizes mass that \(q\) places where \(p\) is small. For a two-component mixture this makes the optimal fit sit on one component once the components are well separated, and for a correlated Gaussian target it makes mean-field variances equal to conditional rather than marginal variances.

  • Does the posterior find the truth? Prior mass, tests and Schwartz's theorem

    Posterior consistency asks whether, under a fixed true distribution, the posterior concentrates on arbitrarily small neighbourhoods of it. Doob’s theorem answers yes for almost every truth under the prior; Schwartz’s theorem gives checkable conditions at a given truth: enough prior mass near it in Kullback–Leibler divergence, and tests that separate it from distant alternatives.

  • Objections to Bayesian inference, and the replies they have drawn

    An essay on the main published lines of criticism of Bayesian statistics: subjectivity of the prior, the impossibility of representing ignorance, frequentist performance, the likelihood principle and stopping rules, model checking, and testing. Each is stated in its strongest form, together with the principal responses and what remains unresolved.

  • Sampling a posterior with a Markov chain: Metropolis–Hastings, Gibbs and what makes them work

    A posterior known only up to its normalizing constant can still be sampled, by running a Markov chain that leaves it invariant. We derive the Metropolis–Hastings and Gibbs kernels from detailed balance, state the convergence results that justify ergodic averages, and look at how the proposal scale controls efficiency.

  • The same p-value at every sample size: Lindley's paradox and the point null

    Hold the z-statistic of a two-sided normal test fixed at 1.96 and let \(n\) grow. The p-value stays at 0.05 while the posterior probability of the point null tends to one. The calculation is elementary; what it says about the two notions of evidence is less so.

  • Priors with infinite mass: when an improper prior still gives a proper posterior

    A flat prior on the real line is not a probability distribution, yet it often produces a perfectly good posterior. The posterior exists exactly when the marginal density \(m(x)\) is finite, and this condition can fail in ordinary models, most notably for the scale parameter of a hierarchical model.

  • Updating beliefs with a loss instead of a likelihood: Gibbs posteriors and the learning rate

    A Gibbs posterior replaces the log-likelihood in Bayes’ formula by a scaled empirical loss, so that the target is defined by the loss rather than by a model. It solves a variational problem that explains its form, but its spread depends on a learning rate \(\omega\) that the formula does not fix. For the median under absolute-error loss the calibrated value is \(\omega = 2f(m)\).

  • Learning the prior from many parallel problems: empirical Bayes for normal means

    With a thousand noisy measurements of a thousand different means, the measurements themselves reveal the distribution the means were drawn from. Estimating that prior and plugging it into Bayes’ rule leads to the James–Stein estimator, which beats the obvious estimator at every parameter value, and, without a parametric prior, to Tweedie’s formula.

  • A prior on probability distributions: three constructions of the Dirichlet process

    The Dirichlet process can be defined through its finite-dimensional Dirichlet marginals, built explicitly by breaking a stick, or generated sequentially by a Pólya urn. Each description makes a different property evident: the moments and conjugacy, the almost sure discreteness, and the logarithmic growth of the number of distinct values.

  • Mixtures without a fixed number of components: Dirichlet process mixtures and collapsed Gibbs sampling

    Placing a Dirichlet process prior on the mixing distribution of a mixture model gives a prior over partitions of the data with an unbounded number of clusters. With a conjugate base measure the mixture parameters integrate out, and a Gibbs sampler on cluster labels needs only counts and posterior predictive densities. The posterior on the number of clusters, however, should be read with care.

  • A chi-square upper bound on the evidence, and the sandwich it makes with the ELBO

    The quantity \(\tfrac1n \log \E_q[(p(y,z)/q(z))^n]\) is an upper bound on \(\log p(y)\) for every \(n > 1\); for \(n = 2\) its gap is \(\tfrac12 \log(1 + \chi^2(\pi\Vert q))\). Minimizing it gives mass-covering approximations and minimum-variance importance proposals, but it is governed by tails and is hard to estimate.

  • Uncertainty about the model itself: Bayesian model averaging for linear regression

    When several models are plausible, Bayes’ theorem applies to the model index as to any other unknown. We derive the closed-form posterior model probabilities under Zellner’s g-prior, average over all \(2^6\) subsets of six predictors, and look at what averaging does to a coefficient’s posterior and to predictive accuracy.

  • Bayesian inference for the normal model with conjugate priors

    Conjugate posteriors for a normal sample in three cases: unknown mean with known variance, unknown variance with known mean, and both unknown, where the marginal posterior of the mean is a Student \(t\) distribution.

  • Maximum a posteriori estimation

    The maximum a posteriori estimate maximizes the log-likelihood plus the log prior density. For a coin with a Beta prior it has a closed form that shrinks the sample proportion towards the prior mode; it also depends on the parameterization, unlike the maximum likelihood estimate.

  • Variational inference and the evidence lower bound

    Variational inference replaces posterior integration by optimization: it maximizes a lower bound on the log evidence over a tractable family of distributions. Two derivations of the bound, the coordinate ascent updates for mean-field families, and a Gaussian example in which the gap can be computed exactly.

Computation

  • Two Hessians for one likelihood: Newton's method, Fisher scoring and IRLS in generalized linear models

    Fisher scoring replaces the observed information in Newton’s method by its expectation. For a generalized linear model this turns each step into a weighted least-squares fit, guarantees an ascent direction, and loses nothing when the link is canonical.

  • Stochastic variational inference: natural-gradient steps from minibatches

    In conditionally conjugate models the natural gradient of the ELBO with respect to the global variational parameter is the coordinate-ascent target minus the current value. Replacing the full-data target by an unbiased minibatch estimate gives stochastic variational inference, whose cost per step does not grow with the number of observations.

  • Slice sampling: uniform draws under the graph of a density

    Sampling from a density is the same as sampling uniformly from the region under its graph. Slice sampling turns this into a two-step Gibbs sampler and makes it practical with Neal’s stepping-out and shrinkage procedures, which preserve the uniform distribution on each slice.

  • Differentiating an expectation through its density: the score-function gradient estimator

    The identity \(\nabla_\theta \E_\theta f(X) = \E_\theta[f(X) \nabla_\theta \log p_\theta(X)]\) turns the gradient of an expectation into an expectation that can be estimated by Monte Carlo. Sufficient conditions for the identity, its failure when the support moves, and exact variances of the estimator with and without baselines, compared with the pathwise estimator.

  • Conditioning out the noise: Rao–Blackwellization of Monte Carlo estimators

    Replacing a simulated quantity by its conditional expectation given part of the simulation never increases variance for independent draws. For a two-block Gibbs sampler the gain can be computed exactly: the asymptotic variance drops by at least the variance of a single draw.

  • Power iteration and the eigengap: what sets the rate

    Repeated multiplication by a symmetric matrix aligns any generic starting vector with the leading eigenvector at the geometric rate \(\lvert\lambda_2/\lambda_1\rvert^k\), while the Rayleigh quotient converges to the eigenvalue at the squared rate. Proofs, failure modes, and the shifts and subspace variants that improve on the basic method.

  • LU, Cholesky, QR and the SVD: what each factorization buys

    Four matrix factorizations, each reducing a problem with a general matrix to problems with triangular, orthogonal or diagonal factors. Existence and uniqueness proofs, what each one costs, and why least squares should be solved with QR rather than with the normal equations.

  • Just-in-time compilation in JAX: tracing, fusion and where the time goes

    An interpreted array program pays a fixed cost for every operation and a full pass over memory for every intermediate array. Just-in-time compilation in JAX removes both by tracing a function into a small intermediate language and compiling it with XLA. The mechanism explains both the restrictions JAX imposes and the timings measured here for an array computation and a sequential recursion.

  • Hamiltonian Monte Carlo: momentum, the leapfrog map and an exact accept–reject step

    Hamiltonian Monte Carlo moves far through the parameter space by simulating a physical system and corrects the simulation error with a single Metropolis step. The leapfrog map composed with a momentum flip is a volume-preserving involution, which makes the method exact; the Gaussian case, solved in closed form, explains the step-size limit \(\varepsilon < 2\sigma_{\min}\) and how the integration time determines the autocorrelation.

  • One posterior, four computations: grid, Laplace, Metropolis and variational inference

    A two-parameter logistic regression is small enough to compute its posterior exactly on a grid, which makes it a clean test case for the three approximations used in larger problems. The Laplace approximation and a random-walk Metropolis sampler recover it well; mean-field variational inference gets the centre right and the spread wrong, for a reason that can be proved exactly.

  • Expectation–maximization as coordinate ascent on the evidence lower bound

    Both steps of EM maximize the same functional \(F(q,\theta) = \log p(x\mid\theta) - \mathrm{KL}(q\Vert p(\cdot\mid x,\theta))\), one coordinate at a time. Read as mean-field variational inference over latent variables and parameters, EM is the limit in which the variational factor for the parameters collapses to a point mass.

  • What is a correlated sample worth? Effective sample size and the integrated autocorrelation time

    For a stationary sequence the variance of the sample mean behaves like \(\gamma_0\tau/n\), where \(\tau = 1 + 2\sum_k \rho_k\) is the integrated autocorrelation time, so \(n\) correlated draws are worth \(n/\tau\) independent ones. We derive this, compute it exactly for an autoregression, and explain why estimating \(\tau\) requires truncating the sum.

  • Control variates: a regression adjustment inside a Monte Carlo average

    A Monte Carlo average can be corrected by any quantity with known mean that moves together with the integrand. The best linear correction is a least-squares regression, its variance is the residual variance, and estimating the coefficient from the same draws costs nothing to first order.

  • Conjugate gradients: conjugate directions, Krylov subspaces and the square-root speed-up

    For a quadratic with positive definite Hessian \(A\), conjugate gradients minimizes the error over a growing Krylov subspace using one matrix–vector product per step. This optimality gives finite termination, a convergence rate governed by \(\sqrt{\kappa}\) instead of \(\kappa\), and fast convergence for clustered spectra.

  • From two uniforms to two Gaussians: the Box–Muller transform and the polar method

    A pair of independent uniforms can be turned into a pair of independent standard normals by reading them as a radius and an angle. We prove that the Box–Muller transform and Marsaglia’s polar method are exact, and look at what each costs.

  • BFGS as a least-change secant update: derivation, positive definiteness and behaviour on quadratics

    The BFGS method replaces the inverse Hessian in Newton’s method by a matrix updated from gradient differences. The update is the smallest change, in a weighted Frobenius norm, that satisfies the secant equation; it preserves positive definiteness under the curvature condition that a Wolfe line search guarantees, and on quadratics it recovers the inverse Hessian in \(n\) steps.

  • Stochastic EM for Gaussian mixtures

    Stochastic EM replaces the E-step by a draw of the latent labels from their posterior. For Gaussian mixtures this gives a simple Markov chain on the parameter space that escapes slow regions of EM but fluctuates around, and does not converge to, a maximum likelihood estimate.

  • EM for Gaussian mixture models

    The EM algorithm for a finite mixture of multivariate normal distributions: responsibilities in the E-step, closed-form M-step updates with complete proofs, the monotonicity of the likelihood, and the unboundedness that makes the maximum likelihood estimator a local notion.

  • Sampling from the multivariate normal distribution

    Draws from \(N_d(\mu, \Sigma)\) are affine images \(\mu + Az\) of standard normal vectors, for any factor with \(AA^\top = \Sigma\). The Cholesky factor makes this cheap and also gives the log-density; when the precision matrix is the natural input, the correct map is \(\mu + L_Q^{-\top} z\).

  • Inverse transform sampling

    If \(U\) is uniform on \((0,1)\) and \(F^{-1}\) is the generalized inverse of a distribution function \(F\), then \(F^{-1}(U)\) has distribution function \(F\). One proof covers continuous, discrete and mixed distributions alike.

  • The expectation–maximization algorithm

    EM maximizes a likelihood with missing or latent data by alternating between computing the expected complete-data log-likelihood and maximizing it. A derivation through a lower bound, the monotonicity proof, what is and is not known about convergence, and a worked Poisson mixture.

Machine learning

  • How many labelings can a classifier family produce? VC dimension and the Sauer–Shelah lemma

    A class of binary classifiers restricted to \(n\) points produces at most \(2^n\) labelings. The VC dimension is the largest n for which all of them occur, and the Sauer–Shelah lemma shows that beyond it the count grows only polynomially, which is what makes uniform laws of large numbers over the class possible.

  • Inducing points as variational parameters: a lower bound for sparse Gaussian processes

    Keeping the exact Gaussian process model and approximating its posterior through inducing values gives a lower bound on the log marginal likelihood, \(\log\mathcal N(y\mid 0, Q_{ff}+\sigma^2 I) - \mathrm{tr}(K_{ff}-Q_{ff})/(2\sigma^2)\). The inducing inputs become variational parameters: adding them can only tighten the bound, and optimizing them cannot overfit the model.

  • The representer theorem: why kernel ridge regression lives in the span of the data

    Penalized empirical risk minimization over a reproducing kernel Hilbert space, an infinite-dimensional problem, has a solution in the \(n\)-dimensional span of the kernel sections at the data. For squared loss this gives kernel ridge regression in closed form, a primal–dual identity for ordinary ridge regression, and an exact leave-one-out shortcut.

  • Is the random seed a tuning parameter? Seed variability and the optimism of the best seed

    A fixed seed makes a computation reproducible, but the result is still a random function of the seed. Treating the computation as a randomized procedure separates seed variability from sampling variability, explains why averaging over seeds helps, and quantifies how much reporting the best of \(k\) seeds overstates performance.

  • Random Fourier features: approximating a kernel by sampling its spectral measure

    By Bochner’s theorem a continuous shift-invariant positive definite kernel, normalized to equal one on the diagonal, is the characteristic function of a probability measure. Sampling frequencies from that measure gives an explicit finite-dimensional feature map whose inner products are unbiased estimates of the kernel, with error of order \(D^{-1/2}\) in the number of features.

  • Normalizing flows: exact densities through invertible maps

    A normalizing flow models a density as the image of a simple base distribution under an invertible map, so the likelihood is available in closed form through the change-of-variables formula. Triangular Jacobians make that formula cheap; maximum likelihood fits the flow by minimizing a Kullback–Leibler divergence that can be computed on either side of the map.

  • Natural gradient descent: steepest descent measured by the Fisher information

    The direction of steepest descent depends on how step lengths are measured. Measuring them by Kullback–Leibler divergence between the corresponding distributions gives the natural gradient \(F(\theta)^{-1}\nabla L(\theta)\), which does not depend on the parametrization, coincides with Fisher scoring, and for exponential families reduces to a plain gradient in the mean parameters.

  • Linear dimensionality reduction as optimization over orthonormal frames

    Following Cunningham and Ghahramani (2015), many linear dimensionality reduction methods can be written as the maximization of an objective over matrices with orthonormal columns. PCA is the case with an eigenvector solution; a Fisher-type trace-ratio objective is not, and a generic Riemannian gradient method on the Stiefel manifold solves it where the usual eigenvector heuristic does not.

  • Inducing points: Gaussian process prediction through a low-rank summary

    Exact Gaussian process regression costs \(O(n^3)\). Routing all information about the latent function through its values at \(m \ll n\) inducing inputs gives the subset-of-regressors, DTC and FITC approximations, which cost \(O(nm^2)\). They share a predictive mean but differ in how they treat variance, which matters most away from the inducing inputs.

  • Gaussian process regression as conditioning a Gaussian vector

    A Gaussian process prior on a regression function turns every prediction problem into the conditioning of a finite Gaussian vector. The posterior is again a Gaussian process, its mean is a kernel expansion over the data, its covariance does not depend on the responses, and the marginal likelihood balances fit against complexity when choosing hyperparameters.

  • AdaBoost as stagewise minimization of the exponential loss

    AdaBoost’s reweighting rule and its step size \(\alpha_m = \tfrac12\log\{(1-\varepsilon_m)/\varepsilon_m\}\) are exactly a greedy coordinate step on the empirical exponential loss. The same identity gives the exponential training-error bound and a bound on the fraction of small margins.

Quantitative finance