1. Motivation: Sequential experimental design with expensive, gradient-free black-box functions is ubiquitous in computational statistics, biomedical research, and engineering. Evaluation budgets are severely constrained (each CFD simulation takes hours, each clinical trial costs millions), making efficient information extraction from each observation critical.

2. Mathematical setup: Let $f: \mathcal{X} \rightarrow \mathbb{R}$ be an unknown expensive objective function where $\mathcal{X} \subseteq \mathbb{R}^d$ is the design space. At time $t$, we have observations $\mathcal{D}_t = {(x_i, y_i)}_{i=1}^t$ where $y_i = f(x_i) + \epsilon_i$ and $\epsilon_i \sim \mathcal{N}(0, \sigma^2)$. We model $f$ with a Gaussian process prior:

   \[f(x) \sim \mathcal{GP}(\mu_0(x), k(x,x'))\]

   yielding posterior mean $\mu_t(x)$ and variance $\sigma_t^2(x)$. We have $K$ parallel agents and evaluation budget $N$.

   Assumptions:

   1. Evaluations are expensive (minutes to hours per query)
   2. Function evaluations are noisy but unbiased
   3. Agents can operate asynchronously with bounded delays $\delta_t^j \leq \tau_{\max}$
   4. Communication costs scale as $O(\log K)$ or $O(K)$

3. Toy example: When $d=2$, $\mathcal{X} = [0,1]^2$, $K=4$ agents, and $N=20$ budget, this reduces to selecting 20 points in the unit square across 4 parallel workers to minimize cumulative regret while handling asynchronous feedback delays.

4. Formal objective: Minimize cumulative regret over $T$ rounds:

   \[R_T = \sum_{t=1}^T \left( f(x^*) - f(x_t) \right)\]

   where $x^* = \arg\max_{x \in \mathcal{X}} f(x)$ is the global optimum.

ALMAB-DC combines three components in a unified framework:

1. Active Learning: A Gaussian process surrogate selects informative candidates via acquisition functions:

   \[x_{t+1} = \arg\max_{x \in \mathcal{X}} U_t(x; \mu_t(x), \sigma_t^2(x))\]

   where $U_t$ is UCB, expected improvement, or max-variance depending on task.

2. Multi-Armed Bandit Controller: Each candidate is treated as a bandit arm. UCB selection:

   \[a_t = \arg\max_{i \in \mathcal{A}} \left[ \hat{\mu}_i(t) + c\sqrt{\frac{\log t}{n_i(t)}} \right]\]

   Thompson Sampling draws from posterior distributions.

3. Distributed Asynchronous Computing: $K$ agents evaluate configurations in parallel. Empirical mean updates:

   \[\hat{\mu}_{a_t}(t+1) = \hat{\mu}_{a_t}(t) + \frac{1}{n_{a_t}(t)+1}(\bar{r}_t - \hat{\mu}_{a_t}(t))\]

   where $\bar{r}_t = \frac{1}{K} \sum_{j=1}^K r_t^{(j)}$ averages across agents.

   Applied to toy example: With $d=2$, $K=4$, $N=20$: GP posterior over unit square identifies high-uncertainty regions, UCB balances exploration/exploitation across 4 arms, dispatcher assigns evaluations to 4 workers asynchronously, results update GP for next iteration.

Step 1 — Prior work:

• Auer et al. (2002): “Finite-time analysis of the multiarmed bandit problem” — established $O(\sqrt{AT\log T})$ UCB regret bounds for single-agent synchronous bandits
• Joulani et al. (2013): “Online learning under delayed feedback” — analyzed bandit regret under bounded delays but without active learning integration
• Li et al. (2018): various distributed hyperparameter optimization papers — applied bandits to HPO but without GP surrogates or active learning

Step 2 — Delta: This paper integrates active learning (GP-based acquisition), multi-armed bandits (UCB/Thompson Sampling), and distributed computing into a single framework. Prior work addressed subsets but not the full combination.

Step 3 — Theory-specific assessment:

• Main contribution is architectural/systems: the unified framework, not new theory
• Regret bounds are standard UCB/Thompson Sampling results from Auer et al./Russo et al., extended heuristically to distributed setting
• No new proof techniques — assembles known lemmas from MAB theory, delayed feedback, and Amdahl’s Law
• No lower bound analysis or optimality claims

The paper provides extensive empirical validation (5 benchmarks, 500 replicates each) but limited theoretical novelty.

Verdict: INCREMENTAL — solid engineering contribution combining existing components, but no fundamental theoretical advances or surprising results.

No substantial proofs are provided. The theoretical analysis relies on standard results:

1. UCB regret bound (from Auer et al. 2002):

   \[R^{\text{UCB}}_T \leq \sum_{i: \Delta_i > 0} \frac{8\log T}{\Delta_i} + \left(1 + \frac{\pi^2}{3}\right)\Delta_i\]

2. Delayed feedback extension (from Joulani et al. 2013):

   \[R^{\text{delay}}_T = O\left(\sqrt{(T + \tau_{\max})K \log T}\right)\]

3. Amdahl’s Law speedup:

   \[S(K) = \frac{T_1}{T_K} = \frac{1}{p + \frac{1-p}{\eta K}}\]

4. AL-MAB regret decomposition (Remark 1, informal):

   \[E[R^{\text{ALMAB}}_T] = E[R^{\text{MAB}}_T] + \varepsilon(q_T)\]

   where $\varepsilon(q_T)$ is surrogate approximation error.

   The paper assembles these known results without rigorous proofs for the combined setting. The “optimal agent count” formula in Eq. (7) is derived via basic calculus but lacks formal validation.

Datasets: 5 benchmarks with 500 replicates each using calibrated surrogate models:

• Case 1: EfficientNet-B0/CIFAR-10 HPO (N=60 evaluations)
• Case 2: CFD airfoil drag minimization (N=50 evaluations)
• Case 3: MuJoCo HalfCheetah RL policy search (N=50 evaluations)
• Case 4: Bayesian dose-response optimization (10 rounds)
• Case 5: Adaptive spatial field estimation (N=30 observations)

Baselines: Grid Search, Random Search, BOHB, Optuna (TPE), ALMAB-DC (Thompson Sampling variant)

Key Results:

• Case 1: 93.4% CIFAR-10 accuracy (+1.7pp over BOHB, +1.1pp over Optuna)
• Case 2: CD = 0.059 drag coefficient (36.9% below Grid Search)
• Case 3: 9,599 expected return (+50% over Grid Search)
• Scaling: 7.5× speedup at K=16 agents consistent with Amdahl’s Law
• All improvements statistically significant (Bonferroni-corrected Mann-Whitney U tests)

Note: Uses surrogate models rather than live experiments for reproducibility and computational tractability.

Limitations:

1. Surrogate evaluation — RESTRICTIVE: All experiments use calibrated surrogate models rather than live systems, limiting real-world validation
2. No theoretical analysis for combined setting — TECHNICAL: Regret bounds assume decomposition between AL and MAB components without rigorous proof
3. Communication model oversimplification — TECHNICAL: $O(\log K)$ or $O(K)$ communication costs ignore network topology and protocol details
4. GP scalability — NATURAL: $O(n^3)$ kernel inversion limits surrogate to moderate dataset sizes

5. Homogeneous agent assumption — RESTRICTIVE: Analysis assumes identical agent capabilities, but real distributed systems are heterogeneous

   Open Problems:

6. Large-scale validation: Demonstrate ALMAB-DC on live expensive experiments (clinical trials, large-scale CFD) rather than surrogates
7. Theoretical unification: Provide rigorous regret bounds for the combined AL+MAB+distributed setting, particularly handling the interaction between surrogate approximation error and bandit regret

Motivation: PAC-Bayes learning rules balance training error against KL divergence to avoid overfitting, but their agnostic consistency properties remain unclear, especially their relationship to Bayesian prediction.

Mathematical setup: Consider i.i.d. samples $S \sim D^m$ from source distribution $D$ over $X \times {0,1}$. For prior distribution $\Pi$ over predictors $h: X \to {0,1}$ and regularization parameter $\lambda > 0$, define the PAC-Bayes learning rule:

\[\hat{Q}_\lambda(S) = \arg \min_Q J_\lambda(Q, S)\]

where:

\[J_\lambda(Q, S) = mH(L_S(Q)) + \lambda \text{KL}(Q \| \Pi)\]

with $H(\cdot)$ the binary entropy function, $L_S(Q) = \mathbb{E}_{h \sim Q}[L_S(h)]$ the expected training error, and constraint $L_S(Q) \leq 1/2$.

Assumptions:

1. Source distribution $D$ is arbitrary (agnostic setting)
2. Prior $\Pi$ is a valid probability measure on predictor space

3. Posteriors $Q$ are absolutely continuous with respect to $\Pi$

   Toy example: When $\lambda = 1$ and $\Pi$ is uniform over two predictors ${h_0, h_1}$ with errors $L^* < 1/2$ and $L’ > L^*$, the rule reduces to selecting the predictor minimizing $mH(L_S(h)) - \log \Pi(h)$.

   Formal objective: Characterize the worst-case limiting error:

   \[L_\infty(\lambda, L^*) = \lim_{m \to \infty} \sup_{(\Pi,D) \in M(L^*)} \mathbb{E}_{S \sim D^m}[L_D(\hat{Q}_\lambda(S))]\]

   where $M(L^*) = {(\Pi, D) : \exists Q^* \text{ with } L_D(Q^*) = L^*, \text{KL}(Q^* | \Pi) \leq 10}$.

The PAC-Bayes learning rule produces a Gibbs posterior with self-consistent temperature:

\[d\hat{Q}_\lambda(h) \propto d\Pi(h) \cdot 2^{-\frac{m}{\lambda}H'(L_S(\hat{Q}_\lambda)) L_S(h)}\]

The method involves:

1. Define regularization function for $\lambda > 0$:

   \[T_\lambda(q) = \min_{0 \leq p \leq 0.5} \lambda \text{KL}(p \| q) + H(p)\]

2. For $0 < \lambda \leq 1$: minimizer is $p^* = 0$, so $T_\lambda(q) = -\lambda \log(1-q)$ and limiting error is:

   \[\ell_\lambda(L^*) = 1 - 2^{-\frac{H(L^*)}{\lambda}}\]

3. For $\lambda > 1$: minimizer is $p^* = \frac{1}{1 + (\frac{1-q}{q})^{\frac{\lambda}{\lambda-1}}}$, giving $T_\lambda(q) = U_\lambda(q)$ where:

   \[U_\lambda(q) = \lambda \text{KL}\left(\frac{1}{1+(\frac{1-q}{q})^{\frac{\lambda}{\lambda-1}}} \| q\right) + H\left(\frac{1}{1+(\frac{1-q}{q})^{\frac{\lambda}{\lambda-1}}}\right)\]

   and limiting error is $\ell_\lambda(L^*) = U_\lambda^{-1}(H(L^*))$.

   Toy example application: For two predictors with uniform prior and $L^* = 0.1$, $\lambda = 0.5$: the limiting error is $\ell_{0.5}(0.1) = 1 - 2^{-H(0.1)/0.5} \approx 0.85$, showing catastrophic overfitting.

Prior work:

1. Zhu and Srebro [2025]: analyzed MDL with discrete priors, showed $\lambda \approx \sqrt{m}$ needed for consistency
2. McAllester [2003], Catoni [2007]: established PAC-Bayes bounds for continuous priors but without agnostic consistency analysis

3. Grünwald and Langford [2004]: showed MDL₁ is inconsistent in agnostic setting

   Delta: This paper extends discrete MDL analysis to continuous PAC-Bayes rules and provides explicit Bayesian interpretations. New contributions include:

4. exact characterization of worst-case limiting error for continuous priors via matching upper/lower bounds
5. proof that $\hat{Q}_\lambda$ corresponds to empirical Bayes with sample-size-dependent prior

6. connection to profile posterior and full Bayesian posterior.

   Theory-specific assessment:

   • Main theorem extends known discrete results to continuous case using similar proof techniques (concentration + entropy concavity)
   • Bounds are tight by construction, matching discrete case exactly when specialized
   • Lower bounds achieved by explicit constructions, not surprising given discrete precedent
   • Bayesian interpretation is new but follows from standard variational analysis

   The theoretical contribution is solid but represents a natural extension rather than breakthrough. The continuous setting requires technical refinements (entropy concavity) but core insights follow established patterns.

   Verdict: INCREMENTAL — clean extension of discrete MDL overfitting analysis to continuous PAC-Bayes setting with useful Bayesian interpretations.

The proofs combine PAC-Bayes concentration with entropy-based analysis:

Upper bounds: Core inequality from PAC-Bayes bound gives with probability $1-\delta$:

\[T_\lambda(L_D(\hat{Q}_\lambda(S))) \leq H(L_D(Q^*)) + \lambda \frac{\log(m+1) + \text{KL}(Q^*\|\Pi)}{m} + O\left(\sqrt{\frac{\log^3 m}{m}}\right)\]

Key step: analyze $T_\lambda^{-1}$ applied to RHS using different regimes for $\lambda$.

Lower bounds: Constructive proof using hard instances:

1. Hypothesis space: Binary sequences $x = x[0]x[1]…$ with predictors $h_i(x) = x[i]$
2. Prior: Universal coding prior $\Pi(h_i) = 1/(i \log^2 i + 10)$

3. Source: $y \sim \text{Ber}(1/2)$, $x[0] = y \oplus \text{Ber}(L^*)$, $x[i] = y \oplus \text{Ber}(L’)$ for $i \geq 1$

   Concentration argument: Show $\hat{Q}_\lambda$ puts zero mass on good predictor $h_0$ by proving existence of competitor $\delta_{h_{\hat{i}}}$ with strictly better objective:

   For $\lambda \leq 1$: Use interpolating predictor with $L_S(h_{\hat{i}}) = 0$

   For $\lambda > 1$: Use empirical risk minimizer among first $k(m) = 2^{mKL(\hat{L}|L’)}$ bad predictors

   Entropy concavity: Critical technical tool - for any posterior $Q$ with mass $q_0$ on $h_0$:

   \[H(L_S(Q)) \geq q_0 H(L_S(h_0))\]

   This enables lower bounding the objective when $Q$ has positive mass on good predictor, leading to contradiction.

Purely theoretical. Empirical validation would require:

1. Implementing PAC-Bayes rule with entropy regularization on synthetic datasets matching the hard instances (binary sequences with coordinate predictors)
2. Comparing limiting errors across different $\lambda$ values and sample sizes
3. Validating the empirical Bayes interpretation by comparing with explicit Bayesian posterior computations
4. Testing consistency claims for $\lambda_m \approx \sqrt{m}$ scaling on real classification datasets

Limitations:

1. RESTRICTIVE: Analysis requires $L_S(Q) \leq 1/2$ constraint which may not hold with arbitrary priors

2. TECHNICAL: Hard instances use discrete hypothesis spaces embedded in continuous setting - unclear if results extend to genuinely continuous spaces like neural networks

3. RESTRICTIVE: Bounds only proven for competitor $Q^*$ with $\text{KL}(Q^*|\Pi) \leq 10$ - larger complexity requires different analysis

4. NATURAL: Agnostic setting assumption is standard but may be pessimistic for many real problems

5. TECHNICAL: Upper bounds for profile posterior $\hat{Q}^{\text{mod}}_\lambda$ remain open - only lower bounds proven

   Open problems:

6. Prove matching upper bounds for profile posterior and full Bayesian posterior to complete consistency characterization
7. Extend analysis to genuinely continuous hypothesis classes (e.g., neural networks) beyond discrete embeddings

1. Motivation: Bayesian model inference on binary model spaces includes graphical models, variable selection, mixture distributions, and decision trees. Traditional methods like reversible jump MCMC and birth-death algorithms suffer from slow exploration, making them impractical for large-scale problems (e.g., 500,000 parameters). Faster algorithms are needed.

2. Mathematical setup: Consider unknown binary vectors $m = (m_1, \ldots, m_k) \in M := {0,1}^k$ where each $m_i \in {0,1}$ and $M$ has $l = 2^k$ elements. Each model $m$ has associated parameters $\theta_m \in \Theta_m$. The joint posterior is:

   \[p(m, \theta_m | y) = \frac{p(y | m, \theta_m)p(\theta_m | m)p(m)}{p(y)}\]

   The marginal posterior is:

   \[p(m | y) = \frac{p(y | m)p(m)}{p(y)}\]

   Key assumptions:

   1. Binary model space structure

   2. Marginal likelihood $p(y

      m)$ is available or can be approximated

   3. Model prior $p(m)$ is specified

3. Toy example: When $k = 2$ with $M = {(0,0), (0,1), (1,0), (1,1)}$ and uniform prior $p(m) = 1/4$, this reduces to comparing 4 models. If posterior probabilities are $\pi(1,1) = 0.01$ and $\pi(0,0) = \pi(0,1) = \pi(1,0) = 0.33$, traditional BD algorithms require many steps to explore properly since they only flip one bit per iteration.

4. Formal objective: Design an MCMC algorithm with limiting distribution equal to the model posterior:

   \[\lim_{s \to \infty} \text{P}(m^{(s)} = m) = p(m | y)\]

The Multiple Jump MCMC (MJ-MCMC) algorithm starts with a birth-death process with rates $Q(m, m’)$ and transforms it into a discrete-time algorithm.

Algorithm steps:

1. Given current state $m^{(s)}$, compute rates $q_i(m^{(s)}) = \min(1, p(m_i

   y)/p(m^{(s)}

   y))$ for all $i = 1, \ldots, k$

2. For user-defined sequence $\varepsilon_s \in (0,1)$, flip each element $m_i^{(s)}$ independently with probability $q_i(m^{(s)})\varepsilon_s$

3. The resulting vector is $m^{(s+1)}$

   Transition probabilities are:

   \[P_{\varepsilon_s}(m, m') = \prod_{i \in H_{m,m'}} q_i(m)\varepsilon_s \prod_{i \notin H_{m,m'}} (1 - q_i(m)\varepsilon_s)\]

   where $H_{m,m’}$ is the set of positions where $m$ and $m’$ differ.

   Toy example application: For the 2D case with $m = (0,0)$ and rates $q_1(0,0) = q_2(0,0) = 0.03$, if $\varepsilon_s = 0.1$, then:

   • Stay at $(0,0)$: probability $(1-0.003)^2 = 0.994$
   • Move to $(1,0)$: probability $0.003 \times 0.997 = 0.003$
   • Move to $(0,1)$: probability $0.997 \times 0.003 = 0.003$
   • Move to $(1,1)$: probability $0.003 \times 0.003 = 0.000009$

Step 1 — Prior work:

• “Bayesian Structure Learning in Sparse Gaussian Graphical Models” (Mohammadi & Wit, 2015): introduced birth-death MCMC for GGMs, but limited to single-edge updates
• “Reversible Jump Markov Chain Monte Carlo Computation and Bayesian Model Determination” (Green, 1995): established trans-dimensional MCMC framework, but suffers from acceptance/rejection steps
• Birth-death algorithms (Jennings & Corcoran, 2016): rejection-free but only allow local moves

Step 2 — Delta: This paper transforms any birth-death process into a discrete-time algorithm that can update multiple (potentially all) model components simultaneously. Key innovation is the transition probability formula (8) that maintains theoretical guarantees while allowing global jumps.

Step 3 — Theory-specific assessment:

• Main theorem is somewhat predictable: scaling birth-death rates by $\varepsilon_s$ and allowing simultaneous updates naturally leads to faster mixing
• Proof technique is routine: uses standard MCMC convergence theory and detailed balance arguments
• Bounds are not particularly tight: $O(\varepsilon)$ convergence rate in homogeneous case is expected
• No known lower bounds are provided for comparison

The algorithm’s practical performance (100-200x speedup) is impressive, but the theoretical contributions are incremental extensions of existing birth-death theory.

Verdict: INCREMENTAL — solid extension of birth-death MCMC with good practical benefits but predictable theoretical results.

The proofs use standard MCMC convergence theory with two main components:

1. Homogeneous case (Theorem 3.1): Uses detailed balance to show the stationary distribution $\pi_\varepsilon$ satisfies:

   \[\pi_\varepsilon(m) P_\varepsilon(m, m') = \pi_\varepsilon(m') P_\varepsilon(m', m)\]

   Key inequality for convergence rate:

   \[\|\pi_\varepsilon - \pi\| = O(\varepsilon)\]

   The waiting time convergence uses:

   \[\varepsilon W_m^\varepsilon \xrightarrow{d} W_m \text{ as } \varepsilon \to 0\]

2. Inhomogeneous case (Theorem 3.2): Applies Robbins-Siegmund theorem for stochastic approximation. Requires:

   \[\lim_{s \to \infty} \varepsilon_s = 0 \text{ and } \sum_{s=1}^\infty \varepsilon_s = \infty\]

   Key technical insight: The algorithm maintains irreducibility because $P_{\varepsilon_s}(m, m’) > 0$ for all $m, m’$ with $\pi(m’) > 0$, enabling global jumps while preserving convergence.

3. Special case (Lemma 3.3): For factorizable posteriors:

   \[p(m | y) = \prod_{i=1}^k p(m_i | y)\]

   The algorithm is exact for any $\varepsilon \in (0,1)$ because element-wise independence matches the algorithm’s independent update structure.

Datasets:

1. Simulated Gaussian graphical models: p=1000 variables, n=400 observations, edge density 0.2%

2. Real mouse gene expression data (GSE15907): p=623 genes, n=653 observations

   Baselines: BD-MPL algorithm (Mohammadi et al., 2024), the current state-of-the-art for large-scale Bayesian structure learning

   Key numbers:

   • 100-200x speedup over BD-MPL on simulated data
   • 10x speedup on real gene expression data
   • Solves 1000-node problems in under 30 seconds
   • AUC-PR scores comparable to baseline across all tested $\varepsilon$ values
   • Correlation >0.988 with BD-MPL edge inclusion probabilities

   Performance tested with homogeneous $\varepsilon \in {0.001, 0.01, 0.3, 0.7, 0.95}$ and inhomogeneous sequences $\varepsilon_s = 0.3/\log_{10}(s+9)$ and $\varepsilon_s = 0.3/(s\log_2(s+1))^{0.4}$.

Limitations:

1. Algorithm can visit low-probability states when $\varepsilon_s$ is large - TECHNICAL (needed for fast mixing but diminishes as $\varepsilon_s \to 0$)
2. Choice of $\varepsilon_s$ sequence presents accuracy-speed tradeoff - TECHNICAL (common in stochastic optimization)
3. Theoretical guarantees only for marginal posterior, not joint posterior - TECHNICAL (extension via Lemma 3.4 is straightforward)
4. Requires closed-form or accurate approximation of Bayes factors - NATURAL (standard requirement in Bayesian model selection)
5. Performance on very dense graphs ($>$50% edge density) not thoroughly tested - RESTRICTIVE (limits applicability to some domains)

6. No convergence rate bounds provided - TECHNICAL (standard limitation in MCMC theory)

   Open problems:

7. Extend framework to continuous model spaces using techniques from Chen (2023)
8. Develop adaptive strategies for automatically choosing optimal $\varepsilon_s$ sequences based on problem characteristics