Multi-objective Bayesian Optimization with Heuristic Objectives for Biomedical and Molecular Data Analysis Workflows
Authors: Alina Selega
Kieran R. Campbell
Venue: Accepted by TMLR ()

Abstract: Many practical applications require optimization of multiple, computationally expensive, and possibly competing objectives that are well-suited for multi-objective Bayesian optimization (MOBO). However, for many types of biomedical data, measures of data analysis workflow success are often heuristic and therefore it is not known a priori which objectives are useful. Thus, MOBO methods that return the full Pareto front may be suboptimal in these cases. Here we propose a novel MOBO method that adaptively updates the scalarization function using properties of the posterior of a multi-output Gaussian process surrogate function. This approach selects useful objectives based on a flexible set of desirable criteria, allowing the functional form of each objective to guide optimization. We demonstrate the qualitative behaviour of our method on toy data and perform proof-of-concept analyses of single-cell RNA sequencing and highly multiplexed imaging datasets for univariate input optimization.

Published in Transactions on Machine Learning Research (03/2023)
Multi-objective Bayesian Optimization with Heuristic Objec- tives for Biomedical and Molecular Data Analysis Workflows
Alina Selega
aselega@lumenfeld.cor
Lumenfeld- Tonenbaum Research Institute Vector Institute
Kieran R. Campbell
Seeranaampiell@hmmmjeelll....
Lumenfeld-Tanenttee Tanenboum Research Institute University of Toronto
Vector Institute
Onlario Institute for Cancer Research
Reviewed on OpenReview: /ttps://openreene net/ forum? id= QspAcsAyis
Abstract
Many practical applications require optimization of multiple, computationally expensive, and possibly competing objectives that are well-suited for multi-objective Bayesian optimization (MOBO). However, for many types of biomedical data, mcasures of data analysis workflow success are often heuristic and therefore it is not known priori which objectives are useful. Thus, MOBO methods that return the full Pareto front may be suboptimal in these cases. Here we propose a novel MОBO method that adaptively updates the scalarization function using properties of the posterior of a multi-output Gaussian process surrogate function. This approach selects useful objectives based on flexible set of desirable criteria, allowing the functional form of each objective guide optimization. We demonstrate the qualitative behaviour of our method on toy data and perform proof-of-concept analyses of single-cell RNA sequencing and highly multiplexed imaging datasets for univariate input optimization.
1
Introduction
The analysis of high-dimensions biological data is often exploratory and unsupervised. For example, gene expression data may be subject to clustering algorithms to find groups representative of mcaningful biological variation. For assays that profile at the patient level, these clusters may represent novel disease subtypes, while for assays at the single-cell level, they may represent. novel cell types.
Despite the importance of these methods, there is no "one-sizo-fitsall" approach to the analysis of such data. Instead, there is a myriad of different possible parameter combinations that govern these workflows and lead to variations in the results and interpretation. For example, in the analysis of single-cell RNA-sequencing (scRNA-seq) technology that quantifies the expression profile of all genes at single-cell resolution a common analysis strategy is to cluster the cells to identify groups with biological significance. However, each workflow for doing so has variations with respect to data normalization, cell filtering strategies, and the choice of clustering algorithm and parameters thereof. Changes to these algorithm and parameter choices produce dramatically different results (Germain al., 2020 Duò et al. 2018) and there is no ground truth available. This motivates important question: how do we optimize these workflows such that the resulting exploratory analysis bost reflects the underlying biology? It is important to note the presence of measurement noise in virlually all biomedical data, which can arise from the technology used for data acquisition or represent underlying biological heterogeneity.
In the adjacent field of supervised machine learning (ML), such optimization over workflows has largely been tackled from the perspective of automated ML (AutoML, He al. (2021)). This comprises a diverse of methods such Bayesian optimization (Snoek et al. 2012) and Neural Architecture Search (Elsken al.

---

Published in Transactions on Machine Learning Research (03/2023)
2018) that attempt to optimize the success of the model with respect to one or more hyperparameter settings. In this context, success is defined as the model accuracy on held out test set, though can also correspond to the marginal likelihood of the data given the model and hyperparameters.
However, in the context of exploratory analysis of genomic data, existing AutoML approaches face three challenges. Firstly, they are almost exclusively unsupervised, meaning there is no notion of accuracy on a test set we may optimize with respect to. Secondly, the majority of methods are not generative probabilistic models (Zappia et al., 2018) so it is impossible to optimize with respect to the marginal or test likelihood. Finally, the objectives used to optimize workflow numerous, conflicting, noisy, due to the underlying noise present in the raw data, and can be highly subjective, due to often being heuristics.
This is demonstrated by attempls to benchmark clustering workflows of scRNAA-eqq data. As said above, there are many parameters that must be set, e.g. which subset of genes and clustering algorithm to use, along with such paramcters as resolution in the case of community detection (Germain et al.) 2020). However, there
is no quantitative way to choose which parameter setting is "best" and so the community turns to number of heuristic objectives to quantify the performance of a workflow. For example, Cui et al. (2021) | attempt to optimize the adjusted Rand index (ARI) with respect to expert annotations and heuristic based around downsampling rare cell types while minimizing runtime. Germain et al. (2020) similarly consider the ARI but also the average silhouette width to maximize cluster purity. Zhang el al. (2019) consider a range of heuristics including agreement with simulated data and robustness to model misspecification.
However, given that these objectives are all heuristic and open to user preference, there no guarantee that all of them are useful and have maxima that align with the meta-objective at hand, which in the above example is the ability to identify biologically relevant population of cells. Conversely, some heuristic objectives may be non-useful they are largely noisy and attribute nothing to the overall optimization problem by not aligning with a meta-objective. For example, in an Imaging Mass Cytometry experiment, which also aims to cluster cells, an antibody that quantifies the expression of a given protein may fail entirely, which would not be identified prior to data analysis. In that case, any objective that included that protein's expression would be irrelevant to the meta-objective, but this would not. be known up front. This motivates the central question we attempt to address: how can we adapt AutoML approaches to optimize unsupervised workflows over multiple hcuristic objectives that are frequently subjective and confficting?
To begin to tackle this question, we introduce MANATEE (Multi-objective bAyesiaN optimizAtion wiTh hluristic objiEctive)). The key idea is that by considering a linear scalarization as a probabilistic weighting over (heuristic) objective inclusion, we may up- or downweight an objective based on desirable or non- desirable properties of its posterior functional form. Consequently, rather than returning the full Pareto front that may include points (parameter values) that maximize potentially non-uscful heuristic objectives, we automatically concentrate on useful region in accordance with the specified properties. In this work, we evaluate method only for univariate input optimization; while we designed our framework to be extensible to multivariate inputs, we do not make any claims with regards to performance in that setting here. The main contributions we present are:
1. Introduce the concept of behaviours B of the posterior functional form of the surrogate objective function f that are desirable if function is useful for overall optimization.
Suggest an example set of such behaviours that may be inferred from the posterior of a multi-output. Gaussian process, used as the surrogate function.
3. Build upon previous MOBO approaches with random scalarizations that compute the distribution of scalarization weights p(A) but instead condition on objective behaviours with p(A(B), inferring which objectives are useful.
Construct a set of example objectives measuring workflow success for real molccular imaging and tran- scriptomic data and show that the proposed procedure compares favourably existing approaches for univariate input optimization.
2

---

Published in Transactions on Machine Learning Research (03/2023)
2 Background
2.1 Bayesian optimization
Bayesian optimization (BO, see Frazier (018 and references therein) attempts to optimize a function g(x) E R for some I E IR that is, in some sense, expensive to evaluate and for which derivative information is not available, precluding gradient-based approaches. Note that while we set the scene in the general case of a multivariate input, in this work we perform experiments with univariate inputs, though our approach can be extended to D > in the future. Applications of BO have become popular in the tuning of ML hyperparaneters (Turner et al. 2021 and indeed entire workflows (uusi et al., 2018) due to the expensive nature of ro-training the models.
At their core, BO approaches propose surrogate function S defined on the same range and domain as g that may be searched efficiently to find points x that either maximize reduce uncertainty about f, or both. This leads to the concept of an acquisition functior acq(x) E IR that may be optimized to find the next æ at which y may be evaluated. While multiple acquisition functions have been proposed, here we focus on the Upper Confidence Bound (UCB) (Auer 2002) defined as:
acquen (z) (x) + (x)
(1)
where (z) and ((t) are the posterior mean and standard deviation of f at after t acquisitions from g while 3e is a hyperparameter that controls the balance between exploration and exploitation. While there are
many possible choices for the surrogate function S, including deep neural networks (Snoek et al. 2015), a
popular choice is a Gaussian process due to its principled handling of uncertainty and capacity to approximate a wide range of functions.
2.2 Gaussian processes
Overview Gaussian processes (GPs) [Williams & Rasmussen 2006) define a framework for performing inference over nonparametric functions. Let m(x) be a mean function and k(x, x²) a positivo-definito covariance function for x, x' e RD. We define f(x) to be a Gaussian process denoted f(x) ~ GP (m(x), k(x, x')) if for any fiite-ddimeniional subset x = [r,,..., xxx]" E NRD the corresponding function outputs f [f( f( follow multivariate Gaussian distribution p(f]x) (O, K), where K is the covariance matrix with entries (K),, = 1z) and we have assumed zero- mero-mean function without loss of generality. The kernel fully specifies the prior over functions, with one popular choice we use throughout the paper being the exponentiated quadratic kernel k(x,x') = exp
(m==")² 12
It is common to model noisy observations y via the likelihood p(ylf), which when laken to be N(f, with noise variance o² leads to the exact. marginalization of
Multi-output Gaussian processes GPs may be extended to model K distinct outputs via the functions {f [Boillla et al. 2007). One construction is to model the full covariance matrix as the Kronecker product between the K x K inter-objective covariance matrix KIO and the data covariance matrix:
cov (fa(z), for (a')) (K (k )e,
(2)
Here the kernel hyperparameter is shared across objectives, though in the following we model objective-speeific observation noises Ek with variances o², as p(u)lfs)) ~ o(fk.o??)
2.3 Multi-objective Bayesian optimization
Multi-objective optimization Multi-objective optimization attempts to simultaneously optimize K objectives 91 ( 9K (I) over 2 E RD, which is common in many real-world settings (Deb) 2014). However, it is rare in practico to be able to optimize all K functions simultancously and instoad is common to attompt to recover the Parelo front. We say a point æ1 Parelo dominated by 32 ifl gx(x)) S ga (x22 Vk = 1,...., K
We will later refer the "sbjeiee acquisttion function".
2 2Commonly referred to tasks, we here refer them as objectives given the application.
3

---

Published in Transactions on Machine Learning Research (03/2023)
and 3 k .....KKK s.t. gx(x1) any other point. The Parelo front is then defined as the set of Pareto optimal points, which intuitively corresponds to the set of equivalently optimal points given no prior preference between objectives.
(x())). A point is said to be Pareto optimal if it is not dominated by
Scalarization functions One popular approach to multi-objective optimization is the use of scolarization functions (see Chugh (2020) for an overview). scalarization function sx(g(z)) parameterized by takes the set of functions g 2 (z)...... 9K (a)] and outputs a single scalar value to be optimized in lieu of g(xr). Roijers et al. (2033) show that if sx monotonically increasing in all go (x) then the resulting optimum x² lies the Pareto front of g.
While many scalarization functions exist, one popular choice is the linear scalarization function (x(())) E (ag)(a)), A 0 Ak. This has the intuitive Explainability interpretation that each A corresponds to the weight of function k, with a larger relative value pulling the optimum of Sx towards the optimum of Ok.
Hypervolume improvement Another multi- objective optimization approach relies on the notion of hypervolume (HV), the volume of the space dom- inated by Pareto front and bounded from below by a reference point, which current work assumes to be known by the practitioner (Daulton et al. 2021 HV is used as a metric to assess the quality of a Pareto front and is sought to maximized in the optimization. Expected HV improvement (EHVI) for set of points be computed using box decomposition algorithms Yang et al. 2019a).
Inter-objective agreement
Multi-objective Bayesian optimization with scalarizations Multi-objective Bayesian optimiza- tion (MOBO) approaches that use scalarizations op- erate under the same conditions as BO, where each evaluation of gx(x) expensive and derivative in- formation is unavailable. An example method is ParEGO (Knowles 2006). which randomly scalarizes objectives with angmented Chebyshev scalarization and uses expected improvement. It was recently ex- tended to gNParEGO (Daulton et al. (220),, which supports parallel and constrained optimization in a noisy sctting. Unlike hypervolume-
Maximum not at boundary
Figure 1: Cartoon illustrating proposed desirable be- haviours of objectives. Explainability captures how much an objective covaries with the parameter x, favouring those with lower observation noise. Inter- objective agreement favours objectives that agree with each other. Maximum not at boundary determines whether the optimum is contained within the user- specified parameter range.
methods which can struggle with > objectives (Balandat et al. 2021 gNParGGO is more suited for such problems.
Paria et al. 2019 propose a MOBO procedure that, rather than maximizing s for a single X, constructs a distribution p()) and minimizes the expected pointwise regret,
R(X) = Ep(2))
max sx(g()) - max (B(x))) #EX
xEX
where is the feature space of x and is the subset of X lying on the Pareto front to computed. The exact region of the Parcto front be considered governed by p()) and the authors provide a bounding box procedure for the user to select p(A), akin in the case of linear scalarization asserting a priori which objectives important. However, to our knowledge, MOBO approach has proposed (()).), inferred
4

---

Published in Transactions on Machine Learning Research (03/2023)
from either the data the posterior over functions, that adaptively up- or downweights objectives based desirable properties.
Multi-objective Bayesian optimization beyond scalarizations For hypervolume- baased methods in the MOBO setting, EHVT has been extended to parallel evaluation of q points, leveraging automatic differentiation and boosting efficiency (Daulton ct al. 2020). As EHVI assumes noisc-frce case and can be affected in noisy seltings, recent work introduced noisy EHVI (NEHVT), which uses ils expectalion under the posterior distribution of the surrogate function values given noisy observations (Daulton et al. 2021) p. NEHVI is more robust to noise than other hypervolume-basedd MOBO methods, is equivalent to EHVI in the noiseless setting, and its parallel formulation (gNEHVT) achieves computational gains and state-of-the-art performance large batch optimization (Daulton et al. 2021
Another approach MOBO utilizes uncertainty, overcoming the issuc of scalability with the number objectives faced by hypervolume- based methods. Predictive entropy search for multi-objective optimization (PESMO) aims to minimize entropy of the posterior distribution over the Pareto set (Hernández-Looatle et al. 2016 Max-value Entropy Search for Multi-objective Optimization (MESMO) employs efficient output space entropy search, improving its computation time over PESMO (Belakaria et al., 2019). Uncertainty-aware Search framework for optimizing Multiple Objectives (USeMO) selects points that maximize a multi-objectivv measure of uncertainty and outperforms existing methods on problems with up to six objectives with faster convergence (Belakaria al. 2020
2.4 Applications of AutoML and Bayesian optimization in molecular biology and genomics
AutoML and BO approaches have previously been successfully applied across multiple problems in molecular biology and genomics. One example is a popular application of BО to generate protein candidates at the sequence level with desirable chemical properties YYnng et al. 2019b). While less common, there are handful of examples of AutoML applications to hyperparameter optimization in genomics. The GenoML project [Makarious et al.) 2021) provides a Python framework centered on open science principles to perform end-to-end AutoML procedures for supervised learning problems in genomics. AutoGeneS (Aliee & Theis 2021) develops a multi-objective optimization framework for the selection of genes for the deconvolution of bulk RNA-sequencing without relying on marker genes and instead optimizing properties of identified cell clusters. However, to knowledge, there bioinformatics and genomics workflows in the absence of well-defined objective functions. In contrast, there are multiple BO techniques that allow user to express a preference between solutions (González et al., 2017/. While these methods could have exciting applications in genomics, we consider an alternative setup, where the user expresses a preference on the functional form of the unsoen objectives but docs not participate during acquisition.
work that tackles the general problem of optimizing
3 Multi-objective Bayesian optimization over heuristic objectives
3.1 Setup
We assume we have access to K noisy, heuristic objectives that at acquisition step t return a measurement yor for input location 2, E X, where X' is a compact subset of R on [a, bj. We introduce surrogate functions fo(x) that model with multi-output GP as described in Section 2.2 with full kernel given by Eq. 2 The choice to fit multi-output GP to data reflects our prior assumption that the heuristic objectives may have a correlated functional form or operate on similar lecallscales. In settings where these assumptions do not hold, using a multi-output GP may be suboptimal. Our framework is applicable to any scalarization function that uses weightings (e.g. Tchobyshev scalarization (Nakayama et al. 2009)) and here we consider a linear scalarization function over objectives Ss (f(r)) (z). Ultimately, we seck to maximize (P(A))ss>(((()))
5

---

Published in Transactions on Machine Learning Research (03/2023)
3.1.1 Acquisition functions
The next point to query F++1 is chosen by maximizing the expectation of the acquisition function. For this propose two approaches:
SA (scalarized acquisition): maximize the expectation of the scalarization of the single-objective acquisition function of each objective By(x))) [s> (a)q(((()))) per Paria et al. (2019);
AS (acquisition of scalarized): maximize the expectation of the single-objective acquisition function of the scalarized objectives Ep(al-) [ac](=)(((=))))
We refer to these as acquisition functions and derive expressions for both in Appendix C While many
choices of single-objective acquisition functions are possible, we use acqucB (Eq. D., the UCB single-objective acquisition function. The SA formulation simplifies to an intuitive interpretation where cach objective's UCB
function value is weighted by the probability of that objective being useful given its behaviours (Appendix C The AS formulation takes into account the multi-objective posterior covariance structure (Appendix C but
has longer computation time that may require approximations when K is large.
3.2 Desirable heuristic objective behaviours
Next, we wish to set p(a|-) to upweight objectives that are inferred as useful based on desirable properties learned from the data. In our framework, we assume that the practitioner has a priorí no preference over objectives, only that some or all may be useful. Instead, they have a preference over the functional form of the objectives, which is expressed via some desirable properties. We begin by considering what properties of given heuristic objective f& (z) may be considered desirable. While many are possible, we suggest three behaviours (Figure
1. Explainability: (x) covaries significantly with æ (i.e. is explained by z). The justification here is that the practitioner has selected heuristic & assuming it will provide insight into the choice of x, so if there is no correlation then it should be downweighted. Given that the data have been scaled to empirical variance 1, the inferred variance of the fitted observation noise o², represents the proportion of variance unexplained by Sa so we define B ==
Inter-objective agreement: fa shares a similar functional form with fwek's + k. with the intuition that it is useful for practitioners to find regions of the input space where multiple heuristics agrec. After fitting the multi-oulput GP, (K19),, ) defines the covariance ween objective & and k' for k & k and (K!9) )x, defines the variance of objective k. We therefore introduce the inter-objective agreement behaviour
Σ max 0,
1
(K10) KK
B
k
K (KIO)
W=1:k
The intuition is that
Vikio (K!0)) represents the correlation between objectives k and k' so B[(0)
(K10
represents the average correlation with other objectives while not penalizing negative correlation worse than correlation.
3. Maximum not at boundary: Within X, fu has maximum that is not at the boundary of z. The useful range of x is specified by the practitioner. Then, if fa is maximized by a boundary value of then either (i) the optimum is outside of the specified range, conflicting with the practitioner's intuition, or (ii) fü is unbounded in x, in which case is not useful for optimization. In the former case, one approach would be to revise the range and repeat the process. Since the derivative of GP is also a GP, we may identify whether a stationary point exists in X by searching for the zeros of the posterior mean derivative ['(x). We therefore define B(a) = hasmax(si, X)), where hasmax returns 1 if fs has a maximum on and otherwise by evaluating the derivatives of the posterior mean of the multi-output GP (derived Appendix D
We will denote set of three above behaviours for given objective k B (B(.) (1) ((2) B

---

Published in Transactions on Machine Learning Research (03/2023)
Acquire data D for K heuristic objectives
2 Fit multi-output Gaussian process surrogate P(
Infer objective behaviours p(B]
X Acquire next point
(i()))
Houristic objective
1) Explainabiiity
2) Inter-task agreement
3) Maximum boundary
@ Update objective weights
fa( Houristic
p (Ag B&) ix p. By Ax (p()
2 Optimize scalarized acquisition function
En, P(Ax BA) [fcg(((>))]
xx (x) Heuristic obiective A
Figure 2: Multi-objective Bayesian optimization with heuristic objectives. A multi-output Gaussian process is fitted (Step 2) to the initial training dataset (Step 1). Objective behaviours are inferred from the posterior form (Step 3) and are used to update the distribution of objective inclusion weights (Step 4). The acquisition function defined as the expectation under the weights distribution is optimized (Step 5) to give the next location sample objectives at (Step 6). Step shows the scalarized acquisition (SA) function approach.
3.3 Incorporating desirable behaviours into scalarization weights
We next consider how to use the set of behaviours B = [B*] to parameterize the scalarization probabilities P(AlB). We assume that d is binary variable Ak that corresponds to whether objective k is useful or otherwise, with p(Ax|Ba given by a Bernoulli distribution. While this construction initially appears restrictive compared existing approaches with random scalarizations where À R.K it has a desirable property that maintains its generality (note that we assume SA be linear scalarization function). Specifically, optimization under our setup returns Pareto optimal points of f (proof presented in Appendix E.
Theorem 3.1. If Ep(A, B.) [Ax] > 0 Vk, then x- arg max, (A lies on the Pareto front of f.
However, how construct p()= = 1|B ( ( directly non-obvious. Instead, we ask how would each objective behaviour appear if we knew that objective was useful (or otherwise)? This allows us to specify conditional distributions over behaviours for a meefulnnnn-seefu objective p( (B( = 1), p( (BB(' = 0) for = 1,2,3 3 and combine with a prior p(Ax = = p((x = 0) to compute the posterior probability that an objective useful given its behaviours:
II, PLB 14 = 1) =1
1(Bb)
Σx=u,1 11, p( ( (p(()pp))))))
With these considerations, suggest distributions for p p(B( Axx); however, emphasize that these are suggestions only and there many possible that would fit the problem.
Explainability: For B(a) the explainability of objective & (i.c. the proportion of variance unexplained by that function), we assume that if that objective is desirable (Ax = 1) then the lower the observation noise, the better and in the non-desirable case (Ax = 0), higher noise is expected. Given the lack of additional assumptions, appeal to the principle of parsimony and propose lincar relationship of the form:
P(B(1)As)) = 2(1 A BB + 2Ax(1 - B( B E [0,
5
otherwise.
7

---

Published in Transactions on Machine Learning Research (03/2023)
Inter-objective agreement: For inler-objective agreement BR our reasoning that identifying solutions where multiple objectives agree can be desirable for practitioner. We thus posit that high inter-objective correlation should be more likely for a desirable objective and vice-versa for a non-desirable one, and again a linear relationship is the most parsimonious:
2 _ + Ax)(1 B if B E [0.1]
B (b6))
otherwise.
Maximum not at boundary: We propose 1>==i Bernoulli(T.) where To, α1 are user-settable iyper-parameters. This means that conditioned on an objective being useful (or otherwise), there a fixed probability of that objective containing a maximum the region. In our experiments, we set #0, क such that useful objective has a maximum in the region with 75% chance, and non-useful one with 25% chance.
p(B (BB ARx
Bernoull((zo) Bernouli(m))
if A5 =
(7)
if A
3.4 MANATEE
Putting these steps logether results in the MANATEE framework, an ilerative MOBO procedure outlined in Figure 2 and Algorithm [ First, the objectives are evaluated at a set of input locations randomly chosen on the paramcter space. Second, the multi-output GP surrogate function with covariance given by Eq. 2 is fitted to all objectives. Then, the objective behaviours are computed from the surrogate function and the distributions over objective weights that indicate whether each objective is useful are updated. Finally, the updated acquisition function optimized (see Section 3.1.1 and Appendix C guiding acquisition of the next point. The procedure is repeated for predetermined mumber of steps. The overall "best" point to be used for downstream analysis may be chosen as that which maximizes the scalarized surrogate function.
Algorithm MANATEE framework
Input: training data D("; input space X on region [a, b]; objectives f(x) = (((())) acquisition func- tion (ff 8x, acq); behaviours B (B("), B(²), B( distributions over behaviours for
objective p( (de); number of iterations for f 1....,
; T do
Fit multi-output Gaussian process p(f(D(t-1)
for k 5- 1,..., K do for i +-++1 1, do Compute objective behaviour _ from the Gaussian process posterior (Section 3.2
Compute p(B(() Xx)) (how likely B(() is given objective is useful/non-usefig using Eq.
6
end for Update p(Ax = 1 [B₂) (probability that objective is useful given its behaviours) using Eq. end for
Update acquisition function A(f(2), sx, acq) using objective weights p(AlB) (Section 3.1.1 Pick next candidate by maximizing acquisition function Ia < arg max A(f(r), $A acq) Acquire data at new location (*x), Vk = l,.... K objectives Update training data with new datapoints D(t) + D((-1) y(E)))))))
end for
3.5 Baselines for experiments
We contrast our method against two baselines and three existing approaches for MOBO: (i) Random acquisition: draw 2, ~ Unif(00 at cach iteration, (ii) Random scalarization: use identical surrogate and acquisition functions as MANATEE-SA
sample z, but draw d ~ Unif(0, 1) rather than conditional on B, (iii) @NEHVI (Daulton et al. 2021 with approximate hypervolume computation to facilitate inference over
8

---

Published in Transactions on Machine Learning Research (03/2023)
A
sptimal Osjective
Offective
Objective
1.
1.0
1.
1.0
9.
0.0
B
Objective
Oéjective 2
Oojective 3
Objective
Objective 5
00x 11bc
0.2
Acquisiion
Acquisition
Acquisiiion
Acquisition
Acquisition step
Figure 3: A 150 random samples of toy data for the objectives, including those on the Pareto front (orange) and otherwise (blue), along with points acquired by MANATEE-SA (red) in an example run. B Inclusion probabilities for each of the objectives a function of acquisition step. Solid line shows the mean and shaded region denotes the 95% confidence interval across runs.
> 5 objectives, (iv) NNarEGO (Daulton et al. 2020), and (v) UScMO (Belakaria et al.; 2020 Since it is not known how noisy real biomedical problems are, supplement our evaluation including NEEVV and gNParEGO, specifically designed for MOBO of noisy objectives, with USeMO, an efficient uncertainty-based state-of-the-art approach for many objectives in the noiseless setting (Belakaria et al. 2020).
In our setup not all objectives are considered useful, so quantifying the hypervolume improvement over all objectives as a commonly used performance measure (Daulton et al., 2020 2021) does not align with the stated goal. Thus, we construct evaluation measures (referred as (neta-objectivvs) that use independent quantify the quality of the acquisitions.
information (e.g. expert labels, see Sections 4.2. 4.3
When a meta-objective h(x;)) is available at every iteration t 1...., T with overall maximum V" = max,ex h(x), to compare among approaches we compute the following metrics: (i) Cu- mulative regret: E (v* - h((x)),, (ii) Full regret: y" max h(x), and (iii) Bayes regret: + E (v" max h(x)), where X1 is the set of x acquired up to time 1. Of these, we place most emphasis on cumulative regret as it quantifies how close each method gets to the optimal solution on average. In contrast, the full and Bayes regrct quantify how close the "best" acquired point gets to y" as mcasured by the max over h of all points acquired so far; however, since the meta-objective our problem setup (and only used for method comparison), it is impossible to quantify max h(x) in practice outside of benchmarking exercises.
in general inaccessible for
4 Experiments
4.1 Toy data experiment
We begin by remonstrating the overall problem setup on toy data on an input space [0, 1]. We consider
5 objectives overall - 3 that act as the useful objectives with maxima around that of a meta-objective at 1/4 given by sin 2T.., max(O, sin 2Tzz), and sin (T(x= 0.05) and 2 that disagree and act as the non-useful objectives given by 2r and -2r. Each objective is augmented with noise (Appendix B.I. Note that on real data we do not know a priorí which objectives are usefu[]] Further, the meta-objective is not specified it may be linear, non-linear, and not necessarily a function of the heuristic objectives it simply needs a maximum at r 1
Qtherwise only useful objectives would he included and standard MORO procedures applied.
9

---

Published in Transactions on Machine Learning Research (03/2023)
Samples from each of these functions can be seen in Figure 3A (blue points). The overall Pareto front (orange points) spans almost the entire region including samples at the very right where one of the non-useful linear objective functions has its maximum. However, when applied to this toy problem, MANATEE quickly begins acquiring samples around the joint maxima of the three useful objective functions (red points). Indeed, tracing the inclusion probabilities p(>x = 1[Bx) across the iterations (Figure 33) demonstrates how MANATEE learns to upweight objectives 1-3 while downweighting 4-5. This demonstrates than when we do not know priori which objectives trust, we may still recover a region of high utility when the Pareto front spans
the full space of conflicting objectives.
4.2 Imaging Mass Cytometry cofactor selection
We next apply MANATEE the selection of cofactors for Imaging Mass Cytometry (IMC) data, new technology that mcasure the expression of up to 40 proteins at subcellular resolution in tissue sections Giesen et al.
2014 In the analysis of mass cytometry data, a cofactor is frequently used to normalize
the data (Ray & Pyne 2012 Wagner et al. 2019) via the transformation ÿ = sinh (w/e). However, to our knowledge no systematic approach exists to set the cofactor and it is typically left as a user-specified parameter.
Table 1: Results for IMC cofactor optimization experiment. CR: cumulative regret, FR: full regret; BR: Bayes regret. M-SA: MANATEE with scalarized acquisition, M-AS: MANATEE with acquisition of scalarized function, RA: random acquisition, RS: random scalarization. ARI: adjusted Rand index, NMI: normalized mutual information. Values are mean (s.d.).
ARI
NMI
Method
CR
FR
BR
CR
FR
BR.
M-SA
0.017(0.005)
0.003(0.003)
0.007(0.006)
0.019(0.009)
0.002(0.005)
0.008(0.009)
M-AS
0.006(0.010) 0.021(0.010)
0.011(0.011)
0.025(0.017)
0.007(0.016)
0.015(0.017)
RA
0.045(0.001)
0.024(0.013)
0.031(0.008)
0.065(0.003)
0.026(0.013)
0.036(0.009)
RS
(0020()
0.003().004)
0.008(0.005)
0.025(0.006)
0.003(0.006)
0.008(0.007)
NEHVI
0.021(0.011) 0.011(0.007) 0.042(0.004)
0.011(0.007) 0.061(0.007)
0.026(0.017)
NNarEGO
0.002(0.003) (0000(
0.017(0.012)
0.049(0.009)
0.005(0.003)
(0.277
USeMO
0.010(0.009) 0.043(0.003)
0.016(0.009)
0.062(0.005)
0.012(0.011)
0.018(0.013)
Here, we consider the standard workflow where (i) the expression data normalized with a given cofactor and (ii) the data is clustered using standard methods with the "best" cofactor being the one that leads to the most biologically relevant cellular populations Given that this problem in general has notion of "test accuracy" with respect to which we could optimize the cofactor, we instead suggest a number of heuristic objectives based around maximizing the correlation of cluster-specific mean expression of known protein marker combinations. For example, the proteins CD19 and CD20 are highly expressed in B lymphocyte cells and lowly expressed in all others. Therefore, a clustering correctly separates B cells from others, the correlation between the mean CD19 and CD20 expression in each cluster should be high as the proteins should either be co-expressed or both not expressed (at the origin), as demonstrated in Appendix F.1 We can apply this logic to range of cell type markers to construct set of heuristic objectives (Appendix
F
To quantify the ability of each clustering to uncover biologically relevant populations, we use expert annotated cell types from Jackson et al. (2020) and assess cluster overlap with the adjusted Rand index (ARI) and normalized mutual information (NMI), which for this experiment form the overall meta-objectives 0 in line
with prior benchmarking efforts of single-cell clustering (Qi et al. 2020 Kiselev et al.) 2017). Note that this is in general unavailable for the analysis of newly generated data and we would only have access to the correlation (heuristic) objectives.
4 All paramcters of the clustering procedure are held constant across colactors to allow for fair comparison. The fact that we can easily specify meta-objectives highlights the ubiquity of the "multiple heuristic objective" issue in bioinformaties.
10

---

Published in Transactions on Machine Learning Research (03/2023)
The results comparing MANATEE to the alternative methods are shown in Table I and the optimization trajectories across acquisitions are shown in Supplementary Figure 4 On the metric of cumulative regret, which as above, is most relevant for the problem setup at hand, MANATEE-SA outperforms the alternative approaches. On full and Bayes regret, MANATEE performs comparably with the basclines. On cumulative regret, NEHVT and USeMO are comparable to random acquisition, suggesting that consistently acquiring close-to-optimal solutions over > 5 noisy objectives challenging even for approximate hypervolume computation and that noisy setting is challenging for methods assuming noise-free observations. Interestingly, find that random scalarization exhibits strong performance on several measures, which may be understood by the fact that the scalarized objective Σ& fa naturally places high weight on regions where many objectives agree, mirmicking similar scenario to our inter-objective agreement criterion. Samples from the objectives along with methods' acquisitions are shown in Supplementary Figure 5 MANATEE samples low cofactor values corresponding to regions where meta-objectives' values are high, while other methods sample throughout the paramcter range. Examining the inclusion probabilitics shows that MANATEE learns to downweight the /CD2 co-expression objective, which is maximized for high cofactor values (Supplementary Figures 5 and 6 We also show each regret metric function of the acquisition step to demonstrate the convergence rate of each method (Appendix A.5 Supplementary Figure 10a). MANATEE-SA reaches low regret values faster than or comparable to other methods.
We further performed ablation experiments of each behaviour and found that no single behaviour drives the performance (Appendix A..).. We also performed cross-valt validation on data splits to demonstrate that MANATEE does not overfit to a given dataset (Appendix A.3)
4.3 Single-cell RNA-seq highly variable gene selection
Single-cell RNA--eqqeencing (scRNA-soq, see Hwang et al. 2018 for overview) quantifies whole- transcriptome gene expression at single-cell resolution. A key step in the analysis of the resulting data is selection of a set of highly variable genes (HVGs) for downstream analysis, typically taken as the "top x%% (Yip et al., 2019). but there are no systematic or quantitative recommendations for selecting this proportion (Luecken & Theis 2019). Therefore, we apply MANATEE to this problem following a clustering workflow similar to the IMC experiment, but by varying the proportion of HVGs used for the analysis and koeping all other clustering parameters fixed. We again propose a number of co-cxpression basod hcuristies (Appendix F.3 and augment. these with measures of cluster purily (mean silhouette width, Calinski and Harabasz score,
score) previously used in scRNA-seq analysis (Germain et al. 2020).
Table 2: Results for scRNA-seq HVG selection oplimization experiment. CR: cumulative regret, FR: full regret; BR: Bayes regret. M-SA: MANATEE with scalarized acquisition, M-AS: MANATEE with acquisition of scalarized function, RAA: random acquisition, RS: random scalarization. ARI: adjusted Rand index, NMI: normalized mutual information. Values are mean (s.d.).
ARI
NMI
Method
CR
FR
BR
CR
FR
BR
M-SA
0.126(0.019)
0.056(0.010)
0.064(0.012)
0.125(0.026)
0.022(0.010)
0.031))
M-AS
0.127(0.026)
0.058(0.015)
0.067(0.015)
0.124(0.036)
0.023(0.013)
0.034(0.017)
RA
0.192(0.020)
(.059(0.009)
0.067(0.013)
0.199(0.025)
0.025(0.009)
0.040(0.015)
RS
0.140(0.018)
0.049(0.008)
0.059(0.010)
0.130(0.021)
0.014(0.007)
0.026(0.012)
NEHVI
0.186(0.038)
0.050(0.008)
0.092(0.045)
0.023(0.009) (..007)
(.073))
PNrEEGO
0.161(0.039)
0.055(0.009)
0.091(0.048)
0.158(0.049)
0.023(0.009)
0.070(0.064)
USeMO
0.057(0.013) (..03)
0.071(0.014)
0.027(0.014) (.230()
0.042(0.019)
For these workflows, no general ground truth clustering or cell types are available. However, a new technology called CITE-seq can simultaneously quantify both the RNA and surface protein expression at single-cell level (Stocckius et al.] 2017). Given that cell types are traditionally defined by surface protein expression (Oostrum et al.) [2019), we use a clustering of the surface protein expression alone as the ground truth following existing work (Liu et al. 2021 The concordance with this clustering acts as the meta-objective
11

---

Published in Transactions on Machine Learning Research (03/2023)
in this experiment, which we benchmark the proposed approaches against. We supply each method with the heuristic objectives above and benchmark the gene proportion acquisitions by contrasting the resulting clusterings with the surface protein-derived ground truth using ARI and NMI as metrics. Once again, these represent only two possible choices of meta-objective and there are many more we could design, highlighting the prevalence of heuristic objectives in the field.
The results are shown in Table 2 the optimization trajectories in Supplementary Figure 7 and regret metrics (Appendix A at each acquisition step in Supplementary Figure 10p. As above, our main focus is on cumulative regret since when deploying in a real-world scenario, we would not have access to the meta-objective. MANATEE performs favourably on cumulative regret compared to the other approaches, though has higher full and Bayes regret. MANATEE learns to strongly upweight the Davies-Bouldin objective, which agrees with the meta-objectives and is less noisy, and downweight the CD3E/ /CDE/CD4 and CD3D FFBD/CDD /CDSA objectives, which are maximized in a different region to the meta-objectives (Supplementary Figures and 9 This allows our method to acquire proportions of HVGs close to the upper bound of the parameter range where the mela-objectives are maximized, unlike gNParEGO and USeMO that don't reach this region (Supplementary Figure Overall, this demonstrates our method to be a promising approach to tackle hyperparameter optimization on real, noisy datasets and achieve competitive performance compared to existing basclinos and sfateoo-thhe-r methods.
Finally, we quantified the effects of different model initializations on the parameter values selected by our method by reanalyzing the out puts of the above experiments, which contained multiple repetitions across different random seeds, controlling the imitial training set and model initialization. In the IMC cofactor selection experiment, MANATEF-SA and MANATER-AS most frequently selected the same cofactor value in 99% and 90% of all repetitions (Supplementary Figure I.).. In the HVGs selection experiment, the interquartile range of percentage values most frequently selected by MANATEE-SA in all repetitions was 5% and that for MANATEE-AS was 8%, while the considered range of optimization was between 1% and 50% of genes (Supplementary Figure 11p). All implementation details are provided in Appendix A.1 and the experimental setup is described in Appendix A.2
5 Discussion
A common theme here is the subjectivity of parameter setting in biological data analysis workflows. Setting these often involves no hcuristic objectives at all, simply relying on itcrative data exploration to find a parameter combination that "works". Even when heuristic objectives are involved -ucc as in the benchmarking analyses of scRNA-seq workflows the precise choice of which objectives to include is fundamentally subjective too.
It important to note that our proposed approach does not remove subjectivity from the analysis. Many important steps, including the chosen behaviours B and their conditional inclusion distributions p(B]A) are set by the user. Therefore, it abstracts the subjectivity by a level, changing the question from "which objectives should I use to benchmark my method?" to "what would the behaviour of good objective function de?". Given that no link is assumed between the specified heuristic objectives and the true meta-objective and that the choice of desirable objective behaviours is given as example only, we make no optimality claims about the ability to explore the Pareto front.
In our sctup, weighting of the acquisition functions may effectively exclude some objectives and lead to potential under-exploration. This may be mitigated by starting the objective weighting process after some set number of steps. Furthermore, the inter-objective agreement behaviour may be potentially limiting if a practitioner is interested in jointly suboptimal solutions that still perform reasonably well all objectives, some of which are competing. If all objectives had the maximum at the boundary, this behavior would no longer be relevant to the optimization; in such cases, the user could be alerted and asked to re-define which behaviours they want to use. As each behaviour's weight depends on the specification of ())),, there is no guarantce that these would have a similar scale, though we found no one behaviour driving the performance on the problems we considered here. Here we have assumed that all tasks are quantifiable for all observations, though in some Bayesian optimization settings this doesn't 't hold KKause & Ong 2011].
12

---

Published in Transactions on Machine Learning Research (03/2023)
Our framework can be used with other scalarization functions that include weighting À which may be more suitable in some scenarios, but here we performed experiments using linear scalarization. We note that if the true meta-objective is not linear in the observed objectives then a linear scalarization will be suboptimal. However, in the context we considered, the true meta-objective may represent complex notions such as the ability to uncover biologically meaningful results, which, while likely nonlinear, is impossible to quantity or specify. Therefore, the linear scalarization represents a trade-off between having access to an intuitive weight A for each function while not necessarily being theoretically optimal. In addition, MANATEE performed competitively with state-of-the-art methods in our experiments with real data analysis pipelines, which represent cases with unknown true meta-objective.
As our general approach is applicable to any multi-objectivee optimization scenario (while tailored to biomedical analyses), we acknowledge that it could be used in highly diverse applications. We note that these may include ethically dubious bioinformatics analyses such as those pertaining to genetic testing of embryos. We strongly caution against any such use without thorough ethical review process.
There several extensions that would serve future steps. We have only considered optimizing E E, but future work can use our method to optimize multiple pipeline parameters as our work generalizes to 2 E R for D > 1 (Appendix Dr. There is much current research in BO methods over both contimuous and categorical domains (Ru et al. 2020), which may better suit. the parameter space of seRNA-seq analysis pipelines (Germain et al.] 2020 A lot of research in BO centres on the incorporation of user input and expert opinions to guide optimization (Häse et al.) 2021 Abdolshah et al. 2019). While we have explicitly considered the opposite problem where a priori it not known which objectives should be upweighted there could be situations where both approaches could be integrated. For example, expert may provide ratings for the results of each scRNA-seq clustering during optimization. In such settings, these ratings could be integrated into our proposed framework by updating the distributions p()[B, Q) over θ such that they confer high weights to functions of expert ratings. Finally, we welcome future work evaluating gains in discovery and accuracy of biological results and computation time arising from using our mothod to choose values for pipeline parameters.
Code availability
Our code is available at
http://mmmmm2mmmmmataeeeeeee www.mmaaaaaaaaaeeee manatee_paper _apper__ppp_____ wttppppppaalllllaapaapaapaa a_nallsasss
Acknowledgments
We thank Dr. Gavia Gray for her advice on PyTorch optimizers and the reviewers for their thoughtful comments and suggestions. This work was supported by a Canadian Institutes of Health Research project grant number PJT175277 to K.R.C, the Canadian Statistical Sciences Institute Ontario Top-up for Postdoctoral Fellows in Data Science to A.S., the Hold'em for Life Oncology Fellowship to A.S., and the Vector Institute Postgraduate Affiliate Award to A.S. We acknowledge the support of the Natural Sciences and Engineering Research Council of Canada (NSERC) RGPIN-2020-04008. to K.R.C. This research was undertaken, in part, thanks to funding from the Canada Research Chairs Program to K.R.C.
References
Majid Abdolshah, Alistair Shilton, Santu Rana, Sunil Gupta, and Svetha Venkatosh. Multi-objective bayesian optimisation with preferences over objectives. In Hanna M. Wallach, Hugo Larochelle, Alina Beygelzimer, Florence d' Aalcé--Buc, Emily B. Fox, and Roman Garnett (eds.), Advances in Neural Information Processing Systems 32: Annual Conference on Neural Information Processing Systems 2019, NewrIPS 2019, December 8-14, 2019, Vancouver, BC, Canada, pp. 12214 12224, 2019. URL https; //proceedings..euue http://eiw..iiitea//aaaaaaaaaann htw.::o02so---btttttttt...tt
Hananeh Alice and Fabian J Theis. Autogenes: Automatic gene selection using multi-objectivve optimization for RNA-seq deconvolution. Cell Systems, 2021.
13

---

Published in Transactions on Machine Learning Research (03/2023)
Peter Auer. Using confidence bounds for explotaaiion-cclloaaiio trade-offs. Journal of Machine Learning Research, 3(Nov):397-473 2002.
Maximilian Balandat, Brian Karrer, Daniel R. Jiang, Samuel Daulton, Benjamin Letham, Andrew Gordon Wilson, and Eytan Bakshy. Botorch: A framework for efficient monte-carlo bayesian optimization. In Hugo Larochelle, Mare' Aurelio Ranzato, Raia Hadsell, Maria-Florina Balcan, and Hsuan-Tien Lin (eds.), Advances Neural Information Processing Systems 83: Anmual Conference on Neural Information Processing Systems 2020, NewrlPS 2020, December 6- 12, 2020, virtual, 2020. URL httpps:/pproceedingg.. www..ius...cpppppp.... Abstroataabbbtttt.....c.cc html
Maximilian Balandat, Brian Karrer, Daniel R. Jiang, Samuel Daulton, Benjamin Letham, Andrew Gordon Wilson, and Eytan Bakshy, Noisy, parallel, multi-objective BO in BoTorch with qEHVI, qNEHVI, and P.raEGG. httpp:://btoore: http://///....iiaaaaaaaaaaaooeia wwwpe/ittttiiteeeeeeeeeeeeeeeeee boj 2021. Accessed: 2022-01-26.
Syrine Belakaria, Aryan Deshwal, and Janardhan Rao Doppa. Max-value entropy scarch for multi-objective bayesian optimization. In Hanna M. Wallach, Hugo Larochelle, Alina Beygelzimer, Florence d' Alchhé Buc, Emily B. Fox, and Roman Garnett (eds.), Advances in Neural Information Processing Systems 32: Annuul Conference Neural Information Processing Systems 2019, NewrIPS 2019, December 8-14, 2019, Vancouver, BC, Canada, pp. 7823 7833, 2019. URL https:://ppeet htttp::/prroceedings neuripa. cc/paper/2019/ hash/ http://2225222220057744488668880
html
Syrine Belakaria, Aryan Deshwal, Nitthilan Kamnappan Jayakodi, and Janardhan Rao Doppa. Uncertainty- aware search framework for multi-objective bayesian optimization. In The Thirty-Fourth AAAI Conference on Artificial Intelligence, AAAI 2020, The Thirly-Second Innovalive Applications of Artificial Intelligence Conference, IAAI 2020, The Tenth AAAI Symposium on Educational Advances in Artificial Intelligence, EAAL 2020, New York, NY, USA, February 7-12, 2020, pp. 10044 10052. AAAI Press, 2020. URL https://aaas http///aaaa.rr//jj//./ ottg//////dde.... phcedexx.pp///aee Att.AAAAAAaiitteeeieaeaaa
Lukas Biewald. Experiment tracking with Weights and Biases, 2020. URL https: htttp://www.wandb.co wwwdw.anddoooo.... com/
Edwin V. Bonilla, Kian Ming Adam Chai, and Christopher K. I. Williams. Multi-task gaussian process prediction. In John C. Platt, Daphne Koller, Yoram Singer, and Sam T. Roweis (eds.), Advances in Neural Information Processing Systems 20, Proceedings of the Twenty-First Annual Conference on Neural Information Processing Syslems, Vancouver, Brilish Columbia, Canada, December 3-6, 2007, pp. 153 160. Curran Associates, Inc., 2007. URL https://proceediing. http::gg.....r......ppppppeeee
hash//
Abstract. html
Richard P Brent. Algorithms for minimization without derivatives. Courier Corporation, 2013.
Tinkle Chugh. Scalarizing functions in Bayesian multiobjective optimization. In 2020 IEEE Congress on Evolutionary Computation (CEC), pp. 1-8. IEEE, 2020.
Yaxuan Cui, Shaoqiang Zhang, Ying Liang, Xiangyun Wang, Thomas N Ferraro, and Yong Chen. Consensus clustering of single-cell RNA-seq data by enhancing network affinity. Briefings in Bioinformatics, 2021.
Samuel Daulton, Maximilian Balandat, and Eytan Bakshy. Differentiable expected hypervolume im- provement for parallel multi-objective bayesian optimization. In Hugo Larochelle, Marc' Aurelio Ran- zato, Raia Hadsell, Maria-Florina Balcan, and Hsuan-Tien Lin (eds.), Advances in Neural Information Processing Systems 33: Annual Conference Newral Information Processing Systems 2020, NewrIPS 2020, December 6-12, 2020, virtual, 2020. URL httpp::/pproceediiei neurips. httpripps.cc/pppe///200.... /hash/ html
Samuel Daulton, Maximilian Balandat, and Eytan Bakshy. Parallel Bayesian optimization of multiple noisy objectives with expected hypervolume improvement. ArXiv preprint, abs/2105.081155 2021. URL
Kalyanmoy Deb. Multi-objective optimization. In Search methodologies, pp. 403 449. Springer, 2014.
14

---

Published in Transactions on Machine Learning Research (03/2023)
Kalyanmoy Deb, Amrit Pratap, Sameer Agarwal, and TAMT Meyarivan. A fast and clitist multiobjective genetic algorithm: Nsga-ii. IEEE transactions on evolutionary compulation, 6(2):182 197, 2002.
Angelo Duò, Mark D Robinson, and Charlotte Soneson. A systematic performance evaluation of clustering methods for single-cell RNA-seq data. F1000Research, 7, 2018.
P.T. Eendebak and A.R. Vazquez. Oapackage: A python package for generation and analysis of orthogonal arrays, optimal designs and conference designs. Journal of Open Source Software, 2019.
Thomas Elsken, Jan Hendrik Metzen, and Frank Hutter. Neural architecture search: A survey, arxiv e-prints, page. ArXiv preprint, abs/ abs/1808.05377, 2018. URL httppss//aaaxvv..orgaabbb///b//a 0037 05377
Peter I Frazier. A tutorial on Bayesian optimization. ArXiv preprint, abs/1807.022 2018. URL https
7/arxiv. www. abs// 1807.
Nicoló Fusi, Rishit Sheth, and Melih Elibol. Probabilistic matrix factorization for automated machine learning- In Samy Bengio, Hanna M. Wallach, Hugo Larochelle, Kristen Grauman, Nicolò Cesa- Bianchi, and Roman Garnett (cds.), Advances in Neural Information Processing Systems 31: An- wual Conference on Neural Information Processing Systems 2018, NewrIPS 2018, December 3-8, 2018, Montréal, Canada, pp. 3352-3361, 2018. URL https: www.gis...rr...... www..ps..ccpaape..... cc/pappr/2018//aate
Awt.ssattbbstt....a.a..ccaaa.... html
Jacob R. Gardner, Geoff Pleiss, Kilian Q. Weinberger, David Bindel, and Andrew Gordon Wilson Gpytorch: Blackbox matrix-matrix gaussian process inference with GPU acceleration. In Samy Bengio, Hanna M. Wallach, Hugo Larochelle, Kristen Grauman, Nicolò Cesa-Bianchi, and Roman Garnett (eds.), Advances in Neural Information Processing Systems 31: Annual Conference on Neural Information Processing Systems 2018, NewrIPS 2018, December 3-8, 2018, Montréal, Canada, pp. 7587-7597, 2018. URL https html
//proceedings.oo www.nnnnii...mm wwtp:/pwwppppppppppeee/e/e/pceea http:///aaa.a..ccooecceeeccaaaaa
Pierre-Luc Germain, Anthony Sonrel, and Mark D Robinson. pipocomp, a general framework for the evaluation of computational pipelines, reveals performant single cell RNA-seq preprocessing tools. Genome Biology, 21(1):1-28, 2020.
Charlotte Giesen, Hao AO Wang, Denis Schapiro, Nevena Zivanovie, Andrea Jacobs, Bodo Hattendorf, Peter J Schüffler, Daniel Grolimund, Joachim M Buhmann, Simone Brandt, et al. Highly multiplexed imaging of tumor tissues with subcellular resolution by mass cytometry. Noture Methods, 11(4):417-42) 2014.
Javier González, Zhenwen Dai, Andreas C. Damianou, and Neil D. Lawrence. Preferential bayesian opti- mization. Tn Doina Precup and Yee Whye Teh (eds.), Proceedings of the 34th International Conference Machine Learning, 1CML 2017, Sydney, NSW, Australia, 6-11 August 2017, volume 70 of Proceedings of Machine Learning Research, pp. 1282 1291. PMLR, 2017. URL httpp://proceedings. wwweeediinnmmmi......mmmmm.meeee p.l..pees//70// gonzalez17a. html
Yuhan Hao, Stephanie Hao, Erica Andersen-Nissen William M. Mauck I11, Shiwei Zheng, Andrew Butler, Maddic J. Lee, Aaron J. Wilk, Charlotte Darby, Michael Zagar, Paul Hoffman, Marlon Stocckius, Efthymia Papalexi, Eleni Mimitou, Jaison Jain, Avi Srivastava, Tim Stuart, Lamar B. Fleming, Bertrand Yeung, Angela J. Rogers, Juliana M. McElrath, Catherine A. Blish, Raphael Gottardo, Peter Smibert, and Rahul Satija. Intograted analysis of multimodal single-cell data. Cell, 2021. doi: 10.1016/j.cell. 2021..0..
Florian Häse, Matteo Aldeghi, Riley J Hickman, Lofc M Roch, and Alán Aspuru Guzik. Gryffin An algorithm for Bayesian optimization of categorical variables informed by expert knowledge. Applied Physics Reviews, 8(3):031406, 2021.
Xin He, Kaiyong Zhao, and Xiaowen Chu. AutoML: A survey of the -tate-of-the-art. Knowledge-Baod Based Baseed Systems, 212:106622, 2021.
15

---

Published in Transactions on Machine Learning Research (03/2023)
Daniel Hornández-Lobato, José Miguel Hernández-Lobato, Amar Shah, and Ryan P. Adams. Predictive entropy search for multi-objective bayesian optimization. In Maria-Florina Balcan and Kilian Q Weinberger (eds.), Proceedings of the 33nd International Conference on Machine Learning, 1CML 2016, New York City, NY, USA, June 19-24, 2016, volume 48 of JMLR Workshop and Conference Proceedings, pp. 1492-1501. JMLR.org, 2016. URL http:///rooceddings mlr. prees// wwwpss//e/eeeeeeeeee w48/hhrnandez--eee 1ow.lezz-oobataaa66o.. html
Byungjin Hwang, Ji Hyun Lee, and Duhee Bang, Single-cell RNA sequencing technologies and bioinformatios pipelines. Experimental 8 molecular medicine, 50(8))1-14, 2018.
Hartland W Jackson, Jana R Fischer, Vito RT Zanotelli, H Raza Ali, Robert Mechera, Savas D) Soysal, Holger Moch, Simone Muenst, Zsuzsanna Varga, Walter P Weber, et al. The single-cell pathology landscape of breast cancer. Nature, 378(7796):615-66)) 620, 2020.
Vladimir Yu Kiselev, Kristina Kirschner, Michael T Schaub, Tallulah Andrews, Andrew Yiu, Tamir Chandra, Kedar N Natarajan, Wolf Reik, Mauricio Barahona, Anthony R. Green, et al. Sc3: consensus clustering of single-cell RNA-seq data. Nature Methods, 14(5):483-4835 2017.
Joshua Knowles. Parego: A hybrid algorithu with on-line landscape approximation for expensive multiobjec- tive optimization problems. IEEE Transactions on Evolutionary Computation, 10(1):50 66, 2006.
Andreas Krause and Cheng Soon Ong. Contextual gaussian process bandit optimization. In John Shawe-Taye Taylor, Richard S. Zemel, Peter L. Bartlett, Fernando C. N. Pereira, and Kilian Q. Wein- berger (eds.), Advances in Neural Information Processing Systems 24: 25th Annual Conference on Neural Information Processing Systems 2011. Proceedings of a meeting held 12-14 December 2011, Granada, Spain, pp. 2447 2455, 2011. URL httpps:/ httpp://proceedings wwwrips.e..peeeee cc/paper//011/hahh 13f 35111ffc688777666611112222221110
Abstract...t..t........ html
Xuan Liu, Sara JC Gosline, Lance Pflieger, Pierre Wallet, Archana Iyer, Justin Guinney, Andrea H Bild, and Jefsrey Chang. stowledge-based classification fine-grained immune cell types in single-cell RNA-seq data. Briefings in Bioinformatics, 2bbb)bbbab 2021.
Malte Luccken and Fabian J Theis. Current best practices in single-cell RNA-soq analysis: a tutorial. Moleculor systems biology, 15(6):88746 2019.
Wesley J Maddox, Maximilian Balandat, Andrew G Wilson, and Eytan Bakshy. Bayesian optimization with high-dimensional outputs. Advances in Neural Information Processing Systems, 34:19274-19287, 2021.
Mary B Makarious, Hampton L Leonard, Dan Vitale, Hirotaka Iwaki, David Saffo, Lana Sargent, Anant Dadu, Eduardo Salmerón Castaño, John F Carter, Melina Maleknia, et al. GenoML: Automated machine lcarning for genomics. ArXiv preprint, htss/203..0221,, 2021. URL
03221
Davis J. McCarthy, Kieran R. Campbell, Aaron T. L. Lun, and Quin F. Willis. Scater: pre-processing, quality control, normalisation and visualisation of single-cell RNA-seq data in R. Bioinformatics, 33:1179-1188 2017. doi: 10.1093/ 10.1093/bioinffe bioinformatics/btv/// /tt777.
Andrew McHutchon. Diffcrentiating Gaussian Processes. Cambridge (ed.), 2013.
Jonas Mockus, Vytautas Tiesis, and Antanas Zilinskas. The application of Bayesian methods for seeking the extremum. Towards global optimization, 2(117-129):2, 1978.
Hirotaka Nakayama, Yeboon Yun, and Min Yoon. Sequential approzimate multiobjective optimization using computational intelligence. Springer Science & Business Media, 2009.
Mare Oostrum, Maik Müller, Fabian Klein, Roland Bruderer, Hui Zhang, Patrick GA Podrioli, Lukas Reiter, Panagiotis Tsapogas, Antonius Rolink, Bernd Wollscheid, et al. Classification of mouse cell types using surfaceome proteotype maps. Nature Communications, 10(1):1-9 2019.
16

---

Published in Transactions on Machine Learning Research (03/2023)
Biswajit Paria, Kirthevasan Kandasamy, and Barnabás Póczos. A flexible framework for multi-objective bayesian optimization using random scalarizations. In Amir Globerson and Ricardo Silva (eds.), Proceedings of the Thirty-Fifth Conference on Uncertainty in Artificial Intelligence, UAI 2019, Tel Avie, Isruel, July 22-25 2019, volume 115 of Proceedings of Machine Learning Research, pp. 766-776. AUAI Press, 2019. URL http://proceeednggs mmr.eingsmsrr....mmm www...rreettt.. https//vv1//pparriia
.tmll
Adam Paszke, Sam Gross, Francisco Massa, Adam Lerer, James Bradbury, Gregory Chanan, Trevor Killeen, Zeming Lin, Natalia Gimelshein, Luca Antiga, Alban Desmaison, Andreas Köpf, Edward Yang, Zachary DeVito, Martin Raison, Alykhan Tejani, Sasank Chilamkurthy, Benoit Steiner, Lu Fang, Junjic Bai, and Soumith Chintala. Pytorch: An imperative style, high-performance deep lcarning li- brary. In Hanna M. Wallach, Hugo Larochelle, Alina Beygelzimer, Florence d' d'Alché- Buc, Emily B. Fox, and Roman Garnett (eds.), Advances in Neural Information Processing Systems 32: Annual Confer- ence on Neural Information Processing Systems 2019, NewrIPS 2019, December 8-14, 2019, Vancou- ver, BC, Canada, pp. 8024-8035, 2019. URL https://pproceedingg www.ggs..urrs.... 0ebea888ee1119222f88877777777772 Awsoo--tbsttttttaa html
F. Pedregosa, G. Varoquaux, A. Gramfort, V. Michel, B. Thirion, O. Grisel, M. Blondel, P. Prettenhofer, R. Weiss, V. Dubourg, J. Vanderplas, A. Passos, D. Cournapeau, M. Brucher, M. Perrot, and E. Duchesnay. Scikit-learn: Machine learning in Python. Journal of Machine Learning Research, 12:2825 2830, 2011.
Ren Qi, Anjun Ma, Qin Ma, and Quan Zou. Clustering and classification methods for single-cell RNA- sequencing data. Briefings in Bioinformatics, 21(4):1196 1208, 2020.
Surajit Ray and Saumyadipta Pyne. A computational framework to emulate the human perspective in flow cytometric data analysis. PloS one, 7(5)633569,, 2012.
Diederik M Roijers, Peter Vamplew, Shimon Whiteson, and Richard Dazeley. A survey of multi-objective sequential decision-making. Journal of Artificial Intelligence Research, 48:67-113, 2013.
Bin Xin Ru, Ahsan S. Alvi, Vu Nguyen, Michael A. Osborne, and Stephen J. Roberts. Bayesian optimisation over multiple continuous and categorical inputs. In Proceedings of the 97th International Conference on Machine Learning, ICML 2020, 13-18 July 2020, Virtual Event, volume 119 of Proceedings of Machine PMLR, 2020. URL http://prroceedings mwr. prosss/111//r200
Learning Research, pp. html
Jasper Snock, Hugo Larochelle, and Ryan P. Adams. Practical bayesian optimization of machine learning algorithms. In Peter L. Bartlett, Fernando C. N. Pereira, Christopher J. C. Burges, Téon Bottou, and Kilian Q, Weinberger (eds.), Advances in Neural Information Processing Systems 25: 26th Annual Conference on Neurul Information Processing Systems 2012. Proceedings of a meeting held December 3-6, 2012, Lake Tahoe, Nevada, United States, pp. 2960 2968, 2012. URL httpps://proceedingg www.oediinnnnuuuuss... neurine..roomm cc/ hbw.:seccdd-bttrt..hh.....
Jasper Snoek, Oren Rippel, Kevin Swersky, Ryan Kiros, Nadathur Satish, Narayanan Sundaram, Md. Mostofa Ali Patwary, Prabhat, and Ryan P. Adams. Scalable bayesian optimization using deep neural networks. In Francis R. Bach and David M. Blei (eds.), Proceedings of the 32nd International Conference Machine Learning, ICML 2015, Lille, France, 6-11 July 2015, volume 37 of JMLR Workshop and Conference Proceedings, pp. 2171 2180. JMrg 2015. URL httt::://prcceedii
mlr ptt.s//33//v3/ http//333///oe11.../ snoek15..... html
Marlon Stoeckius, Christoph Hafemeister, William Stephenson, Brian Houek-Loomis, Pratip K Chattopadhyay, Harold Swerdlow, Rahul Satija, and Peter Smibert. Simultaneous epitope and transcriptome measurement in single cells. Nature Methods, 14(9):865 868, 2017.
William R Thompson. On the likelihood that one unknown probability exceeds another in view of the evidence of two samples. Biometrika, 25(3-4):285 294, 1933.
Ryan Turner, David Eriksson, Michael McCourt, Juha Killi, Eero Laaksonen, Zhen Xu, and Isabelle Guyon. Bayesian optimization is superior to random scarch for machine lcarning hyperparameter tuning: Analysis of the black-box optimization challenge 2020. ArXiv preprint, abs/2104.10201, 2021. URL https: /tttss///arxv.orrr//////0/0...1 wrw./b////0..11////222
17

---

Published in Transactions on Machine Learning Research (03/2023)
Johanna Wagner, Maria Anna Rapsomaniki, Stéphane Chevrier, Tobias Anzeneder, Claus Langwieder, August
Dykgers, Martin Rees, Annette Ramaswamy, Simone Muenst, Savas Deniz Soysal, et al. A single-cell atlas of the tumor and immune ecosystem of human breast cancer. Cell, 177(5):1330 1345, 2019.
Christopher K Williams and Carl Edward Rasmussen. Gaussian processes for machine learning, volume 2. MIT press Cambridge, MA, 2006.
F Alexander Wolf, Philiip Angerer, and Fabian J Theis. Scanpy: large-scale single-cell gene expression data analysis. Genome Biology, 19(1):1-5, 2018.
Kaifeng Yang, Michael Emmerich,, André Deutz, and Thomas Bäck. Efficient computation of expected hypervolume improvement using box decomposition algorithms. Journal of Global Optimization, 75(1): 3 34, 2019a.
Kevin Yang, Zachary Wu, and Frances H Arnold. Machine-learninng-uiiee directed evolution for protein engineering. Nature Methods, 16(8):687-69) 694, 2019b.
Shun H Yip, Pak Chung Sham, and Junwen Wang. Evaluation of tools for highly variable gene discovery from single-cell RNA-seq data. Briefings in Bioinformatics, 20(4):1583- 1589, 2019.
Luke Zappia, Belinda. Phipson, and Alicia Oshlack. Exploring the single-cell RNA-seq analysis landscape with the scRNA-tools database. PLoS Computational Biology, 14(66)41006245, 2018. doi: 10.1371 /..33a/journal.
Allen W Zhang, Ciara O'Flanagan, Elizabeth A Chavez, Jamie LP Lim, Nicholas Ceglia, Andrew McPherson, Matt Wiens, Pascale Walters, Tim Chan, Brittany Hewitson, et al. Probabilistic cell-type assignment of single-cell RNA-seq for lumor microenvironment profiling. Nature Melhods, 16(10):1007-100 1015, 2019.
Luisa M Zintgraf, Timon V Kanters, Diederik M Roijers, Frans Oliehoek, and Philipp Beau. Quality assessment of MORL algorithms: A utility-based approach. In Benelearn 2015: proceedings of the 24th anmual machine learning conference of Belgium and the Netherlands, 2015.
A Implementation details
A.1 Method hyperparameters
Hyperparameters for all methods are summarized in Supplementary Table 3
MANATEE MANATEE (with both AS and SA acquisition functions) was implemented with PyTorch v. 1.9.0 P Paszke et al. 2019) and gpytorch v. 1.6.0 Garrdner et al.) 2018) for the Gaussian process model and inference. Optimization was performed with the LBFGS optimizer. At every acquisition step, the model was initialized and fit to to the current training set 5 times and the model with the highest log-likelihood was kept. If fitting failed, the process would be re-tried a maximum of 20 times before halting. Optimization of the acquisition function was initialized with the maximum of 100 random samples. Above implementation details also apply the random scalarization (RS) bascline. Optimization of the SA and AS acquisition functions performed with the LBFGS optimizer with line search. For MANATEE, maxima of the posterior mean were identified by computing the first derivative of the posterior mean, finding its zeros using Brent's method Brent 2013) implemented by scipy- owt.mize.........mm brenttq with default parameters, and computing the second derivative at those locations. A candidate was declared a maximum not at boundary if its second derivative was less than -10 and if the candidate was at least 0.01 units away from the range extrema.
( (Balandat et al. 2020j. The development version was used to facilitale usage of qNEHVI with the ProneckerMultiTTaskeee model (Maddox et al. 2021). Implementation closely followed the tutorial on multi-objective Bayesian optimization (Balandat
NEEVV The NEHVI approach was implemented with botorch
al. (0021. Batch size was set to 1. rit_gpytorch_me model was called with max max_retries to 20. Other paraneters were set following the tutorial. Reference point was set to the minimum of the initial acquired points. To accommodate > 5 objectives, used approximate hypervolume
18

---

Published in Transactions on Machine Learning Research (03/2023)
Supplementary Table 3: Method hyperparameters used in the experiments. RS: random scalarization.
Method
Parameter
Value
Explanation
p(Ax= =)
0.5
Prior over binary A is set as Bernoulli(O55)
T1, To f
0.75, 0.25
Bernoulli hyperparameters for p A4 =i)
MANATEE
822
<
-10
Upper bound to call a max
Min distance
0.01
Distance from max to extrema
0.1
Kernel lengthscale constraint
o
1
Kernel variance
o Line search function bolfn
0.01
Observation noise variance constraint
LBFGs optimizer arg
UCB 3
0.125 log (2t + 1)
Set as in Paria et al. (0019)
MANATEE GP lils
5
Model inits and lits at each acquisition
RS
GP fit re-tries
Max re-tries to fit model each acquisition
Acquisition samples
100
Initial samples from acquisition function
MC_SAMPLES
128
QMC sampler arg, set as in Balandat et al. (0021)
max retries
20
botorch. fit_gpytorch __e___mdelll____ arg
batch range
[0,
QMC sampler arg, set as in Balandat et al. 2021
GP model
KroneckerHhultiTeskkPP
Set as in Balandat et al. (0021
NEEVV
BATCH_ SIZE
1
Points acquire
NaarEGO
NUM RESTARTS
20
Optimization re-starts, as
Balandat al.
2021
RAW SAMPLES
1024
Acquisition samples, set in
Balandat et al.
2021
batch_ Bachh_imitt
5
optimize_ _cqf arg, set as in
Balandat et al.
2021
maxiter
200
optimize_ace _acf arg, set as in
Balandat et al.
2021
reference point alpha
sample minimum
Lower HV bound
:NEHVI
10
approximate partitioning level
reference point
105
Set as in github. hatt.//////eaaaaa////See//
acquisition function
TS
Thompson sampling [Thompson) 1933
setch_sszze
1
Points to acquire
USeMO
d
1
Input dimensionality
12
beta
log
Set as in github.
0.15
algorithm
NSGA-II Deb et al. (2002) Set as in github. howpu//ooooo/ellaaaaaa//s//s//
function evalutions
2500
Set as in github. com//bellaaarria//SeS000
computation by setting the alpha parameter according to the heuristic based on the number of objectives as proposed in (Daulton et al. 2020).
The @NParEGO approach was implemented with botorch v. 3.6.I.ee337+44ffa2222 Balandat ct al. 2020). Implementation closely followed the tutorial (Balandat al.) 2021). Batch size was set 1. fit_ lrrt__pytorcc__cee model was called with max_retrii set to 20. Other parameters were following the tutorial.
@NParEGO
19

---

Published in Transactions on Machine Learning Research (03/2023)
USeMO USeMO was implemented using the code deposited in the tttppttp//gggttuu.ooo/ooootatttt repository, with main.py py adapted to work with pipelines for the IMC and scRNA-seq experiments. USeMO was run with the Thompson sampling (TS) acquisition function (Thompson 1933) as USeMO-TS and UScMO-EI (EI, expected improvement (Mockus et al.) 1978)) were shown to outperform cxisting methods Belakaria et al, 2020) and TS didn't. require hyperrarrmeeerr to set like the exploralion/exploiatiiin trade-off hyperparameter in the implementation of the EI acquisition function. All other hyperparameters were left as in main.py. py.
A.2 Experimental setup
Supplementary Table 4: Parameters of the experimental procedure.
Experiment Parameter
Value
Number of initial points Number of acquisitions
Toy
30
Replicates
100
Emin
Tmar
Number of objectives
Number of initial points Number of acquisitions 35
IMC
Replicates
98
Zmin
Emar
100
Number of objectives
Number of initial points 5
scRNA-seq
Number of acquisitions 36
Replicates
100
Tmin
0.01
Tmax Number of objectives
0.5
Parameters for all experimental procedures are summarized in Supplementary Table
Toy experiment The Pareto front on toy objectives was computed with the OApackage (Eendebak & Vazquez 2019). The initial dataset contained 5 training points at random locations and MANATEF-SA performed 30 acquisition steps. The experiment was repeated 100 times. The range the optimized paramcter & was between 0 and Toy objectives are described Appendix B.I.
IMC experiment At each acquisition step, data was normalized with the acquired cofactor value and clustered. Mean marker expression in each cluster was computed on the data normalized with cofactor=1. The co-expression objective values were computed a Pearson correlation between the mean expression of marker pair across clusters. The experiment was repeated 98 times. The range of the optimized cofactor value was between and 100. The overall meta-objective maximum for ARI and NMI, used to compute regrets, was set as the maximum of ARI/NMI values computed from all acquisitions by all methods (MANATEE-SA, MANATEE-AS, RA and RS baselines, (NEHVI, gNParEGO, USeMO), including behaviour ablation methods (MANATEE-SA with leave-one-out behaviour). qNEHVI returned an error and completed fewer than the total of 35 acquisitions on 49 runs; NNarEGO returned an error and completed fewer than the total of acquisitions on 80 runs, other methods completed all acquisitions on all runs. The objectives are described in Appendix F.2
20

---

Published in Transactions on Machine Learning Research (03/2023)
M-SA
RA
QNEHVI
M-AS
qNParEGO
USeMO
RS
60
50
Cofactor
40
20
10
0
5
25
30
Acquisition step
Supplementary Figure 4: Optimized cofactor value as a function of acquisition step for all methods. M-SA: MANATEE with scalarized acquisition, M-AS: MANATEE with acquisition of scalarized function, RA: random acquisition, RS: random scalarization. Solid line shows the and shaded region denoles the 95% confidence interval.
scRNA-see experiment At each acquisition step, data subsetted to the top highly variable genes according to the acquired proportion value and clustered. The co-expression objective values were computed as Pearson correlation between the mean expression of a marker pair across clusters. The experiment was repeated 100 times. The range of the optimized highly variable gene proportion was between 0.01 and 0.5. Unsupervised cluster purity metries were computed on the PCA transform of the normalized scR.AA-sq data. The overall meta-objective maximum for ARI and NMI, used to compute regrets, was as the maximum of ARI/NMI values computed from all acquisitions by all methods (MANATEE-SA, MANATEE-AS, RA and RS baselines, NNEVII, &NParEGO, USeMO), including behaviour ablation methods (MANATEE-SA with leave-one-out behaviour). ęNEHVI returned an error and completed fewer than the total of 36 acquisitions on 64 runs, &NOarEGO returned an crror and complcted fower than the total of 36 acquisitions on 62 runs, other methods completed all acquisitions on all runs. The objectives are described in Appendix F.3
Experimental results were tracked with Weights & Biases (Biewald 2020 Reported regrets include acquisi- tions from incomplete runs for şNEHVI and ParBBGO
21

---

Published in Transactions on Machine Learning Research (03/2023)
CD45_CD20_
CD45_ _D68_
CD45 CD3
Vimentin Fibronectn.
CD19_C020_+
panCytolenati.com _Kt
#[Cobetori
0.94
0.00
0.50
Re
#tColastor
2.50
2.56
#IColbutori
#iCoebbort
0.
0.05
0.7
2.04
8.2
0.92
0.80
50
100
10 Cotactor
100
50
100
53 Colactor
150
10 Coractor
103
13 Cotaccor
103
53 Colactor
150
Coractor
Cofactor
(a) MANATEE-SA
Farelo spoial
non-ostimal
NANATEE AS asquied
ECadherin ECoehhettnnnnCCnoleeet
CD45 CD20
CD45, CD45CD66
CD68
CD45 _03
Vimentin Aiborectin. +
CD19 0_,0 CD19_C020_
panCytobenatii CKT
0.20
0.90
MCodertor
2.00
0.44
2.96
2.04
822
0.92
0.92
100
300
130
100
53
100
Cotactor
Cotactor
Cofactor
Cofactor
(b) MANATEE-AS
Forats mpt
non-op0mmal
(NEHVI acquired
Caadherin un panCytokeratinn
CD455 Co0._
co4S_con__
Vimente bronectin... +
CD19 _c020_ CD19_cn20__
panCytobenttl CKS
2.98
riCoftuceer
REcdertigg
5:86
&
0.03
63
0.92
0.90
2.92
677
0.2
Coftctoo
Cofacbor
Cofactor
Cotador
Coflcttor
Cofactor
Cofactor
NEHVV
Poroto optimal
no--ptinaal
ghFor ECO acquined
panEyttberaaiin
Vimenti
panCytabenati_com
0.06
TcCrfhecor
0.90
(Crfecterr
5.40
TCCofeetod
RCCreeetr 0,3
0.85
5.44
0.80
0.2
1000
100
50 Cofactor
100
100
50
130
100
53 Cofactor
100
Cottector
Cottctoo
Coraooo
Coftecto
Cotactor
(d) ParEGG
Pareto outimal
nonapdmmal
USeNO
ECatherin in pantyytleeaaiin
CD45, CD220
CD4S cc...088..
CD4S CD3.
Vinwnir F bronectin
CD19 ca2ccccc..
panCytobareiin CKS,
0.56
0.00
(CCrfeettrd
0.54
0.83
0.7
-22
0.3
0.92
0.9
9:2
50
1000
50
200
50 Cornchor
205
30
100
50
100
50
100
32
100
UOeMO
0.1
6,32
6,30
0.28
8.33
0.26
100
50
103
(f) Meta-objectives
Supplementary Figure 5: 150 random samples from 7 co-expression objectives and meta-objectives (ARI, NMT) in the TMC cofactor selection experiment. Orange indicates samples on the Pareto front. blue indicates Pareto dominated samples, and red indicates points acquired by each method. The acquisitions shown for each method are from the run with the highest average ARI valuc. Pareto optimality of random samples was computed with the OApackage (Eendebak & Vazquez 2019).
22

---

Published in Transactions on Machine Learning Research (03/2023)
CCothent percyleeettt
CDR 0000,
CDO CEAR
CDRS
wererte Naranedr
COLR
powcyskenttecc
64
Actutatlin stop
Acsuls ton sap
Acqu Sll strslap
Acquiaiion AR2
Acquaalor step
Accui Mon 2002
Acqusstonn step
(a) MANATEE-SA
CCuhhelt jantysloratt,
CD4S, CD20,
COHS,
cos5,cto,
Vinentin, Rarenetti
CDER
pancyskorinn.coo
04
Arquisttin trop
Neuktiin.rg vop
Argidtenndo
Arcuitatoo daog
Arqucttin Sto
(b) MANATEE-AS
Supplementary Figure 6: Inclusion probabilities for MANATEE-SA and MANATEE-AS for each of the objectives function of acquisition step in the IMC cofactor selection experiment. Solid line shows the mean and shaded region denotes the 95% confidence interval across runs.
M-SA
RA
QNEHVI
M-AS
qNParEGO
USeMO
RS
0.35
0.30
% HVGs
0.25
0.20
10
15 Acquisition step
35
Supplementary Figure 7: Optimized percentage of highly variable genes as a function of acquisition step for methods. M-SA: MANATEE with scalarizod acquisition, M-AS: MANATEE with acquisition scalarizod function, RA: random acquisition, RS: random scalarization. Solid line shows the mean and shaded region denotes the 95% confidence interval.
A.3 Cross-validation experiments
In both experiments, subsampled and processed data (as described in Appendix B was divided in 70/30% train/test splits in 5-fold cross-valide valsdaaai,,,ion, where the parameter was optimized with MANATEE-SA only on train and cumulative regrets were computed in train and test. To ensure sufficient data, CITE-seq data was subsampled to 3000 cells in the cross-validation experiments. For the highly variable gene selection experiment, the same set of highly variable genes (computed according to the proportion value acquired train at each step) was used when computing regrets test. For both experiments, the initial datasct contained 5 training points at random locations and MANATEE-SA performed 10 acquisitions. Each cross-validation experiment was repeated 100 times. When computing regrets, the overall meta-objective maximum for ARI and NMI was set as the maximum of ARI/NMI values computed from the acquisitions by MANATEF-SA in cach run and each fold. For each run, cumulative regret was averaged over folds. Supplementary Table shows cross-validann validatiio cumulative regrets averaged over runs train and test for both experiments.
23

---

Published in Transactions on Machine Learning Research (03/2023)
optims
2.4
511
252
CH
CB
HUMAN CD3E, HUMAN CDA
HUWAN, CDSD HUWAN, COSA,
-25
11% H4000
0.2
(1% MANGO
-3.0
-75
9.6
0.0
100
-125
-10.0
0.4
-150
HUMAN, SYPRC_HUNAN.COOE
NUMAN_COD3 WWI.CCCCOOHAAA HUMALUNGALL 1
WWWAN_CDDDMMUOON CUNAL
MM HUMAN, M54AL
RS AVGaI
-12
08
To HVCs
% HVG
MANATEE-SA
2.4
51
HUMAN CDDE HUMAN CD4
HUWAN COSD HUMAN, COSA.
-25
71% 10/00
1.2
-5.
11% HV/
-75
2.0
-10.0
-125
-150
0.4
IVVG
HUMAN, WWW.CPPPPPMMMCC MUMAAAANNEE,L
HUMAN_ CWwAANCCO_HAACAA BWT.CCCOO_UAAAA HUMAUNS4AL
HUMAN_CC69_HUUUN NUMAN_M64AL
-12
-0.8
To HVCS
HVGG
(b) MANATEE-AS
2.4
CN
HUMAN CDDE BWW__CCD__HHMMAAA H.M
HUMAN, COOD, JUWAN CORAL
-25
11% 10000
0.2
-5.0
-7.
0.0
-10.
-125
HUMAA PTPRC JRIMAN [D10
HUMAN, co19 HUMAA NSSA3
HUMAN (DOD JAMAN CD6B
HUMAW cooa NUWAN MS-K1
18% HANS0
3.8
FIS HWOG
FIS HW060
2.7
2.6
% 1M
STNGS
(c) qNEHVI
Supplementary Figure 8: 150 random samples from 9 objectives and meta-objectives (ARI, NMI) in the % HVGs selection experiment. (cont.)
A.4 Behaviour ablation experiments
In leave-one-out behaviour experiments, MANATEE-SA was run without one of each behaviours the time. Supplementary Table T shows results for the IMC cofactor selection experiment, with each row indicating
24

---

Published in Transactions on Machine Learning Research (03/2023)
2.4
51
252
HUUMN CDOE, HUMAN CDA
HUWAN, COSD HUWAN, COSA,
-25
(1% H4000
0.2
(1% MANGO
-3.0
-75
0.0
100
-125
-10.0
0.4
-150
HUMAN, WWW._COD._HHMA 1
HUMAN_CO66 M54AL.
8%% AWS51
-06
-12
-08
to HVCs
HVG
(d) q NarFGGO
2.4
511
250
Cal
CB
HUMAN CDDE HUMAN CDA
HUWAN COSD HUWAN, COSA
-15
71% 10/00
1.2
-5.0
10 2000
-75
0.0
-12.0
-125
-150
HUMAN, WWW.CPPPPHMMMC JUNAAMAAAANN,L
HUMAN_ WWW.__C19__UPPP 1
HUNAN_CODD
HUMAN_CO6H_HUU MS4AL
2.8
To HVCs
HVG
USeMO
A11
104
(f) Meta-objectives
Supplementary Figure 8: 150 random samples from objectives and meta-objectives (ARI, NMI) in the % HVGs selection experiment. Orange indicates samples on the Parcto front, blue indicates Pareto dominated samples, and red indicates points acquired by each method. The acquisitions shown for each method are from the run with the highest average ARI value. Pareto optimality of random samples was computed with the OApackage [Eeddebak & Vazquez 2019
regrets for MANATEE-SA without said behaviour. Supplementary Table 7 similarly shows results for the scRNA-seq highly variable gene selection experiment.
A.5 Computing regret curves
Regret curves for cumulative regret (CR), full regret (FR), and Bayes regret (BR) (denoted
CR(I"),F((1)) BR(()) were computed function of the acquisition step 1' = ===,,.... 1,..., T as follows:
25

---

Published in Transactions on Machine Learning Research (03/2023)
HUMAN LONANLLEEENNNNCCLL
HUMAN__CO3, HWW.HHHAAAACO COBA,
0.8
026
10
30
KUMAN (WWP)PPPCCCCAA HOGGN HUWAN CD19 HUMAN M
HUMAN WWW.A_COOOHHUMMMMCCOMMM CDER
HU6A M54AL
Acqusiliion step
1.0
1163
1.0
22
30
10
20
32
30
26
30
0
13
20
32
Asquatttion step
Acquss lton step
Acquiattion stap
Acquarion step
(a) MANATEE-SA
1961
0.0
024
Acqukition stop
HUMAR PTRRC H HUMAN HUMAN,
HUMARLCB3O
HUNGALCD6S, HUMAN
0.8
22
20 Acquii Bonn step
32
Acquaition step
Acquuattion step
Acquastoon step
(b) MANATEE-AS
Supplementary Figure 9: Inclusion probabilities for MANATEF-SA and MANATEF-AS for each of the objectives function of acquisition step in the % HVGs selection experiment. Solid line shows the mean and shaded region denotes the 95% confidence interval across runs.
Supplementary Table 5: 5-fold cross-validation mean cumulative regret on train and test data splits. IMC denotes the IMC cofactor optimization experiment, scRNA-seq denotes the highly variable gene proportion optimization experiment. ARI: adjusted Rand index, NMI: normalized mutual information.
IMC ARI train
0.015 + 0.003
IMC ARI test
0.017 0.003
IMC NM1 train
0.021 0.005
IMC NMI test
0.023 = 0.005
scRNA-seq ART train
0.095 + 0.015
scRNA-seg ARI test
0.089 + 0.015
sc scRNA-seq NMI train
0.126 + 0.021
scR.AA-ee NMI test
0.126 = 0.020
CR((') 1 J Σ (0' M(x))). V' 1,.,,T
FR(t') = 3²* max h(x), W 1,.,T x x
BR(t)= = E
(v'
max h(r)), Vt' = I...T T
x x
1=1
26

---

Published in Transactions on Machine Learning Research (03/2023)
Supplementary Table 6: Results for the behaviour ablation experiments for IMC cofactor optimization. CR: cumulative regret, FR: full regret; BR: Bayes regret. ARI: adjusted Rand index, NMI: normalized mutual information. Values are mean (s.d.).
ARI
NM1
Ablated
CR
FR
BR
CR
FR
BR
Explainability
0.007(0.004) 0.016(0.004) 0.002(0.003)
0.017(0.006) (.002(0.005) 0.007(0.007)
Inter-obj agreement 0.018(0.005) 0.003(0.005) 0.007(0.006) 0.020(0.009) 0.003(0.008) 0.008(0.010)
Max not at boundary 0.018(0.007) 0.004(0.006) 0.008(0.007) 0.020(0.012) 0.004(0.011) 0.009(0.012)
Supplementary Table 7: Results for the behaviour ablation experiments for scRNA-seq highly variable gene sclection. CR: cumulative regret, FR: full regret; BR: Bayes regret. ARI: adjusted Rand index, NMI: normalized mutual information. Values are mean (s.d.).
ARI
NMI
Ablated
CR
FR
BR
CR
FR
BR
Explainability
0.126(0.023) 0.055(0.012) (.061() 0.121(0.032) (.021( 0.030(0.012)
Inter-obj agreement (0.188) (.055( (.061)) 0.126().025) 0.020(0.010) 0.030(0.014)
Max not at boundary 0.120(0.022) 0.055(0.012) 0.063(0.013) 0.113(0.033) 0.020(0.010) 0.030(0.013)
Note that these cumulattie quantitics up to step 1' for t' = =......... . T are normalized by the number of acquisition steps so far, leading to most curves decreasing with acquisition step. However, these quantities are not guarantced decrease with acquisition step for all methods as they are computed w.r.t. the meta-objectives (ARI, NMI), while all methods make acquisitions by having access only to heuristic objectives, which do not necessarily agree with the meta-objectives in our setup. The heuristic objectives and meta-objectives are also noisy functions owing to being computed using real noisy genomics data.
Regret curves for all methods for the IMC colactor selection experiment are shown in Supplementary Figure 10 and those for the % HVGs selection experiment are shown in Supplementary Figure 10 b.
Cumulative regret
Full regret
Bayes regret
Cumulative regret
Full regret
Bayos regret
0.250
0.25
0.25
0.04
0.04
0.225
0.20
0.20
0.03
0.00
0.03
0.375
0.200
0.15
0.15
0.159
0.10
0.02
0.00
0.01
0.325
0.05
0.07
0.05
0.5
0.25
0.20
0.05
0.04
0.04
INN
020
0.15
IWN
0.03
0.35
0.10
0.02
0.00
0.5
20
20
20
20
20
Acquisiion step
Acquisition step
Accuisiton step
Acquistion step
Acquisition step
Acquistion Step
(a) IMC colactor selection
(b) % HVGs selection
Supplementary Figure 10: Normalized cumulative CR, BR, FR for ART and NMT meta-objectives for all methods in both experiments. Solid line shows the mean and shaded region denotes the 95% confidence interval. CR: cumulative regret, FR: full regret; BR: Bayes regret. M-SA: MANATEE with scalarized acquisition, M-AS: MANATEE with acquisition of scalarized function, RA: random acquisition, RS: random scalarization. ARI: adjusted Rand index, NMI: normalized mutual information.
27

---

Published in Transactions on Machine Learning Research (03/2023)
M-SA
Method
Method
M-AS
M-AS
80
100
0.0
0.1
0.2
0.3
0.4
0.5
Cofactor
% HVGs
(a) Cofactor values selected in runs of the IMC cofac- tor selection experiment
(b) Percentage values selected in runs of the % HVGs selection experiment
Supplementary Figure 11: Boxplots showing modes of the acquisitions made in each run by two MANATEE variants for a) IMC cofactor selection experiment and b) % HVGs selection experiment. The mode of the acquisitions made in a run represents the parameter value selected by a method in that run. Modes were computed on the acquisition values rounded decimal places. Each run (repetition) of an experiment had different random seed controlling the initial training and initialization of the Gaussian process model. M-SA: MANATEE-SA, M-AS: MANATEE-AS.
A.6 Computing best achieved hypervolume
For each experiment, the hypervolume was computed with the compute compute_hypers hypervolume ( function of the DominatedPartitionnig BoTTorh class. At each acquisition step, we computed the dominated hypervolume of the dataset acquired so far w.r.t. the reference point. The acquired dataset included the initial points and the acquisitions made up to the current step. All acquisitions ions were scaled to [0, 1] to account for different. ranges of the objectives. To scale the Calinski and Harabasz score and the Bouiden score, we computed score) as the
their upper and lower bounds (note that we define our objective as the negative maximum and minimum acquired value across all methods in all runs, respectively. We defined the lower bound of the Calinski and Harabasz score and the upper bound of the negative Davies-Bouldin For all iesetttooobbaedd objectives and the silhouette width, defined the lower and upper bounds as -1 and 1, respectively. For each experiment, we defined the reference point as the lower bound of each objective. Specifically, for the correlation-based objectives and the silhouette width the lower bound was -1, for the Calinski and Harabasz score the lower bound was 0, and for the Davics-Bouldia score the lower bound was set as the minimum of all acquisitions.
O.
The curves showing the best achieved hypervolume at each iteration for all methods in both experiments are shown in Supplementary Figure 12 The best achieved hypervolume at each iteration was computed the maximum hypervolume at that iteration across all runs. For the cofactor selection experiment, the maxima were computed across 98 runs for MANATEE-S3 RS, RA, USeMO; 96 runs for MANATEE-AR AS; 49 succossfully terminated runs for (NEHVI; and 17 succossfully terminated runs for ęNParEGO. For the highly variable gene selection experiment, the maxima were computed across 100 runs for RA, USeMO; 99 runs for MANATEF-SA, MANATEF-AS, RS; 36 successfully terminated runs for NEHVT; and 37 successfully terminated runs for PNO.EGGO
28

---

Published in Transactions on Machine Learning Research (03/2023)
M-5A
RA
QNEHVI
M-SA
RA
qNEHVI
M-AS
qNParEGO
USeMO
M-AS
qNParEGO
USeMO
RS
RS
78.0
1000
77.5
Hypervolumt volume
77.0
Hypervolume
950
76.5
76.0
850
75.5
75.0
800
74.5
750
74.0
15
25
15
30
Acquisition step
Acquisition step
(a) Best achieved hypervolume in the IMC cofactor
(b) Best achieved hypervolume in the % HVGs selec- experiment
selection experiment
Supplementary Figure 12: Best achieved hypervolume at each acquisition step for a) IMC cofactor selection experiment and b) % HVGs selection experiment. Dominated hypervolume of the acquired dataset far computed for each acquisition step and the maximum selected across all runs. Each run (repetition) of an experiment had a different random seed controlling the initial training set and initialization of the Gaussian process model. M-SA: MANATEE-SA M-AS: MANATEE AS, RA: random acquisition, RS: random scalarization.
B
Data processing
B.1 Toy data
The 5 objectives in the loy experiment had the following functional forms defined on x [0, 1]:
y1 (ar) = max(0. sin 2 + 6
(2(x) sin (x 0.05) +, (O,, 0.3²)
3. Va(x) = sin 2πx + ((,
y1(2) = 2x + 6 N(0,00.8))8)
y5(2) 2x e, (0.00.1
B.2
IMC data
IMC data and expert annotated ground truth clustering used in the experiments come from JJackson al. (2020). Data was randomly subsampled to 5000 cells and the heavy metal markers that weren't 't conjugated to antibodies were removed. Data clustered with the scikit-learn (Pedregosa et al. 2011 k-mcans algorithm with k 10.
B.3 CITE-seq data
CITE-scq data used in these experiments come from Stocckius al. (2017) retrieved using the Singlecel1MultiMoaaal Bioconductor R package with data version 1.0.0. Cell surface antibody expres- sion normalized using the logNormCounts from the scuttle package (McCarthy al. 2017] and
29

---

Published in Transactions on Machine Learning Research (03/2023)
clustered using Seurat 4.1.0 (Hao et al. 2021 with top 10 principal components input and reso- lution parameter set to 0.8. Intra-cellular single-cell RNA-seq data was filtered for genes with at least. 100 reads and further processed with scanpy (Wolf et al. 2018). Data was randomly subsampled to 1000 cells (except for cross-validation experiment, where it was subsampled to 3000 cells) and normalized using pp- www.moooeeeea..eee tttalttellllllleeeoaaaleiilellla with target_ sum=1e4, pP- pp. logip, and pp- pp. scale. scanpy was used select. highly
variable genes, compute the neighbourhood graph with 10 neighbours and top 40 principal components, compute the PCA decomposition with default arguments, and compute Leiden clustering with resolution parameter set to 0.8.
C Derivations of acquisition functions
C.1 Expectation of scalarization of the single-objective acquisition function of objectives (SA)
We define the SA acquisition function ()((1)) (sx(aq9ccn((())) and derive the following expression:
(AAS (a(u(()) + (Bo( Bp(A1B) Σ Au(a(())+ +VBO(()))
EP((1) (ux(() Box(x))
K
Σ p(Ax = ([Ba)(((()) + (zx
For the SA acquisition function, the optimized expression simplifies to the sum over K objectives of each objective's acqUCB (fk (x)) function value weighted by the probability of that objective being useful conditioned its behaviours p(),- 1[Bx). While the SA formulation only takes into account the posterior variance of each fR, &o = 1.... K, fitting a multi-output Gaussian process to data is still desirable as it ensures that the posterior form of f reflects our assumptions about the heuristic objectives, namely, that they should be maximized by similar values of x.
C.2 Expectation of single-objective acquisition function of the scalarized objectives (AS)
We define the AS acquisition function as (A())) [acquen(s>)) (sx(f(x)))) and wish to derive the following expression:
E EA(B
acQUCB
Σ
Aafllco))
First, notice that,
(Nffi(s), (
N X X"AE(=)),
where ((x) is the posterior mean and E(x) is the posterior covariance of f evaluated at some x. Then, their sum distributed as:
Σ
Aama(1)) Σ )aja(t))> ~ N Σ
(F +2 Σ AeAZ(((()
&
ISKCJKK
Now, we use this to derive:
30

---

Published in Transactions on Machine Learning Research (03/2023)
Σ Aatalaa) ]
k
(p((1)) acqUCB
B)
Σ (x VB.
Xx) +2
E
V
ISijjSK
Σ p(h 1Bx))u(()) + VB. B)AB))
()) +2
Σ (1))
ISiJSSK
For the AS acquisition function, the UCB is applicd to the weighted sum of fk. h = l,..., expression contains the sum the posterior means weighted by p(), = 1BBa) as before, but the variance term takes into account the posterior covariance of
J....,K, K.
The optimizod
The proposed acquisition functions can differ c.g. within regions of search space, where those objectives that have negalive posterior covariance with others but still have a large weight due to other behaviours, have large posterior variance. In these cases, the SA formulation will assign larger value to those regions than the AS formulation.
C.2.1 Computing the variance term
We compute the expectation of the variance term w.r.t. p(\ by exhaustively considering all K-dimensional binary vectors X", I,..., 2 which is feasible our experiments with and objectives. Specifically, we evaluate the expectation as:
(P(
x +2
Σ ((g(
ISKUSK
Σ
B(1BB)
Σ
(XE) (xn(x) +2 +2
Σ N1AEZ1((0)
V
x
ISKJJKK
Σ X
Σ (AP) (x) +2
Σ X))))(())
ISKGSK
In problems with large K where this approach becomes infeasible, the expectation can be approximated with S Monte Carlo samples of Xi ~ p(>n [B)) V and I,..., as follows:
E PA(15)
A Ex) +2 Σ Ax))En(())
ISiCjK
(A2)² E +2
Σ {11)Z1(())
ISKJEK
Derivatives of a multi-output Gaussian process
The posterior mean of single-output Gaussian process with N noisy observations y at a new location x, is given by (Williams & Rasmusssn 2006):
((x,)) (K(X,X) + o2IN)-1
31

---

Published in Transactions on Machine Learning Research (03/2023)
D.1 First derivative
For the exponentiated quadratic kernel function, the first derivative of the posterior mean (which is also the mean of the distribution over derivatives of the GP posterior functions) is (McHutchon 2013
8f. ór.
8K(I,, X)
a -- (K(I,, @a)
dz
(K(X, X) + o,.
X
[r, "1,,.., ]N 1T
For D-dimensional input x; diagonal matrix A collects lengthscales l3. 1...... D on the diagonal, X an N D matrix and O represents element-wise multiplication. The resulting derivative of. with each element. corresponding to the d-th input dimension.
or. is D-dimensional,
D.2 Second derivative
The second derivative of the posterior mean is given by
f- 8K(,,, X)
8(x,)²
Sx
Sx,
The second derivative of k(z,, xx)) is D x D matrix given by (Mclluccoon 2013
82 (-A-1 + A (x. 21)(xx 2:)(x=-11)) A ') x(2,,
(x,
6(x,)²
We stack these to compute a X X matrix (K(z.. X) and multiply by a to compute the second o(x.)² ya
derivative of the posterior mean. The resulting derivative 8 is a D x D matrix with (i,j)-th element f
8[2,}
(the second derivative w.r.t. dimensions and j x.). corresponding to ((x.):
D.3 One-dimensional input case
For one-dimensional input and lengthscale 1, the second derivative of the posterior mean w.r.t. I, simplifies to:
f²f.
(- K(nn X" o X P'ToOK(r.... K(x.,))
8x?
D.4 Multi-output Gaussian process
For a multi-output Gaussian process with M objectives, we arrange observations as MN x array: = [y1,...., NN1 ,,. WN We also augment the auxiliary matrix X to become an MN matrix: X [r. .., EN .. II...., x. =N
The multi-output. kernel is defined Kmulti K10 & K and the additive noise term D IN, where D is the diagonal matrix with task-sp observation noises.
The first derivative of the posterior mean for multi-output Gaussian process is D-dimensional for each of the M tasks. The second derivative returns a x D matrix for each of the M tasks. In the computation of the second derivative, the last two dimensions of 8 (guulii (x..)) are flattened similarly to y before multiplication by a.
8(2.))
We note that the desired gradients can also be accessed via the automatic differentiation engine instead of using the analytically derived expressions we have done here.
32

---

Published in Transactions on Machine Learning Research (03/2023)
E Proof of theorem
Theorem 3.1 Tf (AA(A > 0 vk, then x² = arg max, B) lies the Pareto front of f.
Proof: Let dR E {0,1} be binary variable for k ,,..,, with da p(>x (Ba) given by Bernoulli distribution parameterized by function of B₁. Let us denote its expectation 인h = (Ax))B [A > 0 Ak
and denote À = [AI,..., Ax]. We further assume Ag is independent from d conditioned on its behaviours B₃ vk, k': k + k'. We also denote the set of all behaviours = [Bx] Then, for a linear scalarization function 8x and by linearily of expectation:
Emnsacoo)
E BAA ([(x)) B)(AB)
fn(x)) EMAAss(AN
E(((x) B dR BI.....x))) [dc]
Σ>e AR
Σ
p(A Ax BI,...... Bxx)
(A......... K)W
(8)
Σ A P(A [Bb)
AnE(O.1)
IM( (1 (L(s(A=1((B)) +0 p(A, = O(B)))
IMOAA 1I
As the 0 Vk, Ex fx())o is monotonically increasing in all fx(x.. If any scalarization function 8사 is monotonically increasing in all coordinates fz( x² arg max, sx (f(x)) lies on the Pareto front of f for f K Roijers et al. 2013 Zintgraf et al. 2015 Paria et al. 2019). Thus, the solution to max, 15 s (f()) max, Σx fA (z) lies on the Parcto front of f. In other words, maximizing expected 8x (f()) under (A)B) with our construction returns Pareto optimal solutions.
F Cluster mean co-expression as a heuristic
F.1
Overview
The cluster mean co-expression heuristic is demonstrated in Figure 13 After clustering the single-cell data, consider the expression of two proteins which should be markers for given cell type, i.c. they should either both be co-expressed or not expressed. An example of this is shown in Figure 13 (right). Consequently, the correlation in the cluster means is high. Conversely, if the clustering does not capture the cell types well, the correlation in the cluster means will be low (Figure 13 left). The opposite logic applies two proteins should be mutually exclusively expressed: the correlation of cluster mean expression should be minimized.
F.2
IMC experiment
The protein pairs used to construct co-cxpression heuristic objectives are listed in Supplementary Table
F.3
scRNA-seq experiment
The gene pairs used to construct co-expression heuristic objectives listed in Supplementary Table
33

---

Published in Transactions on Machine Learning Research (03/2023)
Cellular populations mixed
Cellular populations separated
Cluster
Expression protein 2
0.8
1
:
2
0.4
3
4
0.0
0.0
0.4
0.8
0.0
0.4
0.8
Expression protein
Supplementary Figure 13: Mean expression of two example proteins after two sets of clustering. Left: the cluster means poorly separate into double-positive and double-negatin ive populations as would be expected if the two proteins are markers for the same cell type. Left: the ideal situation, where clusters only co-express both proteins simultaneously or not at all.
Supplementary Table 8: Co-expression of marker pairs used as objectives in the IMC cofactor selection experiment.
Protein pair
Co-expression direction Cell type
E-Cadherin, pan- pan-Cytokeratin
Epithelial
CD45, CD20 CD45, CD68 CD45, CD3 Vimentin, Fibronectin
B cell
Myeloid
T cell
Stromal cell
CD19, CD20 pan- Man-Cytokeratin, Cytokeratin
B cell
Basal epithelial
Supplementary Table 9: Co-expression of marker pairs used as objectives the scRNA-seq highly variable gene selection experiment.
Gene pair CD3E, CD4 CD3D, CD8A PTPRC, CD68
Co-expression direction Cell type(s)
Regulatory T cell
Cytotoxic T cell
Myeloid
CD19, MS4A1
+
cell
CD3D, CD68
T/myeloid
CD68, MS4A1
Myeloid/B
34


Reviews:
[{"review_id": "XOLsxpwojF", "reviewer": "Reviewer_ZVWs", "summary": "", "questions": "", "limitations": "", "rating": null, "confidence": null, "soundness": null, "presentation": null, "contribution": null, "strengths": "", "weaknesses": "", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "Full disclosure, I have reviewed this paper before and find that most of both my requested changes and that of the other reviewers have been made. Good work overall.\n\n*Strengths:*\n\n- Maximizing for predictive confidence (high explainability) and not moving towards the boundary in multi-objective problems is certainly a good idea.\n\n- Overall, I think I now better understand the subtle idea of a meta-objective that we’re slowly moving the optimization towards. This seems generally applicable across many areas of biological design campaigns, especially when our proxy tasks are defective in some way.\n\n- Experimental tasks are realistic and seem grounded in biological motivation.\n\n- Botorch code is provided in a common framework.\n\n- Thank you for providing the algorithm box.\n\n*Weaknesses:*\n\n- Only single dimensional problems are considered in the experiments. All things considered this isn’t that bad of an issue.\n\n- The regret curves are only somewhat demonstrating superiority over observations (e.g. over the course of a design campaign), which is probably just most reflective of the noisiness of the tasks.\n\n- This is somewhat of a fundamental weakness that can probably be addressed in the discussion: Inter-objective agreement seems like it might limit the practitioner’s need to make tradeoffs between competing objectives in favor of jointly suboptimal but reasonable performing objective values.\n", "was_revised": false, "final_justification": "", "rebuttal": "", "claims_and_evidence": "Yes", "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": "", "paper_topic_and_main_contributions": "", "excitement": null, "reproducibility": null, "ethical_concerns": "", "main_review": "", "summary_of_the_paper": "", "summary_of_the_review": "", "strength_and_weaknesses": "", "limitations_and_societal_impact": "", "review_text": "", "summary_of_contributions": "This paper proposes MANATEE, an acquisition function for multi-objective Bayesian optimization, that infers positive qualities of different objective functions in terms of practitioner preferences, rather than being solely reliant on the hypervolume or other scalarizations of the several objective functions. In the experiments here, MANATEE focuses in on points with low predictive noise (what they call explainability), correlation between objectives, and non-boundary points. \n\nThe authors perform experiments with synthetic data and two biology tasks – mass cytometry cofactor detection and single cell gene selection. The method performs reasonably well across the objectives they measure across both cumulative regret (Fig. 10) and the objectives they were looking for (Tables 1 and 2). The biology tasks seem well connected to real problems despite being one dimensional and are designed around the concept of a meta-objectives, which softly measure optimization goals.\n", "claims_and_evidence_yn": "Yes", "claims_explanation": "", "audience_yn": "Yes", "audience_explanation": "", "requested_changes": "Requested changes:\n\n(medium) Can you doublecheck the cumulative regret plots in Fig. 10b? I’m not entirely sure they are going to be leveling off, but should continue upwards, or at least asymptoting?\n\n(minor) Fig 10a: label these as cumulative regret, full regret, bayes regret.\n\n(minor) Fig 5,8: maybe only display the pareto front or show the other points with a lower alpha value\n\n(minor) fig 7: these should probably be all on the same plot to maximize differences\n\n(medium) I would also suggest that the authors plot best achived hypervolume across iterations for the tasks as that seems to be a common concern.\n\nAppendix d: I still think that this can be all achieved with autograd but it’s really a minor point.\n", "broader_impact_concerns": "n/a", "additional_comments": "", "extra_scores": {}, "extra_text": {}}, {"review_id": "EurF13oeZ3", "reviewer": "Reviewer_QD2q", "summary": "", "questions": "", "limitations": "", "rating": null, "confidence": null, "soundness": null, "presentation": null, "contribution": null, "strengths": "", "weaknesses": "", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "N/A", "was_revised": false, "final_justification": "", "rebuttal": "", "claims_and_evidence": "Yes", "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": "", "paper_topic_and_main_contributions": "", "excitement": null, "reproducibility": null, "ethical_concerns": "", "main_review": "", "summary_of_the_paper": "", "summary_of_the_review": "", "strength_and_weaknesses": "", "limitations_and_societal_impact": "", "review_text": "", "summary_of_contributions": "I have reviewed the previous version of this submission.\n\nThe current paper has addressed the outstanding concerns from last reviewing round satisfactorily. I recommend accepting the paper.", "claims_and_evidence_yn": "Yes", "claims_explanation": "", "audience_yn": "Yes", "audience_explanation": "", "requested_changes": "N/A", "broader_impact_concerns": "N/A", "additional_comments": "", "extra_scores": {}, "extra_text": {}}]

Author rebuttal:

Meta-review: This resubmission was only reviewed by two reviewers because we were unable to reach one of the original reviewers.

Both reviewers find that the resubmission addresses all reviewer comments raised for the original submission and is publishable.
Decision: Accept (poster)