Variational Search Distributions
Authors: Daniel M. Steinberg
Rafael Oliveira
Cheng Soon Ong
Edwin V. Bonilla
Venue: ICLR 2025 Poster (2025)
Keywords: black box optimization
Bayesian optimization
variational inference
generative models
level set estimation
biological sequence design
protein engineering
fine tuning
Primary area: probabilistic methods (Bayesian methods, variational inference, sampling, UQ, etc.)
TL;DR: A variational inference framework for learning conditional generative models of (biological) sequences for black-box optimization and more.

Abstract: We develop VSD, a method for conditioning a generative model of discrete, combinatorial designs on a rare desired class by efficiently evaluating a black-box (e.g. experiment, simulation) in a batch sequential manner. We call this task active generation; we formalize active generation's requirements and desiderata, and formulate a solution via variational inference. VSD uses off-the-shelf gradient based optimization routines, can learn powerful generative models for desirable designs, and can take advantage of scalable predictive models. We derive asymptotic convergence rates for learning the true conditional generative distribution of designs with certain configurations of our method. After illustrating the generative model on images, we empirically demonstrate that VSD can outperform existing baseline methods on a set of real sequence-design problems in various protein and DNA/RNA engineering tasks.

Published as a conference paper at ICLR 2025
VARIATIONAL SEARCH DISTRIBUTIONS
Daniel M. Steinberg, Rafael Oliveira, Cheng Soon Ong & Edwin V. Bonilla Data61, CSIRO, Australia
steinberr,, rafsel. wtw.el..ddsanttoeeolleeeee iveiiir cheng snon.ong, edwin..boniLLa)@aatalll
ABSTRACT
We develop variational search distributions (VSD), a method for conditioning a generative model of discrete, combinatorial designs rare desired class by efficiently evaluating a black-box (e.g. experiment, simulation) in a batch sequen- tial manner. We call this task active generation; we formalize active generation's 's requirements and desiderata, and formulate a solution via variational inference. VSD uses off-the-sheif gradient based optimization routines, can learn powerful generative models for desirable designs, and can take advantage of scalable predic- tive models. We derive asymptotic convergence rates for learning the true condi- tional generative distribution of designs with certain configurations of our method. After illustrating the generative model on images, we empirically demonstrate that VSD can outperform existing baseline methods set of real sequence-design problems in various protein and DNA/RNA engineering tasks.
INTRODUCTION
We consider a variant of the active search problem (Garnett et al., 2012; Jiang et al., 2017; Vanchi- nathan et al., 2015), where we wish find members (designs) of rare desired class batch sequential manner with fixed black-box evaluation (e.g. experiment) budget. We call sequential active learning of generative model of these designs active generation. Examples of rare designs are compounds that could be useful pharmaceutical drugs, or highly active enzymes for catalyzing chemical reactions. We assume the design space is discrete or partially discrete, high-dimensional and practically innumerable. For example, the number of possible configurations of single protein is (100) (see, e.g., Sarkisyan et al., 2016). Learning generative model of these designs allows us to circumvent the need for traversing the whole search space.
We are interested in this active generation objective for variety of reasons. We may wish to study the properties of the "fitness landscape" " (Papkou et al., 2023) to gain better scientific understanding of a phenomenon such natural evolution. Or, may not be able to completely specify the constraints and objectives a task, but we would like to characterize the space of, and generate new feasible designs. For example, want enzymes that degrade plasties in industrial setting, but we may not yet know the exact conditions (c.g. temperature, pH), some of which may be anti- correlated with enzyme catalytic activity. Alternatively, if we know these multiple objectives and constraints, we may only want to generate designs from a Pareto set.
Assuming we can take advantage of a prior distribution over designs, we formulate the search prob- lem as inferring the posterior distribution over rare, desirable designs. Importantly, this posterior can be used for generating new designs. Specifically, we use (black-box) variational inference (VI) (Ran- ganath et al., 2014 4)), and so refer to our method as variational search distributions (VSD). Our major contributions are: (1) we formulate the batch active generation objective over a (practically) innu- merable discrete design space, (2) we present variational inference algorithm, VSD, which solves this objective, (3) we show that VSD performs well theoretically and empirically, and (4) we connect active generation to other recent advances in black-box optimization (BBO) of discrete sequences that use generative models. VSD uses off-the-shelf gradient based optimization routines, is able to learn powerful generative models, and can take advantage of scalable predictive models. In our ex- periments we show that VSD can outperform existing baseline methods on a set of real applications. Finally, evaluate our approach on the related sequential BBO problem, where we want to find the globally optimal design for a specific objective and show competitive performance when compared with state-of-the of-the-art methods, e.g., based latent space optimization (LSO) (Gruver et al., 2023).

---

Published as a conference paper at ICLR 2025
2
METHOD
In this section we formalize our problem and describe its requirements and desiderata. We also develop proposed solution, based on variational inference, which will refer to variational search distributions (VSD).
2.1 THE PROBLEM OF ACTIVE GENERATION
We are given design space X', which can be discrete or mixed discrete-continuous and high dimen- sional, and where for each instance that we choose x E X, we measure some corresponding property of interest (so-called fitness) y R. For example, in our motivating application of DNA/RNA or protein sequences, X' = DM where V'ss is the sequence vocabulary (e.g.. amino acid labels, VI = 20) and M is the length of the sequence. However, we do limit the application of our method se- quences. Using this framing, real world experiment (e.g., measuring the activity of an enzyme) can be modeled as an unknown relationship,
y y= f((x) +
(1)
for some black-box function (the experiment), f., and measurement error e E IR, distributed accord- ing to p(e) with [e] O. Instead of modeling the whole space, X. we are only interested in a set of events which we choose based on fitness, S C X. In particular for active generation we wish to learn generative model, qxx., that only returns samples x(*) E S by efficiently query- ing the black-box function in Equation 1. We assume that S are events in high dimensional space, and that we have access to a prior belief, p(x), which helps narrow in on this subset of X. We are given an initial dataset, DN == ((y):xn)) of Xn E S. Given p(x) and DN we aim to generate batches of unique candidates, {xM [x for black-box (experimental) evaluation in a series of rounds, E (1...., T}, where and we desire Xbt S. After each round DN is augmented with the black-box results of the batch,
which may contain only a few instances
(1100))
ie. DN + NN U {(XMr (he) As we shall see later, our solution allows us to satisfy the follow- ing requirements and additional desiderata for active generation.
Requirements & Desiderata. Active generation requirements (R) and other desiderata (D).
(R1) Rare feasible designs, S. events in X that need be identified
(RS) Generative models, x(*) q(x), that task- specifie for rare, feasible designs
(D1) Guaranteed convergence for certain choices of priors, variational distributions and predic- tive models
(R2) Sequential non-myopic candidate generation, x E for sequential black-box evaluation
(R3) Discrete search over (combinatorially) large
(D2) Gradient based optimization strutegies for can- didate searching
design spaces, e.g. x E Vil
(R4) Batch generation of to ((1000) diverse can- (D3) Scalable predictive models that enable high- didate designs per round
throughput evaluationlexperiment
Like active search (Garnett et al., 2012) in our case we are interested in the solution space of the super level-set, SSLS {x f.(x) 7) for threshold IR (e.g., wild-type fitness). As we only have access to noisy measurements, our task is to estimate the super level-set distribution, p(x|y > T), using active generation. Estimating this distribution over SSLS is computationally and statistically challenging and, therefore, we cast this as variational inference problem. We also consider the case of BBO for which SgBO == argmax, f((x). and we show that we can accommodate this in our variational framework by iteratively raising per round, l. We visualize the properties and models involved active generation as applied to a continuous "fitness landscape" in Figure 1.
2.2 VARIATIONAL SEARCH DISTRIBUTIONS
We cast the estimation of p(xly > r) sequential optimization problem. A suitable objective for a round, t, is to minimize a divergence,
Φ = argmin Dip(xly > [x((()))))
(2)
2

---

Published as a conference paper at ICLR 2025
(a) argmax f-(x)
(b) fx f((x) +}
(c) p(x)
(d) p(x(v >> T)
Figure 1: Fitness landscape properties and models. (a) A noise-less fitness landseape, f-(x), and the maximum fitness design, Snno (x*). the white 'x'. (b) The super level-set, SSLS, of all fit designs as the white hatched area. (c) Prior belief p(x). (d) The density/mass function of the super level-set, p(xly 7), blue contours. Our goal is to sequentially estimate a generative model for the distribution of the super level-set (d), We assume a noisy relationship between f. and so the super level-set will not have a hard boundary, and p(xly T) will be defined over all X.
where (xx()) is parameterized distribution from which sample candidate designs Xbc= (R5), and which we aim to match to p(xly y > r). The difficulty is that we cannot directly evaluate or empirically sample from p(xly > T). However, if we consider the reverse Kullback-Leibler (KL) divergence,
((x@)) log > x1x) p(x)
argmin [ ( +) = argmin E (x40)
log
(3)
where we have expanded p(xly > r) using Bayes rule and dropped the constant term p(y T), note that we no longer require evaluation of p(xly > т) directly. We recognize the right-hand side of Equation 3 as the well known (negative) variational evidence lower bound (ELBO),
(4) For this we assume access to prior distribution over the space of designs, p(x), that may be in- formed from the data at hand (or pre-trained) Henceforth, as we will develop a sequential algorithm, we will denote this prior as p(x(Do). We note the relationship between log pvy > T|x) and the prob- ability of improvement (PI) acquisition function from Bayesian optimization (BO) (Kushner, 1964),
(o) := Egtxie)]logg(( > t x)] - [xxx,x(([))]l((x)]
log ply > r|x) = log p(u r|x, Dw) log Ep(x)x [1 r]] oogori(x DN, r). (5) Here 1 {false, true)
{0, 1} is the indicator function and p(y)x, is typically esti- mated using the posterior predictive distribution of a Gaussian process (GP) given data, DN- So p(y > r|x, DN) V ((UN
/)/oN (x)), where (-( cumulative standard normal distribu- tion function, and (x)).. o.y (x) are the posterior predictive mean and variance, respectively, of the GP. We refer this estimation strategy as GP-PI. and rewrite the ELBO accordingly,
(6) We refer to the method that maximizes the objective in Equation 6 as variational search distributions (VSD), since we using the variational posterior distribution a means of searching the space of fit designs, satisfying (RI), (R2) and (R4). It is well known that when the true posterior is member of the variational family indexed by the above variational inference procedure has the potential to recover the exact posterior distribution. To recommend candidates for black- bbx evaluation we sample a of designs from our search distribution each round,
(d, DN) Ex(x)o) (og DN, - D (a(x|)))[l((())))
B (x(())) b=1
(xbt) )
where & argmax T,
We discuss the relationship between VSD and BO in Appendix G. In general, because of the dis- crete combinatorial nature of our problem, we cannot use the ce-parameterization trick (Kingma & Welling, 2014) to estimate the gradients of the BLBO straightforwardly. Instead, we use the score function gradient, also known as REINFORCE (Williams, 1992; Mohamed et al., 2020) with standard gradient descent methods (D2) such as Adam (Kingma & Ba, 2014),
VrLnno (,, Dx) By (llg a (x, DN, r) log (x(()) logp(x (D))) V log (x(())]
3

---

Published as a conference paper at ICLR 2025
Here use Monte-Carlo sampling to approximate the expectation with suitable variance duction scheme, such as control variates (Mohamed et al., 2020). We find that the exponentially smoothed average of the ELBO works well in practice, and is the same strategy employed in Daulton et al. (2022). VSD implements black-box variational inference (Ranganath et al., 2014) for parame- ter estimation, and despite the high-dimensional nature of X, we find we only need (1000) samples to estimate the required expectations for ELBO optimization on problems with M = =(100) satis- rying (R3). Note that Equation 6 8 do not involve any data (Dx) directly, only indirectly through the acquisition function. Hence the scalability of VSD is dependent on the complexity of training the underlying estimator of p(y) x, DN).
2.3 CLASS PROBABILITY ESTIMATION
So far our method indirectly computes PI by transforming the predictions of a GP surrogate model, p(y|x, D.). as Equation 5. Instead we may choose to follow the reasoning used by Bayesian optimization by density-ratio estimation (BORE) in Tiao et al. (2021); Oliveira et al. (2022); Song al. (2022), and directly estimate the quantity we care about, p(y > T T|x, Dv). This can be ac- complished using class probability estimation (CPE)
[ > 7] {0, 1 so ply > r|x, Dx) = p(= = 1|x, DN) s (x)),, where Te : X - [0, 1]. We can estimate the class probabilities using proper scoring rule (Gneiting & Raftery, 2007) such as log-loss,
the labels
(C (N)
1 N'
Σm In log (xx) + (1 zn. log(1 (x.))))
N
where D [Z2n Xn xn) J.y
The VSD objective and gradient estimator using CPE then become, (o, 0) = By(xjoy To [x) - DKI. (xxl((((()))) [x())
(10)
V V 8) = En(x|o))(( [log no (x) - log (x(() + log V log [x]]))
(11)
into which we plug 07 argmin, (C ). We refer to this strategy as CPE-PI. Using a CPE enables the use of more scalable estimators than GP-PI, satisfying our desideratum (D3). This is crucial if we choose to run more than a few rounds of experiments with B = ((100)). Since VSD is a black-box method we may choose to use CPEs that are non-differentiable, such as decision tree ensembles. The complete VSD algorithm is given in Algorithm and depicted in Figure B.I. We have allowed for a threshold function, t. = f+v y DN ).)), that can be used to modify the threshold each round. For example, an empirical quantile function 7. Qy(7)) where " E (0, as in Tiao et al. (2021). Or constant 7 for estimating a constant distribution of the super level-set.
Algorithm 1 VSD optimization loop with CPE-PI.
Require: Threshold function fr and n/ dataset DN. black-box f., prior p(x(Do). CPE To (x)), variational family g(x (o), budget T and B.
1: function FITMODELS(DNN1))
2: D - {(zn, x.) ) where zn = 1 > 7]
3: 0 argming (CPF (0, )
4:
6" + argmax, CELBO
5:
return 6-,0"
6:
for round E (l,.... .......,TT T] do
7:
To F ff y E DNJ.
8:
0; - FITMODELS(DN,1))
9:
{xo1] B
/b=1 +- (xx
10:
{ybe) B +- (ff(Xx))) Eb Ib=
11:
DN f DN U (xxx, )b
12:
T. + f+ e D.t.
13: 0-,0- - FITMODELS(D
14: return 0²,0"
2.4 THEORETICAL ANALYSTS
We show that VSD sampling distributions converge to a target distribution that characterizes the level set given by t, satisfying (D1) in two general settings. We first derive results assuming is drawn
4

---

Published as a conference paper at ICLR 2025
from Gaussian process, i... f. ~ GP(O,,) k), with residive-eemidefiiite covariance (or kemel) function k : x X' + R (Appendix E), using GP-PI as the CPE for VSD. These results are then extended to probabilistic classifiers based on wide neural networks (NNs) (Appendix F) by means of the neural tangent kernel (NTK) for the given architecture (Jacot et al., 2018). For the analysis, set B = 1 and N = f, though having B > 1 should improve the rates by multiplicative factor. Theorem 2.1. Under mild assumptions (E.1 to E.5), the variational distribution of VSD equipped with GP-PI converges to the level-set distribution in probability at the following rate:
DKLlp(xly> > Te, D:)llp(xly>> > Ts, fo)] Op(t-1/2).
(12)
This result is based on showing that the GP posterior variance vanishes at an optimal rate of O(t in our setting (Lemma E.5). We also analyze the rate at which VSD finds feasible designs, or "hits", compared to an oracle with full knowledge of f.. After T rounds, the number of hits found by VSD is HT Σi= 1y > Tt- 1]. where 1/t follows Equation and x ~ p(xly > T--1; 1). The number of hits, H7, from agent that fully knows f. is the same but for generating conditioned on f. with X, ~ p(x > TC-1, Using this definitoon and Theorem 2.1, we prove the following.
Corollary 2.1. Under the settings in Theorem 2.1, also have that:
Bi|HT - 17 E O(VT).
(13)
B[Hr] is related to the empirical recall measure (18) up to the normalization constant, but it does not account for repeated hits, which are treated as false discoveries (false positives) under recall. Lastly, for NN-based CPEs, we obtain convergence rates dependent on the spectrum of the NTK (Proposition F.2), which we instantiate for infinitely wide ReLU networks below. For the full results and proofs, please see Appendix E for the GP-based analysis and Appendix F for the NTK results.
Corollary 2.2. Let Th be modeled via a fully connected ReL.U network. Then, under assumptions on identifiability and sampling (F.1 to F.6), in the infinite-width limit, VSD with CPE-PI achieves:
DxI.[p(xly Te, De)llp(xly Te, f.)] E Op ((2CaFF3)
(14)
This result finally indicates that, when equipped with flexible NN-based CPBs, VSD is also capable of recovering the target distribution for arbitrary sequence lengths in combinatorial problems.
3 RELATED WORK
We will consider related work firstly in terms of methods that have similar components to VSD, then secondly in terms of related problems to specification of active generation. VSD can be viewed as one of many methods that makes use of the variational bound (Staines & Barber, 2013),
max f.(x) 2 max [fx( x
(15)
The maximum is always greater than or equal to the expected value of a random variable. This bound is useful for black-box optimization (BBO) of f., and becomes tight if q(x()) - 8(x").
See Appendix G for more detail and VSD's relation to BO. Other well known methods that make use of this bound are natural evolution strategies (NES) (Wierstra al., 2014), variational opti- mization (VO) (Staines & Barber, 2013; Bird et al., 2018), estimation of distribution algorithms (EDA) (Larrañaga & Lozano, 2001; Brookes et al., 2020), and Bayesian optimization with proba- bilistic reparametrization (BOPR) (Daulton et al., 2022). For learning the parameters of the vari- ational distribution, é, they variously make use of maximum likelihood estimation or the score function gradient estimator (REINFORCE) (Williams, 1992). Algorithms that explicitly modify Equation 15 to stop the collapse of (x to a point mass for batch design include design by adaptive sampling (DbAS) (Brookes & Listgarten, 2018) and conditioning by adaptive sampling (CbAS) (Brookes et al., 2019). They use fixed samples x*) from g(x ( for approximating the expectation, and then optimize @ using a weighted mexiiuu--ikeelhoon or variational style proce- dure. Though DbAS and CbAS were formulated for offline noo-seqqunttal)) tasks, they have often been used a sequential setting (Ren al., 2022). We can take a unifying view of algorithms that use surrogate model for f. by recognizing the general gradient estimator,
Ey(xler) ( [w(x) V log {x([))]
(16)
5

---

Published as a conference paper at ICLR 2025
Method
w(x)
o'
Fixed x(*) ~ g(((o') ) per round?
VSD
log To. (x) + log p(xx Du) log (x((())
No (REINFORCE)
DbAS Ie. (x)p(x Da) //(x///(///////
CbAS
(1-1
Yes (importance Monte Carlo)
TT- (x)
Yes (Monte Carlo)
BORE
TT-- (x)
No (REINFORCE)
BOPR
(x, Dx)
No (REINFORCE)
Table I:
How related methods can be adapted from Equation 16. VSD, CbAS and DbAS may also use cumulative distribution representation of (x)(x, DN, 7) in place of To- (x)(
where we give each component in Table 1. For experiments BORE has been adapted to discrete X' by using the score function gradient estimator, which we denote by BORE", while CbAS and DbAS have been adapted to use CPE their original derivations use PI acquisition function.
A number of finite horizon methods have been applied to biological sequence BBO tasks, such as Amortized BO (Swersky et al., 2020), GFlowNets (Jain al., 2022), and the reinforcement learning based DynaPPO (Angermueller et al., 2019). LSO-like methods (Gómez-Bombaree Tripp et al.. 2020; Stanton et al., 2022; Gruver et al., 2023) tackle optimization of sequences by encoding them into a continuous latent space within which candidate optimization or generation takes place. Selected candidates are decoded back into sequences before black box evaluation; see González-Duque et al. (2024) for comprehensive survey. VSD does not require a latent space nor an encoder, and as such can be seen as an amortized variant of probabilistic reparameterisatioon methods (Daulton al., 2022) or continuous relaxations (Michael et al., 2024). Heuristic stochastic search methods such as AdaLead (Sinai et al., 2020) and proximal exploration (PEX) (Ren et al., 2022) have also demonstrated strong empirical performance on these tasks. We compare the prop- erties of the most relevant methods to our problem in Table
et al., 2018;
In contrast to just finding the maximum using BBO, active generation considers another problem generating samples from a rare set of feasible solutions. Generation methods that estimate the super level-set distribution, p(xly 7), include CbAS, which optimizes the forward KL divergence, [xxly(xxly > (x][( using importance weighted cross entropy estimation (Rubinstein, 1999). Batch-BORE (Oliveira al., 2022) also optimizes the reverse KL divergence and uses CPE, but with Stein variational inference (Liu & Wang, 2016) for diverse batch candidates (with a continuous relaxation for discrete variables). There a rich literature on the related task of active learning and BO for level set estimation (LSE) (Bryan et al., 2005; Gotovos et al., 2013; Bogunovic et al., 2016; Zhang et al., 2023a). However, we focus on learning generative model of discrete space.
For active generation VSD, CbAS and DbAS all use an acquisition function defined in the original domain, X, to weight gradients (see Equation 16) for learning a conditional generative model, from which Xbc are sampled. An alternative is to use guided generation, that is to train an unconditional generative model, and then have a discriminative model guide (condition) the samples from the unconditional model at test time. This plug-and-play of a discriminative model has shown promise for controlled image and text generation of pre-trained models (Nguyen et al., 2017: Dathathri et al., 2020; Li et al., 2022; Zhang al., 2023b). LaMBO (Stanton et al., 2022) and LaMBO-2 (Gruver et al., 2023) take a guided generation approach to solve the active generation problem. LaMBO uses an (unconditional) masked language model auto-encoder, and then optimizes sampling from its latent space using an acquisition function as guide. LaMBO-2 takes a similar approach, but uses a diffusion process the unconditional model, and modifies Langevin sampling de-noising process with an acquisition function guide.
4 EXPERIMENTS
Firstly we test our method, VSD, on its ability to generate complex, structured designs, x, in a single round by training it to generate subset of handwritten digits from flattened MNIST images (LeCun et al., 1998) in Sec. 1. We then compare VSD on two sequence design tasks against existing baseline methods. The first of these tasks (Sec. 4.2) is to generate many unique, fit sequences as possible using the datasets DHFR (Papkou et al., 2023), TrpB (Johnston et al., 2024) and TFBIND8 (Barrera et al., 2016). These datasets contain near complete evaluations of X, and to our knowledge DHFR and TrpB are novel in the machine learning literature. The second (Sec. 4.3) is a more traditional black-box optimization task of finding the maximum of unknown function; using

---

Published as a conference paper at ICLR 2025
Butch [Nue] (R4)
Generative 9(x/9) (R5)
Gradient descent (D2)
Goneral acq arg//r fi.
Rare E S (RI)
Sequential (R3)
Discrete X (R3)
Guaranteed (DI)
Scalable (D3)
Amortization
Method
(Daulton 2022)
BORE (Tiao al., 2021)
Barch BORE (Oliveira al., 2022)
DbAS (Brookes Listgarien, 2018)
X
ChAS (Brookes al., 2019)
X
Amortized RO (Swersky al., 2020)
GetowNets (Jain al., 2022)
DynaPPO (Angermaeller
2019)
AdaLead (Sinai 2020)
PEX (Ren al., 2022)
GGS (Kirjner 2024)
1,50 e.g. (Tripp al., 2020)
LaMBO (Stanton 2023)
LaMBO-2 (Grvver al., 2023)
VSD (ours)
Table 2: Feature table of competing methods: / has feature, X does not have feature, -arttallyy has feature, requires only simple modification. We follow Swersky et al. (2020) in their definition of amortization referring to the ability to use g(x((}_1) for warm-starting the optimization of Φz-
datasets AAV (Bryant et al., 2021), GFP (Sarkisyan et al., 2016) and the biologically inspired Ehrlich
functions (Stanton al., 2024). The corresponding datasets involve | {4, 20), 4 S M < 237 and 65, 000 < | < 20237 We discuss the settings and properties of these datasets in greater detail in Appendix C. For the biological sequence experiments we run predetermined number of experimental rounds, T = 10 or 32 for the Ehrlich functions. We set the batch size to B = 128, and use five different seeds for random initialization. We compare against DbAS (Brookes Listgarten, 2018), CbAS (Brookes et al., 2019), AdaLead (Sinai et al., 2020), and PEX (Ren et al., 2022) all of which we have adapted use CPE, BORE (Tiao et al., 2021) which we have adapted to use the score function gradient estimator, and naïve baseline that uses random samples from the prior, p(x(Do). To reduce confounding, all methods share the same surrogate model, acquisition functions, priors and variational distributions where possible. We compare against LaMBO-2 (Gruver et al., 2023) on the Ehrlich functions, uses its own surrogate and generative models.
4.1 CONDITIONAL GENERATION OF HANDWRITTEN DIGITS
Our motivating application for VSD is to model the space of fit DNA and protein sequences, which are sringrreressenationss of complex 3-dimensional structures. In this experiment we aim to demon- strate, by analogy, that VSD can generate sequences that represent 2-dimensional structures. For this task, we have chosen to 'unroll' (reverse the order of every odd row, and flatten) down-scaled (14 x 14 pixel, 8-bit) MNIST (LeCun et al., 1998) images into sequences, x, where A 196 and IV = 8. We then train long short-term memory (LSTM) recurrent neural network (RNN) and decoder-only causal transformer generative models on the entire MNIST training set by maximum likelihood (ML). These generative distributions are used as the prior models, p(x for VSD and we detail their form in Appendix C.3. The task then to use VSD one round estimate the pos-
9 8
3
0
0
5
0
6
B
4
3
3
3
3
5
3
3
3
3
3
5
S
6
5
0
/
3
7
>
7
9
8
0
7
6
0
5
8
3
5
3
5
5
5
3
3
8
5
5
S
7
4
6
2
2
3
7
4
2
9
8
4
5
3
3
5
3
9
6
5
8
2
6
4
7
6
4
9
3
4
5
5
3
3
r
5
3
5
5
8
9
3
c
5
3
9
9
0
9
4
6
5
0
3
6
5
7
9
3
3
3
3
3
3
5
3
3
3
3
9
3
5
3
6
r
8
7
6
0
r
2
9
3
5
3
3
3
3
3
3
5
6
6
1p
8
4
7
6
3
6
1
6
b
4
4
S
3
3
5
3
3
5
3
3
5
$
3
3
3
9
8
4
4
4
9
? 0
9
0
9
0
$
3
5
5
3
3
5
5
3
3
5
(a) LST Prior
(b) Transformer Prior
(c) LSTT Posterior
(d) Transformer Posterior
Figure 2: (a) and (b) are samples from the LSTM and transformer priors, respectively. (c) and (d) show samples from the L.T and transformer VSD variational distributions respectively. We also report the samples mean scores according to the CPE probabilities.
7

---

Published as a conference paper at ICLR 2025
terior p(xly (3. 5)) using CPE trained on labels 2n 1g [3.5]]. We convolutional architecture for the CPE given in Appendix C.4, and it achieves a test balanced accuracy score of 99%. We parameterize the variational distributions, ((x()), in the same way as the priors, and initialize these distribution parameters from the prior distribution parameters. During training with ELBO the prior distribution parameters are locked, and we run training for 5000 iterations. This is exactly lines 8 and 9 in Algorithm 1. Samples are visualized from the resulting variational distri- butions with the corresponding priors in Figure We see that the prior LSTM and transformer are able to generate convincing digits once the sampled sequences are 're-rolled', and that VSD is able to effectively refine these distributions, even though it does not have access any data directly only scores from the CPE. Both the LSTM and transformer yield qualitatively similar results, and have similar mean scores from the CPE.
4.2 FITNESS LANDSCAPES
In this setting we wish to find fit sequences SSLS, we fix over all rounds. We only consider the combinatorially (near) complete datasets to avoid any pathological behavior from relying on machine learning oracles (Surana et al., 2024). Results presented in Figure 3. The primary measures by which we compare methods are precision, recall and performance,
Precision, min (t,,tt,, |S]
ybr T] [xbr & 1,
(17)
B
Recall, min (rin(TB,SS)
[ybrr > 7] 1 [xb $ X2 x
(18)
Performance, =
Σ> yor [xor & X2 1--1
(19)
bm
Here 2g C X' is the set of experimentally queried sequences by the bth batch member of the rth round, including the initial training set. These measures are comparable among probabilistic and non probabilistic methods. Precision and recall measure the ability a method to efficiently explore S, where min [S]] is the size of the selected set at round t (bounded by the number of good solutions), and min [B [S]] is the number of positive elements possible in the experimental budget. Performance measures the cumulative fitness of the unique batch members, but unlike Jain et al. (2022) we do not normalize this measure.
For exact experimental settings we refer the reader Appendix C.I. We set T to be that of the wild-type sequences in the DHFR and TrpB datasets, and use - = 0.75 for TFBIND8. We find that a uniform prior over sequences, and a mean field variational distribution (Equation 22) are adequate for these experiments, as is simple MLP for the CPE. Results are presented in Figure 3. VSD is the best performing method by most of the measures. We have found the AdaLead and PEX evolutionary-search based methods to be effective on lower-dimensional problems (TFBIND8 being the lowest here), however we consistently observe their performance degrading as the dimension of the problem increases. We suspect this is direct consequence of their random mutation strategies being suited to exploration to the learned generative models employed by VSD, CbAS, and DbAS. Our modified version of BORE (which is just the expected log-likeliin likeelioooo component of Equation 10) performs badly in all cases, and this is direct consequence of its variational distribution collapsing non-batch setting this behavior is not problematic, but shows the importance of the KI. divergence of VSD in this batch setting. We replicate these experiments in Appendix D.1 using GP-PI, also backed by our guarantees. In all cases VSD's results remain similar or improve slightly, whereas the other methods results remain similar or degrade. We report on batch diversity scores Appendix D.3.
low dimensions, but less efficient in higher dimensions compared
point mass. In
4.3 BLACK-BOX OPTIMIZATION
In this experiment we use VSD on the related task of BBO for finding SBBO- We set 7, adaptively by specifying it as an empirical quantile Q', of the observed target values at round t,
(1=
(20)
8

---

Published as a conference paper at ICLR 2025
CAA
Pounal
1200
Sacul
Toand
Roard
Board
Addeed
(a) DHFR
(b) TrpB
(c) TEBINDS
Figure 3: Fitness landscape results. Precision (Equation 17), recall (Equation 18) and performance (Equation 19) higher is better - for the combinatorially (near) complete datasets, DHFR, TrpB and TFBIND8. The random method implemented by drawing samples uniformly.
where Pt-1 percentile from the previous round, and 17 E [0, 1] is an annealing parameter for that trades off exploration and exploitation. Performance is measured by simple regret It, which quantifies the fitness gap between the globally optimal design and the best design found,
B.t max (36s T==1;=11
(21)
Here y" is the fitness value of the globally optimal sequence x". We use the higher dimensional AAV (w==-19544) GFP 411=2.12) and Ehrlich functions (1"=1) datasets/benchmarks to show that VSD can scale to higher dimensional problems. X' AAV and GFP completely intractable fully explore experimentally. and so as the ground-true -tuthh black-box. function. We use the CNN-based oracles from Kirjner et al. (2024) for these experiments. However, we note here that some oracles used these experiments do not predict well caatof-diitributtio (Surana et al., 2024), which limits their real-world applicability. The Ehrlich functions (Stanton al., 2024) are challenging biologically inspired closed-form simulations that cover all X. We compare against genetic algorithm (GA), CbAS and LaMBO-2 (Gruver et al., 2023) for sequences of length M (15, 32, 64) using the POLI and POLI-BASELINES benchmarks and baselines software (González- Duque et al., 2024). For these experiments we use CNNs for the CPEs all experimental settings are in Appendix C.2.
predictive oracle trained on all the original experimental data
The results are summarized in Figure 4 and 5. Batch diversity scores for these experiments are pre- sented in Appendix D.3, and for HOLO Ehrlich function implementation results see Appendix D.2. VSD is among the leading methods for all experiments. VSD takes better advantage of the more complex variational distributions than CbAS and DbAS since it can sample from the adapted vari- ational distribution while learning it. We can see that AdaLead, PEX and often BORE all perform worse than random for reasons previously mentioned. Simple regret can drop below zero for AAV & GFP since oracle is used the black box function. but the global maximizer is taken from the experimental data. VSD outperforms CbAS on the Ehrlich function benchmarks, and competitive with LaMBO-2. We also present ablation study in Appendix D.4.

---

Published as a conference paper at ICLR 2025
Sinple rapars
Giople ooper
GFP
Gingee
RourJ
Tounal
AAV
Strple
Sarple
Round
Tound
Found
(a) Independent
(b) LSTM
(c) Transformer
Figure 4: AAV GFP BBO results. Simple regret (Equation 21) lower better on GFP and AAV with independent and auto-regressive variational distributions. The PEX and AdaLeea results are replicated between the plots, since they are unaffected by choice of variational distribution.
CXP
Smople rogect
Groppe oogge
C:5
C2
13A9 BN
Bound
Round
Rourd
(a) M 15
(b) M 32
(c) M 64
Figure Ehrlich function (POLI implementation) BBO results. VSD and CbAS with different variational distributions: mean field (MF), LSTM and transformer (TFM), compared against genetic algorithm (GA) and LaMBO-2 baselines.
CONCLUSION
We have presented the problem of active generation sequentially learning generative model for designs rare class by efficiently evaluating a black-box function and method for ef- ficiently generating samples which call variational search distributions (VSD). Underpinned by variational inference, VSD satisfies critical requirements and important desiderata, and we show that VSD converges asymptotically to the true level-set distribution at the same rate as a Monte-Carlo estimator with full knowledge of the true distribution. We showcased the benefits of our method empirically on set of combinatorially complete and high dimensional sequential-dest problems and show that it can effectively learn powerful generative models of fit designs. There is close connection between active generation and black-box optimization, and with the advent of pow- erful generative models we hope that our explicit framing of generation of fit sequences will lead to further study of this connection. Finally, our framework can be generalized to more complex applica- tion scenarios, involving learning generative models over Pareto sets, Spareto in multi-objective set- ting, or other challenging combinatorial optimization problems (Bengio et al., 2021). such as graph structures (Annadani et al., 2023), and mixed investigating as future work directions. For the code implementing the models and experiments in this paper, please https: //gi thub. com cs ro- funml /va riat ona sea rch.
biological
variables. All of which are worth
10

---

Published as a conference paper at ICLR 2025
ACKNOWLEDGMENTS
This work funded by the CSIRO Science Digital and Advanced Engineering Biology Future Sci- ence Platforms, and was supported by resources and expertise provided by CSIRO IMT Scientific Computing. We would like to thank the anonymous reviewers, and especially reviewer xamp for the constructive feedback and advice, which greatly increased the relevance and quality of this work.
REFERENCES
Yasin Abdasi-Yadkori. Online Learning for Linearly Parametrized Control Problems. Phd, Univer- sity of Alberta, 2012.
Christof Angermueller, David Dohan, David Belanger, Ramya Deshpande, Kevin Murphy, and Lucy Colwell. Model-basee reinforcement learning for biological sequence design. In International conference on learning representations (ICLR), 2019.
Yashas Annadani, Nick Pawlowski, Joel Jennings, Stefan Bauer, Cheng Zhang, and Wenbo Gong. BayesDAG: Gradient-based posterior inference for causal discovery. Thirty-seventh Confer- ence on Neural Information Processing Systems (NeurIPS), 2023.
Maximilian Balandat, Brian Karrer, Daniel R. Jiang, Samuel Daulton, Benjamin Letham, An- drew Gordon Wilson, and Eytan Bakshy. BoTorch: A Framework for Efficient Monte-Carlo Bayesian Optimization. In Advances in Neural Information Processing Systems 33, 2020.
Luis A Barrera, Anastasia Vedenko, Jesse V Kurland, Julia M Rogers, Stephen S Gisselbrecht, Elizabeth J Rossin, Jaie Woodard, Luca Mariani, Kian Hong Kock, Sachi Inukai, et al. Survey of variation in human transcription factors reveals prevalent dna binding changes. Science, 351 2016.
(6280):1450-11))
David G. T. Barrett and Benoit Dherin. Implicit gradient regularization. In International Conference on Learning Representations (ICLR). OpenReview, 2021.
Heinz Bauer. Probability theory and elements of measure theory. Academic Press, 2nd edition, 1981.
Yoshua Bengio. Andrea Lodi, and Antoine Prouvost. Machine learning for combinatorial opti- mization: A methodological tour d d'horizon. European Journal of Operational Research, 290(2): 405-421, 2021. ISSN 0377-2217.
Thomas Bird, Julius Kunze, and David Barber. Stochastic variational optimization. arXiv preprint arXiv:1809.0483 04855, 2018
Ilija Bogunovic, Jonathan Scarlett, Andreas Krause, and Volkan Cevher. Truncated variance reduc- tion: A unified approach to bayesian optimization and level-set estimation. Advances in neural information processing systems, 29, 2016.
Stéphane Boucheron, Gábor Lugosi, and Pascal Massart. Concentration inequalities: A Nonasymp- totic Theory of Independence. Oxford University Press, 2013.
David Brookes, Hahnbeom Park, and Jennifer Listgarten. Conditioning by adaptive sampling for robust design. In International conference on machine learning, pp. 773-782, PMLR, 2019.
David Brookes, Akosua Busia, Clara Fannjiang, Kevin Murphy, and Jennifer I Listgarten. A view of estimation of distribution algorithms through the lens of expectationmmxxiiizztoon In Proceed- ings of the 2020 Genetic and Evolutionary Computation Conference Companion, pp. 189-190, 2020.
David H Brookes and Jennifer Listgarten. 2018.
Design by adaptive sampling.
arXiv preprint
Brent Bryan, Robert C Nichol, Christopher R Genovese, Jeff Schneider, Christopher Miller, and Larry Wasserman. Active learning for identifying function threshold boundaries. Advances in neural information processing systems, 18, 2005.

---

Published as a conference paper at ICLR 2025
Drew H Bryant, Ali Bashir, Sam Sinai, Nina K Jain, Pierce J Ogden, Patrick F Riley, George M Church, Lucy Colwell, and Erie D Kelsic. Deep diversification of an aav capsid protein by machine learning. Nature Biotechnology. 39(6):691-696, 2021.
Lin Chen and Sheng Xu. Deep neural tangent kernel and laplace kernel have the same RKHS. In International Conference on Learning Representations, 2021.
Sayak Ray Chowdhury and Aditya Gopalan, On Kernelized Multi-armed Bandits. In Proceedings of the 34th International Conference on Machine Learning (ICML), Sydney, Australia, 2017.
Sumanth Dathathri, Andrea Madotto, Janice Lan, Jane Hung, Eric Frank, Piero Molino, Jason Yosin- ski, and Rosanne Liu. Plug and play language models: A simple approach to controlled text generation. In International Conference on Learning Representations, 2020.
Samuel Daulton, Xingchen Wan, David Eriksson, Maximilian Balandat, Michael A Osborne, and Eytan Bakshy. Bayesian optimization over discrete and mixed spaces via probabilistic reparam- eterization. Advances in Neural Information Processing Systems (NeurIPS), 35:12760- 2022.
Lester E. Dubins and David A. Freedman. A sharper form of the Borel-Cantelli lemma and the strong law. The Annals of Mathematical Statistics, 36(3):800-807 1965.
Audrey Durand, Odalric-Ambrym Maillard, and Joelle Pineau. Streaming kernel regression with provably adaptive mean, variance, and regularization. Journal of Machine Learning Research, 19, 2018.
Rick Durrett. Probability: theory and examples. Cambridge University Press, New York, NY, 5th edition, 2019.
Ronen Eldan, Dan Mikulincer, and Tselil Schramm. Non-asymptotic approximations of neural networks by Gaussian processes. In 34th Annual Conference on Learning Theory (COLT), volume 134 of Proceedings of Machine Learning Research, 2021.
Henry E. Fleming. Equivalence of regularization and truncated iteration in the solution of ill-posed image reconstruction problems. Linear Algebra and Its Applications, 130(C): 1990.
E. Portfa-Portugués. Notes for Nonparametric Statistics. Bookdown, 2024. URL https://
bookdown. org/ /tt///////pp////////////// /NPP-CNM/ Version 6.9.1. ISBN 978- -79-8-----2953711
Roman Garnett. Bayesian optimization. Cambridge University Press, 2023.
Roman Garnett, Yamuna Krishnamurthy, Xuehan Xiong, Jeff G. Schneider, and Richard P. Mann. Bayesian optimal active search and surveying. Proceedings of the 29th International Confer- ence on Machine Learning, ICML 2012, Edinburgh, Scotland, UK, June 26 - July 1, 2012. icml.cc / Omnipress, 2012.
Tilmann Gneiting and Adrian E Raftery. Strictly proper scoring rules, prediction, and estimation. Journal of the American statistical Association, 102(477):359-3777 2007.
Rafael Gómez-Bombarede Jennifer N Wei, David Duvenaud, José Miguel Hernández-Lobato, Benjamín Sánchez-Lengeling, Dennis Sheberla, Jorge Aguilera-Iparraguine Timothy D Hirzel, Ryan P Adams, and Alán Aspuru-Guzik. Automatic chemical design using data-driven contin- uous representation of molecules. ACS central science, 4(2):268 276, 2018.
Miguel González-Duque, Richard Michael, Simon Bartels, Yevgen Zainchkovskyy, Søren Hauberg, and Wouter Boomsma. A survey and benchmark of high-dimensional Bayesian optimization of discrete sequences. arXiv preprint arXiv:2406.04473 04739, 2024.
Miguel González-Duque, Simon Bartels, and Richard Michael. poli: libary of discrete sequence objectives, January 2024.
URL. https:://gt http//gttbtt.tt////// com/
poli.
Alkis Gotovos, Nathalic Casati, Gregory Hitz, and Andreas Krause. Active learning for level set estimation. In Proceedings of the Twenty-Thing international joint conference on Artificial Intel- ligence, pp. 1350,, 2013.
12

---

Published as a conference paper at ICLR 2025
Nate Gruver, Samuel Stanton, Nathan Frey, Tim GJ Rudner, Isidro Hotzel, Julien Lafrance-Vanasse, Arvind Rajpal, Kyunghyun Cho, and Andrew G Wilson. Protein design with guided discrete diffusion. Advances in neural information processing systems (NeurIPS), 36, 2023.
Kurt Homik, Maxwell Stinchcombe, and Halbert White. Multilayer feedforward networks uni- versal approximators. Neural Networks, 2(5):359-366, 1989.
Jiri Hron, Yasaman Bahri, Jascha Sohl-Dickstein, and Roman Novak. Infinite attention: NNGP and NTK for deep attention networks. In Proceedings of the 37th International Conference on Machine Learning (ICML). PMLR 108, 2020.
Arthur Jacot, Franck Gabriel, and Clément Hongler. Neural tangent kernel: Convergence and gen- eralization in neural networks. In Advances Neural Information Processing Systems, Montreal, Canada, 2018
Moksh Jain, Emmanuel Bengio, Alex Harnandez-GGrria, Jarrid Rector-Brooks, Bonaventure FP Dossou, Chanakya Ajit Ekbote, Jie Fu, Tianyu Zhang. Michael Kilgour, Dinghuai Zhang, et al. Biological sequence design with gflownets. In International Conference on Machine Learning, pp. 9786-9801. PMLR, 2022.
Shali Jiang, Gustavo Malkomes, Geoff Converse, Alyssa Shofner, Benjamin Moseley, and Roman Garnett. Efficient nonmyopic active search. In International Conference on Machine Learning, pp. 1714-1723. PMLR, 2017.
Kadina E. Johnston, Patrick J. Almhjell, Ella J. Watkins-Dulaney, Grace Liu, Nicholas J. Porter, Jason Yang, and Frances H. Arnold. combinatorially complete epistatic fitness landscape in an enzyme active site. Proceedings of the National Academy of Sciences, 121(32):e240044911, 2024. doi: 2...0733/paas24404991222
Diederik Kingma and Jimmy Ba. Adam: method for stochastic optimization, arXiv preprint arXiv:1412.6958 2014.
Diederik P. Kingma and Max Welling. Auto-encoding variational bayes. In Yoshua Bengio and Yann LeCun (eds.), 2nd International Conference on Learning Representations, ICLR 2014, Banff. AB, Canada, April 14-16, 2014, Conference Track Proceedings, 2014.
Andrew Kirjner, Jason Yim, Raman Samusevich, Shahar Bracha, Tommi S Jaakkola, Regina Barzi- lay, and lla Fiete. Improving protein optimization with smoothed fitness landscapes. In The Twelfth International Conference Learning Representations, 2024.
Harold J Kushner. A new method of locating the maximum point of an arbitrary multipeak curve in the presence of noise. Journal of Basic Engineering, 1964.
Pedro Larrañaga and Jose A Lozano. Estimation of distribution algorithms: A new tool for evolu- fionary computation, volume 2. Springer Science & Business Media, 2001.
Yann LeCun, Léon Bottou, Yoshua Bengio, and Patrick Haffner. Gradient-based learning applied to document recognition. Proceedings of the IEEE, 86(1) (6(1))2227--224 1998.
Jachoon Lee, Lechao Xiao, Samuel S. Schoenholz, Yasaman Bahri, Roman Novak, Jascha Sohl- Dickstein, and Jeffrey Pennington. Wide neural networks of any depth evolve linear models under gradient descent. In 33rd Conference on Neural Information Processing Systems (NeurIPS 2019), Vancouver, Canada, 2019.
Xiang Li, John Thickstun, Ishaan Guirajani, Percy S Liang, and Tatsunori B Hashimoto. Diffusion- Im improves controllable text generation. Advances in Neural Information Processing Systems, 35:4328-4333 2022.
Chaoyue Liu, Libin Zhu, and Mikhail Belkin. On the lincarity of large non-linear models: When and why the tangent kemel is constant. In Advances in Neural Information Processing Systems, 2020.
Qiang Liu and Dilin Wang. Stein variational gradient descent: general purpose bayesian inference algorithm. Advances neural information processing systems, 29, 2016.
13

---

Published as a conference paper at ICLR 2025
Odalric-Ambrym Maillard. Self- normaliaaiiino techniques for streaming confident regression. HAL, <hal-01349, 2016.
Richard Michael, Simon Bartels, Miguel Sonzález-Duque, Yevgen Zainchkovskyy. Jes Frellsen, Søren Hauberg, and Wouter Boomsma. A continuous relaxation for discrete Bayesian optimiza- tion. arXiv preprini urXiv:2404. /7452, 2024.
Shakir Mohamed, Mihaela Rosca, Michael Figurnov, and Andriy Mnih. Monte carlo gradient esti- mation in machine learning. Journal of Machine Learning Research, 21(132):1-62, 2020.
Mattes Mollenhauer, Ingmar Schuster, Stefan Klus, and Christof Schütte. Singular value decompo- sition of operators on reproducing kernel hilbert spaces. In Oliver Junge, Oliver Schütze, Gary Froyland, Sina Ober-Blöbaun and Kathrin Padberg- Gehle (eds.), Advances in Dynamics, Opti- mization and Computation, pp. 109-131, Cham, 2020. Springer International Publishing.
Michael Murray, Hui Jin, Benjamin Bowman, and Guido Montufar. Characterizing the spectrum of the NTK via a power series expansion. In International Conference on Learning Representations (ICLR), Kigali, Rwanda, 2023.
Mojmír Mutný and Andreas Krause. Efficient high dimensional Bayesian optimization with addi- tivity and quadrature Fourier features. In Conference on Neurul Information Processing Systems (NeurIPS), Montreal, Canada, 2018
Anh Nguyen, Jeff Clune, Yoshua Bengio, Alexey Dosovitskiy, and Jason Yosinski. Plug & play generative networks: Conditional iterative generation of images in latent space. In Proceedings of the IEEE conference on computer vision and pattern recognition, pp. 4467-4477,200
Rafael Oliveira, Louis Tiao, and Fabio T Ramos. Batch bayesian optimisation via density-ratio estimation with guarantees. Advances in Neural Information Processing Systems (NeurIPS), 35: 29816-29829, 2022.
Andrei Papkou, Lucia Garcia-Pastor, José Antonio Escudero, and Andreas Wagner. rugged yet easily navigable fitness landscape. Science, 382(6673))cadh3866)) 2023.
Mary Phuong and Marcus Hutter.
Formal
algorithms for transformers.
arXiv preprint
arXiv:2207.09233, 2022.
Rajesh Ranganath, Sean Gerrish, and David Blei. Black box variational inference. In Artificial intelligence and statistics, pp. 814-822. PMLR, 2014.
Carl E. Rasmussen and Christopher K. I. Williams. Gaussian Processes for Machine Learning. The MIT Press, Cambridge, MA, 2006.
Zhizhou Ren, Jiahan Li, Fan Ding, Yuan Zhou, Jianzhu Ma, and Jian Peng. Proximal exploration for model-guided protein sequence design. In International Conference on Machine Learning, pp. 18520-18536, PMLR, 2022.
Reuven Rubinstein. The cross-entropy method for combinatorial and continuous optimization. Methodology and computing in applied probability, 1:127-190, 1999.
Karen S Sarkisyan, Dmitry A Bolotin, Margarita V Meer, Dinara R Usmanova, Alexander S Mishin, George V Sharonov, Dmitry N Ivankov, Nina Bozhanova, Mikhail Baranov, Onuralp Soyle- mez, et al. Local fitness landscape of the green fluorescent protein. Nature, 533(7603):397 2016.
Bernhard Schölkopf and Alexander J. Smola. Learning with Kernels: Support Vector Machines, Regularization, Optimization, and Beyond. The MIT Press, 2001.
Shayle R Searle. Matrix algebra useful for statistics. John Wiley & Sons, 1982.
Sam Sinai, Richard Wang. Alexander Whatley, Stewart Slocum, Elina Locane, and Eric D Kelsic. Adalead: A simple and robust adaptive greedy search algorithm for sequence design. arXiv preprint arXiv:2010.02144: 2020.
14

---

Published as a conference paper at ICLR 2025
Samuel L. Smith, Benoit Dherin, David G. T. Barrett, and Soham De. On the origin of implicit regularization stochastic gradient descent. In Proceedings of the 10th International Conference on Learning Representations (ICLR 2021). OpenReview, 2021.
Jiaming Song, Lantao Yu, Willie Neiswanger, and Stefano Ermon. A general recipe for likelihood- free bayesian optimization. In International Conference on Machine Learning, pp. 80384-20404 PMLR, 2022.
Daniel Soudry, Elad Hoffer, Mor Shpigel Nacson, Suria Gunasekar, and Nathan Srebro. The implicit bias of gradient descent separable data. Journal of Machine Learning Research, 19, 2018.
Niranjan Srinivas, Andreas Krause, Sham Kakade, and Matthias Seeger. Gaussian process optimiza- tion in the bandit setting: no regret and experimental design. In Proceedings of the 27th Interna- tional Conference on International Conference on Machine Learning, ICML' 10, pp. 1015-1022, Madison, WI, USA, 2010. Omnipress. ISBN 978 605589077.
Joe Staines and David Barber. Optimization by variational bounding. In ESANN, 2013.
Samuel Stanton, Wesley Maddox, Nate Gruver, Phillip Maffettone, Emily Delaney, Peyton Green- side, and Andrew Gordon Wilson. Accelerating Bayesian optimization for biological sequence design with denoising autoencoders. In International Conference on Machine Learning, pp- PMLR, 2022.
Samuel Stanton, Robert Alberstein, Nathan Frey, Andrew Watkins, and Kyunghyun Cho. Closed-form test functions for biophysical sequence optimization algorithms. arXiv preprint arXiv:2407.00236, 2024.
Shikha Surana, Nathan Grinsztajn, Timothy Atkinson, Paul Duckworth, and Thomas Barrett. Overconfident oracles: Limitations of in silico sequence design benchmarking. In ICML 2024 Al for Science Workshop, 2024.
Kevin Swersky, Yulia Rubanova, David Dohan, and Kevin Murphy. Amortized bayesian optimiza- tion discrete spaces. In Conference Uncertainty in Artificial Intelligence, pp. 769-778. PMLR, 2020.
Shion Takeno, Yu Inatsu, Masayuki Karasuyama, and Ichiro Takeuchi. Posterior sampling-base based Bayesian optimization with tighter Bayesian regret bounds. Proceedings of the 41 st Inter- national Conference on Machine Learning (ICML 2024), volume 235, Vienna, Austria, 2024. PMLR.
Louis C Tiaa,, Aaron Klein, Matthias W Seeger, Edwin V Bonilla, Cedric Archambeau, and Fabio Ramos. Bore: Bayesian optimization by density-ratio estimation. In International Conference on Machine Learming, pp. 10289-10300, PMLR, 2021.
Brandon Trabucco, Xinyang Geng, Aviral Kumar, and Sergey Levine. Design-bench: Benchmarks for data-driven offline model-based optimization. CoRR, abs/2202.084500 2022.
Austin Tripp, Erik Daxberger, and José Miguel Hernández-Loobtoo. Sample-efficient optimization in the latent space of deep generative models via weighted retraining. Advances in Neural Infor- mation Processing Systems, 33:1 1259-11272. 2020.
Sattar Vakili and Julia Olkhovskaya. Kernelized reinforcement learning with order optimal regret bounds. In Advances in Neural Information Processing Systems (NeurIPS), volume 36, 2023.
Sattar Vakili, Kia Khezeli, and Victor Picheny. On information gain and regret bounds in Gaus- sian process bandits. In Arindam Banerjee and Kenji Fukumizu (eds.), Proceedings of The 24th International Conference on Artificial Intelligence and Statistics, volume 130 of Proceedings of Machine Learning Research, pp. 82-90. PMLR, 2021.
Hastagiri P Vanchinathan, Andreas Marfurt, Charles-Antoine Robelin, Donald Kossmann, and An- dreas Krause. Discovering valuable items from massive data. In Proceedings of the 21th ACM SIGKDD International Conference on Knowledge Discovery and Data Mining, pp. 1204, 2015.
15

---

Published as a conference paper at ICLR 2025
Daan Wierstra, Tom Schaul, Tobias Glasmachers, Yi Sun, Jan Peters, and Jürgen Schmidhuber. Natural evolution strategies. The Journal of Machine Learning Research, 15(1):949-995 980, 2014.
Ronald J Williams. Simple statistical gradient-followngg algorithms for connectionist reinforcement learning, Machine learning, 8:229-256, 1992.
James T Wilson, Riccardo Moriconi, Frank Hutter, and Mare Peter Deisenroth. The reparameteri- zation trick for acquisition functions. arXiv preprint arXiv: 17Xiv:1712.00424, 2017.
Yuan Yao, Lorenzo Rosasco, and Andrea Caponnetto. On early stopping in gradient descent learn- ing. Constructive Approximation, 26(2):289-315, 2007.
Fengxue Zhang, Jialin Song, James C Bowden, Alexander Ladd, Yisong Yue, Thomas Desautels, and Yuxin Chen. Learning regions of interest for Bayesian optimization with adaptive level-set estimation. In Andreas Krause, Emma Brunskill, Kyunghyun Cho, Barbara Engelhardt, Sivan Sabato, and Jonathan Scarlett (eds.). Proceedings of the 40th International Conference on Ma- chine Learning, volume 202 of Proceedings of Machine Learning Research, pp. 41579-41595. PMLR, 23-29 Jul 2023a.
Hanqing Zhang, Haolin Song, Shaoyu Li, Ming Zhou, and Dawei Song. survey of controllable text generation using wransfomer-basee pre-trained language models. ACM Computing Surveys, 56(3):1-37, 2023b.
16

---

Published as a conference paper at ICLR 2025
A
ACRONYMS
ACRONYMS
BBO black-box optimization. 1. 5, 8, 10, 16, 19-21, 23
BO Bayesian optimization. 3,5, 6, 35
BOPR Bayesian optimization with probabilistic reparametrization. 5-7,
BORE Bayesian optimization by density-ratio estimation. 4, 6-9, 16-18, 20, 22
CbAS conditioning by adaptive sampling. 5-10, 16-18, 18, 20, 22
CPE class probability estimation. 4-9, 16-19, 21, 22, 35
DbAS design by adaptive sampling. 5-9, 16-18, 22
EDA estimation of distribution algorithms.
ELBO evidence lower bound. 3, 4.
ES evolution strategies.
GA genetic algorithm. 9
GP Gaussian process. 3-5, 8, 18, 20, 39
KL Kullback-Leibler. 3, 6, 8
LSE level set estimation.
LSO latent space optimization. 6,7
LSTM long short-term memory. 7, 8, 10, 17, 18, 20-23
ML maximum likelihood. 7,
NES natural evolution strategies.
NTK neural tangent kernel.
PEX proximal exploration. 6-10, 17, 18,
PI probability of improvement, 3-6, 8, 18, 20, 29, 35
RNN recurrent neural network. 7, 18
UCB upper confidence bound.
VI variational inference.
VO variational optimization. 5
VSD variational search distributions. 1-10, 16-18, 20-23, 29, 31, 35, 36
17

---

Published as a conference paper at ICLR 2025
B DEPICTION OF ACTIVE GENERATION
See Figure B.1 for graphical depictions of active generation as implemented by VSD, compared to batch BO. Active generation using VSD follows Algorithm I sequentially approximate p(xly > T). It uses samples from the current learned approximation of this distribution, g(x [x]) for proposing candidates to evaluate each round. Unmodified batch BO in the discrete setting without any specialization requires a list of candidates, from which a batch of candidates are selected per round using surrogate model with a batch acquisition function, e.g. see Wilson et al. (2017). The surrogate model's hyper-parametens, 0, are estimated by minimizing negative log marginal likeli- hood, LNLML Dx). Mcchanistically, active generation learns a generative model of valuable can- didates to circumvent the requirement of having to select candidates from a list. This is especially important for searching over the space of sequences, X, where enumerating feasible candidates is often intractable. Furthermore, active generation naturally lends itself to large batch sizes, while explicit batch optimization is often computationally intensive and is limited to (110) candidates per batch Wilson al. (2017). As mentioned in section alternatives active generation for spe- cializing batch BO to the discrete domain also include LSO (Tripp et al., 2020) and amortized BO (Swersky al., 2020) among others.
Medd tolsing
Medi bandin
ErestA PL)
PA
Fly
Unaadded
Libelal
Lanciel
Papeiren 5000
Dxperiment Aodt
(a) VSD Active Generation
(b) Batch BO
Figure B.I: Depictions of (a) active generation as implemented by VSD, and (b) batch Bayesian optimization as applied to discrete sequences. Please see the text for a discussion of the differences between these approaches.
C EXPERIMENTAL DETAILS
We use three well established datasets; green fluorescent protein (GFP) from Aequorea Victo- ria (Sarkisyan et al., 2016), an adeno-associated virus (AAV) Bryant et al. (2021) I); and DNA binding activity to human transcription factor (TFBIND8) (Trabucco et al., 2022; Barrera et al., 2016). These datasets have been used variously by Brookes & Listgarten (2018); Brookes et al. (2019); Angermueller al. (2019); Kirjner al. (2024); Jain et al. (2022) among others. The GFP task is to maximize fluorescence, this protein consists of 238 amino acids, of which 237 can mutate. The AAV task to maximize the genetic payload that can be delivered, and the associated protein has 28 amino acids, all of which can mutate. A complete combinatorial assessment is infeasible for these tasks, and so use the convolution neural network oracle presented in Kirjner al. (2024) as in-silico ground truth. TFBIND8 contains complete combinatorial assessment of the effect of changing 8 nucleotides on binding to human transcription factor SIX6 REF R1 (Barrera al., 2016). The dataset we use contains all 65536 sequences prepared by Trabucco et al. (2022).
We also use two novel datasets from recent works that experimentally assess the (near) complete combinatorial space of short sequences. The first dataset measures the antibiotic resistance of Es- cherichia coli metabolic gene folA, which encodes dihydrofolate reductase (DHFR) (Papkou et al., 2023), Only sub-sequence of this gene is varied (9 nucleic acids which encode 3 amino acids), and so a near-complete (99.7%) combinatorial scan is available. For variants that have no fitness (resistance) data available, give score of -1. The next dataset is near-complete combinatorial scan of four interacting amino acid residues near the active site of the enzyme tryptophan synthase (TrpB) (Johnston et al., 2024), with 159,129 unique sequences and fitness values, we use -0.2 for
18

---

Published as a conference paper at ICLR 2025
the missing fitness values (we do not the authors' imputed values). These residues explicitly shown to exhibit epistasis - or non non-additive effects on catalytic function which makes navigating this landscape more interesting challenge from an optimization perspective.
Finally, we use the recently proposed Ehrlich functions (Stanton et al., 2024) benchmark. These functions are challenging closed form biological analogues, specifically designed to test BBO meth- ods on high dimensional sequence design tasks without having to resort to physical experimenta- tion or machine learning oracles. We use the POLI and COLI-BASELINS software package for the benchmark and baselines (González-Duque et al., 2024), and test on both the original HOLO imple- mentation (Stanton et al., 2024) as well as the native POLI implementation of these functions.
The properties of these datasets and benchmarks are presented in Table C.I.
Dataset
M
TFBIND8
65,536
65,536
TrpB
20
159,129
160,000
DHFR
9
261,333
262,144
AAV
20
28
42,340
2028
GFP
20
237
51,715
20
Ehrlich-11 15
20
15
20
20
Ehrlich-33 32
20
32
200a
20
Ehrlich-64
20
64
2084
2064
Table C.I: Alphabet size, sequence length, and number of available sequences for each of the datasets use in this work.
We optimize VSD, CbAS, DbAS and BORE for minimum of 3000 iterations each round (5000 for all experiments but the Ehrlich functions) using Adam (Kingma & Ba, 2014). When we use a CPE, AdaLead's x parameter set to 0.5 since the CPE already incorporates the appropriate threshold.
C.I FITNESS LANDSCAPES SETTINGS
For the DHFR and TrpB experiments we set maximum fitness in the training dataset to be that of the wild type, and 7 to be slightly below the wild type fitness value (so we have 10 positive examples to train the CPE with). We use a randomly selected Nusin = 2000 below the wild-type fitness to initially train the CPE, we also explicitly include the wild-type. The thresholds and wild- type fitness values are; DHRF: T = -0.1, V/w O, TrpB: T = 0.35, V/w = 0.409. We follow the same procedure for the TFBIND8 experiment, however, there is no notion of wild-type sequence in this data, and so we set T = 0.75, and Ytrain max = 0.85. We use a uniform prior over sequences, p(x) = - Categ(Em.) 1 |V]- since these are relatively small search spaces, and the sub- sequences of nucleic/amino acids have been specifically selected for their task. Similarly, we find that relatively simple independent (mean-field) variational distributions of the form in Equation 22 and MLP based CPEs work best for these experiments (details in Sec. C.4).
C.2 BLACK-BOX OPTIMIZATION SETTINGS
We follow Kirjner et al. (2024) in the experimental settings for the AAV and GFP datasets, but we modify the maximum fitness training point and training dataset sizes make them more amenable to sequential optimization setting. The initial percentiles, schedule, and max training fitmess values are; AAV: Po 0.8. 17 = 0.7, ymax 5. GFP: Po 0.8, 17 0.7 ymax = 1.9. We aim for pr 0.99. The edit distance between x" and the fittest sequence in the CPE training data is 8 for GFP, and 13 for AAV. We again use random Nirain = 2000 for training the CPEs, which this case are CNNs architecture specifics are in Sec. C.4.
For the Ehrlich function experiment, we use sequence lengths of M = (15,32, 64] with motifs for the shorter sequence lengths, and 8 motifs for M 64. All use motif length of 4 and quantization of 4. B = 128, T = 32 and only 128 random samples of the function are used for D these resampled for each seed. As before, 5 different random seeds are used for these trials, and for VSD we use an the same scheduling function for т, as in Equation 20, with Po = 0.5 and n 0.87 (so Pr 0.99). The lower initial percentile is used since the training dataset is much
19

---

Published as a conference paper at ICLR 2025
smaller than in the other experiments, and find allowing for more exploration initially improves VSD's performance.
In these higher dimensional settings, find that performance of the methods heavily relies using an informed prior (in the case VSD and CbAS), or initial variational distribution (in the case of DbAS and BORE). To this end, we follow Brookes al. (2019) and fit the initial variational distribution to the CPE training sequences (regardless fitness), but we use maximum likelihood. For the more complex variational distributions (LSTM and transformer), to over-fit - so we implement early stopping and data augmentation techniques. Then for VSD and CbAS we copy this distribution and fix its parameters for the remainder of the experiment for use as a prior. We also use this prior for the Random method, but AdaLead and PEX use alternative generative heuristies. For these experiments we use the simple independent variational distribution and the same LSTM and causal decoder-only transformer models from Sec. 4.1.
have be careful
C.3 VARIATIONAL. DISTRIBUTIONS
In this section we summarize the main variational distribution architectures considered for VSD, BORE, CbAS and DbAS, and the sampling distributions for the Random baseline method. Some- what surprisingly, we find that we often obtain good results for the biological sequence experiments using simple independent (or mean-field) variational distribution, especially in lower dimensional settings,
M (x(()) = II m=1
(22)
where E V and 6m E R Vl, However, this simple mean-field distribution was not capable of generating convincing handwritten digits or, in some cases, higher-dimensional sequences. We have also tested variety of transition variational distributions,
M x1- = II m=1
Categ (zem /stts//assssasasooo (NN (X
(23)
where NNm (x;-1,,) 6) is the mh vector output of a neural network that takes a sequence from the previous round, X;-1, as input. We have implemented multiple neural net encoder/decoder archi- tectures for NN Ø). but we did not consider architectures of the form NN,,. (6) since the variational distribution in Equation 22 can always learn a Φm = NN (d'). We found that none of these transition architectures significantly outperformed the mean-field distribution (Equation 22) when it was initialized well (e.g. fit to the CPE training sequences), see Sec. D.4 for results. We also implemented regreseegees variational distributions of the form.
M
g(x () (solegg(x/sofmmxx())) II ( [1 v1m)) where,
(24)
6m- (1m)))
q(mm 111m-1, 21m)) O
(soft/// (DTransformerr( (11mm 1:mm111
For LSTM RNN and decoder-only transformer with causal mask, for the latter Phuong & Hutter (2022, Algorithm 10 & Algorithm 14) for maximum likelihood training and sampling imple- mentation details respectively. We list the configurations of the LSTM and transformer variational distributions in Table C.2. We use additive positional encoding for all of these models. When using these models for priors or initialization of variational distributions, we find that over-fitting can be an issue. To circumvent this, we use early stopping for larger training datasets, or data augmentation techniques for smaller training datasets (as the case of the Ehrlich functions).
C.4 CLASS PROBABILITY ESTIMATOR ARCHITECTURES
For the fitness landscape experiments the smaller combinatorially complete datasets use two-hidden layer MLP, with an input embedding layer. The architecture is given in Figure C.2 (a). For the larger dimensional AAV and GFP datasets and Ehrlich function benchmark, use the
20

---

Published as a conference paper at ICLR 2025
Configuration
Digits
AAV
GFP
Ehrlich 15
Ehrlich 32
Ehrlich 64
LSTM
Layers
5
4
4
3
3
Network size
128
32
32
32
64
Embedding size
4
10
10
10
10
10
Transformer
Layers
1
2
2
Network Size
256
64
64
32
64
128
Attention heads
8
2
2
2
3
Embedding size
32
20
20
20
30
Table C.2: LSTM and transformer network configuration.
convolutional architecture given in Figure C.2 (b). On all but the Ehrlich benchmark, five fold cross validation was used to select the hyper parameters before the CPEs are trained on the whole training set for use in the subsequent experimental rounds. For the Ehrlich benchmark do not validation to select the CPE hyper parameters but we do use an additive ensemble of 10 randomly initialized CNNs for the CPE following LaMBO-2. Model updates performed by retraining the whole query set.
D ADDITIONAL EXPERIMENTAL RESULTS
D.1 FITNESS LANDSCAPES GAUSSIAN PROCESS PROBABILITY OF IMPROVEMENT
Here we present additional fitness landscape experimental results, where we have used GP as surrogate model for p(y|x, x, Dv) in conjunction with complementary Normal CDF as the PI acquisition function. This is one of the main frameworks supported by our theoretical analysis. VSD, DbAS, CbAS and BORE can make use of the GP-PI acquisition function, and so BORE is BOPR (Daulton et al., 2022) in this instance since we are not using CPE. PEX and AdaLead only use the GP surrogate, as per their original formulation. The GP uses simple categorical kernel with automatic relevance determination from Balandat et al. (2020),
1 M
Fm mmm]
k(x, x') exp
(25)
1m
where and 4m are hyper-parameters controlling scale and length-scale respectively. See Figure D.3 for the results.
D.2 EHRLICH FUNCTION HOLO RESULTS
See Figure D.4 for BBO results on the original HOLО Ehrlich function implementation (Stanton al., 2024). We present additional diversity scores for these and the POLI implementation in Sec. D.3.
D.3 DIVERSITY SCORES
The diversity of batches of candidates is common thing to report in the literature, and to that end we present the diversity of our results here. We have taken the definition pair-wise diversity from (Jain et al., 2022) as,
Diversity, = B(B - X.ED Ex.
Σ Σ
Lev(x;, x;),
(26)
where Lev X' + No is the Levenshtein distance. caution the reader as to the interpretation of these results however, as more diverse batches often do not lead to better performance, precision, recall or simple regret (as can be seen from the Random method results). Though insufficient diver- sity also explain poor performance, experiment are presented in Figure D.5, and black-box optimization for AAV & GFP in Figure 6.6 and Ehrlich functions in Figure D.7.
in the case of BORE. Results for the fitness landscape
21

---

Published as a conference paper at ICLR 2025
Sequential
Sequential(
Embedding num num_embeddiinge n enbedding_dimes
Embedding
num_ _um_enbeddingeaa, umbeddinggdim-11
Drepout (p=C.2),
Bropout (p-0.2), Convld
Flatten (),
LeakyReIU(), Linear
in_channels=10, out_channels-16, hernel_size=3 7,
in_features-8 + M, out_features=32
) LeakyReLU MaxPoolld(
LeakyReLU(), Linear
kernel_size=2 4, stride=2 4,
fe_ fut_features=1
Convid
in_channels=16, out_chhamees==16,
kernel size=7,
LeakyReLU MaxPoolld
kernel_size=2 stride=2 or 4,
Flatten LazyLinear( out_features=128
TeakyRelu li ninear
in iettureers1
out fat_uretu
(a) MLP architecture
(b) CNN architecture
Figure C.2: CPE architectures used for the experiments in PyTorch syntax. A = V = M, GFP uses max pooling kernel size and stride of 4, all other datasets and benchmarks 2. The Ehrlich function benchmark uses and ensemble 10 randomly initialized CNNs that are additively combined. The Ehrlich-15 functions use kernel size of 3, all other BBO experiments use kernel size 7. LaMBO-2 uses the same kernel size as our CNNs for the Ehrlich functions.
D.4 ABLATIONS VARIATIONAL AND PRIOR DISTRIBUTIONS
In Figure D.8 we present ablation results for VSD using different priors and variational distributions. We use the BBO experimental datasets for this task as they are higher-dimensional and so more sensitive to these design choices. We test the following prior and variational posterior distributions:
IU
Independent categorical variational posterior distribution of the form in Equation 22, and uniform prior distribution, p(x) II Categ(
-
Independent categorical prior and variational posterior of the form in Equation 22, The prior is fit using ML on the initial CPE training data.
LSTM
LSTM prior and variational posterior of the form Equation 24. The prior is fit using ML on the initial CPE training data.
DTFM Decoder-only causal transformer prior and variational posterior of the form Equation 24. The prior fit using MI. on the initial CPE training data.
22

---

Published as a conference paper at ICLR 2025
2012
Raard
1200
Paculh Lauagee1
Roand
Roard
Board
Addeod
(a) DHFR
(b) TrpB
(c) TEBINDS
Figure D.3: Fitness landscape results using GP-PI. Precision (Equation 17), recall (Equation 18) and performance (Equation 19) higher is better for the combinatorially (near) complete datasets, DHFR and TrpB and TFBIND8. The random method is implemented by drawing samples uni- formly.
Strole COgAS
Sincin PHUNE
Sinoto Oounn
62
LANOC
6H4
OMULTS
Roand
Sound
Rourd
(a) M 15
(b) M 32
(c) = 64
Figure D.4: Ehrlich function (HOLO implementation) BBO results. VSD and CbAS with different variational distributions; mean field (MF), LSTM and transformer (TFM), compared against genetic algorithm (GA) and LaMBO-2 baselines.
TAE
Independent categorical prior and a cansition-style auto-encoder variational posterior the two-hidden layer MLPs for the encoder and decoder. The
form Equation 23, where prior is fit using ML on the initial CPE training data.
TCNN
Independent categorical prior and a transition-style convolutional auto-encoder variational posterior of the form Equation 23, where we use a convolutional encoder, and transpose convolutional decoder. The prior is fit using ML on the initial CPE training data.
We use the informed (hoomedd-indeendenn priors with the transition variational distributions since they are somewhat counter-in intuuttive as priors themselves.
23

---

Published as a conference paper at ICLR 2025
Draf
RoLnd
(a) DHFR
(b) TrpB
(c) TFBIND8
Figure D.5: Fitness landscape diversity results. Higher is more diverse, as defined by Equation 26.
50
GFP
Diverddy
Round
Rouns
Routc
AAV
Degrott
Deerdity
Round
Rrund
Round
(a) Independent
(b) LSTM
(c) Transformer
Figure D.6: Black-box optimization results for diversity on GFP and AAV with independent and auto-regressive variational distributions. Higher is more diverse, defined by Equation 26. The PEX and AdaLead results are replicated between the plots, since they are unaffected by choice of variational distribution.
From Figure D.8 we can see that while using an uninformative prior works in the lower-dimensionall fitness landscape experiments, using an informative prior is crucial for these higher dimensional problems. We found a similar result when using this uninformative prior with CbAS, or using a uniform initialization with DbAS and BORE. The methods not able to make any significant progress within the experimental budget given. The independent and transition variational distri- butions achieve similar performance, whereas the auto-regressive models generally outperform all others. This is because of the LSTM and transformer's superior generalization performance when generating sequences measured both when training the priors (on held-out sequences) and during VSD adaptation.
E THEORETICAL ANALYSIS FOR GP-BASED CPES
In this section, we present theoretical results concerning VSD and its estimates when equipped with Gaussian process regression models (Rasmussen & Williams, 2006). We show that VSD sam- pling distributions converge to target distribution that characterizes the level set given by The approximation eгror mainly depends on the predictive uncertainty of the probabilistic model with respect to the true underlying function f.. For the analysis, will assume that f. is drawn from Gaussian process, i.c., f. ~ (0, k), with positive-semideffiit covariance (or kernel) function k X x X R. In this case, can show that the predictive uncertainty of the model converges
24

---

Published as a conference paper at ICLR 2025
POLI
Hourc
HOLO
Doanty
Plantc
Roanc
Roanc
(a) M 15
(b) M 32
(c) M 64
Figure D.7: Black-box optimization results for diversity on the POLJ and HOLO implementations of the Ehrlich functions. Higher is more diverse, defined by Equation 26.
(in probability) to zero the number of observations grows. From this result, prove asymptotic convergence guarantees for VSD equipped with GP-PI-based CPEs. These results form the basis for analysis of CPEs based on neural networks (Appendix F).
E. GAUSSIAN PROCESS POSTERIOR
Let f. ~ GP(O, k) be a zero-mean Gaussian process with positive-semiddefnite covariance function
A X' x X + R. Assume that are given a set Dx = {(x;, M)) of N > observations fl(x;) + es. where (0, o?) and x; E X. The GP posterior predictive distribution at any X' is then given by (Rasmussen Williams, 2006):
f(x))DD ~ (X))) (x)))
(27)
mN(x) = kN(x) (KN +
(28)
(x x') k(x,x') - -x())) (x(xx)))) (KN o -NN(xx
(29)
oN (x) kxN(
(30)
where kx(x) = [bx x, ) E IRN KN = [x(x,, x)x,,, NNN and yN = [u: E V
Batch size. In the following, we will assume a batch of size B = 1 to keep the proofs simple. With this assumption, every iteration > 1, we have N = observations available in the dataset. We would, however, like to emphasize that sampling a batch of multiple observations, instead of a single observation, per iteration should only improve the convergence rates by a constant (batch- size-dependent) multiplicative factor. Therefore, our results remain valid as upper bound for the convergence rates VSD in the batch setting.
E.2 BACKGROUND
We will consider underlying probability space (2, 21, IP), where Ω is the sample space, 2 denotes the o-algebra of events, and IP is a probability measure. For any event E M, we have that P[AA] [0. 1] quantifies the probability of that event. For events involving random variable, e.g., X
(2, 2) + (R, Bz), where BR denotes the Borel o-algebra of the real line with its usual topology, will let:
P(x 0] = Pl[w E x(w) 0)]
(31)
25

---

Published as a conference paper at ICLR 2025
14
2.5
Simple regiet
VSD,DTPM
X50, 1STV
VAD,DTPM
VSD
VSD
VSDJUU
-2
VSD. 7AE
V50, TOW
Simple
1.0
VSDJ
VSOJU
-
-
VSD 7A
V50, V50,7CN TOWN
0.0
2
10
Round
Round
200
Diversity
DTM
125
VSD DTM
VSDJ
V50_7NE
2.0
V50.7M
30
Round
(a) AAV
(b) GFP
Figure D.8: Ablation results for the AAV and GFP BBO experiments. VSD is trialed with different prior and variational posterior combinations, "T" indicates simple independent informed prior and posterior, "IU" is the same but with uniforn prior, "LSTM" and "DTFM" the LSTM and decoder only transformer prior and posteriors, "TCNN" and "TAE" are transition convolutional encoder-decoder and auto-encoder posteriors, with informed independent priors. See text for details.
We will also use conditional expectations, i.c., given a o-sub-algebra of 21, the conditional expec- tation G[ is (6-meaurrable random variable such that:
VAEG, f ElX S xP
(32)
We will denote by (3 an increasing filtration on 2. For instance, we could set 3t as the a- algebra generated by the random variables in the algorithm (i.e., the candidates, target observations, etc.) at time t. For more details on the measure-theoretic definition of probability, we refer the reader classic textbooks in the area (e.g. Bauer, 1981; Durrett, 2019)
We will use the following well known notation for asymptotic convergence results. For given strictly positive function g + R. we define Og((t) as the of functions asymptotically bounded by (up to a constant factor) as:
[h((t]
Og((t) = h: +
lim sup t
(33)
g())
and for convergence in probability we use its stochastic counterpart:
Oz(g(t) P:N (2, + (R, Bx)
lim lim sup IP
Ip(t)
C
0
(34)
C-+00
t1000
g(1)
which is equivalent to:
|p g(t)
0, DC E (0)
C
Se,
T2T
(35)
26

---

Published as a conference paper at ICLR 2025
for some Tg N. For almost convergence, may also say that sequence of random variables Pe. E N, is almost surely ((g(t)) if Plp: E O(g(t)) = 1. A deeper overview on these notations and their properties can be found in Carcia-Portugués (2024). In addition, as it is common in the multi-armed bandits literature, we use the variants O and Op to denote asymptotic rates which are valid up to logarithmic factors.
E.3 AUXILIARY RESULTS
We start with a few technical results which will form the basis for our derivations. The following re- cursive relations allow us to derive convergence rates for the variance of GP posterior by analyzing how much it reduces per iteration.
Lemma E.1 (Chowdhury & Gopalan (2017, Appendix F)). The posterior mean and covariance
functions of a Gaussian process given 2 observations obey the following recursive identities:
x(x,,,)
Me(x)= 14 (x)
(yt pe-1(x))
(36)
0 oi (xt)
h (x. X1 k4 (Xx, x')
1x x) 1x,(xx,))
(37)
o² +
1 Xt)
of(x) r2-1(x)
(38)
o? ? (xt)
for x, x' E X.
We will also make use of the following version of the second Borel-Cantelli lemma adapted from Durrett (2019, Thr. 4.5.5) and its original statement in Dubins & Freedman (1965).
Lemma E.2 (Second Borel-Cantelli lemma). Let {A be a sequence of events where A, E St. for all t N, and let Xt : w + 1/ E for E St. Then the following holds with probability 1:
Σtsxe Xt
lim T-000
L < 00;
(39)
EL P|Ac 3t-1]
assuming P[A1/3o) > 0. In addition, lim Σ P[A = co, then 1.
The next result provides us with upper bound the posterior variance of a Gaussian process which is valid for any covariance function.
Lemma E.3. Let k : X X + IR be any posttive-semmidefnite kernel on X, and let is XXXXX x X be a kernel defined as:
+ R
xx(x,),, 0,
== x'
x(x,x')=
(40)
x X',
for x, E X. Göven any set of observations {xi, V/ v. for 1, denote by of the predictive variance of a GP model with prior covariance given by k, and let ōة denote the predictive variance of a GP model configured with k as prior covariance function, where both models are given the same set of observations. Then the following holds for all t 2 0:
0200(x (x)
o2(x) S 0?(x) xx()))
V/x E EX,
(41)
+² Nx(
where Nz(x) denotes the number of observations at x, and (x) = xx)) : k(), for X
Proof. It is not hard to show that ie defines valid ever k positive semidefinite. We will then focus on proving the main statement by an induction argument. The proof that the statement holds for the base case at /. = 0 is trivial given the definition:
covariance function when-
xx
xxxxxx T(xx, x) = (x),
V X.
(42)
Now assume that, for given t > 0, it holds that of (x) S (x),, for all x X. We will then check the inequality remains valid at + 1. By Lemma F, we have that:
1x, x(+1)
[x) = ox() of (x++1) oz
(43)
27

---

Published as a conference paper at ICLR 2025
For any x E X' such that x + X++1, we know that k x 2 0, so that (again by Lemma E.1):
by (x, x(+1) 1) S 1:2 (x, Xx+1) 0,
(44)
which shows that:
Vx f X1+1, (x) S oz S <x,(x) (x))
(45)
At X4+1, can rewrite f+1 (x+111++1)) as:
(x
(x++1) 07(x++1) + o?
(46)
We then check the difference:
ox (x
( 0 (x++1)
(x(+1) + t?? (x(x++1) + t
of (x (6) + 02) - o. -2 (xH +2)
(x, + (0) + +
(47)
(of ox(of((x+11) (x++1))
+ + (0) (x++) +
20,
since ox S (x(+1 by our assumption for time I. Therefore, have shown that:
(x) S (x) = o - o (x). Vx E X'
(48)
From the conclusion above and the base case, the inequality in the main result follows by induction.
Now we derive an explicit form for 87. Note that this case corresponds to an independent Gaussian model, i.e., f((x) II f(xx)) whenever x # x', for f. ~ GP(O, k). For any 1 2 1, this model's predictive variance at any x E X is given by:
(x). o.
x +
(x)
(49)
a² + +x)11(x))
(x,)
+ o²
x x
Looking the reciprocal, have that:
1xx x
2 1,
Vx E X.
(50)
(x
0i (x,)
o²
Therefore, every observation at x is simply adding factor of ---² to o, (x). Unwrapping this recursion leads us to:
Vt 1
11xx= Vx X. o?
1
(51)
ox)
(x)
The result in Lemma E.3 then follows as the reciprocal of the above, which concludes the proof.
Lemma E.4. Let f. (O k) for given ls X' X - R, where o supxxx k(x,x) oo, and | < co. Then f. is almost surely bounded, and: sup X
S a <xxv2logg2
(52)
Proof. The result follows by application of a concentration incquality for the maximum of finite collection of sub-Gaussian random variables (Boucheron et al., 2013, Sec. 2.5). Note that (.f((x) bxxx is collection of | Gaussian, and therefore lub-Gaussian, random variables with Gaussian parameter given by o& 2 (x)). for all X. Applying the maximal inequality for a finite collection Sub-Gaussian random variables (Boucheron et al., 2013, Thr. 2.5), we have that:
E
max f-(x) S Saxv log | < oo. XEX
(53)
By symmetry, we know that fxxx) is also sub-Gaussian with the same parameter, so that the bound remains valid for maxxex fx(x.. As a consequence, the expected value of the maximum of [fx(x) is upper bounded by the same constant. On finite set, the maximum and the supremum coincide. As the expected value the supremum is finite, the supremum must be almost surely finite by Markov's inequality, and therefore f. is almost surely bounded.
28

---

Published as a conference paper at ICLR 2025
E.4 ASYMPTOTIC CONVERGENCE
The main assumption we will be working with in this section is the following. Assumption E.1. The functionn a sample from a Gaussian process fa ~ GP(O, k), where k: X x X + R is bounded seniiie--emmdeffiitee kernel X.
The next result allows us to derive a convergence rate for the posterior variance of a GP as a function of the sampling probabílities. This result might also be useful by itself for other sampling problems involving GP-based approximations.
Lemma E.5. Let {x h22 be a sequence of X-valued random variables adapted to the filtration (31h221- For a given X, assume that the following holds:
T -
JL E VT T,
1 BT 00,
(54)
for a some sequence of lower bounds & {B (IN. Then, under Assumption E.I, given observations at {x:)i=> the following holds with probability 1:
(x) E O(
(55)
In addition, if B, -> oo, then lim Bio S o.
Proof. At any iteration i, the posterior variance of of GP model is upper bounded by worst case assumption of no correlation between observations (see Lemma E.3), In this case, we have that:
(x)
(x) - o²) =
(56)
a? + (x)
where (x) kxxx x)
(x x). and N; (x) S t denotes the total number of observations
taken at x as of iteration L. Without loss of generality, assume that o (x) = 1.
The only random variable to be bounded in Equation 56 is No. Let Xt 1[x, x). that:
N
X1
Σ 1x
x], t 121. 1.
(57)
We now apply the second Borel-Cantelli lemma (Lemmaa E.2) to N,. Namely, let Ñ, denote the sum of conditional expectations of {xi) given available data, i.c.:
N,
E 8:-1]
Σ [[ x| 3
Σ P(x: = x 8:-1].
(58)
By Lemma E.2, know that the following holds for some L R:
N LC N
lim 10000
(59)
Hence, N₁ is asymptotically equivalent to N₁. Applying this fact G2, we have that:
lim o
lim (x)
t 0² + N₁
B
lim
1-000 o² LN:
lim B,02
(60)
(-700 o² + LB,
< lim min LBr,
T, 19000
< 60,
29

---

Published as a conference paper at ICLR 2025
which holds with probability 1. Lastly, note that, if Be co, then L = 1 by Lemma E.2, and the last limit above becomes o²- The main result then follows by an application of Lemma E.3 and the definition of the big-O notation (see Equation 33).
M
We assume finite search space, which the case for spaces of discrete sequences of bounded length. However, we conjecture that our results can be extended to continuous or mixed discrete- continuous search spaces via a discretization argument under further assumptions on the kernel (e.g., ensuring that f. is Lipschitz continuous, as in Srinivas et al. (2010)). Assumption E.2. The search space X' is finite, | < oc.
We assume that our family of variational distributions is rich enough to be able to represent the PI-based distribution p(xly > Tr. Dz), which is the optimum of our variational objective when the optimal classifier is given by GP-PI. Although this assumption could be seen as strong, note that, due Gaussian noise, the classification probability p'y > x, D:) should be reasonably smooth function of x, which facilitates the approximation of the resulting posterior by generative model.
Assumption E.3. For every 2 O, p(xly > Tr, D:) is member of the variational family, i.e.:
TD ( ( (x > Tr, D.)] = =.
(61)
The next assumption is a technical one to ensure that the thresholds will not diverge to infinity. Assumption E.4. The sequence of thresholds almost surely bounded.2
sup
T S T 00.
(62)
(EN
We can now state our main result regarding the GP-based approximations learned by VSD.
Theorem E.1. Let assumptions E.1 to E.4 hold. Then the following holds with probability for VSD equipped with GP-PI:
df(x) O(
(63)
at every x E X such that p(x) > O.
Proof. Let (x(x) == ply T2 x, Dz). For any given x E X' where p(x) > 0, by Assumption E.2, we have that the next candidate will be sampled according to:
Vt 2 0, P(x++1
3e] p(xly Tt, D.)
lx(x)p(x)
(64)
(ex))
2 (x)x(x)),
where we used the fact that (xp()) have derive a lower bound on (x) apply Lemma E.5 and derive a convergence rate.
[x((x)] (,p())(((x))
< since (xx) < 1, for x E As psp(x) > 0, we only
A lower bound on (x(x) is given by:
Vt O, (((x) V/ (x) (x) +
|elloo + (
V/
(65)
Oe
where (-(-) denotes the cumulative distribution function of standard normal random variable, and |o denotes the essential supremum of a function under P (the probability measure of the underlying abstract probability space). Therefore, if lim;-soc ll < co, will have that lim;- (x() > O, and the sum in Lemma E.5 will diverge.
By Jensen's inequality for conditional expectations, have that:
11 0, || |B[
B.]1l
B[|lA||oo 3.).
(66)
Recall that for convergent sequences lim and lim sup coincide. We do require he known, only finite.
30

---

Published as a conference paper at ICLR 2025
As fl implies that:
3r]] bllbo]
< oo (cf. Lemma E.4), application of Markov's inequality
lim Ff Fr) 2 a| lim B[llf = 0.
(67)
a-+00
a->00 a
Furthermore, my := ||l 3:] also defines a non-negative martingale, and by the martingale
convergence theorem (Durrett, 2019. Thr. 4.2.11), lim M, = mx = EIIl ft 8xel is well
defined and [|][|] 8x = fflos < co. Again, by Markov's s inequality, for any a > 0,
we have that:
[limm-- [Plllool E [lim fo 3ell [I[I|]
lim |l 10o 2 fllc)
(68)
alE[1| f. [oo]
Therefore, for any a > O and any given x E X, with probability at least 1 a. the following holds:
lim Ppxx: = x 8 2 p(x) lim (x-((x) +
t000
||+t- 110 +
2 p(x)
lim v 1000
dc
(69)
alF [T[[1[Salloo) + Tx
2 p(x)
=: (a) > O.
Hence, for any ea E bz do (a)), there is N. N. such that P(x, =x 8 b Ea 0, for all / 2 Na. As result, Ex= P(x which asymptotically diverges at a rate proportional to t. By Lemma E.5 and the definition of the big-O notation, for any x E X, we then have that:
11] 2 (b - Eg) Na), for all 1 > Nu.
IP sim |Lo (x)] S o < c0
2 1 - 1
Va > 0,
(70)
1+000
Taking the limit as a oo, we can finally conclude that:
I sup (x) < 00 = 11
(71)
t-
i.e, of is almost surely O(L-)), which concludes the proof.
Remark F.1. The convergence rate in Theorem I is optimal and cannot be further improved. As
shown by previous works in the online learning literature (Muiný & Krause, 2018; Takeno et al., 2024), a lower bound on the GP variance at each iteration t Z 1 is given by (x) Z (o +t) (assuming k(x,x) the same point x E X (see Takeno et al., 2024, Lem. 4.2). Therefore, the lower and upper bounds on the asymptotic convergence rates for the GP variance differ by only up multiplicative constant.
1), which is the case when every observation in the dataset was collected at
The result Theorem E. now allows us derive a convergence rate for VSD's approximations to the level-set distributions. To do so, however, we will require the following mild assumption, which satisfied by any prior distribution which has support on the entire domain X.
Assumption E.5. The prior distribution is such that p(x) > O, for all x E X.
Theorem 2.1. Under mild assumptions (E.1 1 to E.5), the variational distribution of VSD equipped with GP-PI converges to the level-set distribution in probubility at the following rate:
DKL(P(xly> > Te, De)llp(xll > Tr, f.)] E Op(t 112).
(12)
Proof. We first prove an upper bound for the KL divergence in terms the PI approximation error. We then derive bound for this term and apply Theorem E. to obtain convergence rate.
KI. bound formulation. Let (x) == p④y > Te Tx|x, D.) and E (x) = p④y r.|x, f-), for x X.
From the definition of the KL divergence, we have that: DKL (p(xly > Te, D.)llp(xlp) > n. f.)] Ep(x [lo (log Epixlp>>. [log (()) - log f (x)
T4, D,) logp(x > Te, f-]]
+ log Ep(x) ( (x)] - log Ep(x) (x))
E7 (x)) }
(x(x)
Ep(x) (6] (x)] (72) Ep(()) [e(())
log
+ log
31

---

Published as a conference paper at ICLR 2025
For logarithms, we know that log(1 + a) S for all a > -1, which shows that:
(x(
1 (x(x) - (x(x))
(((x) (x)
log
log
(73)
Eg (x)
E (x)
(x)
pp(x) (x)]
1+ E Bp(x) 17 (x) [((x]]
Bx(x) 1ez (x) (x(x)
log
log
(74)
pp(x)) [x)]
1x(x) [x)
pp(x) le
Combining the above into Equation 72 yields:
[((x)) - E7 (xx)
E Bx(x) [ey (x) (x(x
DKL(p(xly > Tr, Dx)llp(xly > Te. f.) & (.
(x(x)
[exx)](((x)) (75)
The denominator in the expression above is such that:
lifilloo tz A1 0, (x)) = ply > Ti xx f.) = V f((x) Te V
VX Ex. (76)
O.
de
By Lemma E.4, know that [| follso] < co, which implies that f-llo < o0| by Markov's inequality. Next, we derive a bound for the approximation error term.
Error bound. We now derive upper bound for the difference (x, then show that it asymptotically vanishes. Applying Taylor's 's theorem to L we can bound AG, as function of the approximation error between the mean or and the true function as:
4xx) E7 (x) and
(x) [x + o?
f((x T2
Vt 2 0,
Ab((x)
y
V
Oe
(x)
f((x - Tt
2
(x +
Or
(xHo(x) fx(x)) (x) + o2 T2
+f -.))
(77)
22T
Te [o +
(x.H((x) f((x) V (x 02 + [x)
V
(xl - f((x) + o (x) + mml)
Vx E X,
02
(C) /2 1, and we used the fact that o, S vof(x) + of S 5x(xx) + ar to
since suPreR obtain the last two inequalities.
de
Convergence rate. derive convergence rate, given any x X and / 2 0, we have that:
(x)) | [.] a. [.B(((e(x) f(x)) 3.) + (E( 3.) + +md)
(78)
We know that [x]] 3r] is almost surely bounded, and by Jensen's inequality, it also holds that: (l] - fa(x)) 3r] (x))..
(79)
Applying Theorem 1. we then have that:
A E Op(t-1/2) 1/2).
(80)
Since || & B[llf e Op(1), we also have that:
E (xxx) E O0((1).
(81)
Lastly, we know that (x) E O+(1) by Equation 76 and the observation that |ffAllo E Oz(1). The main result then follows by combining the rates above into Equation 75.
32

---

Published as a conference paper at ICLR 2025
E.5 PERFORMANCE ANALYSIS
At every iteration / 1, VSD samples from (an approximation to) the target p(xly Tz-1, DI and obtains an observation Y: pyy) x:). A positive hit consists an event yt > Tt-1, where - is computed based on the data available in D or constant. Therefore, probability of positive hit for given realization f. as:
compute the
Plu:
Dr-1; [C]
Pp D. (p(u m-11x,
(82)
Then the expected number of hits Hr after T 2 1 iterations given by:
Birr fl]
Eptxbrr I.D, w)p(y Tx-1 f-)]. Σ
(83)
We will compare this quantity with the expected number of hits H7 obtained by sampling distri- bution with full knowledge of the objective function f::
T
E H f
Σ
Ep ppll > Te-1/x, f.)
(84)
The next result allows us bound the difference between these two quantities. Corollary 2.1. Under the settings Theorem 2. we also have that:
B[ - H7ll E O(VT).
(13)
Proof. For all T 2 1, we have that:
H - H = E
Ex(xp 1.D, s)lp(y > n-11x, - Ep(xpsrr wollllu > T--1/x, fi)l
E
Σ
ply T-11x, (p(xlv De-1) p(xly > f.))
XEX
ex(x)
(x))
(85)
E
ply >
(Ep(x)l((xx))
XEX
[x)(
E
p(x)
min (E [(((x')] E (07
XEX
since p④y f.) < 1. for all > 1. As both || and |l in Op(1),
(Eptx) (et (x']]} is lower bounded by some constant. As (xx(x) Op(1-1/2), for T large enough and some C > O, we then have that:
1 22CVT Ve
B H
(86)
which follows by application of the Euler-Maclaurin formula, since J 1 dl = 2/T 2 and Ji
the remainder term asymptotically vanishes.
If the oracle achieves E(H;) = T, the error bound in Corollary 2.1 suggests an increasing rate of positive hits by VSD as [E[lT]] > 1 1T-1 for some constant C. > and large enough T. Therefore, VSD should asymptotically achieve a full rate of positive hit per iteration in the single-point batch setting consider. Note, however, that the results above do not discount for repeated samples, though should still indicate that VSD achieves a high discovery rate the course of its execution.
Remark E.2.
33

---

Published as a conference paper at ICLR 2025
F VSD WITH NEURAL NETWORK CPES
In this section, we consider VSD with class probability estimators that are not based on GP regres- sion, which the case for the previous section, while specifically focusing on neural network models. We will, however, show that with a kernel-based formulation we are able to capture the classification models based on neural networks which use. This is possible by analyzing the be- havior of affiite-widtth neural networks (Jacot et al., 2018; Lee et al., 2019), whose approximation error with respect to the finite-width model can be bounded (Liu et al., 2020; Eldan al., 2021).
Although our classifiers are learned by minimizing the cross-entropy (CE) loss, we can connect their approximations with theoretical results from the infinite-width neural network (NN) literature, which are mostly based on the mean squared enror (MSE) loss. Recall that, given dataset Di-
{(x,, zn) with binary labels 2n E {0, 1}. the cross-entropy loss for a probabilistic classifier TT. : X 0. parameterized by is given by":
1 N
(0 Dīy)
In log (x + (1 zn) log(1 - =x(xx)))
(87)
n=1
The MSE loss for the same model corresponds to:
N
LM ( = N
Σ (Zn - Te TT(xx))) (xx))²
(88)
The following result establishes a connection between the two loss functions.
Proposition F1. Given a binary classification datasel Di of size N 2 1, the following holds for the cross-entropy and the mean-square error losses:
CCrE(O, D 2 LMsE(0, DR), N
(89)
Proof. Applying the basic logarithmic inequality log(1 + a) S a, for all > -1, to the cross- entropy loss definition yields:
1
Core(O, :=
nn log To(xx) (1 zn. log(1 (x()))
N
N
Zn (x)) - -1)- 2n (x))
(90)
N
2zn (x(n) Zn -To(xn)
n=1
2n 22, + Tx(()))
Now note that 2n = 22 for 2, E 10. and (x,))) 2 T )2 as Tx(x)) E [0, 1]. for all
NJ. NJ. Making these substitutions Equation 90, we obtain:
N
(O 1
2 2 n=1
+ Tx(x))² 2xx LMSE (o)
(91)
which concludes the proof.
The result in Proposition F.1 suggests that minimizing the cross-entropy loss will lead us to minimize the MSE loss as well, since the latter upper bounded by the former. This result provides us with theoretical justification to derive convergence results based on the MSE loss, which has been better analyzed in the NN literature (Jacot et al., 2018; Lee et al., 2019), as a proxy to establish convergence guarantees for the CE-based VSD setting.
We implicitly assume that 0 < (x(xx) < 1, for n (1..... N). so that the CE loss is well defined. This assumption can, however, be relaxed when dealing with the MSE loss, which remains well defined otherwise.
34

---

Published as a conference paper at ICLR 2025
F.1 LINEAR APPROXIMATION. VIA THE NEURAL TANGENT KERNEL
For this analysis, we will follow a frequentist setting. Namely, let n." denote the unknown true clas- sifier, i.e., π (x) == ply > - fc), for x E X. We assume that x" is unknown, fixed element of a reproducing kernel Hilbert space (RKHS) associated with given kernel (Schölkopf & Smola, 2001). In the case of infinite-width neural networks, we know that under certain assumptions the NN trained via gradient descent under the MSE loss will asymptotically converge to a kernel ridge regression solution whose kernel is given by the neural tangent kernel (NTK, Jacot et al., 2018). This asymptotic solution is equivalent to the posterior mean of Gaussian process that assumes no ob- servation noise. However, for a finite number of training steps & < oo, the literature has shown that gradient-based training provides a form of implicit regularization, which we use to ensure robustness to label noise. Moreover, although our analysis will be based on the NTK, the approximation error between the infinite-widd and the finite-width NN vanishes with the square root of the network width for most popular NN architectures (Liu et al., 2020). Therefore, we can assume that these approximation guarantees will remain useful for wide-enough, finite-width NN models.
Reproducing kernel Hilbert spaces. The RKHS F& associated with a positive- semidfii-e ker- nel k X' x X + R is a Hilbert space of functions over X' with an inner product (+,) and corresponding norm ||-lI& == (rr1), such that, for every π E Fk. the reproducing property (x) = (T, k(., x))x holds for all
(Schölkopf & Smola, 2001).
Implicit regularization. Several results in the literature have shown that training overparameter- ized neural networks via gradient descent provides a form of implicit regularization on the learned model (Fleming, 1990; Yao et al., 2007; Soudry et al., 2018; Barrett & Dherin, 2021), with some of the same behavior extending to the stochastic gradient setting (Smith et al., 2021). In earlier works, Fleming (1990) showed direct equivalence between an early stopped gradient-descent linear model and the solution of a regularized least-squares problem with a penalty on the parameters vector Eu- clidean norm. In the NTK regime, the network output predictions at iteration s e N of gradient descent are given by (Lee et al., 2019):
n (x) (x) + (x) KN' e PoKx') - (o) X X X
(92)
where To represents the network's initialization, V > 0 denotes the learning rate, the kernel k corresponds to the NTK associated with the given architecture, and the data is represented by Xx := {x C X' and nn == [ E {0, 11 Rearranging terms and performing basic algebraic manipulations, the equation above can be shown to be equivalent to:
(x) = (x) + (x) (KN Σx) (ZN mo(X))), x E X,
(93)
where En == (N - 1)= corresponds to a data-dependent regularization matrix. The above is equivalent to the solution of the regularized least-squares problem, we show below.
Lemma F.1. Assume ks : X x X + R is positive definite and = 0. Then Equation 93 solves the following regularized least- squares problem:
n E argmin TEFi
Σ (T (x() - 21)² + || n 1
(94)
where RN x 1 N. nN (x((x),.... .., ,(xx)),.,,((xx))) s((xN),, and xp(x) k(,,x) E FR
represents the kernel's canonical feature map. for x E X.
Proof. The least-squares loss can be rewritten as:
N (N(T) = Σ - z:)² + +|T = N zNl2 + ||R 1/2 7
(95)
where |1-112 denotes the Euclidean norm of a vector. Taking the functional gradient with respect to r E Fk and equating it to zero, we have that an optimal solution á satisfies:
Vln (+) = 2 (
-zN)= + y क = 0.
(96)
35

---

Published as a conference paper at ICLR 2025
The minimum-noon norm solution is then given by:
मे = ( Dy -
(97)
where A denotes the Moore-Penrose pseudo-inverse of an operator A : FR + Fk.
Let &N U V represent the singular value decomposition of & in the RKHS (Mollen- hauer et al., 2020), where UN := [u1,.... "y]; VN == [Vi,..., Vyl, with [u, C Fx and {v: C R denoting the left and right singular vectors, respectively, and Sy NNN corre- sponds to the diagonal matrix of singular values of ΦN There are N non-zero singular values, since k is assumed to be positive definite, and we have the correspondence KN = 2 f = VxAN V with AN = S representing the diagonal matrix of eigenvalues of Kv. which is full-rank for positive-defintt kernel with distinct entries Xy = {x)] C X. Applying the SVD derive the pseudo-inverse in Equation 97 then yields:
D N & PN(C V 1 NZN
(UNANU + UNSN (E VNANV V UUNSS
UNSN
.UN(AN + AN(E USAN
USNN 1) ") SS NNN
(98)
VU N V VN + ANN(C -PAAN 1) NZN
@NN(KN
- -1)-1)-
PN ΣN) ZN,
which coneludes the proof.
For our analysis, we will assume that the classifier network is zero initialized with O, noting that the least-squares problem can always be solved for the residuals (x(x) and then have mo added back to the solution. We refer the reader to Lee et al. (2019) for further discussion on the effect of the network initialization.
Approximation for finite-width networks. For fully connected, convolutional or residual net- works equipped with smooth activation functions (e.g., sigmoid or tanh), Liu et al. (2020) showed where
that the approximation error between the linear model and the finite-width NN is m denotes the minimum layer width, and the Ô notation corresponds to the O-notation with loga- rithmic factors suppressed. NTK results for other activation functions and different neural network architectures, such as multi-head attention (Hron et al., 2020), are also available in the literature.
F.2 ASSUMPTIONS
In the following, we present a series of mild technical assumptions needed for our theoreti- cal analysis of NN-based CPEs. For this amalysis, we mainly assume that the true classifier T Tx(x) = pyy > x, f.) is a fixed, though unknown, element of the RKHS F& given by a bounded NTK which is formalized by the following two assumptions. As in the GP case, assume batches of size B = 1 to simplify the analysis.
Assumption F.1. There x²" E F₁ such that:
" +x(x) = ply > S.), E X.
(99)
For a rich enough RKHS, such assumption is mild, especially given that most popular NN architec- tures offer universal approximation guarantees (Hornik al., 1989).
Assumption F.2. The NTK k corresponding the network architecture in wa is positive definite and bounded in X.
We will also assume that the threshold is fixed to simplify the analysis. However, our results should asymptotically hold for time-varying thresholds as long as the limit lime-soo Ti exists.
Assumption F.3. The threshold is fixed, i.e., To IR, for all 2 1.
36

---

Published as a conference paper at ICLR 2025
The following assumption on label noise should always hold for Bernoulli random variables (Boucheron et al., 2013). Any upper bound on the sub-Gaussian parameter should suffice for the analysis (c.g., o for Bernoulli variables).
Assumption F.4. For all N and all E X, label noise Ç p(y|x, x, f.), is re-sub-Gaussian
1 7] xx" xx(xx, with
Ta'0?
Va E
E (eexp(a() I exp
(100)
for some og 2 O.
The next assumption ensures a sufficient amount of sampling is asymptotically achieved over the domain X, which we still assume is finite.
Assumption F5. For any t the variational family is such that sampling probabilities bounded away from 0, i.e.:
3& b>00: 0: Vt T
(x()))) b, vx E X.
(101)
The assumption above only imposes mild constraints the generative models q (6), that prob- abilities for all candidates x E X are never exactly O, though still allowed to be arbitrarily small.
Assumption F.6. The learning rate Vt at each round t is such that:
0 O<U.SuS Amax K.)
VEEN,
(102)
for V. > where )mux(() denotes the maximum eigenvalue of matrix, and K, denotes the NTK matrix evaluated at the training points available at iteration t 2 1.
This last assumption ensures that a gradient descent algorithm is convergent (Fleming, 1990), though it to bound the spectrum of the implicit regularization matrix E, after finite number of training steps 8 < oo (a.k.a. early stopping), which is needed for our results. We highlight that, under mild assumptions on the data distribution, (M,, E O(1) w.r.t. the number of data points (Murray et al., 2023), so that the bound in Assumption F.6 will not vanish. Therefore, such assumption is easily satisfied by maintaining sufficiently small learning rate.
F.3 APPROXIMATION ERROR FOR NN-BASED CPES
Similar to the GP-PI setting, we will assume batch size of 1, so that we can simply use the iteration index é > for our estimators. We recall that convergence rates for the batch setting should only be affected by a batch-size-depeee dependent multiplicative factor, preserving big-O convergence rates. We start by defining the following proxy variance:
2 1, (x) k(x, x) (x) (K: + E.) (x),), XEX,
(103)
where E, the implicit regularization matrix in Equation 93 due early stopping. The proxy variance is then equivalent to a GP posterior variance under the assumption of heteroscedastic (i.e., input dependent) Gaussian noise with covariance matrix given by E,. Given its similarities, we have that if enough sampling is asymptotically guaranteed, we can apply the same convergence results available for the GP-Pl-based CPE, i.e., of E O((")) almost surely.
Lemma F.2. Let assumptions F2, E.5 and E.6 hold. Then the following almost surely holds for the proxy variance:
O(
(104)
Proof. We first observe that à? (103) upper bounded by the posterior predictive variance of GP model assuming i.i.d. Gaussian noise (cf. Sec. E. with variance p² satisfying:
Amax (E.), VI VIEN,
(105)
37

---

Published as a conference paper at ICLR 2025
which is such that:
Σx K -I
(i(Ussk,))+ 1 8Vt
(106)
since eA > A, for any Hermitian matrix A, where 2 denotes the Loewner partial ordering in the space of positive-seiidefinte matrices, i.e., A > B if and only if A is positive semidefinite. Noting that the sum of sampling probabilities at any point x E X' diverges as - oo by Assump- tion F.5, the result then follows by applying Lemma E.5 to the GP predictive variance upper bound with noise variance set to p² == (sv.)- (Assumption F.6).
Lemma F.3. Let assumptions F.I to F.6 hold. Then, given any E (0, 1], the following holds with probability least -- for the approximation error between tee and n²:
vI 2 1, [xx(x) T (x)] < (x), x XEX,
(107)
where 3,(8) = ||~0||~ oc 2p log (eg((dttI p- (K, and p == 8V.
Proof. The result above is direct application of Theorem 3.5 in Maillard (2016) which provides an upper confidence bound on the kernelized least-squares regressor approximation error (another version of the same result is also available in Durand al. (2018, Thr. 1)).
Let 4 == x² (xx) denote the label noise in observation i, for {l...., t}. Expanding the definition of ft (93) with To 0, given any E N and x E X, we can decompose the approximation error as:
|m- (x) (x)]
|m² (x) 1x(x) T (K, E |2y
= n (x) k((x) " (K, + E.) ( + C.)1
(108)
< n (x) - (x(x) + E.) + (x) (K, + E,)
where applied the triangle inequality obtain the last line. Analyzing the two terms the right-hand side, by the reproducing property, we now have for the first term:
(x'(x) kx (K, + Σ.) -; = KTT, (I- P. + E.) L (x)))al
KTT,, (I+ D.E, Φ. (x())al
S ||="||<|| [I ( $.
(x))x
(109)
||="|la p(x) (I E. 1b -xx((x)
S || ((x]ll//(xx) E,E, $ -l(((x)
= =|+(")l0=(xx),
where the second equality follows by an application of Woodbury's identity, the first inequality is due to Mauchh-Schwaate the second inequality is due to the fact that A -2 3 A whenever A E I. and the last line follows from the definition of o?. For the remaining in term (108), we have that:
(x) TK, + -t,l = ((((x), (K, + 5.)=tt())
(()(()), (I Σ, +b - 1 - S.)+1
I((I + 24 1/2 p((), (I P,E, -.)--112),7.7 (c)a
S Vx(x)T (I
E, 1 )
$ $. -1 E, Sellk
1(110)
d((x)) c. x 4. [I P.E, Φ
16 - 'S:,
where we applied the identity T (BB + -A-- (I+B A A BB) A which holds for
an invertible matrix A (Searle, 1982), to obtain the second equality, and the upper bound follows by the Cauchy-Schwan inequality. For the norm of the noise-dependent term, will apply
38

---

Published as a conference paper at ICLR 2025
concentration bound by Abbasi- Yadkori (2012), which first requires a few transformations towards non-time- varying regularization factor. Applying the SVD o, V, (Mollenhauer al., 2020), in the proof of Lemma F.1, we have that Σ, V,D,V, where D; A -I and A, S. which leads to:
[I .,
Sell -1 4,-, 5
"UD = =
(111)
V. " D.) D.)
where we applied the identity (I AB) - A(I + BA) -1 (Searle, 1982). For the eigenvalues of Er, have the following lower bound:
D. A -I)
M e -
Age uisA (U
(112)
1 enne.
V/8
1
eI,
V4S
where the first inequality due to -- Hence, setting = e_ sv. we have that:
A
A, and the last inequality holds by Assumption F.6.
(I +,
pc, V.S((++)))
s vpl S.V.
(113)
=p-ic $ (pI + -
[IloII + +
By Abbasi- Yadkori (2012, Cor. 3.6), given any 5 (0, 1), then have that the following holds with probability at least
ddtt(I+ p D-1K,)
2 1, [[ +, T)
1
S < 20 log
(114)
Finally, combining the bounds above into Equation 108 leads to the result in Lemma F.3.
For the next result, we need to define the following quantity:
ET
max C
log det(I + 'K(tr)),
(115)
2
where K(Xr) = [k(x, x') xxx''X E R [T|x|TT Note that Êr corresponds to the maximum information gain of GP model (Srinivas et al., 2010) with covariance function given by the NTK, assuming Gaussian observation noise with variance given by Then & is mainly dependent on the eigenvalue decay of the kernel under its spectral decomposition (Vakili al., 2021). For the spectrum of the NTK, a few results are available in the literature (Murray et al., 2023).
Proposition F.2. Let assumptions F. to F.6 hold. Then, given E (0, 1}, the following holds with probability at least 8 for VSD equipped with wide enough NN-based CPE model ā,.
DRI|p(xly > Tt. (x)llp(xly Ti, f)) S) (V
(116)
39

---

Published as a conference paper at ICLR 2025
Proof. The result follows by applying the same steps in the proof of Theorem 2.1. We note that (xx = π (x) > 0, due to observation noise, so that E (xxx= E Op(1). Similarly, Lemma F3 implies that [x) n" (x) s Be (x)) with probability at least 1 - å simultancously over all x E X, so that ratio-dependent terms in Theorem 2.1 should remain bounded in probability. The upper bound in the result then follows by noticing that in our case (x)( (x()))x(x)) with
high probability, where oN O0 by Lemma F.2, and 3e(6) (vVE)) by Lemma F.3 and the definition of & in Equation 115.
The result above tells us that VSD equipped with an NN-based CPE can recover similar asymp- totic convergence guarantee to the one we derived for the GP-PI case, depending on the choice of NN architecture and more specifically on the spectrum of its associated NTK. In the case of a fully connected multi-layer ReLU network, for example, Chen & Xu (2021) showed tween the RKHS of the ReLU NTK and that of the Laplace kernel ((x, x') = exp( exp(-Cllx - x'I))- As the latter is equivalent to a Matérn kernel with smoothness parameter set to 0.5 (Rasmussen & Williams, 2006), the corresponding information gain bound & E (t ), where here denotes the dimensionality of the domain X' (Vakili & Olkhovskaya, 2023). In the case of discrete sequences of length M, the dimensionality of is determined by M. Hence, we have proven Corollary 2.2.*
equivalence be-
Corollary 2.2. Let Te be modeled via a fully connected ReLU network. Then, under assumptions on identifibiiiiy and sampling (F.1 to F6), in the infinite-width limit, VSD with CPE-PI achieves:
DKL(P(xlV > Tú, Do)llp(xly > Tu. f.). E Öp (TाTTT))
(14)
Similar steps can be applied to derive convergence guarantees for VSD with other neural network architectures based on the eigenspectrum of their NTK (Murray et al., 2023) and following the recipe in, e.g., Vakili et al. (2021) or Srinivas et al. (2010).
G VSD AS A BLACK-BOX OPTIMIZATION LOWER BOUND
A natural question to ask is how VSD relates to the BO objective for probability of improve- ment (Garntt, 2023, Ch.7),
argmax log oprg, DN,
(117)
x
Firstly, we can see that the expected log-likelihood of term of Equation 6 lower-bounds this quantity.
Proposition G.1. For a parametric model, g(x(o), given & E + S RD and q E P: X' +b - 10, 1]. maxlog pr DN, 2 max Baxxooxxtl [log ( DN,
(118)
and the bound becomes tight as q(x()))
5(x²) a Dirac delta function at the maximizer x
Taking the argmax of the RHS will result in the variational distribution collapsing to delta distri- bution at x² for an appropriate choice of (x)))). The intuition for Equation 118 is that the expected value of random variable is always less than equal its maximum. The proof of this is in Daulton et al. (2022); Staines & Barber (2013). Extending this lower bound, we can show the following.
Proposition G.2. For divergence D P(X) x P(X) + [0, oo), and a prior po E P(X).
max log (x, T) 2 max (x [log o (x, DN, r)] D[g
(119)
We can see that this bound trivially true given the range of divergences, and this covers VSD as a special case. However, this bound is tight if and only if Pa concentrates as a Dirac delta at x; with an appropriate choice of q(x()). In any case, the lower bound remains valid for any choice of informative prior Po or even a uninformed prior, which allows us to maintain the framework flexible to incorporate existing prior information whenever that is available.
*Here Öp suppresses logarithmic factors, as in O. and holds in probability.
40


Reviews:
[{"review_id": "Ha686hWof3", "reviewer": "Reviewer_xamp", "summary": "This paper casts sequential black-box optimization as a variational inference (i.e. amortized optimization) problem, and uses this perspective to unify a collection of different black-box optimization algorithms under a common theoretical framework and presents some proof of concept results on easy sequence optimization tasks.\n\n#### UPDATE - 12/02/2024 ####\nAfter extended discussion and an exceptionally thorough response from the authors, my concerns have been addressed and I am recommending acceptance. I believe this paper presents a very nice conceptual view of the topic of active generation and the empirical results raise interesting questions. I encourage my fellow reviewers to review the updated manuscript and re-evaluate their scores. I have left my original review unaltered for any interested readers.", "questions": "The following questions are sincere:\n\n- Who is the audience for this paper? \n\n- What questions is this paper answering?\n\n- What does the variational inference framing get us in the end? Access to a set of tools for theoretical analysis?", "limitations": "", "rating": 8, "confidence": 5, "soundness": 4, "presentation": 3, "contribution": 3, "strengths": "This paper demonstrates a clarity of thought and composition that is commendable, I particularly enjoyed the related work section.\n\nLikewise I do not have any major concerns regarding the technical soundness of the results presented.\n\nAs a good conceptual introduction to the topic, I think this draft could be useful to researchers new to the topic with some revisions.", "weaknesses": "I have two general impressions of this paper. \n\nFirst, it seems like the authors have not really chosen a direction for the paper. There are at least three different directions here, A) a unifying view of sequential black box optimization algorithms, B) a practical algorithm for sequential BBO, and C) theoretical analysis of convergence rates of a particular sequential BBO algorithm under strong assumptions. I would suggest you pick no more than two directions, preferably one. I actually think this particular subfield could really benefit from a more holistic perspective of the work that has been done, as I constantly see minor variations of these algorithms in my social media feed and review stack with no apparent awareness of the relationships between them. From what I can tell from this draft, it seems that A and C likely play more to your strengths.\n\n\nSecond, the authors seem blissfully unaware of a substantial body of work on this topic. To be quite candid, the paper reads like it was written circa September 2021. This is not mere rhetoric. The most recent baseline the authors consider was published at ICML 2021. It is also odd that two of the baselines you did include, DbAS and CbAS, are not even designed for the sequential setting. As a very active researcher in this exact area, I struggle to understand who this paper is for and how the authors pictured their place in the broader dialogue on this topic. I am sure you worked very hard on this paper and I commend your effort, but I honestly believe the best advice I can give you is to talk to more people working on this topic, preferably from outside your immediate academic circle. While it is difficult to hear this feedback, one of the functions of peer review is to reveal \"unknown unknowns\". I want to be sure this review is constructive, so I will provide some key references if you are serious about diving into this topic. You should also consider making use of tools like [Connected Papers](https://www.connectedpapers.com/) to improve your literature review process and avoid this situation in the future. \n\nYou can start with [A survey and benchmark of high-dimensional Bayesian optimization of discrete sequences](https://arxiv.org/abs/2406.04739). This work is the most up-to-date complete survey on the topic I have seen, and the benchmarking rigor is notably good. This paper is associated with two repositories, [poli](https://github.com/MachineLearningLifeScience/poli) and [poli-baselines](https://github.com/MachineLearningLifeScience/poli-baselines). The former contains a suite of test functions that are much more up to date than the combinatorially complete landscapes considered in this paper, and the latter contains a suite of baseline solvers. You may even want to consider contributing your method as a solver to poli-baselines at some point.\n\nSome key axes of variation to consider: \n\nHow is the optimization problem solved? Most fall into one of three categories, directed evolution (which you seem to be familiar with based on your inclusion of AdaLead and PEX), generative search with explicit guidance, e.g. [2, 3, 4, 5, 6], and generative search with implicit guidance [7, 8], which can also be seen as a kind of amortized search. I could cite more papers but I believe I have made my point. Algorithms also differ in their handling of constraints, and their approach to managing the feedback covariate shift induced by online active data collection by an agent. \n\nIn particular I will draw your attention to [a tutorial for LaMBO-2](https://github.com/prescient-design/cortex/blob/main/tutorials/4_guided_diffusion.ipynb) if you want to start considering more up to date baselines, however I would recommend using the solver interface provided in poli-baselines for actual experiments. You may also be interested in Ehrlich functions if you would like a convenient test function that is much more difficult to solve than small combinatorially complete landscapes but still easy to work with [9]. Ehrlich functions are available in [a small standalone package](https://github.com/prescient-design/holo-bench) or [as part of the poli package](https://machinelearninglifescience.github.io/poli-docs/using_poli/objective_repository/ehrlich_functions.html).\n\nWhile I'm sure this is not the outcome you hoped for, science is a dialogue, and good science requires awareness of what is happening outside your academic niche. Hopefully my feedback is clear and actionable enough to benefit this work and your progression as a scientist.\n\nReferences\n\n- [1] González-Duque, M., Michael, R., Bartels, S., Zainchkovskyy, Y., Hauberg, S., & Boomsma, W. (2024). A survey and benchmark of high-dimensional Bayesian optimization of discrete sequences. arXiv preprint arXiv:2406.04739.\n- [2] Tripp, A., Daxberger, E., & Hernández-Lobato, J. M. (2020). Sample-efficient optimization in the latent space of deep generative models via weighted retraining. Advances in Neural Information Processing Systems, 33, 11259-11272.\n- [3] Stanton, S., Maddox, W., Gruver, N., Maffettone, P., Delaney, E., Greenside, P., & Wilson, A. G. (2022, June). Accelerating bayesian optimization for biological sequence design with denoising autoencoders. In International Conference on Machine Learning (pp. 20459-20478). PMLR.\n- [4] Gruver, N., Stanton, S., Frey, N., Rudner, T. G., Hotzel, I., Lafrance-Vanasse, J., ... & Wilson, A. G. (2023). Protein design with guided discrete diffusion. Advances in neural information processing systems, 36.\n- [5] Maus, N., Jones, H., Moore, J., Kusner, M. J., Bradshaw, J., & Gardner, J. (2022). Local latent space bayesian optimization over structured inputs. Advances in neural information processing systems, 35, 34505-34518.\n- [6] Maus, N., Wu, K., Eriksson, D., & Gardner, J. (2022). Discovering many diverse solutions with bayesian optimization. arXiv preprint arXiv:2210.10953.\n- [7] Tagasovska, N., Gligorijević, V., Cho, K., & Loukas, A. (2024). Implicitly Guided Design with PropEn: Match your Data to Follow the Gradient. arXiv preprint arXiv:2405.18075.\n- [8] Chen, A., Stanton, S. D., Alberstein, R. G., Watkins, A. M., Bonneau, R., Gligorijevi, V., ... & Frey, N. C. (2024). LLMs are Highly-Constrained Biophysical Sequence Optimizers. arXiv preprint arXiv:2410.22296.\n- [9] Stanton, S., Alberstein, R., Frey, N., Watkins, A., & Cho, K. (2024). Closed-Form Test Functions for Biophysical Sequence Optimization Algorithms. arXiv preprint arXiv:2407.00236.", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "", "was_revised": false, "final_justification": "", "rebuttal": "", "claims_and_evidence": "", "methods_and_evaluation_criteria": "", "theoretical_claims": "", "experimental_designs_or_analyses": "", "relation_to_broader_scientific_literature": "", "essential_references_not_discussed": "", "other_strengths_and_weaknesses": "", "other_comments_or_suggestions": "", "questions_for_authors": "", "summary_of_paper": "", "summary_of_recommendation": "", "technical_quality": "", "clarity_of_presentation": "", "potential_impact": "", "robotics_focus": "", "quality_of_limitations": "", "questions_for_rebuttal": "", "reasons_to_accept": "", "reasons_to_reject": "", "questions_to_authors": "", "extra_scores": {}}, {"review_id": "MUV7SmyOso", "reviewer": "Reviewer_1vHR", "summary": "The authors propose a black box varoiational inference approachfor discrete designs generation. The authors derive asymptotic convergence rates for learning the true conditional generative distribution of designs. Compelling results on high dimensional sequence-design problems are demonstrated.", "questions": "* What's 'x' in the title of Figure 1?\n* What are the limitations of this approach?\n* How is diversity within a batch enforced? \n* The reverse KLD is known to result in mode collapse. Why wasn't this an issue?\n* Which variation reduction method did you use for the gradient estimator?", "limitations": "", "rating": 6, "confidence": 3, "soundness": 3, "presentation": 2, "contribution": 3, "strengths": "* The problem is important as it has applications in pharmaceutical drugs/enzyme design.\n* The paper paper is well written and the method is sound\n* Experimental results on high dimensional datasets demonstrate superiority of the approach", "weaknesses": "* The method lacks novelty, it's based on putting together blocks that have already been proposed in the litterature\n* The paper clarity can be improved with an overview plot of the method", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "", "was_revised": false, "final_justification": "", "rebuttal": "", "claims_and_evidence": "", "methods_and_evaluation_criteria": "", "theoretical_claims": "", "experimental_designs_or_analyses": "", "relation_to_broader_scientific_literature": "", "essential_references_not_discussed": "", "other_strengths_and_weaknesses": "", "other_comments_or_suggestions": "", "questions_for_authors": "", "summary_of_paper": "", "summary_of_recommendation": "", "technical_quality": "", "clarity_of_presentation": "", "potential_impact": "", "robotics_focus": "", "quality_of_limitations": "", "questions_for_rebuttal": "", "reasons_to_accept": "", "reasons_to_reject": "", "questions_to_authors": "", "extra_scores": {}}, {"review_id": "gPHpFf2S41", "reviewer": "Reviewer_z76B", "summary": "The paper develops the variational search distribution method to solve the active search problem in biological design. VSD estimates the super level-set distribution in a sequential manner by generating batches of data points for sequential experiments. Empirical results on optimizing protein fitness in several datasets showcase the effectiveness of VSD.", "questions": "- How robust are the results to the selection of the threshold $\\tau$ and the batch size $B$?\n- While the reviewer is not familiar with the field, could the authors give some intuitions about the difference between VSD and active learning approaches like Bayesian optimization, and why VSD is better?", "limitations": "", "rating": 6, "confidence": 2, "soundness": 3, "presentation": 2, "contribution": 3, "strengths": "- The paper formulates the batch active search problem in the variational inference framework and provides theoretical guarantees to the learned distribution based on the sequentially attained data.\n- Experimental results on real-world biological datasets demonstrate the practical use of the algorithm and its effectiveness to solve the problem.", "weaknesses": "- The precision of VSD and most other methods is decreasing with more rounds in TrpB and TFBIND8 datasets while the recall values are in general low. However, an ideal method should achieve a better estimation of the ground truth super level-set distribution as more samples are collected. This may be due to the initial training set size being too large or the fitness landscape being easy to model. How do the models perform with a smaller initial training set size?\n- How is VSD compared with the simple and commonly used directed evolution method?", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "", "was_revised": false, "final_justification": "", "rebuttal": "", "claims_and_evidence": "", "methods_and_evaluation_criteria": "", "theoretical_claims": "", "experimental_designs_or_analyses": "", "relation_to_broader_scientific_literature": "", "essential_references_not_discussed": "", "other_strengths_and_weaknesses": "", "other_comments_or_suggestions": "", "questions_for_authors": "", "summary_of_paper": "", "summary_of_recommendation": "", "technical_quality": "", "clarity_of_presentation": "", "potential_impact": "", "robotics_focus": "", "quality_of_limitations": "", "questions_for_rebuttal": "", "reasons_to_accept": "", "reasons_to_reject": "", "questions_to_authors": "", "extra_scores": {}}, {"review_id": "ciTPqvfntP", "reviewer": "Reviewer_9cf6", "summary": "The authors present a novel variational method for learning to sample from rarely observed events, aiming to minimize a distance between the distribution of interest, namely $p(x∣y>t)$, and its parametric variational counterpart $q(x|\\phi)$. The problem is reformulated to leverage the “whole dataset,” not just rarely observed events, and is expressed as Equation (5), which comprises two terms: $log p (y>t∣x)$ and the negative KL divergence between$q(x|\\phi)$ and $p(x)$. The authors' final proposal is to estimate $p(y>t∣x)$ using a parametric function instead of a simple PI estimate. The variational distribution is optimized by a REINFORCE gradient estimator.", "questions": "1.Since your algorithm heavily relies on another model ($p(x | D)$), I would be highly interested in better understanding the influence of a good prior on your variational distribution.\n2. Regarding the GFP experiments, do you sample already existing sequences ? What is the influence of the relative poor performance of the oracle on ood data on the interpretation of the results  ?\n3. How can you explain that only a very simple prior such as a mean field performs on average better ?  It seems quite logical for GFP for instance where a wild type exists, however it is less intuitive for datasets without wild type.\n\nTypo: the recall and precision have the same expression.", "limitations": "", "rating": 6, "confidence": 3, "soundness": 4, "presentation": 3, "contribution": 3, "strengths": "The paper is clear, well-written, and aligns with well-established benchmarks in the field, such as CBAS (Brookes et al.). \nThe model is supported by convergence analysis and an extensive set of well-handled experiments.", "weaknesses": "While the model description is clear, the model comprise a parametric distribution $p(x|D_0)$ which might be the biggest model shortcoming originating from the model own formulation. \n\nIts major impact is that it reweights the gradient estimates of $q(x|\\phi)$. Intuitively, how would that compare simply to the iterative strategy of Cbas ?", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "", "was_revised": false, "final_justification": "", "rebuttal": "", "claims_and_evidence": "", "methods_and_evaluation_criteria": "", "theoretical_claims": "", "experimental_designs_or_analyses": "", "relation_to_broader_scientific_literature": "", "essential_references_not_discussed": "", "other_strengths_and_weaknesses": "", "other_comments_or_suggestions": "", "questions_for_authors": "", "summary_of_paper": "", "summary_of_recommendation": "", "technical_quality": "", "clarity_of_presentation": "", "potential_impact": "", "robotics_focus": "", "quality_of_limitations": "", "questions_for_rebuttal": "", "reasons_to_accept": "", "reasons_to_reject": "", "questions_to_authors": "", "extra_scores": {}}]

Author rebuttal:

Decision: Accept (Poster)