張 丁 吳盈鋒 何迎花 劉 科 汪啟勝1, 何建華
1(中國科學(xué)院上海應(yīng)用物理研究所 上海201800)2(中國科學(xué)院大學(xué) 北京100049)3(中國科學(xué)院上海高等研究院 上海201210)4(武漢大學(xué)高等研究院 武漢430072)
生物大分子晶體學(xué)是在原子分辨率水平上解析生物大分子結(jié)構(gòu)最主要的技術(shù)手段之一。相比于常規(guī)實驗室的X射線光源,同步輻射的高通量性使得之前絕大多數(shù)很難或者不能進(jìn)行的實驗都能成功開展下去,基于同步輻射的生物大分子晶體學(xué)也受益頗多。除了光源帶來的諸多優(yōu)越性外,隨著實驗硬件、軟件的突破以及實驗方法的創(chuàng)新和發(fā)展,同步輻射生物大分子晶體學(xué)也飛快發(fā)展成熟[1]。但隨著生命科學(xué)領(lǐng)域的迅猛發(fā)展,海量的新蛋白質(zhì)結(jié)構(gòu)的分析測定需要依靠同步輻射生物大分子晶體學(xué)的實驗方法,對同步輻射生物大分子晶體學(xué)光束線站的需求日益增加[2],光束線必須要能夠極其可靠高效地進(jìn)行高通量實驗來滿足這種需求,重要的手段就是實現(xiàn)光束線的自動化運(yùn)行,建立自動化程度高的集成系統(tǒng)已經(jīng)成為幾乎所有同步輻射生物大分子晶體學(xué)實驗站發(fā)展的重要工作[3]。
一個完整的生物大分子晶體學(xué)實驗步驟相對繁復(fù),包括復(fù)雜的光束線優(yōu)化、樣品安裝和準(zhǔn)直、實驗條件確定、樣品衍射情況檢查、實驗策略的制定優(yōu)化、數(shù)據(jù)采集以及后續(xù)的數(shù)據(jù)處理與分析。作為整個實驗流程的第一個環(huán)節(jié),實驗前光束線上光學(xué)元件的配置優(yōu)化對于用戶實驗是十分重要的。在此之前,已經(jīng)有一些實驗室提出了光束線自動優(yōu)化的方案,Pugliese等[4]將模糊邏輯應(yīng)用到光束線自動優(yōu)化上;Olivier等[5]使用波前分析法作為優(yōu)化策略,但這些方案都不能保證全局優(yōu)化,直到2015年杜永華[6]將遺傳算法成功應(yīng)用到新加坡光源光束線全局自動優(yōu)化上,隨后他又采用差分進(jìn)化算法,并使用LabVIEW開發(fā)出光束線自動優(yōu)化程序AI-BL1.0[7]。上海光源的衍射線站借鑒新加坡光源的成功經(jīng)驗,將該程序應(yīng)用到上海光源基于EPICS系統(tǒng)(Experiment Physics and Industrial Control System)的控制平臺上[8]。上述案例中的優(yōu)化目標(biāo)都是樣品位置的光通量,然而光束性能指標(biāo)除了光通量之外,還包括光斑尺寸、光束位置等,因此光束線的優(yōu)化問題是一個典型的多目標(biāo)優(yōu)化問題。非支配排序遺傳算法II(Nondominated Sorting Genetic Algorithm II,NSGA-II)是一種建立在Pareto最優(yōu)解理論上的多目標(biāo)全局優(yōu)化算法,非常適用于光束線復(fù)雜系統(tǒng)的全局優(yōu)化問題上。同時上海光源晶體學(xué)光束線站的運(yùn)動控制和狀態(tài)數(shù)據(jù)獲取都是基于EPICS實現(xiàn)的,Python和EPICS的接口已經(jīng)被各光源和光束線廣泛應(yīng)用,因此利用Python開發(fā)優(yōu)化程序,不僅能使該程序方便應(yīng)用到晶體學(xué)光束線的自動優(yōu)化上,同時也能使程序有效地推廣和應(yīng)用到其他線站的優(yōu)化上。
本文提出一種基于NSGA-II的多目標(biāo)優(yōu)化策略,并將其應(yīng)用于上海光源P2生物防護(hù)蛋白質(zhì)晶體學(xué)線站(BL10U2)[9]的自動優(yōu)化上,用Python實現(xiàn)的優(yōu)化程序在BL10U2上完成了初步測試,并達(dá)到了預(yù)期目標(biāo)。
同步輻射光的性能指標(biāo)包括光子能量、光通量、光斑大小、光束位置等。良好的光束質(zhì)量是獲得可靠實驗數(shù)據(jù)的基本保證,也是提高線站運(yùn)行效率的前提條件。光束線從電子儲存環(huán)引出同步輻射光,再按實驗要求對其進(jìn)行調(diào)制,并將其傳輸?shù)綄嶒炚具M(jìn)行實驗[10]。狹縫、單色器、聚焦鏡等光束線設(shè)備對同步光進(jìn)行準(zhǔn)直、單色化、聚焦和傳輸,將特定性能的同步光傳輸?shù)綄嶒炚?。光束線上每個部件的姿態(tài)都會影響實驗站樣品處的同步光性能,因此,光束線優(yōu)化問題是一個典型的多目標(biāo)全局優(yōu)化問題。
實際問題中大都具有多個目標(biāo)需要同時滿足,即在同一個問題模型中同時存在幾個非線性目標(biāo),這些目標(biāo)函數(shù)同時需要進(jìn)行優(yōu)化處理,并且通常是互相沖突的,此類問題被稱為多目標(biāo)優(yōu)化問題(Multi-objective Optimization Problem,MOP)[11]。多目標(biāo)優(yōu)化問題的數(shù)學(xué)描述由決策向量、目標(biāo)函數(shù)、約束條件組成,一般表述如下:
式中:x為決策向量;y表示目標(biāo)函數(shù),也叫目標(biāo)向量;X為決策空間,由決策向量構(gòu)成;Y是目標(biāo)空間,由目標(biāo)向量構(gòu)成;e(x)≤0是約束條件,表示決策向量的可行取值范圍。
MOP由多個目標(biāo)函數(shù)組成,它會產(chǎn)生多個解,這些解不能像單目標(biāo)優(yōu)化問題那樣直接比較出它們之間的優(yōu)劣,被稱為目標(biāo)函數(shù)的非支配解(Nondominated Solutions),也叫做Pareto最優(yōu)解(Pareto optimal solutions)。多個Pareto最優(yōu)解組成Pareto最優(yōu)解集,求得Pareto最優(yōu)解集就是MOP要達(dá)到的最終目的,決策者根據(jù)實際情況和自身需要從解集中選取合適的Pareto最優(yōu)解。以下為多目標(biāo)優(yōu)化的重要定義:
定義1:Pareto支配
對于最小化多目標(biāo)優(yōu)化問題,對于n個目標(biāo)函數(shù)fi(x),i=1,…,n,任意給定兩個決策向量xa、xb,如果有以下兩個條件成立,則稱xa支配xb。
1)對于?i∈1,2,…,n,都有fi(xa)≤fi(xb)成立。
2)?i∈1,2,…,n,有fi(xa)<fi(xb)成立。
定義2:Pareto最優(yōu)解
如果一個解x*不被其他任意一個解支配,則被稱之Pareto最優(yōu)解。
定義3:Pareto最優(yōu)解集
多目標(biāo)優(yōu)化問題中所有Pareto最優(yōu)解組成的解集稱為Pareto最優(yōu)解集(Pareto Set)。
定義4:Pareto前沿
Pareto Set中每個解對應(yīng)的目標(biāo)向量組成的集合稱之為Pareto前沿(Pareto Front),簡稱為PF。
光束線的優(yōu)化主要是通過改變決定光學(xué)元件姿態(tài)的運(yùn)動電機(jī)位置,來調(diào)節(jié)光束性能(光通量、光斑位置等)。對照多目標(biāo)優(yōu)化問題數(shù)學(xué)模型來看,一組電機(jī)位置對應(yīng)著一個決策向量,該組電機(jī)位置決定的光束性能則是對應(yīng)的目標(biāo)向量,約束條件則是各個電機(jī)的行程范圍。光束線優(yōu)化問題可以抽象成尋找多個目標(biāo)函數(shù)的極值問題,多個目標(biāo)函數(shù)無法兼顧,最終優(yōu)化得到的結(jié)果是這些目標(biāo)函數(shù)的Pareto最優(yōu)解集,也就是若干組非支配的電機(jī)位置,決策者根據(jù)實際需要選擇合適的電機(jī)位置。
1994年,Deb等[12]在遺傳算法的基礎(chǔ)上對其進(jìn)行了改進(jìn),從而提出了NSGA,它是一種基于Pareto最優(yōu)解的多目標(biāo)優(yōu)化算法,該算法對種群中的個體進(jìn)行分層,處在同層的非支配個體共享一個虛擬適應(yīng)度值,并采用共享函數(shù)法來保證個體多樣性。盡管NSGA突破了傳統(tǒng)的遺傳算法,可以求解任意多的目標(biāo)函數(shù)問題,并且在實際生產(chǎn)生活中針對于多目標(biāo)求解問題得到了廣泛的應(yīng)用,但是它還是存在一些問題,主要體現(xiàn)在支配排序的計算復(fù)雜度高、缺少精英機(jī)制、需要指定共享參數(shù)這三個方面。2000年,Deb又在NSGA的基礎(chǔ)上提出了一種帶有精英機(jī)制的快速非支配排序遺傳算法(NSGA-II)[13],該算法采用了快速非支配排序算法和擁擠比較算子,計算復(fù)雜度比NSGA極大地降低;通過計算擁擠距離和定義擁擠偏序,使得解的分布性更佳,保持了種群的多樣性;引入了精英保留策略,有效防止了精英個體的丟失。
NSGA-II的關(guān)鍵點在于快速非支配排序、擁擠度和擁擠算子以及精英策略,以下是它們簡單的介紹[14]。
2.2.1快速非支配排序
找出初始種群中的Pareto最優(yōu)解集,將它們的Pareto等級設(shè)置為1,然后將它們從種群中刪除,在余下種群中Pareto最優(yōu)解集,將它們的Pareto等級設(shè)置為2,以此類推,通過支配關(guān)系對整個種群進(jìn)行排序,得到所有個體的Pareto等級。
2.2.2擁擠度和擁擠算子
遺傳算法有自動收斂的性質(zhì),設(shè)置擁擠度nd能夠使同一Pareto等級中的個體相互間分開,從而保證了種群的多樣性。先令同一非支配層里所有個體的nd為0,然后對每一個目標(biāo)函數(shù)fm,根據(jù)目標(biāo)函數(shù)值的大小對個體進(jìn)行排序,排序后兩個邊界個體的擁擠度nd為無窮大,其余個體的擁擠度:
累加個體在所有目標(biāo)函數(shù)的擁擠度得到該個體最終的擁擠度。
利用每個個體的Pareto等級和擁擠度可以得到其擁擠比較算子,概括來說就是:兩個不同Pareto等級的個體,Pareto等級更低的個體更優(yōu);兩個相同Pareto等級的個體,擁擠度更大的個體更優(yōu)。
2.2.3精英策略
NSGA-II采用精英策略將優(yōu)良個體遺傳到下一代,具體操作為:在父代種群Pt基礎(chǔ)上通過遺傳算子來產(chǎn)生子代種群Ct,將父子種群混合組成新種群Rt,對Rt進(jìn)行非支配排序,按Pareto等級由低到高依次將非支配集F1,F(xiàn)2,…,F(xiàn)m放入新的父代種群Pt+1中,直到Pt+1的大小超出原有種群大小,按擁擠度從小到大從Pt+1中剔除Fm中的個體直到Pt+1的大小和原有種群大小相等。
算法流程如圖1所示。NSGA-II具體過程如下:
圖1 NSGA-II流程圖Fig.1 Flowchart of NSGA-II
1)初始化種群。隨機(jī)產(chǎn)生規(guī)模大小為Np的初始父代種群:
式中:rand(0,1)為0到1的隨機(jī)數(shù);n為決策向量里的變量個數(shù);uj、lj分別為變量的上、下限。
2)快速非支配排序。
3)計算擁擠度。
4)交叉、變異得到子代種群。
5)混合父代種群和子代種群,確定混合種群中個體的Pareto等級和擁擠度。
6)采用精英策略,選擇混合種群優(yōu)秀的個體生成新的父代種群。
7)重復(fù)步驟4)~6),直到達(dá)到設(shè)定的進(jìn)化代數(shù),即可獲得Pareto最優(yōu)解集。
NSGA-II中種群個體(決策向量)的Pareto等級和擁擠度都是由其對應(yīng)的目標(biāo)向量決定的,因此,在算法實現(xiàn)中,必須保證決策向量和目標(biāo)向量的是一一對應(yīng)的,也就是要保證一組電機(jī)位置和由其決定的光束狀態(tài)信息的對應(yīng)關(guān)系。上海光源晶體學(xué)光束線站的運(yùn)動控制和光束狀態(tài)信息采集都是基于EPICS系統(tǒng)[15]實現(xiàn)的,該系統(tǒng)的核心是IOC軟件,IOC包括分布式運(yùn)行數(shù)據(jù)庫(Distributed Runtime Database)和輸入輸出支持模塊。EPICS利用IOC運(yùn)行數(shù)據(jù)庫中的記錄完成對輸入和輸出信號量的控制,記錄值的標(biāo)識稱為過程變量(Process Variable,PV),記錄值在數(shù)據(jù)交換中充當(dāng)基本單元,在客戶端與服務(wù)器端之間建立通信連接??蛻舳伺c服務(wù)器端模型建立通道訪問協(xié)議(基于TCP/IP協(xié)議),并提供相應(yīng)的應(yīng)用接口函數(shù)庫,以供數(shù)據(jù)讀寫、連接監(jiān)控和訪問監(jiān)控等服務(wù)。
電機(jī)的運(yùn)動控制采用EPICS中的Motor模塊[16],一個Motor記錄代表一個電機(jī)實體,記錄中的域?qū)?yīng)著電機(jī)屬性,其中域VAL為控制電機(jī)運(yùn)動的位置,域RBV為電機(jī)位置的回讀值。光通量和光束位置的獲取和處理則是由QBPM[17]和數(shù)字信號處理器組合完成。QBPM探測電極上獲取到光電流信號,然后模擬信號被轉(zhuǎn)換為數(shù)字信號后提供給IOC分析,從而得到最終光通量和光束位置。因此,電機(jī)運(yùn)動、電機(jī)位置獲取以及光束狀態(tài)信息的獲取都可以通過操作PV來實現(xiàn)。
NSGA-II的實現(xiàn)由Python完成,優(yōu)化程序通過PyEpics接口與EPICS連接,PyEpics不僅有對EPICS自帶接口函數(shù)的封裝,還有專門的電機(jī)控制模塊和PV操作模塊,使用者可以自由選擇最方便的模塊進(jìn)行操作。優(yōu)化程序接入EPICS后,可以控制電機(jī)運(yùn)動,獲取電機(jī)位置回讀值和光束信息,并在此基礎(chǔ)上執(zhí)行算法流程。自動優(yōu)化的軟硬件結(jié)構(gòu)如圖2所示。
圖2 光束線自動優(yōu)化軟硬件結(jié)構(gòu)Fig.2 Software and hardware structure of automatic beamline optimization
程序開始后,在設(shè)定的電機(jī)行程范圍內(nèi)隨機(jī)生成初始父代種群,根據(jù)種群個體使各個電機(jī)運(yùn)動到相應(yīng)位置,然后獲取此時的光通量和光束位置作為該個體的目標(biāo)向量,確定種群個體的Pareto等級和擁擠度,父代種群通過遺傳操作得到子代種群,根據(jù)子代種群信息控制電機(jī)運(yùn)動和獲取相應(yīng)光束狀態(tài)信息,混合父子種群并確定該混合種群中個體的Pareto等級和擁擠度,隨后依據(jù)精英策略,在混合種群中選擇個體生成新的父代種群,循環(huán)更新父代種群,直到達(dá)到設(shè)定的進(jìn)化代數(shù),迭代結(jié)束,從得到的Pareto最優(yōu)解集選出所需要的最終結(jié)果。圖3為優(yōu)化程序流程圖。
圖3 自動優(yōu)化程序流程圖Fig.3 Flowchart of automatic optimization program
本程序部署在上海光源P2生物防護(hù)蛋白質(zhì)晶體學(xué)光束線(BL10U2)。該線站具有P2生物安全防護(hù)功能,可滿足我國對傳染性病毒晶體結(jié)構(gòu)研究的需要,同時它也具有與其他的生物大分子晶體學(xué)光束線相同的功能,是上海光源里具有代表性的生物大分子晶體學(xué)光束線。光束線主要設(shè)備包括白光狹縫(Slit1)、水平偏轉(zhuǎn)鏡(M1)、雙平晶液氮冷卻單色器(DCM)、水平預(yù)聚焦鏡(M2)、次級光源狹縫(Slit2)、垂直聚焦鏡(M3)以及水平聚焦鏡(M4);除了以上光學(xué)設(shè)備,還包括束線調(diào)試及診斷必不可少的熒光靶、位置靈敏探測器及光強(qiáng)檢測等設(shè)備。布局如圖4所示。
圖4 上海光源BL10U2布局示意圖Fig.4 Layout of the optical equipment of BL10U2
由于P2線站還處于調(diào)試階段,測試優(yōu)化目標(biāo)選為次級光源點的光通量和光束位置,水平預(yù)聚焦鏡M2為被優(yōu)化設(shè)備,水平預(yù)聚焦鏡采用橢圓柱面反射鏡,對光束進(jìn)行水平方向進(jìn)行聚焦,形成水平方向的一級聚焦光斑。該鏡子涉及到的運(yùn)動方向和對應(yīng)控制電機(jī)如表1所示。次級光源點的光束性能狀態(tài)則是由水平預(yù)聚焦鏡后的QBPM獲取。
表1 水平預(yù)聚焦鏡M2的運(yùn)動方向及相應(yīng)控制電機(jī)Table 1 The movement directions and corresponding control motors of M2
先以不同種群大?。?5、20、25)進(jìn)行了3次測試,進(jìn)化代數(shù)設(shè)置為50,3個種群找到的Pareto最優(yōu)解的個數(shù)分別為9、13、14。種群數(shù)越大,能找到的Pareto解也越多,但后兩組的差距很小,可知當(dāng)種群規(guī)模到達(dá)20后,繼續(xù)增大種群規(guī)模對Pareto最優(yōu)解集的影響較小。在此基礎(chǔ)上,以種群大小為20、進(jìn)化代數(shù)為80增加了1次測試,得到Pareto解的個數(shù)為20。4次測試的交叉概率為0.9,變異概率為0.2,測試結(jié)果列于表2中,4組測試的Pareto前沿如圖5所示。
表2 測試條件及結(jié)果Table 2 Test conditions and results of the experiment
測試結(jié)果表明:1)優(yōu)化的目標(biāo)是同時尋找光通量的極大值和光束實際位置到參考位置偏差的極小值。在圖5中,橫坐標(biāo)以QBPM總探測電流來代表光通量,縱坐標(biāo)是光束位置偏差,從圖5中可以看出,兩個優(yōu)化目標(biāo)之間存在著制約關(guān)系。本次測試采取折中的方法,程序預(yù)先設(shè)定的優(yōu)化目標(biāo)以光束位置為標(biāo)準(zhǔn),將目標(biāo)位置定為距參考位置10μm的位置,程序運(yùn)行結(jié)束后會自動從Pareto最優(yōu)解集中選取出離目標(biāo)位置最近的位置和相應(yīng)的光通量作為優(yōu)化結(jié)果,4次測試的優(yōu)化結(jié)果如表3所示,4次測試的光通量都在1.8×10-6A左右的強(qiáng)度;2)在此之前上海光源晶體學(xué)光束線的調(diào)光工作主要是由線站工程師手動完成的,每次只能對單個運(yùn)動電機(jī)進(jìn)行優(yōu)化,即對每個電機(jī)在其行程范圍內(nèi)依次調(diào)整位置,直到光束狀態(tài)滿足實驗要求,這種調(diào)光模式效率非常低,整個優(yōu)化過程少則幾小時,多則數(shù)日。本文開發(fā)的優(yōu)化程序能同時對多個電機(jī)進(jìn)行全局優(yōu)化,所有測試的最長優(yōu)化時間在30 min左右,相較于手動調(diào)節(jié),優(yōu)化效率大大提高;3)圖5中,當(dāng)進(jìn)化代數(shù)為50時,種群大小為15得到的Pareto最優(yōu)解個數(shù)較少且散亂分布;增大種群規(guī)模到20后,數(shù)據(jù)點增多,但其分布有所起伏,不太均勻;繼續(xù)增大種群規(guī)模到25,數(shù)據(jù)點雖然沒有增加太多,但其分布更加均勻;保持種群規(guī)模20不變,增加進(jìn)化代數(shù)到80,得到數(shù)據(jù)點最多,且其分布也十分均勻和緊密。增大種群規(guī)模和進(jìn)化代數(shù)可以得到數(shù)量更多和質(zhì)量更優(yōu)的Pareto最優(yōu)解集,也能更好地滿足多樣的實驗需求。
表3 測試優(yōu)化結(jié)果Table 3 Optimization results of tests
圖5 4組測試的Pareto前沿(a)種群大小15,進(jìn)化代數(shù)50,(b)種群大小20,進(jìn)化代數(shù)50,(c)種群大小25,進(jìn)化代數(shù)50,(d)種群大小20,進(jìn)化代數(shù)80Fig.5 Pareto front of tests (a)Population size 15,generations 50,(b)Population size 20,generations 50,(c)Population size 25,generations 50,(d)Population size 20,generations 80
本文設(shè)計并開發(fā)了基于NSGA-II的光束線多目標(biāo)自動優(yōu)化程序,進(jìn)一步完善了上海光源光束線自動優(yōu)化方案。基于Python實現(xiàn)的優(yōu)化程序不僅能適用于上海光源生物大分子晶體學(xué)線站,還具備向上海光源其他線站推廣的潛能。自動優(yōu)化程序結(jié)合上海光源P2生物防護(hù)蛋白質(zhì)晶體學(xué)線站軟硬件條件,以次級光源點的光通量和光束位置為優(yōu)化目標(biāo)完成了部署與測試。測試結(jié)果表明:該優(yōu)化程序能準(zhǔn)確找到預(yù)定的優(yōu)化目標(biāo)。測試耗時均在30 min以內(nèi),相較于之前的單目標(biāo)優(yōu)化方案,基于NSGA-II的自動優(yōu)化程序不僅保證了自動優(yōu)化效率,還進(jìn)一步提高了光束線智能優(yōu)化的能力。
測試結(jié)果還表明:為了獲得更好的Pareto最優(yōu)解集,則需要花費(fèi)更多的優(yōu)化時間。后續(xù)為了實現(xiàn)真正的全局優(yōu)化,待優(yōu)化的光學(xué)元件會增多,相應(yīng)電機(jī)也會增多,電機(jī)運(yùn)動時間將會大大增加,為了獲得更準(zhǔn)確的光束信息,采樣時間也會增加,所以為了確保優(yōu)化效果,相關(guān)硬件性能的提高也是非常有必要的。同時測試時P2線站還處于調(diào)試階段,優(yōu)化程序得到了初步成功的結(jié)果,但仍需在P2線站穩(wěn)定運(yùn)行狀態(tài)下做進(jìn)一步的驗證。