Metropolis-Hastings Sampler for Hyperparameter Estimation

• Metropolis-Hastings (MH) sampling for hyperparameter estimation is a core MCMC method that samples complex posterior distributions in Bayesian inference.
• A key technique involves adaptively adjusting the MH proposal scale using the Robbins-Monro process to target optimal acceptance probabilities for improved efficiency.
• This adaptive scaling method offers robust, hands-off tuning applicable to both univariate and multivariate random-walk proposals in various Bayesian models and automated software.

The Metropolis-Hastings (MH) sampler is a core algorithm within Markov Chain Monte Carlo (MCMC) methods for performing Bayesian inference and, in particular, hyperparameter estimation. In the context of hyperparameter estimation, one often faces the task of efficiently sampling from a posterior distribution over hyperparameters that may be high-dimensional, multi-modal, or defined implicitly through complex likelihoods. The foundational work of Garthwaite, Fan, and Sisson presents a robust, theoretically backed, and practically validated adaptive scheme for tuning the scaling parameter of random-walk MH samplers using the Robbins-Monro (RM) process. This approach contributes automatic, fast, and stable tuning of the proposal scale, directly impacting efficiency and reliability in hyperparameter inference tasks.

1. Adaptive Scaling via Robbins-Monro Stochastic Approximation

The algorithm adapts the proposal scale $\sigma$ in a random-walk MH sampler to target a desired overall acceptance probability (OAP), denoted $p^*$. The Robbins-Monro process is employed as a sequential root-finding scheme under stochastic observations to implicitly solve for the optimal scale $\sigma^*$ such that the acceptance rate $p(\sigma^*) = p^*$. At each iteration:

• If the proposed move is accepted, $\sigma_i$ is increased.
• If rejected, $\sigma_i$ is decreased. The update is defined by: $$\sigma_{i+1} = \begin{cases} \sigma_i + \dfrac{c (1 - p^*)}{i}, & \text{if accepted} \ \sigma_i - \dfrac{c p^*}{i}, & \text{if rejected} \end{cases}$$ where $c > 0$ is a problem-dependent steplength constant. The adaptation step-size is proportional to $1/i$, ensuring diminishing adaptation and ergodicity of the chain.

2. Optimal Steplength Constant Estimation

The effectiveness of the Robbins-Monro search is controlled by the selection of $c$, optimally set to: $$c^* = -\frac{1}{\left. \frac{d p(\sigma)}{d \sigma} \right|_{\sigma = \sigma^*}}$$ Direct computation is generally intractable because $p(\sigma)$ is typically unknown and difficult to differentiate analytically. The paper introduces robust empirical estimators for $c^*$:

• Univariate normal proposals:

$c^* \approx \frac{\sigma^*}{p^* (1-p^*)}$

• Multivariate normal proposals (dimension $m$):

$$c^* \approx \sigma^* \left\{ \left(1 - \frac{1}{m^*}\right) \frac{(2\pi)^{1/2} e^{\alpha^2/2}}{2\alpha} + \frac{1}{m^* p^* (1-p^*)} \right\}$$

where $\alpha = -\Phi^{-1}(p^*/2)$ and typically $m^* = m$, with adjustments for heavy-tailed cases. These estimators are empirically shown to be stable across a wide array of distributions, including heavy-tailed and multimodal targets.

3. Proposal Distributions and Application Range

The method is specialized for random-walk proposals:

• Univariate: $N(x, \sigma^2)$
• Multivariate: $N(\mathbf{x}, \sigma^2 \Sigma)$

Where $\Sigma$ is an estimated (or fixed) covariance matrix, ideally proportional to the posterior’s covariance. The approach is robust to the target distribution, including non-Gaussian and bounded distributions, heavy tails, and mixtures, provided the proposal is normal. For heavy-tailed targets (e.g., low-degree $t$, Cauchy), only $m^*$ (in multivariate settings) may require adjustment.

4. Algorithmic Implementation

Univariate case:

1. Initialize $\sigma_1$, typically with $p^* = 0.44$.
2. For $i \geq n_0$ (with $n_0 \approx 5/[p^*(1-p^*)]$), propose $y \sim N(x, \sigma_i^2)$, accept/reject, and update $\sigma$ as above, using $c = \sigma_i / [p^* (1-p^*)]$.
3. Restart adaptation if $\sigma$ increases/decreases by more than a factor of 3, or after a fixed number of iterations.

Multivariate case:

1. Initialize $\sigma_1$, covariance $\Sigma_m$ (identity or based on prior), often with $p^* = 0.234$.
2. For $i > 200$, propose $$\mathbf{y} \sim N(\mathbf{x},\,\sigma_i^2 \Sigma_i)$$, accept/reject, and update $\sigma$ using the estimator for $c$ from (18).
3. Regularly re-estimate empirical covariance and regularize as necessary.

For Metropolis-within-Gibbs or block updating, separate Robbins-Monro searches are conducted per component or block.

5. Demonstrated Performance: Simulation and Real Data

Simulated Data

A wide spectrum of univariate target distributions—Gaussian, $t$, Cauchy, double exponential, mixtures, Gamma, Beta—was evaluated over hundreds of replicates.

• The final estimated $\hat{\sigma}^*$ closely matches theoretical optima.
• Acceptance rates tightly concentrate around the target value even for challenging distributions.

Multivariate experiments (e.g., 50-dimensional $MVN$ targets) showed:

• Rapid and robust convergence to optimal scaling, even with ill-conditioned or misspecified covariance estimates.
• Efficiency and mixing comparable to cases where the true optimal scaling is known.

Practical Application

In a real-data logistic additive mixed model with 306 parameters, the Robbins-Monro adaptation yielded:

• Effective, reproducible scaling for a large and diverse set of parameter blocks.
• Consistently optimal acceptance rates and improved mixing compared to non-adaptive manual tuning.

6. Integration in Adaptive MCMC and Hyperparameter Estimation

The Robbins-Monro scaling approach is modular and highly efficient for integration in:

• Adaptive/blockwise/regionally adaptive MCMC methods, providing local adaptation for different blocks or regions of the posterior.
• Automated MCMC software, enabling hands-off robust scaling over a range of models.
• Hyperparameter estimation, especially for proposal scale and variance hyperparameters in hierarchical models, where optimal exploration is crucial and manual tuning impractical.

The adaptability and theoretical guarantees make the method especially suitable for modern automated Bayesian inference platforms, as well as practical high-dimensional or complex hierarchical problems that require hyperparameter adaptation during sampling. By embedding Robbins-Monro scaling searches, MCMC routines can self-calibrate, minimize pilot tuning, and react quickly to the local geometry of the target distribution.

Summary Table: Adaptive Robbins-Monro Scaling

Aspect	Key Details / Steps
Update rule	$\sigma_{i+1} = \sigma_i + \frac{c(1-p^*)}{i}$ (accept); $-\frac{c p^*}{i}$ (reject)
Steplength estimator (univariate)	$c = \sigma^*/[p^* (1-p^*)]$
Steplength estimator (multivariate)	See equation (18) / (14); depends on dimension, $\alpha = -\Phi^{-1}(p^*/2)$
Proposal types	Random walk $N(x, \sigma^2)$ and $N(\mathbf{x}, \sigma^2 \Sigma)$
Target acceptance rates	$p^* = 0.44$ (univariate), $p^* = 0.234$ (multivariate)
Adaptation criteria	Diminishing ($\propto 1/i$); restart if scale shifts too much or fails to converge
Strengths	Rapid, hands-off adaptation; robust to model, dimension, covariance; minimal tuning
Hyperparameter role	Scales, block variances, and other tuning constants estimated online

Conclusion

Adaptive optimal scaling of the MH algorithm using the Robbins-Monro process provides a robust, practical solution for hyperparameter estimation in Bayesian models. It eliminates manual scale adjustments, ensures near-optimal sampler efficiency, and is universally applicable across univariate and multivariate settings—integration into modern MCMC-based Bayesian inference toolkits fundamentally enhances both usability and computational performance for a broad range of statistical and machine learning applications (1006.3690).