基于改進(jìn)粒子群優(yōu)化算法的牙齒正畸路徑規(guī)劃方法

徐曉強(qiáng)，秦品樂(lè)，曾建朝

（中北大學(xué)大數(shù)據(jù)學(xué)院，太原 030051）

（*通信作者電子郵箱zjc@nuc.edu.cn）

0 引言

正畸（Orthodontics）又稱牙齒矯正，是指利用各種牙齒矯治器對(duì)患者牙齒的位置進(jìn)行重新排列，使其恢復(fù)到與正常牙齒形態(tài)相對(duì)一致的過(guò)程［1］。牙齒矯治器包括傳統(tǒng)有托槽矯治器、無(wú)托槽隱形矯治器等。近年來(lái)，隨著無(wú)托槽隱形矯治器相較于前者具有更舒適、美觀、節(jié)約矯治時(shí)間、容易清潔等特點(diǎn)［2］，并且患者通過(guò)無(wú)托槽隱形矯治系統(tǒng)可以提前看到矯治完畢后的牙齒具體排列情況。因此，越來(lái)越多的患者選擇無(wú)托槽隱形矯治器進(jìn)行牙齒矯正，這使得運(yùn)用計(jì)算機(jī)技術(shù)對(duì)按數(shù)學(xué)模型排列的虛擬牙齒正畸技術(shù)成為研究的熱點(diǎn)。復(fù)雜的牙齒正畸路徑規(guī)劃作為上述研究熱點(diǎn)的一個(gè)重要步驟，既要保證矯治路徑具有可行性，同時(shí)也要保證矯治結(jié)果具有優(yōu)化性。專業(yè)地說(shuō)，牙齒正畸路徑規(guī)劃是保證患者的全部牙齒矯治過(guò)程在符合正畸學(xué)的前提下，使多顆牙齒在經(jīng)過(guò)數(shù)個(gè)矯治階段后達(dá)到最終理想位置，并且在該過(guò)程中相鄰牙齒間不發(fā)生碰撞，因此需要尋找一條盡可能減少患者矯治時(shí)間、移動(dòng)距離、旋轉(zhuǎn)角度等因素的最優(yōu)路徑作為最終牙齒正畸移動(dòng)路徑。文獻(xiàn)［3］提取了人工牙三維三角網(wǎng)格數(shù)據(jù)的特征點(diǎn)建立了牙齒選擇系統(tǒng)，通過(guò)計(jì)算機(jī)數(shù)控技術(shù)完成了全口義齒的虛擬排牙方法。文獻(xiàn)［4］提出單顆牙齒自動(dòng)分割方法優(yōu)化了傳統(tǒng)方法，通過(guò)搜索分割路徑、圖像形態(tài)，擬合牙弓線，匹配最優(yōu)路徑確定牙縫邊界路徑，細(xì)化封閉分割輪廓完成整個(gè)牙齒的分割。此方法對(duì)于嚴(yán)重畸形的牙齒具有較好的正畸效果。文獻(xiàn)［5］提出了一種基于數(shù)字圖像與計(jì)算機(jī)斷層掃描模擬的虛擬近似真實(shí)正畸牙齒移動(dòng)，并通過(guò)給牙齒特定位置施加力來(lái)減少平移及旋轉(zhuǎn)次數(shù)。文獻(xiàn)［6］根據(jù)壓力和張力與特定的信號(hào)傳導(dǎo)相關(guān)因子，建立了局部梯度來(lái)重塑牙周，以更快地移動(dòng)牙齒，減少治療時(shí)間和依賴時(shí)間的不良后果。

隨著對(duì)遺傳算法、神經(jīng)網(wǎng)絡(luò)算法等啟發(fā)式算法的深入研究，發(fā)現(xiàn)上述算法大多存在計(jì)算量大、操作復(fù)雜、算法運(yùn)行時(shí)間過(guò)長(zhǎng)等問(wèn)題，越來(lái)越多的研究者提出運(yùn)用粒子群優(yōu)化（Particle Swarm Optimization，PSO）算法來(lái)解決路徑規(guī)劃的相關(guān)問(wèn)題。PSO 算法是Eberhart等［7］通過(guò)模擬鳥(niǎo)類行為規(guī)律，提出的一種全局優(yōu)化算法。PSO 算法因其結(jié)構(gòu)簡(jiǎn)單、效率較高等優(yōu)點(diǎn)已經(jīng)廣泛應(yīng)用于各類優(yōu)化問(wèn)題［8-9］，近年來(lái)對(duì)改進(jìn)粒子群算法在路徑規(guī)劃領(lǐng)域的研究陸續(xù)被提出。文獻(xiàn)［10］提出了斥力場(chǎng)下粒子群，利用柵格法進(jìn)行初步規(guī)劃以及根據(jù)不同的障礙物所構(gòu)建數(shù)學(xué)模型、適應(yīng)度函數(shù)，從而得到了移動(dòng)機(jī)器人最優(yōu)路徑。文獻(xiàn)［11］在光柵化的基礎(chǔ)上構(gòu)建了包括隨機(jī)地表模型、山地地形模型的數(shù)學(xué)模型，并分析了路徑威脅因素，引入遺傳算法的思想改進(jìn)了粒子群算法最終得到了低威脅的干擾機(jī)路徑。文獻(xiàn)［12］提出了一種基于高斯參數(shù)更新規(guī)則和粒子重新初始化的改進(jìn)粒子群算法，從縮短路徑長(zhǎng)度及提高路徑光滑度兩方面優(yōu)化了路徑規(guī)劃問(wèn)題。文獻(xiàn)［13］通過(guò)個(gè)體粒子進(jìn)化速度與群體離散度來(lái)調(diào)整粒子群的慣性權(quán)重，以此避免其過(guò)早收斂，并利用自然選擇方法改進(jìn)量子行為粒子群的優(yōu)化算法，將這兩種優(yōu)化的算法引入移動(dòng)機(jī)器人路徑規(guī)劃，并通過(guò)實(shí)驗(yàn)證明方法的可行性。文獻(xiàn)［14］以ER3A-C60 機(jī)器人作為研究對(duì)象，建立了運(yùn)動(dòng)學(xué)和動(dòng)力學(xué)模型，用參數(shù)合適的PSO算法求得了最小綜合誤差得到了一條最優(yōu)軌跡。

受上述改進(jìn)算法的啟發(fā)，針對(duì)虛擬口腔正畸治療系統(tǒng)中牙齒移動(dòng)路徑規(guī)劃問(wèn)題，本文提出了一種基于正態(tài)分布的簡(jiǎn)化均值粒子群的牙齒正畸路徑規(guī)劃方法。首先建立了單顆牙齒及整體牙齒的數(shù)學(xué)模型，并將牙齒正畸路徑規(guī)劃問(wèn)題轉(zhuǎn)化為帶約束的優(yōu)化問(wèn)題；其次，引入正態(tài)分布及均值粒子群，結(jié)合簡(jiǎn)化粒子群的位置項(xiàng)，提出了一種基于正態(tài)分布的簡(jiǎn)化均值粒子群優(yōu)化（Simplified Mean Particle Swarm Optimization based on Normal distribution，NSMPSO）算法；最后，構(gòu)造了包含平移路徑長(zhǎng)度、旋轉(zhuǎn)角度、碰撞檢測(cè)以及牙齒在單階段的移動(dòng)量、旋轉(zhuǎn)量的高安全性的適應(yīng)度函數(shù)，實(shí)現(xiàn)了牙齒正畸移動(dòng)路徑的規(guī)劃。在三大基準(zhǔn)測(cè)試函數(shù)上將NSMPSO 算法與其他三種粒子群算法進(jìn)行對(duì)比，結(jié)果表明，本文改進(jìn)的算法在收斂速度和收斂精度上均有明顯的提高，通過(guò)Matlab 的仿真實(shí)驗(yàn)表明了利用該數(shù)學(xué)模型和改進(jìn)算法求得的最優(yōu)路徑安全可靠，可以為醫(yī)生輔助診斷提供有效的幫助。

1 建立數(shù)學(xué)模型

1.1 單顆牙齒及整體牙齒的數(shù)學(xué)模型

本文以分割后的牙齒為研究對(duì)象。首先，采用有向包圍盒（Oriented Bounding Box，OBB）對(duì)分割后的單顆牙齒建立包圍盒，以包圍盒的中心點(diǎn)作為牙齒的特征點(diǎn)建立局部坐標(biāo)系，X、Y、Z軸的方向根據(jù)單顆牙齒包圍盒的長(zhǎng)寬高來(lái)確定。圖1為尖牙、側(cè)切牙、后磨牙在三維數(shù)學(xué)模型下單顆牙齒包圍盒的局部坐標(biāo)系。然后，以全部牙齒為目標(biāo)，選取兩顆后磨牙、兩顆尖牙的特征點(diǎn)所確定的平面作為基準(zhǔn)面來(lái)建立牙齒整體的世界坐標(biāo)系。圖2 為整體牙齒在三維數(shù)學(xué)模型下的世界坐標(biāo)系。最后，用每顆牙齒的包圍盒及中心點(diǎn)代替牙齒進(jìn)行整個(gè)矯正過(guò)程的展示，可以更直觀、形象地了解每個(gè)階段內(nèi)的牙齒位置的變化情況。在某一矯正階段內(nèi)，由于全部牙齒的移動(dòng)和旋轉(zhuǎn)是同時(shí)進(jìn)行的，則可將該階段的牙齒移動(dòng)路徑離散化表示為一系列路徑點(diǎn)的連線。

圖1 單顆牙齒的三維坐標(biāo)系Fig.1 Three-dimensional coordinate system of single tooth

圖2 整體牙齒數(shù)學(xué)模型Fig.2 Mathematical model of whole teeth

1.2 牙齒移動(dòng)路徑規(guī)劃的數(shù)學(xué)建模

將n顆牙齒的移動(dòng)路徑規(guī)劃分為m個(gè)階段，m由每個(gè)階段可移動(dòng)、旋轉(zhuǎn)的最大值及醫(yī)生經(jīng)驗(yàn)共同決定。牙齒在矯治過(guò)程中，在滿足單階段內(nèi)患者可以承受的牙齒最大移動(dòng)量、旋轉(zhuǎn)量并且無(wú)碰撞發(fā)生的前提下，確定每個(gè)階段的所有牙齒的位置信息，最終規(guī)劃出一條牙齒從初始位置到最終位置的最短路徑，即為可執(zhí)行的牙齒移動(dòng)最優(yōu)路徑。因此可將牙齒正畸路徑規(guī)劃問(wèn)題等價(jià)為帶約束的優(yōu)化問(wèn)題：其中：n為牙齒個(gè)數(shù)，m為矯治階段數(shù)。(xij，yij，zij)表示第i顆牙齒在第j個(gè)矯治階段結(jié)束后的位置坐標(biāo)，目標(biāo)函數(shù)之一f1表示牙齒平移路徑長(zhǎng)度，是通過(guò)計(jì)算空間中兩點(diǎn)間的歐氏距離構(gòu)建的。計(jì)算患者所有牙齒在全部矯治階段內(nèi)的移動(dòng)量之和，值越小表示牙齒在矯治過(guò)程中需要移動(dòng)的量越小，從而減少了患者的矯治時(shí)間。目標(biāo)函數(shù)f2表示牙齒旋轉(zhuǎn)角度，(αij，βij，γij)為第i顆牙齒第j個(gè)階段時(shí)包圍盒的三邊與世界坐標(biāo)系的三軸之間的夾角。單顆牙齒包圍盒的局部坐標(biāo)軸X1、Y1、Z1在全部牙齒建立的世界坐標(biāo)軸X2、Y2、Z2的投影分別為X1X2、Y1Y2和Z1Z2，組成的三個(gè)夾角即為α、β、γ，計(jì)算患者全部牙齒在整個(gè)矯治階段內(nèi)的旋轉(zhuǎn)量之和即為目標(biāo)函數(shù)二。根據(jù)正畸學(xué)及醫(yī)生臨床經(jīng)驗(yàn)可知，每顆牙齒在某個(gè)矯治階段內(nèi)移動(dòng)范圍不能大于0.5 mm、旋轉(zhuǎn)角度不能超過(guò)3°。考慮到每顆牙齒在規(guī)定范圍內(nèi)進(jìn)行移動(dòng)或者旋轉(zhuǎn)之后，可能出現(xiàn)某顆牙齒與左右相鄰的牙齒發(fā)生碰撞擠壓導(dǎo)致矯治失敗的情況，因此在每個(gè)矯治階段檢測(cè)牙齒間是否發(fā)生碰撞至關(guān)重要。綜上，用g1、g2、g3表示在對(duì)目標(biāo)函數(shù)f1、f2進(jìn)行優(yōu)化時(shí)的約束條件。約束條件g1是限制每顆牙齒在前后兩個(gè)階段的移動(dòng)量小于每個(gè)階段單顆牙齒能夠承受的最大移動(dòng)量；約束條件g2是限制每顆牙齒在前后兩個(gè)階段的旋轉(zhuǎn)量小于人類能夠承受的每個(gè)階段的最大旋轉(zhuǎn)量；約束條件g3是檢測(cè)矯治過(guò)程中相鄰牙齒是否碰撞，式中ci表示判斷第i顆牙齒與其左右相鄰的牙齒是否發(fā)生碰撞。如果兩顆牙齒沒(méi)有發(fā)生碰撞，則ci等于0；否則ci不為0。因此使得約束條件g3中的fc等于0，則滿足牙齒在矯治過(guò)程中無(wú)碰撞發(fā)生。

2 用NSMPSO算法實(shí)現(xiàn)牙齒正畸路徑規(guī)劃

2.1 粒子群算法

在D維搜索空間中，隨機(jī)初始化z個(gè)粒子，這些粒子可視為優(yōu)化問(wèn)題的可行解，并根據(jù)預(yù)先設(shè)定的適應(yīng)度函數(shù)來(lái)評(píng)價(jià)粒子所在位置的好壞。粒子i的當(dāng)前位置為Xi=(xi1，xi2，…，xiD)，當(dāng)前速度為Vi=(vi1，vi2，…，viD)。在每次迭代中，pi=(pi1，pi2，…，pid，…，piD)表示粒子迄今為止經(jīng)歷的最優(yōu)位置，pg=(pg1，pg2，…，pgd，…，pgD)表示整個(gè)種群迄今為止經(jīng)歷的最優(yōu)位置。粒子的速度和位置的更新公式如下：

其中：下標(biāo)i=1，2，…，z，下標(biāo)d=1，2，…，D，t是當(dāng)前迭代次數(shù)，w是慣性權(quán)重。c1和c2是學(xué)習(xí)因子，r1、r2～U（0，1）的隨機(jī)數(shù)。為了使粒子在搜索空間內(nèi)充分搜索，在迭代過(guò)程中限定粒子速度向量上下限為［vmin，vmax］，位置向量上下限為［xmin，xmax］。

2.2 基于正態(tài)分布的簡(jiǎn)化均值粒子群算法

在D維搜索空間中，粒子均按照一定的速度和方向運(yùn)動(dòng)，因此忽略粒子在運(yùn)動(dòng)過(guò)程中的體積與質(zhì)量，可將其看作一個(gè)點(diǎn)。在連續(xù)解空間過(guò)程中，粒子逐漸向全局最優(yōu)解附近靠攏，可將多維可行解空間看作服從多元正態(tài)分布。在真實(shí)的迭代過(guò)程中，粒子在可行解空間的分布不完全意義上服從正態(tài)分布，相比之下仍存在一定的偏離，但鑒于這一部分的偏差的存在對(duì)整體實(shí)驗(yàn)結(jié)果的影響不大，可忽略不計(jì)，并且絕大多數(shù)粒子的分布是服從正態(tài)分布的，所以引入正態(tài)分布作為真實(shí)粒子在搜索空間分布的近似值是合理可行的。因此，本文引入正態(tài)分布和均值粒子群［15］的思想，結(jié)合簡(jiǎn)化粒子群［16］中的位置項(xiàng)，提出了一種基于正態(tài)分布的簡(jiǎn)化均值粒子群（NSMPSO）算法。改進(jìn)后粒子的位置公式如式（4）所示：

其中w為慣性權(quán)重。式（4）將簡(jiǎn)化粒子群的位置項(xiàng)擴(kuò)展為四項(xiàng)：第一項(xiàng)表示向粒子當(dāng)前位置學(xué)習(xí)項(xiàng)。第二、三項(xiàng)中用粒子個(gè)體最優(yōu)位置、種群全局最優(yōu)位置的線性組合取代粒子的個(gè)體最優(yōu)位置和全局最優(yōu)位置，用慣性權(quán)重調(diào)整影響大小。第四項(xiàng)為正態(tài)分布項(xiàng)，其中，μ表示粒子當(dāng)前位置、粒子個(gè)體最優(yōu)位置及全局最優(yōu)位置三者的平均值，根據(jù)式（5）計(jì)算；σ表示粒子當(dāng)前位置、粒子個(gè)體最優(yōu)位置及全局最優(yōu)位置三者的標(biāo)準(zhǔn)差，根據(jù)式（6）計(jì)算；l1和l2是［0，1］范圍內(nèi)均勻產(chǎn)生的隨機(jī)數(shù)，K是由指數(shù)函數(shù)和三角函數(shù)組成的常數(shù)項(xiàng)，用來(lái)控制標(biāo)準(zhǔn)差的影響程度，根據(jù)式（7）計(jì)算。為了避免出現(xiàn)標(biāo)準(zhǔn)差為0的情況，規(guī)定當(dāng)標(biāo)準(zhǔn)差閾值t小于1E-05時(shí)，標(biāo)準(zhǔn)偏差的值取平均值的一半。綜上所述，本文提出的NSMPSO 算法增強(qiáng)了種群的多樣性，加入正態(tài)分布項(xiàng)進(jìn)一步改善了PSO 算法在解決高維問(wèn)題時(shí)求解精度不高、容易陷入局部最優(yōu)等缺點(diǎn)。

2.3 構(gòu)造適應(yīng)度函數(shù)

構(gòu)造一個(gè)合理的、安全性高的適應(yīng)度函數(shù)不僅可以作為衡量?jī)?yōu)化得出的牙齒正畸移動(dòng)路徑優(yōu)劣程度的標(biāo)準(zhǔn)，同時(shí)也是檢測(cè)本文改進(jìn)的粒子群算法在收斂速度及收斂精度方面是否有所提高的判斷準(zhǔn)則。根據(jù)牙齒正畸移動(dòng)路徑規(guī)劃問(wèn)題的特點(diǎn)，本文對(duì)它的評(píng)價(jià)取決于5 個(gè)方面，分別是平移路徑長(zhǎng)度、旋轉(zhuǎn)角度、是否碰撞以及牙齒在單階段的移動(dòng)、旋轉(zhuǎn)量。

平移路徑長(zhǎng)度評(píng)價(jià)函數(shù)對(duì)應(yīng)目標(biāo)函數(shù)f1，旋轉(zhuǎn)角度評(píng)價(jià)函數(shù)對(duì)應(yīng)目標(biāo)函數(shù)f2，用式（1）中的g1、g2表示對(duì)牙齒在矯治過(guò)程中某一階段移動(dòng)距離和旋轉(zhuǎn)角度的限制約束。本文對(duì)約束條件g3通過(guò)添加常數(shù)項(xiàng)進(jìn)一步改進(jìn)，函數(shù)f3作為評(píng)價(jià)是否發(fā)生碰撞的標(biāo)準(zhǔn)，定義如下：

其中，S表示給定的一個(gè)較大常數(shù)項(xiàng)，其作用是當(dāng)相鄰的牙齒發(fā)生碰撞時(shí)，較大的常數(shù)S可以使得整個(gè)f3的值明顯增大。如果兩顆牙齒沒(méi)有發(fā)生碰撞，則f3=0。

綜上，結(jié)合f1、f2、f3、g1、g2定義本文算法的適應(yīng)度函數(shù)為：

其中，i=1，2，…，z；δ、η、λ、μ、θ分別為路徑長(zhǎng)度、旋轉(zhuǎn)角度、是否碰撞以及牙齒在單階段內(nèi)的移動(dòng)、旋轉(zhuǎn)量在適應(yīng)度函數(shù)中的權(quán)重，且為大于0 的實(shí)數(shù)。式（9）中f1和f2是基于目標(biāo)函數(shù)1 和2 構(gòu)建的，設(shè)置其權(quán)重略大于f3、g1、g2的權(quán)重。分析適應(yīng)度函數(shù)Fi可知，函數(shù)值越小說(shuō)明得到的牙齒正畸路徑越優(yōu)。

2.4 NSMPSO算法在路徑規(guī)劃中的實(shí)現(xiàn)

2.4.1 解的編碼問(wèn)題

PSO 算法采用實(shí)數(shù)編碼，一個(gè)粒子代表著所有牙齒在全部階段內(nèi)的矯治結(jié)果，也就是說(shuō)每個(gè)粒子對(duì)應(yīng)著一個(gè)可行解。由于每個(gè)粒子的位置由每顆牙齒在每個(gè)矯治階段中平移所確定的坐標(biāo)(xij，yij，zij)和旋轉(zhuǎn)所確定的(αij，βij，γij)組成。牙齒個(gè)數(shù)為n，矯治階段數(shù)為m，因此粒子的位置是6×n×m維變量。這樣，粒子的編碼結(jié)構(gòu)如下所示：

2.4.2 算法實(shí)現(xiàn)的具體步驟

步驟1 獲取患者牙齒數(shù)據(jù)，提取牙齒特征點(diǎn)作為初始位置坐標(biāo)，本文采用β函數(shù)來(lái)擬合理想牙弓曲線，從而確定牙齒矯正結(jié)束后的最終位置坐標(biāo)［17］，公式如下：

其中：W是患者兩個(gè)磨牙特征點(diǎn)之間的距離，H是患者中切牙特征點(diǎn)到兩顆磨牙特征點(diǎn)之間連線的最短距離，e是將患者尖牙的特征點(diǎn)坐標(biāo)代入式中所確定。

步驟2 按照2.4.1節(jié)對(duì)粒子進(jìn)行編碼，初始化種群并設(shè)定相關(guān)參數(shù)，其中包括患者牙齒個(gè)數(shù)n，需要矯治階段數(shù)m，種群大小z，每個(gè)粒子維數(shù)D，最大迭代次數(shù)Tmax。初始化粒子的位置，公式如下：

其中，r為［0，1］區(qū)間的隨機(jī)數(shù)。

步驟3 根據(jù)式（9）計(jì)算每個(gè)粒子的適應(yīng)度值，并更新粒子的個(gè)體最優(yōu)值和全局最優(yōu)值。

步驟4 根據(jù)式（5）計(jì)算粒子當(dāng)前位置、粒子個(gè)體最優(yōu)位置及全局最優(yōu)位置的平均值，根據(jù)式（6）計(jì)算三者的標(biāo)準(zhǔn)差，根據(jù)式（7）計(jì)算常數(shù)項(xiàng)K。添加正態(tài)分布項(xiàng)，對(duì)種群中的每個(gè)粒子根據(jù)式（4）進(jìn)行位置更新。

步驟5 判斷算法迭代次數(shù)是否達(dá)到所設(shè)定的最大迭代次數(shù)，若是，停止迭代，輸出牙齒正畸移動(dòng)最優(yōu)路徑；否則，返回步驟3）繼續(xù)執(zhí)行。

算法流程如圖3所示。

圖3 NSMPSO算法流程Fig.3 Flowchart of NSMPSO algorithm

3 實(shí)驗(yàn)環(huán)境及分析

3.1 測(cè)試函數(shù)

為了驗(yàn)證本文提出的NSMPSO 算法在收斂精度及收斂速度方面上的提升，將NSMPSO 算法與基本粒子群優(yōu)化（PSO）、均值粒子群優(yōu)化（Mean Particle Swarm Optimization，MPSO）和動(dòng)態(tài)調(diào)整慣性權(quán)重的簡(jiǎn)化均值粒子群優(yōu)化（Simplified Mean Particle Swarm Optimization with Dynamic adjustment of inertia weight，DSMPSO）算法進(jìn)行對(duì)比。本實(shí)驗(yàn)選取三大基準(zhǔn)函數(shù)作為測(cè)試函數(shù)，其中：Sphere 函數(shù)是單峰函數(shù)，主要用于檢驗(yàn)算法的收斂速度和求解精度；Ackley 函數(shù)和Griewank 函數(shù)是多峰函數(shù)，主要用于檢驗(yàn)算法的全局搜索能力。測(cè)試函數(shù)描述如下。

1）Sphere測(cè)試函數(shù)：

其中，搜索區(qū)間為［-100，100］，當(dāng)f1(x)=0時(shí)，Sphere 測(cè)試函數(shù)取得全局最優(yōu)值。

2）Ackley測(cè)試函數(shù)：

其中，搜索區(qū)間為［-32，32］，當(dāng)f2(x)=0時(shí)，Ackley 測(cè)試函數(shù)取得全局最優(yōu)值。

3）Griewank測(cè)試函數(shù)：

其中，搜索區(qū)間為［-600，600］，當(dāng)f3(x)=0時(shí)，Griewank 測(cè)試函數(shù)取得全局最優(yōu)值。

3.2 參數(shù)設(shè)置

在實(shí)驗(yàn)中，測(cè)試函數(shù)的維數(shù)D=30，最大迭代次數(shù)genmax=150，種群粒子數(shù)m=30。對(duì)于PSO 算法和MPSO 算法，w=0.8、c1=c2=2。對(duì)于DSMPSO算法，c1=c2=2，wmax=0.9、wmin=0.4、f1(x)==0.1、a=1、b=2。對(duì)于NSMPSO 算法，w?。?.4，0.9］區(qū)間的隨機(jī)數(shù)，c1=c2=1.494。

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

3.3.1 收斂性分析

本實(shí)驗(yàn)為了驗(yàn)證NSMPSO 算法的收斂性，分別與PSO、MPSO、DSMPSO 算法在三大測(cè)試函數(shù)上進(jìn)行對(duì)比分析。四種算法迭代150次的收斂效果如圖4所示。

從圖4中可以看出，對(duì)于單峰Sphere 函數(shù)，基本PSO 算法由于全局搜索能力較弱，陷入局部最優(yōu)導(dǎo)致算法過(guò)早停滯，在迭代150 次后未收斂，而MPSO 算法收斂速度較慢，DSMPSO、NSMPSO 算法均能較快找到全局最優(yōu)解，相比之下，NSMPSO算法適應(yīng)值小于其他三種算法，求解精度相對(duì)較高。對(duì)于多峰Ackley 函數(shù)，基本PSO 算法由于難以處理復(fù)雜的多峰函數(shù)導(dǎo)致算法停滯，MPSO 算法由于跳出局部最優(yōu)能力較弱，經(jīng)過(guò)140 次迭代后適應(yīng)值開(kāi)始趨于穩(wěn)定，相比其他兩種改進(jìn)算法，收斂速度過(guò)慢，而DSMPSO 和NSMPSO 算法均能在50 次內(nèi)到達(dá)全局最優(yōu)值或者達(dá)到最優(yōu)值附近可接受的范圍內(nèi)，分別在40 次、25 次左右開(kāi)始趨于收斂。對(duì)于多峰Griewank 函數(shù)，基本PSO 算法在迭代150 次過(guò)后，整體適應(yīng)值變化不大，沒(méi)有達(dá)到收斂精度。MPSO 算法在解決復(fù)雜多峰函數(shù)時(shí)穩(wěn)定性不強(qiáng)，DSMPSO 算法以及NSMPSO 算法分別在迭代20 次、15 次左右后適應(yīng)值穩(wěn)定收斂，從圖4（c）可知，NSMPSO 算法找到全局最優(yōu)解所需迭代次數(shù)最少，收斂速度最快。

綜上所述，本文提出的基于正態(tài)分布的簡(jiǎn)化均值粒子群算法（NSMPSO）無(wú)論對(duì)于單峰函數(shù)還是對(duì)于多峰函數(shù)來(lái)說(shuō)，在尋找最優(yōu)解過(guò)程中均能有效地避免種群陷入局部最優(yōu)解，并且能夠在較少的迭代次數(shù)內(nèi)達(dá)到穩(wěn)定收斂，比其他PSO 算法具有更快的收斂速度，在求解精度方面也較為出色，全局搜索能力更佳。

圖4 三個(gè)函數(shù)收斂曲線對(duì)比Fig.4 Convergence curve comparison of three functions

3.3.2 有效性分析

為了進(jìn)一步驗(yàn)證本文提出的NSMPSO 算法的有效性，選取某患者下牙頜14 顆牙齒作為矯治對(duì)象，選取牙齒特征點(diǎn)作為初始位置，根據(jù)理想牙弓曲線確定牙齒矯治結(jié)束時(shí)的位置。利用NSMPSO 算法在約束條件g1、g2、g3的約束下對(duì)目標(biāo)函數(shù)f1和f2進(jìn)行優(yōu)化，從而得到每顆牙齒在每個(gè)階段時(shí)的坐標(biāo)位置，然后在Matlab上進(jìn)行仿真實(shí)驗(yàn)。

首先，牙齒在矯治過(guò)程中垂直方向上的移動(dòng)量較小，因此本實(shí)驗(yàn)忽略牙齒在垂直方向上的移動(dòng)，將原本為一條三維空間的牙弓曲線投影到由患者兩顆磨牙、兩顆尖牙以及中切牙所確定的牙頜平面上，從而簡(jiǎn)化成一條二維平面上的曲線。其次，患者在矯治過(guò)程中所有牙齒是同時(shí)進(jìn)行移動(dòng)、旋轉(zhuǎn)的，并且需要確保過(guò)程中牙齒間不能發(fā)生碰撞、擠壓，否則容易導(dǎo)致矯治失敗，這也是將碰撞檢測(cè)作為牙齒正畸路徑規(guī)劃的重要原因之一。本實(shí)驗(yàn)利用OBB 進(jìn)行碰撞檢測(cè)，在二維平面上用包圍框代替牙齒可以更直觀、形象地表示牙齒移動(dòng)的可行性，確保牙齒正畸的每一階段都是合理移動(dòng)的。最后，以理想牙弓曲線作為牙齒正畸的標(biāo)準(zhǔn)，牙齒正畸路徑規(guī)劃的結(jié)果如圖5所示。

圖5 牙齒初始位置與牙齒矯正結(jié)束位置對(duì)比Fig.5 Comparison of initial position of teeth and position of teeth after orthodontics

圖5 表示牙齒初始位置與牙齒矯正結(jié)束位置對(duì)比，其中實(shí)線包圍框表示牙齒初始位置，虛線包圍框表示牙齒矯正結(jié)束位置，包圍框的中心點(diǎn)即為每顆牙齒的特征點(diǎn)。為了方便敘述，從左到右依次對(duì)牙齒進(jìn)行編號(hào)為1～14，可以看出該患者牙齒排列不齊，里出外進(jìn)，在一塊狹小的區(qū)域內(nèi)多個(gè)牙齒相互重疊遮擋。牙齒正畸路徑規(guī)劃過(guò)程如圖6 所示，其中，牙齒從初始位置到矯正結(jié)束一共需要7 步，每一階段矯正結(jié)束時(shí)的位置在圖6中表示，由于1 號(hào)和14 號(hào)第一磨牙、2 號(hào)和13 號(hào)第二磨牙相對(duì)移動(dòng)、旋轉(zhuǎn)較小，因此在矯治第三階段時(shí)第一磨牙、第二磨牙已經(jīng)達(dá)到理想位置，其他牙齒在經(jīng)過(guò)七個(gè)階段的矯治后，牙齒之間排列整齊，沒(méi)有縫隙，10 號(hào)右側(cè)尖牙和4 號(hào)前磨牙、3 號(hào)前磨牙沒(méi)有之前被擠出的跡象，并且所優(yōu)化出來(lái)的結(jié)果在牙齒矯治過(guò)程中，各個(gè)相鄰的包圍框沒(méi)有發(fā)生重疊部分，意味著牙齒之間沒(méi)有發(fā)生碰撞，說(shuō)明本文提出的改進(jìn)粒子群算法NSMPSO應(yīng)用到牙齒正畸路徑規(guī)劃是合理可行的。

圖6 牙齒正畸路徑規(guī)劃結(jié)果Fig.6 Orthodontic path planning results

圖7 表示四種算法迭代150 次的牙齒正畸路徑規(guī)劃收斂效果對(duì)比?？梢钥闯觯捎肞SO 和MPSO 算法求解牙齒正畸移動(dòng)最優(yōu)路徑時(shí)效果不佳，很難找到全局最優(yōu)解。與DSMPSO 相比，本文提出的NSMPSO 算法有更快的牙齒正畸路徑規(guī)劃收斂速度，在迭代45 次左右算法開(kāi)始收斂，表示NSMPSO 算法可以規(guī)劃出一條有效的牙齒正畸移動(dòng)最短無(wú)碰撞路徑。

圖7 牙齒正畸路徑規(guī)劃收斂圖Fig.7 Orthodontic path planning convergence diagram

4 結(jié)語(yǔ)

針對(duì)虛擬口腔正畸治療系統(tǒng)中牙齒移動(dòng)路徑規(guī)劃問(wèn)題，本文提出了一種基于NSMPSO 算法的牙齒正畸路徑規(guī)劃方法。首先根據(jù)牙齒運(yùn)動(dòng)特性，將問(wèn)題轉(zhuǎn)化成帶約束的多目標(biāo)優(yōu)化問(wèn)題；其次，本文在簡(jiǎn)化粒子群算法的基礎(chǔ)上，引入正態(tài)分布及均值粒子群的思想，提出了基于正態(tài)分布的簡(jiǎn)化均值粒子群算法（NSMPSO）；最后，構(gòu)造了包含平移路徑長(zhǎng)度、旋轉(zhuǎn)角度、碰撞檢測(cè)以及牙齒在單階段的移動(dòng)量、旋轉(zhuǎn)量的適應(yīng)度函數(shù)，得到牙齒從開(kāi)始位置到最終位置的最優(yōu)路徑解。在三大測(cè)試函數(shù)上進(jìn)行對(duì)比實(shí)驗(yàn)，驗(yàn)證了提出的NSMPSO 算法與其他粒子群算法相比具有最高的收斂速度及最高的求解精度。同時(shí)在Matlab中進(jìn)行仿真實(shí)驗(yàn)表明，利用該數(shù)學(xué)模型和改進(jìn)算法求得的最優(yōu)路徑安全可靠，可以為醫(yī)生輔助診斷提供有效的幫助，具有一定的實(shí)用價(jià)值。