貝葉斯最佳化與自動化實驗設計

本文作者雷博文，主要基於在2021年發表於 npj Computational Materials 上的文章 Bayesian optimization with adaptive surrogate models for automated experimental design ［10］。

引言

無論是在日常生活之中，還是在學校課堂上，或者在學術研究之中，實驗一直以來都幫助著人們理解各種形形色色現象之間的因果關係。通常來說，實驗就是在一定條件下，透過控制變數產生不同的結果，來驗證某種假設。比如物理可以透過實驗可以探索物體運動的規律，材料學則可以分析某個化合物的性質。但早期的研究更注重透過觀察來探究自然規律，比如古希臘科學家亞里士多德就透過大量的觀察來進行研究，屬於論證科學。再比如，阿基米德，他成功地結合了邏輯推理和觀察實驗這兩種方法，但這只是為了檢驗邏輯推論 ［18］。

隨著時間的推移，人們逐漸開始不滿足於已有的研究模式，於是培根等人開始推崇實驗科學，強調“從實驗中獲得知識”。之後，伽利略也呼籲“科學的真理不應該在古代聖人的蒙著灰塵的書上去找，而應該在實驗中和以實驗為基礎的理論中去找”。他初步開創了“觀察——假設——推理——實驗——得到結論”的研究方法，促進了科學研究的發展。如今，做實驗已經成為無數物理，生物，化學，材料等學科學者經常要進行的工作。

Humen-centric vs Data-centric

但隨著時代的發展，科技的進步，人們嘗試解決的問題越來越困難，維度越來越高。許多時候，每次實驗可能都對應著大量的時間和金錢的消耗，傳統的最直接的列舉法，及嘗試所有的可能情況，開始變得不切實際。所以人們開始對實驗有了更高的要求，實驗設計開始變得越來越重要。比如如果能透過更少的實驗達到同樣的效果，就能節約許多時間和金錢，可能一些之前預算內無法完成的研究變得可能。拿今天我們主要要聊的材料學來說，一個常見的問題是“在眾多可能的材料結構中找到性質最優的那個”。如圖 1 所示，我們想要在

和

中找到具有最優性質的結構。其中

，

和

分別可以在紅色，藍色和黑色區域選擇一種元素，基於材料學知識排除一些不穩定的組合後，還剩下403種可行的結構。於是問題就變成了如何在這403種材料裡找到我們需要的那個。通常情況下，預算是無法支援我們嘗試完全部的403種可能的結構。

除了列舉法，另一種方式是人們根據自己的經驗 （human-centric） 來安排進行實驗的化合物結構的先後順序。這種方式優於簡單粗暴的列舉法，尤其是對於經驗豐富的學者來說。但這可能對進行實驗的人會有較高的要求，同時當面對一些新的問題，以往的經驗也不一定奏效，於是可能還是需要不少的嘗試才能找到目標的結構。

不管是列舉法還是根據經驗的實驗設計，都有一個問題，就是沒有充分利用已有的資訊。假設我們已經進行了

次實驗，並且記錄好了相關的資料。那麼當我們準備進行

次實驗時，前

次的結果是可以對我們這次的實驗選擇進行指導的。比如從直覺上，如果我們發現了某種結構的化合物性質比較好，那麼和它類似的化合物可能這種性質也會比較好。

這時候問題就變成了如何才能更充分利用已知資訊。隨著統計學、機器學習演算法的發展，在對複雜函式進行建模的問題中，計算機相對於人腦的優勢越來越明顯。人們已經可以利用計算機和設計好的模型演算法，更好地從資料中發掘出有價值的資訊。“Taking human out of the loop“ ［14］ 和 ”data-centric“ 的呼聲越來越高，材料發現問題中的自動化實驗設計的概念 ［11］ 也迅速成為一個活躍的研究領域 ［1， 7， 9］。貝葉斯最佳化 （Bayesian optimization，BO） 就是其中一種讓資料來指揮實驗設計（自動化實驗設計）的方法。

淺談貝葉斯最佳化

這裡先介紹一下貝葉斯最佳化 （BO），這裡舉一個簡單的例子，比如我們想要找未知函式

的最大值，於是我們基於已經觀測了的函式值對

進行建模，然後根據訓練的模型找到最可能是最大值點的

然後進行取樣，接著將得到的新的樣本點加入已觀測資料之中。不斷重複這個過程，直到迴圈終止條件得到滿足。如圖 2 所示，藍色曲線是真實的未知函式，橘黃色曲線和區域是根據模型擬合出來的曲線及其置信區間，黑色圓點代表已觀測的資料，紅色曲線對應計算的取樣函式 expected improvement，而灰色三角形則為下一個取樣點。可以看到，每一次灰色三角形都在紅色曲線的最高點，然後將新的資料加入資料集中並更新橘黃色曲線和區域，然後得到新的紅色曲線，然後不斷迴圈這個過程。

BO 是一個常用的最佳化目標函式的流程。當我們想基於最少的資訊，來知道某個函式的最大值點（或者最小值點），BO 是一個很好的選擇。近年來，BO 在實驗設計和超引數的除錯中受到了極大的歡迎。具體來說，BO 一般包括以下步驟：

從演算法中我們可以發現，使用 BO 進行自動化實驗設計，有兩個非常重要的點需要提前設定：（1）基於資料進行訓練的機率模型；（2）決定下一個取樣點的採集函式。機率模型是用來建立自變數（實驗前可知的特徵）和因變數（感興趣的材料性質）的關係。而採集函式則是我們選擇下一個樣本的標準，代表了我們如何定義出現最大值的可能性。

關於採集函式，常用的包括但不限於 expected improvement，probability of improvement， upper confidence bound，以及 Thompson sampling。使用者可以基於自己的需求選擇其中之一。這裡我們以 expected improvement （EI） 為例進行後續分析。EI 會從期望的角度選擇對

函式值提高最多的樣本點，其公式如下：

其中，

是已觀測的函式值中的最大值，

是對基於觀測資料得到的後驗分佈取期望，

。

跳出高斯過程，也許別有洞天

前面已經提到，除了採集函式，BO 另一個關鍵點是其使用的機率模型。假設我們有自變數

以及因變數

，我們想要找到兩者之間的關係，並且當新的

到來時能較好地預測其對應的

。通常我們會假設因變數是關於自變數的函式再加上一個額外的噪音項

。通常一個非常常用的選擇是高斯過程迴歸 （Gaussian Process Regression， GPR），不管是 python 還是 R 中都有成熟的包，所以使用起來比較方便，同時效果也不錯。

並非萬能的高斯過程迴歸

在 GPR 中，我們會給未知函式

加一個高斯過程的先驗（

是均值函式，

是所選的核函式）：

當新的輸入

，其預測後驗分佈為：

常用的均值函式是常數函式，而對於核函式，我們有很多選擇，包括但不限於：

Radial-basis function （RBF） kernel：可用於擬合平穩和各向同性的特徵。

RBF kernels with automatic relevance detection~（ARD）：在 RBF kernel 的基礎上為每個特徵分配不同的尺度引數，從而進行變數選擇。

Dot-product kernels：可用於捕捉非平穩的特徵。

Deep kernels：是一種更加靈活的核函式，但是引數較多，一般需要更多資料進行訓練來達到較好的效果。

雖然 GPR 非常流行並且在許多場景下有不錯的效果，但其並不是萬能的，也有一些缺點。比如平穩以及各向同性的 GPR 面對高維問題時表現會下降，尤其當許多不重要的自變數混在了真實變數之中但建模者並不知道哪些是關鍵變數的時候。可是這種情況在材料學中卻十分常見。所以如何在建模的時候，同時進行變數選擇就變得十分重要了。近年來有一些基於 GP 的模型在嘗試變數選擇，比如前面提到過的 automatic relevance detection ［2］。還有基於混合高斯模型和貝葉斯模型平均的方法 ［15］。但在所有的這些基於 GP 的模型中，他們使用的核函式通常給擬合的函式引入了平滑性和連續性的假設。但在許多材料學的問題中，真實函式可能並不是平滑的或者函式值會有突然的跳躍。除了上述的研究，還有許多關於靈活的非平穩協方差核函式的文獻 ［12］，比如 Deep kernels ［17］，但由於引數比較多，其在資料量較少的情況下不一定能獲得好的表現，而在許多材料學研究中，昂貴的實驗開銷導致資料量不會太大。所以也許跳出並非萬能的 GPR，會有不一樣的風景。

貝葉斯多元自適應迴歸樣條

貝葉斯多元自適應迴歸樣條 （Bayesian multivariate adaptive regression splines， BMARS） ［5］ 就是一個不錯的替代 GPR 的選擇。BMARS 是一個非常靈活的非引數迴歸模型，透過產生乘積樣條基函式對未知函式進行建模並自動識別變數之間的非線性相互作用。模型的具體形式如下：

表示基函式

對應的係數，而

的構造如下：

其中，

，

代表變數的索引，對於每個基函式，變數集

是不重複的。

其模型訓練和引數估計主要基於 reversible jump Metropolis-Hastings algorithm ［6］。其中模型透過隨機使用下列三個動作之一進行更新，包括 （a） 分割節點位置的變化； （b） 建立新的基函式； （c） 去除一個已有的基函式，然後再根據新的基函式，計算相關的接受機率。

除了更靈活的擬合能力，另一方面，我們還可以在 BMARS 構建的基函式中得到關於變數重要性的資訊。具體來說，我們可以透過每個變數被包含在基函式中的次數來判斷其重要性，次數越多則相對應的特徵越重要。

整合學習 & 貝葉斯累加回歸樹

除了模型混合 （Model Mixing），整合學習 （Ensemble Learning） ［13］ 提供了另一種組合模型的方法，並且取得了很好的效果。其主要想法是構造多個弱學習器然後將它們聚合成一個更強的學習器 ［8］。在某些情況下，單個模型本身很難準確地勾畫整個未知的複雜機制。因此，在整合學習框架中使用分而治之的方法是更好的策略，它允許每個模型擬合函式的一小部分。整合學習處理複雜資料的強大效能使其成為 BO 流程中機率模型的絕佳候選，然後其在 BO 的問題中並沒有被充分研究。一個非常合適的貝葉斯整合學習模型便是貝葉斯累加回歸樹 （Bayesian Additive Regression Trees， BART） ［4］，同時具有整合學習和貝葉斯框架的優勢。它是一種基於樹的模型，沒有固有的關於平滑性的假設， 所以在對材料學中經常遇到的非平滑目標函式進行建模時，是一種更靈活的模型。

BART 也屬於非引數迴歸模型，其具體公式如下：

其中

代表二叉迴歸樹，

代表二叉迴歸樹

的葉節點

的均值向量，

則是將

分配給

的函式。

BART 透過給各個二叉迴歸樹新增限制其複雜度的先驗，使得每棵樹都比較小，只擬合未知函式的一小部分，從而實現整合學習的效果。模型的訓練主要基於 ［4］ 提出的一個 backfitting MCMC algorithm，其抽樣的效果很好所以已經被廣泛使用。

同時，BART 也可以進行自動變數選擇。每個自變數的重要性取決於其作為拆分規則的頻率，頻率越高則相應的因素越重要。［3］ 進一步提出了一種基於排列的變數選擇方法，這也是確定因子顯著性的一個很好的選擇。

當實驗設計遇到 BO

透過在 BO 中使用 BMARS 或 BART 作為擬合黑盒函式的機率模型，［10］ 提出了一個高效的自動實驗設計的平臺，旨在減少所需的試驗次數以及尋找和發現最佳材料結構者的總費用。框架如圖 4 所示。

在這個工作流程中，實驗者從一個初始已知資料集開始，樣本大小可以很小（比如兩個）。然後，實驗者基於觀察到的資料集訓練 BMARS 或 BART 並從得到的後驗分佈中進行抽樣。基於這些樣本，計算每個可能的最佳材料結構（還未進行實驗的材料結構）對應的採集函式。接著選擇得分最高的候選者並進行相應的實驗。得到新的結果後，將新資料加入已知資料集。如果此時迴圈停止的準則被滿足，則結束迴圈，否則則基於擴充後的資料更新機率模型並指導下一輪的實驗。

在這個完全自動化的框架中，我們需要提供是初始樣本和迴圈停止標準。起始資料集可以是該專案之前的一些已經獲得的相關資料。如果我們沒有這種資訊，我們可以隨機進行少量實驗來填充資料庫並初始化框架中使用的機率模型。對於迴圈停止標準，材料學研究中常用的是找到滿足要求的材料或者達到了最多能承受的實驗次數（實驗預算）［15］。

第一回合：測試函式上的初步對決

為了對比 BMARS，BART 和 GPR 的表現，［10］ 選擇了2個常用的測試函式（Rosenbrock 函式 和 Rastrigin 函式）作為第一回合的考題，目標是找到這兩個函式的全域性最小值點。如圖 5 （a） 所示，Rosenbrock 函式有一個狹長的拋物線形平坦山谷，定位到這個山谷並不難，但是找到全域性最優點就並不簡單了。而圖 5 （b） 則對應著Rastrigin 函式，我們可以看到函式曲面存在著眾多的區域性最小值點，這使得這個任務變得更加具有挑戰性，因為演算法很容易陷入這些區域性最優值。

具體參與對比的 GPR 模型包括如下核函式：Radial-basis function （GP RBF），RBK kernel with automatic relevance detection （GP RBK ARD） ［2］，非平穩的 dot-product kernel （GP Dot） ［16］，更靈活的 deep kernel （GP DKNet） ［17］。 同時也和貝葉斯模型平均進行了比較 （BMA1 & BMA2） ［15］。

圖 6 展示了各模型具體的結果，曲線代表了各個每次迭代後以觀測值中的最小值，所以曲線下降越快代表著模型表現越好。其中藍色實線和紅色實線分別代表 BMARS 和 BART 的結果。可以發現在第一回合，BMARS 優勢明顯，同時 BART 也展現出和許多基於 GPR 的模型旗鼓相當的表現。

第二回合：MAX phases 中再次碰撞

為了進一步比較各個模型的效果，第二回合選擇在了一個實際的材料發現的問題，如圖 1 所示。這是一個前面簡單提到過的任務，我們想要在

和

中找到性質最優的結構。其中

，

和

分別可以在圖 1 中紅色，藍色和黑色區域選擇一種元素，剔除一些不穩定的組合後，可能具有最優性質的結構共計403種。圖 7 展示了找到最優結構所需的平均試驗次數，次數越少說明表現越好。我們可以看到，BART 和 BMARS 一直是表現的最好的幾個模型。

更多的討論和資料分析結果詳見 ［10］。

結語

隨著資料科學的不斷髮展，越來越多之前 human-centric 的任務可以透過 data-centric 的演算法得到更好的解決，從而促進相應學科的發展。這個過程其實也同時反過來促進了資料科學的發展。此文就是一個這樣的例子，透過更好地從歷史資料中提取有效資訊，來設計更合理高效的實驗計劃。相信不少學科中還有許多這種可以和資料科學很好地結合同時相互促進的任務，它們都在等待著一雙發現問題的眼睛。

參考文獻

［1］ M。 Aldeghi， F。 H ̈ase， R。 J。 Hickman， I。 Tamblyn， and A。 Aspuru-Guzik。 Golem： An algorithmfor robust experiment and process optimization， 2021。

［2］ S。 A。 Aye and P。 Heyns。 An integrated Gaussian process regression for prediction of remaininguseful life of slow speed bearings based on acoustic emission。Mechanical Systems and SignalProcessing， 84：485–498， 2017。

［3］ J。 Bleich， A。 Kapelner， E。 I。 George， and S。 T。 Jensen。 Variable selection for bart： an applicationto gene regulation。The Annals of Applied Statistics， pages 1750–1781， 2014。

［4］ H。 A。 Chipman， E。 I。 George， R。 E。 McCulloch， et al。 Bart： Bayesian additive regression trees。The Annals of Applied Statistics， 4（1）：266–298， 2010。 ［5］ D。 G。 Denison， B。 K。 Mallick， and A。 F。 Smith。 Bayesian mars。Statistics and Computing，8（4）：337–346， 1998。

［6］ P。 J。 Green。 Reversible jump markov chain monte carlo computation and Bayesian model deter-mination。Biometrika， 82（4）：711–732， 1995。

［7］ F。 Hase， M。 Aldeghi， R。 Hickman， L。 Roch， E。 Liles， M。 Christensen， J。 Hein， and A。 Aspuru-Guzik。 Olympus： a benchmarking framework for noisy optimization and experiment planning。Machine Learning： Science and Technology， 2021。

［8］ B。 Krawczyk， L。 L。 Minku， J。 Gama， J。 Stefanowski， and M。 Wo ́zniak。 Ensemble learning for datastream analysis： A survey。Information Fusion， 37：132–156， 2017。

［9］ A。 G。 Kusne， H。 Yu， C。 Wu， H。 Zhang， J。 Hattrick-Simpers， B。 DeCost， S。 Sarker， C。 Oses，C。 Toher， S。 Curtarolo， et al。 On-the-fly closed-loop materials discovery via bayesian activelearning。Nature communications， 11（1）：1–11， 2020。

［10］ B。 Lei， T。 Q。 Kirk， A。 Bhattacharya， D。 Pati， X。 Qian， R。 Arroyave， and B。 K。 Mallick。 Bayesianoptimization with adaptive surrogate models for automated experimental design。npj Computa-tional Materials， 7（194）， 2021。7

［11］ P。 Nikolaev， D。 Hooper， F。 Webber， R。 Rao， K。 Decker， M。 Krein， J。 Poleski， R。 Barto， andB。 Maruyama。 Autonomy in materials research： a case study in carbon nanotube growth。npjComputational Materials， 2（1）：1–6， 2016。

［12］ C。 J。 Paciorek and M。 J。 Schervish。 Nonstationary covariance functions for gaussian processregression。 InNIPS， pages 273–280。 Citeseer， 2003。

［13］ O。 Sagi and L。 Rokach。 Ensemble learning： A survey。Wiley Interdisciplinary Reviews： DataMining and Knowledge Discovery， 8（4）：e1249， 2018。

［14］ B。 Shahriari， K。 Swersky， Z。 Wang， R。 P。 Adams， and N。 De Freitas。 Taking the human out ofthe loop： A review of bayesian optimization。Proceedings of the IEEE， 104（1）：148–175， 2015。

［15］ A。 Talapatra， S。 Boluki， T。 Duong， X。 Qian， E。 Dougherty， and R。 Arr ́oyave。 Autonomousefficient experiment design for materials discovery with bayesian model averaging。Physical ReviewMaterials， 2（11）：113803， 2018。

［16］ C。 K。 Williams and C。 E。 Rasmussen。Gaussian processes for machine learning， volume 2。 MITpress Cambridge， MA， 2006。

［17］ A。 G。 Wilson， Z。 Hu， R。 Salakhutdinov， and E。 P。 Xing。 Deep kernel learning。 InArtificialintelligence and statistics， pages 370–378。 PMLR， 2016。

［18］ Zengjian Lv。Galileo-the pioneer of experimental science。

http：//

html。rhhz。net/kjdb/2017

0818。htm

， 2017。