黃俊華，唐 平，陳松齡，江小平，梁英蓬

(1.廣東工業(yè)大學(xué)自動化學(xué)院，廣州 510006；2.中山大學(xué)附屬第一醫(yī)院口腔科，廣州 510080；3.廣東工業(yè)大學(xué)醫(yī)院，廣州 510006)

1 概述

口腔種植體通過外科手術(shù)植入人體缺牙牙骨中完成缺牙修復(fù)。傳統(tǒng)的口腔種植外科手術(shù)術(shù)前設(shè)計(jì)都是基于X線平片(根尖片和全景片)，臨床操作中常因設(shè)計(jì)估計(jì)不足而發(fā)生牙槽骨意外穿孔等嚴(yán)重并發(fā)癥，而且往往給上部結(jié)構(gòu)修復(fù)帶來困難，導(dǎo)致種植修復(fù)的最終失敗[1]。隨著醫(yī)學(xué)圖像三維可視化等計(jì)算機(jī)輔助技術(shù)的發(fā)展，根據(jù)CT影像數(shù)據(jù)制作的手術(shù)導(dǎo)板增加了手術(shù)定位的準(zhǔn)確性和操作的可控性。國內(nèi)外分別推出VISIT導(dǎo)航系統(tǒng)和 SIM/plant等口腔種植軟件[2-3]。文獻(xiàn)[1]通過重建CT影像以人工交互的方式測量牙槽骨相關(guān)數(shù)據(jù)，但該方法工作步驟繁瑣，容易造成因主觀因素而導(dǎo)致的定位參數(shù)不準(zhǔn)確。

結(jié)合量子計(jì)算和遺傳算法的量子遺傳算法有較好的全局搜索能力，已廣泛地應(yīng)用于解決組合優(yōu)化和多屬性決策規(guī)劃問題。文獻(xiàn)[4]提出用于解決背包問題的遺傳量子算法(Genetic Quantum Algorithm, GQA)，種群個體的多樣性是基于量子比特和量子態(tài)疊加特性表示的。文獻(xiàn)[5]提出一種利用量子染色體相位比較來更新旋轉(zhuǎn)角的量子遺傳算法。文獻(xiàn)[6]提出一種基于量子位 Bloch球面的量子遺傳算法(Bloch Quantum Genetic Algorithm, BQGA)，該算法中種群的染色體有并行的3條基因鏈，優(yōu)化性能優(yōu)于普通量子進(jìn)化算法。

現(xiàn)有的口腔種植系統(tǒng)智能程度不高，而且鮮有運(yùn)用人工智能解決種植手術(shù)的學(xué)術(shù)報(bào)道?；谏鲜龇治觯瑸榱颂岣叻N植手術(shù)的精度以及種植系統(tǒng)的智能處理能力，針對口腔醫(yī)學(xué)螺旋CT數(shù)據(jù)量較大的特點(diǎn)，本文提出一種改進(jìn)的 Bloch量子遺傳算法(Improved Bloch Quantum Genetic Algorithm, IBQGA)，用于完成口腔種植體的智能設(shè)計(jì)。

2 量子遺傳算法

2.1 基本量子遺傳算法

基本量子遺傳算法可以從以下方面進(jìn)行說明：

(1)量子比特的表示

量子比特是量子計(jì)算的基本單位?；贐loch球面的量子比特可以用φ和θ一對相位角度表示，其對應(yīng)的向量表達(dá)式如下：

其中，|·>為量子表示符號。

本文采用基于Bloch球面的量子比特作為BQGA種群個體的基因位。設(shè)pi為基于Bloch球面的染色體，則其染色體編碼為：

其中，i是進(jìn)化種群中第i個染色體；φij、θij為第i個染色體第j號基因的一組相位角度；i=1,2,…,m，j=1,2,…,n，m是種群規(guī)模，n是量子位數(shù)。

(2)量子旋轉(zhuǎn)門

在量子遺傳算法中，對量子比特狀態(tài)實(shí)施遺傳操作是由量子旋轉(zhuǎn)門來實(shí)現(xiàn)的，由酉矩陣與種群個體的量子比特編碼運(yùn)算實(shí)現(xiàn)。對于基于Bloch球面的量子比特，相應(yīng)的量子旋轉(zhuǎn)門是3×3的酉矩陣U：

其中，Δφ和Δθ都是量子旋轉(zhuǎn)門的旋轉(zhuǎn)角度，可以控制變換后的量子比特所處的狀態(tài)；旋轉(zhuǎn)角Δφ、Δθ的正負(fù)和大小控制著進(jìn)化算法的收斂速度，如果幅值過大，會導(dǎo)致種群早熟；若幅值過小，會使收斂速度減慢。

文獻(xiàn)[6]中的 BQGA算法充分利用量子算法的多態(tài)性，但實(shí)際進(jìn)化的是染色體中的某一基因鏈，造成最終優(yōu)化群體中只有小部分優(yōu)化個體。文獻(xiàn)[7]提出一種引入進(jìn)化穩(wěn)定策略的量子競爭決策算法，用以解決旅行商問題。該算法利用競爭力函數(shù)評價(jià)競爭者的競爭力，從而提高算法的全局優(yōu)化能力。

2.2 改進(jìn)算法

針對上述 BQGA存在的問題，結(jié)合本文應(yīng)用，現(xiàn)提出下列3種改進(jìn)方法：

(1)把種群分為不同的進(jìn)化群體，即特征群體，使各特征群體獨(dú)自進(jìn)化。

考慮到不同的植入點(diǎn)對種植方案的差異，故以植入點(diǎn)為分類標(biāo)準(zhǔn)，把種群細(xì)分為不同的特征群體。每個個體在自己所屬的特征群體中進(jìn)化：即與同一特征群體中的其他個體比較，執(zhí)行遺傳操作和交叉操作。

(2)針對不同的特征群體及特征群體的進(jìn)化程度，實(shí)施自適應(yīng)調(diào)整進(jìn)化步長的遺傳操作。

對于量子遺傳算法的遺傳操作，文獻(xiàn)[4-5]通過遍歷查詢表生成量子旋轉(zhuǎn)門角度，此方法適合二進(jìn)制編碼的量子進(jìn)化算法。雖能推廣到一般的量子算法中，但涉及的多路判斷影響算法的效率。在不同特征群體的基礎(chǔ)上，相應(yīng)進(jìn)化特征群體中各個體的植入角度。本文采用如下自適應(yīng)調(diào)整旋轉(zhuǎn)角度的生成方法：

其中，φb和 θb為最優(yōu)解中量子位概率幅所對應(yīng)的相位；φc和θc是當(dāng)前解中量子位概率幅所對應(yīng)的相位；若A=0時，sgn(A)取正負(fù)號均可；若B=0時，sgn(B)取正負(fù)號均可。Δφ0和Δθ0同為種群進(jìn)化的步長。fij是第j代種群某特征群體中個體i的適應(yīng)值。fjmax和fjmin分別為當(dāng)代此群體的最大優(yōu)適應(yīng)值和最小適應(yīng)值。

(3)引入隨機(jī)錯位的交叉操作

在許多遺傳算法的應(yīng)用中容易發(fā)現(xiàn)，當(dāng)種群進(jìn)化到一定程度時，種群中各染色體相似或相同[8]，種群的多樣性減弱，難以突破局部最優(yōu)所帶來的影響。而適當(dāng)?shù)慕徊娌僮髂鼙ＷC遺傳算法良好的搜索性能，有利于保持種群的多樣性?？紤]到經(jīng)典遺傳算法中，二進(jìn)制編碼的染色體執(zhí)行交叉操作的特點(diǎn)是，分別隨機(jī)地選擇兩父代的某一個基因位作為交叉基因位，父代之間的交叉基因位相互交換，生成2個新的子代染色體，從而子代染色體隱含父代染色體的特征信息。針對醫(yī)學(xué)影像數(shù)據(jù)數(shù)據(jù)量大的特點(diǎn)，本文算法引入隨機(jī)的錯位交叉操作。隨機(jī)的錯位交叉操作的算法如下：

(1)隨機(jī)生成一個0～1之間的數(shù)rand。

(2)由父代 V1、V2 分別按下式生成子代 V1’、V2’：

3 種植體定位參數(shù)優(yōu)化算法設(shè)計(jì)

3.1 口腔種植定位問題的描述

設(shè)口腔種植定位問題有一個解集 SP={S1, S2,…,Sn}，n為解的總數(shù)，而每一個解Si都有一個評價(jià)種植方案 Si優(yōu)劣的綜合評價(jià) criticism(Si)。每一個解Si={Positioni, Anglei, Li, Ri}，i=1, 2,…, n；每個種植方案有如下特征屬性：植入點(diǎn)位置 Position，植入方向Angle，種植體長度L以及種植體半徑R。解中的植入點(diǎn)位置是主特征屬性，可根據(jù)各個解的植入點(diǎn)屬性的差異，把解集細(xì)分不同的小集合SPj，SPj即為特征群體?？谇环N植體定位設(shè)計(jì)就是要在解集SP中尋找一個這樣的Si，Si的各項(xiàng)特征屬性使其對應(yīng)的綜合評價(jià)criticism(Si)最好。

口腔種植手術(shù)中種植體定位，由種植體的半徑、長度、植入點(diǎn)與方位角度決定。上述4個種植體信息構(gòu)成具體的種植方案，是口腔種植體定位模型的參數(shù)。假設(shè)已知植入點(diǎn)的坐標(biāo)(x0, y0, z0)，植入方向角度ρ和 ω，種植體的半徑和長度分別為 r和 L；則可建立以植入點(diǎn)(x0, y0, z0)為原點(diǎn)的球面坐標(biāo)系，分別以植入方向角ρ和ω為此球面坐標(biāo)系的2個角度坐標(biāo)，以種植體長度L為此球面坐標(biāo)系的長度坐標(biāo)，口腔種植體定位模型如圖1所示。建立的球面坐標(biāo)系是以種植植入點(diǎn)(x0, y0, z0)為原點(diǎn)的種植方案搜索空間。

圖1 口腔種植體定位模型示意圖

下列數(shù)學(xué)語言可以描述口腔種植體定位模型：

口腔的影像數(shù)據(jù)量龐大且包含噪聲信息，因而不能直接作為搜索空間，還需要考慮定位模型和原數(shù)據(jù)空間的坐標(biāo)轉(zhuǎn)換。故需要提取缺牙區(qū)域的影像數(shù)據(jù)作為原數(shù)據(jù)空間。通過數(shù)字圖像處理的方法去噪，簡化區(qū)域數(shù)據(jù)；對提取的數(shù)據(jù)進(jìn)行坐標(biāo)變換，映射到搜索空間。因?yàn)榭谇粌?nèi)牙列和牙槽骨的分布及形態(tài)具有一定的特征，把口腔牙列分為左上/下和右上/下 4個部分(分別記為I、IV、II和III區(qū))。只需按照下式進(jìn)行適當(dāng)?shù)淖鴺?biāo)變換，口腔中缺牙區(qū)域的數(shù)據(jù)就能映射到相應(yīng)的搜索空間，解決口腔種植體的定位問題：

其中，x、y、z為原數(shù)據(jù)空間中某點(diǎn)的3個直角坐標(biāo)，sx、sy、sz為搜索空間下的坐標(biāo)。Tn為屬于n區(qū)的原數(shù)據(jù)空間映射到搜索空間的變換矩陣。

運(yùn)用坐標(biāo)變換和利用上面的定位模型，把種植方案的植入方向優(yōu)化限制為2個植入角度的組合優(yōu)化，且 2個角度的優(yōu)化范圍限制為 0°～90°，而不是 360°乃至空間中的搜索優(yōu)化，能有效地減少進(jìn)化算法搜索最優(yōu)解的難度，有利于提高問題的優(yōu)化效率。數(shù)據(jù)的提取與映射流程如圖2所示。

圖2 數(shù)據(jù)的提取與映射流程

3.2 量子染色體編碼

一個染色體代表一個種植方案，分別由植入點(diǎn)、種植方向和種植體半徑長度等參數(shù)決定。為了簡化染色體結(jié)構(gòu)，現(xiàn)考慮給定種植體半徑和長度的情況，即每條染色體由植入點(diǎn)編碼 1、植入點(diǎn)編碼 2、角度編碼1以及角度編碼2這4個部分組成。每條染色體的量子比特均采用如式(1)表示的基于 Bloch球面的實(shí)數(shù)相位角度編碼。

3.3 適應(yīng)值函數(shù)

合理的種植體定位需要種植體外圍有足夠密度的牙骨支撐，以及植入方向符合生物力學(xué)、吻合度，此外，還需要考慮牙槽骨的骨密度、形態(tài)、寬度、高度和骨缺損等情況[9-10]。因此，本文算法的適應(yīng)值函數(shù) Fitness分為支撐部分 Fsup和方向部分 Fdir，F(xiàn)sup評價(jià)此種植方案的種植體外圍牙骨密度的大小；Fdir作為種植體對應(yīng)的植入方向是否滿足生物力學(xué)的評價(jià)。最優(yōu)的種植方案有如下特點(diǎn)：適應(yīng)值函數(shù)的支撐部分和方向部分的適應(yīng)值總和最大，即max{Fsup+Fdir}。

Fsup作為評價(jià)種植體周圍的牙骨灰度，假設(shè)植入坐標(biāo)(x0, y0, z0)，角度ρ、ω，半徑r和長度L已知，其數(shù)學(xué)表述如下：

其中，value(x, y, z)表示的是點(diǎn)(x, y, z)對應(yīng)的骨灰度值；X、Y、Z為原數(shù)據(jù)空間的各維度的最大值。

種植方向的適應(yīng)度Fdir，由2個定位方向角度來評價(jià)。牙長軸與種植體植入方向的夾角應(yīng)在一定角度范圍以內(nèi)，為合理的種植方向，符合牙合力傳導(dǎo)的規(guī)律[10]。定位方向的適應(yīng)度可通過引入與植入方向和牙長軸的夾角相關(guān)的懲罰函數(shù)決定。如果植入方向和牙長軸的夾角較大，則使適應(yīng)度值為0；若夾角在口腔醫(yī)學(xué)允許范圍之內(nèi)，可增大對應(yīng)的適應(yīng)度值。

本文算法步驟如下：

Step1 隨機(jī)初始化種群，并建立一個用以記錄種植點(diǎn)、植入方向和最優(yōu)適應(yīng)值的空表 list。初始化迭代截止代數(shù)。

Step2 計(jì)算種群中各個體的適應(yīng)值。然后根據(jù)植入點(diǎn)的不同劃分特征群體，記入 list。其中，各特征群體需記錄最優(yōu)及最小適應(yīng)值。

Step3 根據(jù)list表中當(dāng)前特征群體的最優(yōu)適應(yīng)值以及最小適應(yīng)值，分別按式(2)～式(4)更新群體中個體，直至當(dāng)前特征群體更新完畢。

Step4 重復(fù)Step3，直至全部群體被更新為止。

Step5 更新當(dāng)前代數(shù)，并返回 Step2，直到滿足收斂條件為止。

4 實(shí)驗(yàn)結(jié)果與分析

為了證明本文算法在解決種植體定位問題具有優(yōu)于經(jīng)典優(yōu)化算法的特性，下面將本文算法與 GA、BQGA進(jìn)行比較。生成一個模擬CT醫(yī)學(xué)影像數(shù)據(jù)的三維數(shù)組，作為測試種植體定位算法的實(shí)驗(yàn)數(shù)據(jù)。實(shí)驗(yàn)數(shù)據(jù)是一個隨機(jī)生成的服從正態(tài)分布的大小為10×10×10的三維數(shù)組，僅設(shè)置 3個對角相鄰，且數(shù)組標(biāo)號遞增的數(shù)組元素值為10，已知最優(yōu)的坐標(biāo)及對應(yīng)的最優(yōu)角度(π/4)和理論最優(yōu)值(37)。下面利用GA、BQGA與本文算法測試上述測試數(shù)據(jù)，得出各優(yōu)化算法的綜合比較結(jié)果如表1所示，其中，qq表示平均角度；pp表示角度平均相對誤差；kk表示平均優(yōu)化值。

表1 各優(yōu)化算法的綜合比較結(jié)果

3種算法種群規(guī)模均為 200，截止進(jìn)化代數(shù)均為200。GA中的交叉概率為 0.8，變異概率為 0.3，而BQGA與本文算法的交叉概率均為 0.6，變異概率均為 0.1。進(jìn)化算法的優(yōu)化值相對誤差比較如圖 3所示。

圖3 進(jìn)化算法的優(yōu)化值相對誤差比較

從表1的實(shí)驗(yàn)數(shù)據(jù)中，可以得出：GA的優(yōu)化效率不如 BQGA和本文算法高。在 10次實(shí)驗(yàn)中，GA只有2次優(yōu)化接近最優(yōu)方案，即植入坐標(biāo)優(yōu)化為最優(yōu)的植入坐標(biāo)，總體優(yōu)化值的平均相對誤差是0.416 6；而 BQGA和本文算法均能收斂且優(yōu)化值的平均相對誤差均優(yōu)于GA。在本文實(shí)驗(yàn)的GA中，種群中的坐標(biāo)以及角度的優(yōu)化信息都依賴于交叉以及變異算子。即使增大種群數(shù)量或進(jìn)化代數(shù)，也難以保證GA能搜索到理論的最優(yōu)坐標(biāo)。而角度優(yōu)化僅在坐標(biāo)得以優(yōu)化下討論才有意義，因此，植入點(diǎn)坐標(biāo)的優(yōu)化程度決定了最終優(yōu)化解的優(yōu)劣。在種群進(jìn)化過程中，如果種群中坐標(biāo)及植入角度得到充分的優(yōu)化，最終得到的優(yōu)化解就越接近問題的最優(yōu)解。與 BQGA相比，由于本文算法考慮到種群個體間的差異性，本文算法在角度平均相對誤差及優(yōu)化值的平均相對誤差上均優(yōu)于BQGA算法。

針對中山大學(xué)第一附屬醫(yī)院提供的 CT影像資料，結(jié)合本文研究的種植定體位模型，將本文算法應(yīng)用于種植體方案設(shè)計(jì)優(yōu)化上，完成種植體的智能生成，最終優(yōu)化結(jié)果的矢狀切面圖如圖4所示，其中，牙骨中的直線部分為種植體。

圖4 本文算法優(yōu)化種植體方案的矢狀切面圖

由圖4可見，本文算法能正確生成缺牙區(qū)域的種植體種植方案。但就其植入方向等指標(biāo)的合理性，仍有待提高與改進(jìn)。造成這樣的原因可能是圖像處理的方法及其參數(shù)的選取設(shè)置不當(dāng)造成的，也可能是量子遺傳算法適應(yīng)值函數(shù)的設(shè)置問題。

5 結(jié)束語

本文提出一種改進(jìn)的量子遺傳算法。利用基于Bloch球面的量子基因編碼種群多樣性，引入交叉操作，使該算法的全局搜索能力和種群的多樣性得以保證。實(shí)驗(yàn)結(jié)果表明，該算法比GA和BQGA有較好的搜索能力和收斂性能，能解決種植體的定位問題。今后的研究內(nèi)容主要有：如何合理設(shè)置進(jìn)化算法的參數(shù)與約束，如種植體長度半徑以及適應(yīng)函數(shù)的定型等。

[1]康 璐, 黃遠(yuǎn)亮, 顧立栩, 等.計(jì)算機(jī)輔助設(shè)計(jì)軟件在口腔種植外科的應(yīng)用研究[J].口腔頜面外科雜志, 2008,18(4): 265-268.

[2]吳 婷, 廖文和, 戴 寧, 等.口腔種植導(dǎo)板計(jì)算機(jī)輔助技術(shù)研究[J].生物醫(yī)學(xué)工程學(xué)雜志, 2011, 28(1): 1-6.

[3]王培志, 夏 露, 陳 寧, 等.種植導(dǎo)航模板的計(jì)算機(jī)輔助設(shè)計(jì)和制造[J].中國口腔種植學(xué)雜志, 2010, 15(3):128-130, 133.

[4]Han Kuk-Hyun, Kim Jong-Hwan.Genetic Quantum Algorithm and Its Application to Combinatorial Optimization Problems[C]//Proc.of IEEE Conference on Evolutionary Computation.[S.l.]: IEEE Press, 2000.

[5]李士勇.一種基于相位比較的量子遺傳算法[J].系統(tǒng)工程與電子技術(shù), 2010, 32(10): 2219-2222.

[6]李士勇, 李盼池.量子計(jì)算與量子優(yōu)化算法[M].哈爾濱:哈爾濱工業(yè)大學(xué)出版社, 2009.

[7]劉 勇, 馬 良, 寧愛兵.量子競爭決策算法及其在旅行商問題中的應(yīng)用[J].計(jì)算機(jī)應(yīng)用研究, 2010, 27(2):586-589.

[8]陳國龍, 陳火旺, 郭文忠, 等.基于隨機(jī)錯位算術(shù)交叉的遺傳算法及其應(yīng)用[J].模式識別與人工智能, 2004,17(2): 250-255.

[9]周 苗, 陳松齡, 陳建靈, 等.螺旋CT在牙種植術(shù)前評估和設(shè)計(jì)中的應(yīng)用[J].實(shí)用醫(yī)技雜志, 2005, 12(6):1562-1563.

[10]孫婷婷, 潘 瑾, 崔明軍, 等.基于CT影像的牙種植模板相關(guān)的頜骨解剖學(xué)研究[J].口腔頜面外科雜志, 2007,17(3): 252-255.