楊龍成,陸繼慶,劉冀成
(成都信息工程學(xué)院,四川成都610225)
有研究表明,經(jīng)顱磁刺激(Transcranial magnetic stimulation,TMS)對人體安全,無痛,無損傷,重復(fù)性好,且療效明顯[1-5]。但由于磁刺激儀的磁刺激線圈對磁場的匯聚能力差、刺激重復(fù)性差、定位不準(zhǔn),且磁刺激技術(shù)本身的缺陷影響了經(jīng)顱磁刺激的進一步發(fā)展。研究者們經(jīng)過研究,設(shè)計出各種不同線圈陣列[6-7],而線圈的各特性參數(shù)決定刺激的準(zhǔn)確性和刺激深度[8],因此經(jīng)顱磁刺激過程中頭部內(nèi)電磁場分布的準(zhǔn)確模擬受線圈參數(shù)的影響,所以要提高電磁場模擬的精度和實現(xiàn)磁場更好的聚焦性,就需要將線圈參數(shù)進行全面優(yōu)化,逐步改善經(jīng)顱磁刺激系統(tǒng)磁聚焦性能。本文設(shè)計新型的簡單經(jīng)顱磁刺激線圈,對影響磁場分布的線圈可調(diào)參數(shù)進行分析,用高級語言編寫的優(yōu)化算法及其與CST軟件的外部通信接口,完成對經(jīng)顱磁刺激的優(yōu)化求解,達到改善磁刺激磁場的聚焦效果。
模型由兩個邊長相等的正方形子線圈組成,分別水平放置在上下兩個不同平面,如圖1。
圖1 兩方框形線圈CST結(jié)構(gòu)圖Fig.1 Two square coils of CST structure graph
兩子線圈電流有4種方案,分兩種情況討論線圈邊長和線圈重疊間距對磁場分布影響。
(1)線圈注入方向相同的1A電流
圖2是幾組線圈參數(shù)下磁感應(yīng)強度分布圖。圖a與b、c與d、e與f兩兩各為一組,線圈邊長依次為100 cm、160 cm、160 cm,線圈重疊間距依次為0.9×100 cm、0.9×160 cm、0.75×160 cm,兩線圈所在平面垂直距離均為3 cm。
由圖2分析線圈邊長和線圈重疊間距對磁場分布的影響為:子線圈邊長越大越好,子線圈重疊間距對磁場的分布影響不大。
圖2 平面線圈模型生成的B場分布Fig.2 Magnetic field distribution established by plane coil model
(2)線圈注入方向相反的1A電流
圖3是不同線圈參數(shù)組合的磁感應(yīng)強度分布圖。圖a與b、c與d中線圈的邊長為80 cm,前組線圈重疊間距0.9×80 cm、后組0.4×80 cm,兩線圈平面垂直距離均為3 cm。
由圖3分析線圈邊長和線圈重疊間距對磁場分布的影響為:線圈邊長越小,磁場分布越好,線圈重疊間距對磁場分布影響也不大;隨著線圈Z方向距離的增加,磁場出現(xiàn)中心凹陷情況。
上述研究表明,子線圈注入反向電流時,磁場分布較好。
通過以上模型的討論和研究,在圖1所示的兩方框形線圈模型基礎(chǔ)上加入第3個線圈,線圈的形狀可為長方形,或其他模型。如圖4,3個子線圈所在平面兩兩間的間距相同。
圖3 平面線圈模型生成的B場分布Fig.3 Magnetic field distribution established by plane coil model
圖4 三方框形線圈CST模型Fig.4 Three square coils of CST model
圖4線圈模型中共有12種電流方案,根據(jù)第三個子線圈是否注入電流,分兩種情況研究磁場。其中有6種模型的聚焦性能比對比線圈好。圖5列出了這6種電流方案,圖a、b、c、d為新加入線圈注入電流的情況,e和f為無電流情況。
圖5 線圈電流組態(tài)Fig.5 The coil current configuration
圖6給出兩組線圈參數(shù)模型的磁場強度分布圖。每兩個線圈所在平面的垂直距離為3 cm。圖a、b中線圈邊長80 cm,左邊線圈+1 A電流,右邊線圈-1 A電流,重疊距離為0.4×80 cm,第三個線圈邊長為32 cm的正方形,電流+5 A;圖c、d中邊長60 cm,左邊線圈-1 A,右邊線圈1 A,重疊距離0.1×60 cm,第三個線圈的長60 cm,寬為42 cm,電流為+1 A。
圖6 平面線圈模型生成的B場分布Fig.6 Magnetic field distribution established by plane coil model
從圖中可知這兩組參數(shù)的磁場分布聚焦都較好,磁場的聚焦性能與上方加入的第三個線圈的電流密切相關(guān),電流越大聚焦越好。相比于其他4種模型此模型的磁場分布整體要好。
從遺傳算法[9-13]和粒子群算法[14-17]的大量研究中可知兩種算法搜索機制的異同之處和各自的優(yōu)缺點,都屬于全局搜索算法。結(jié)合粒子群算法中粒子是基本單位且無需編碼,遺傳算法交叉和變異操作對最優(yōu)解有破壞影響,但粒子群算法可對最優(yōu)解保存,粒子群算法主要處理連續(xù)問題的這些特性實現(xiàn)本文的混合優(yōu)化算法(PSO-GA)。此混合算法是在粒子群算法的每一次迭代后期加入遺傳算法的交叉與變異[18]操作。
混合算法中的變異是高斯變異。高斯表達式如公式(1)所示,式子中g(shù)bestfitness表示粒子全局最優(yōu)解、pbestfitness[i]表示粒子局部最優(yōu)解的適應(yīng)度值、fitness[i]表示粒子適應(yīng)度值。粒子替換公式為公式(2)和(3)所示,公式(2)為迭代次數(shù)在1/2最大次數(shù)前的變異粒子替換公式,公式(3)為迭代次數(shù)大于或等于1/2最大迭代次數(shù)的變異粒子替換公式,xir表示第i個粒子第r維的位置,gbestir表示全局最優(yōu)解第r維位置。
通過人們公認(rèn)的純數(shù)學(xué)測試函數(shù)[19-20]的有效評估,得出PSO-GA在多峰搜索中存在明顯的優(yōu)勢,其性能更適合磁場的優(yōu)化。
線圈模型的優(yōu)化對象主要是線圈的參數(shù)特性,如線圈的半徑,注入線圈的電流及線圈的空間位置等。通常我們根據(jù)計算線圈的先驗知識確定線圈參數(shù)的解空間,將此解空間作為優(yōu)化算法初始種群的范圍及計算精度的依據(jù),然后優(yōu)化算法依據(jù)其制定的初始種群產(chǎn)生方法產(chǎn)生初始種群,CST與高級語言的數(shù)據(jù)通信接口將初始種群所表示的線圈參數(shù)輸入CST,CST軟件根據(jù)所輸入的參數(shù)建立模型,計算模型及進行結(jié)果處理,最后通過數(shù)據(jù)接口將優(yōu)化算法所需的場結(jié)果輸出,用于優(yōu)化算法的適應(yīng)度計算,此過程就為CST軟件與優(yōu)化算法的單次通信,而算法種群個體通常不止一個,所以計算一代種群需進行這樣的通信多次。圖7給出混合優(yōu)化算法與CST通信的聯(lián)合求解整體步驟流程:
圖7 線圈模型優(yōu)化流程圖Fig.7 Optimization coil model for processing
上述提出的線圈陣列模型中3方框形線圈模型結(jié)構(gòu)簡單,磁場分布較好,所以為了得到用于TMS的激勵線圈,同時驗證CST與優(yōu)化算法接口的可用性,對此模型進行了優(yōu)化。
適應(yīng)度函數(shù)反映了目標(biāo)函數(shù)的特性,線圈優(yōu)化中算法的每一個個體都代表一組線圈參數(shù),同時對應(yīng)一種線圈模型,即對應(yīng)一種場分布。本文對磁感應(yīng)強度進行了歸一化處理,防止出現(xiàn)磁場值和磁場分布同比例增長帶來的計算或?qū)Ρ鹊腻e誤。根據(jù)大量的計算經(jīng)驗,本文使用的場值適應(yīng)度函數(shù)如下公式(4)所示。
線圈優(yōu)化計算中,CST計算空間設(shè)置為100 cm×100 cm×100 cm,磁感應(yīng)強度取值平面為Z=-4平面。公式(4)中BarvN表示輸出平面上,輸出點B場歸一化值的平均值,Barvn表示計算平面中心向外延伸4 cm×4 cm方形范圍內(nèi)B場輸出值的歸一化平均值,N表示從CST輸出的B場取值點總數(shù),N0.98表示B場輸出平面磁感應(yīng)強度歸一化值大于等于0.98的點,同理可得N0.9,N0.8。計算適應(yīng)度算法中利用N0.98與計算平面中心周圍10 cm×10 cm范圍內(nèi),歸一化值大于或等于0.98點的總數(shù)比較,將聚焦點限制在計算平面中心周圍6 cm×6 cm平面內(nèi)且去除了多峰情況。此適應(yīng)度函數(shù)反映了實際應(yīng)用所需磁場分布的特性,將磁場問題轉(zhuǎn)換為數(shù)學(xué)模型,用于TMS優(yōu)化計算。
優(yōu)化線圈的匝數(shù)為100匝,電流取+1 A,0 A和-1 A三個值。優(yōu)化參數(shù)如表1所示。
表1 3方框形線圈優(yōu)化參數(shù)表Table1 3 square coils optimization parameter table
按照上述優(yōu)化所得參數(shù),在CST軟件中建立模型,如圖8所示,計算空間與優(yōu)化B場取值空間一致20 cm×20 cm×20 cm,其他計算條件與優(yōu)化計算設(shè)置一致。計算Z=-4平面歸一化B場幅值分布和二維等高線圖(見圖9)。
為了對比優(yōu)化效果,在單線圈尺寸和刺激強度盡量接近的情況下選取邊長為8 cm對比線圈,注入電流大小為1 A。圖10為對比單線圈的歸一化B場幅值分布和二維等高線圖。從圖11可同時看出兩種模型刺激強度和聚焦程度的對比。
圖8 a優(yōu)化模型b單線圈模型Fig.8 a Optimization model b Single coil model
圖9 優(yōu)化模型Fig.9 Optimization model
圖10 單線圈模型Fig.10 Single coil model
圖11 a優(yōu)化模型B場分布b單線圈B場分布Fig.11 Optimization model of B distribution Single coil model of B distribution
上述對比結(jié)果顯示,優(yōu)化模型無論在聚焦性還是刺激強度上都優(yōu)于對比單線圈模型,因此該模型可用于TMS中。
結(jié)合遺傳算法與粒子群算法的基本特性設(shè)計了混合優(yōu)化算法;用高級語言編程控制CST軟件,實現(xiàn)CST軟件與高級語言程序的數(shù)據(jù)通信接口及聯(lián)合CST軟件的TMS優(yōu)化求解;對當(dāng)前新型的線圈陣列模型中影響磁場分布的參數(shù)進行了討論,用TMS優(yōu)化算法優(yōu)化了其中磁場分布集中且結(jié)構(gòu)簡單的模型,優(yōu)化結(jié)果對比顯示,3方框形線圈的刺激強度和聚焦程度都有不同程度的提高,得到用于TMS的優(yōu)化激勵線圈模型,為全面優(yōu)化激勵線圈的空間結(jié)構(gòu)研究提供了參考依據(jù)與實用價值。
References)
[1] Angela Ziluk,Azra Premji,Aimee J.Nelson.Functional connectivity from area 5 to primary motor cortex via paired-pulse transcranial magnetic stimulation[J].Neuroscience Letters,2010,484:81 -85.
[2] 喬清理,王明時,田心.經(jīng)顱磁刺激的原理方法和應(yīng)用[J].中國生物醫(yī)學(xué)工程學(xué)報,2004,21(4):259 -265.
[3] 王小明,龍存國,吳壁華,張國元,李秋茹,趙小瓊,鄭霞清.反復(fù)經(jīng)顱磁刺激安全性的實驗研究[J].中國臨床康復(fù),2003,7(13):1896.
[4] 吳原,巖永書朋,矢野直次.磁刺激對脊髓前角細(xì)胞興奮性的影響[J].臨床神經(jīng)電生理學(xué)雜志,2003,12(1):3.
[5] Jacob Jolij,Victor A.F.Lamme.Transcranial magnetic stimulation-induced‘visual echoes’are generated in early visual cortex[J].Neuroscience Letters,2010,484:178 -181.
[6] 劉冀成,黃卡瑪,胡雅毅,華偉.功能磁刺激線圈陣列設(shè)計與場分布計算[J].航天醫(yī)學(xué)與醫(yī)學(xué)工程,2004,17(5):365 -369.
[7] 劉冀成,黃卡瑪,郭蘭婷,張弘,胡雅毅.基于遺傳算法的經(jīng)顱磁刺激線圈陣列設(shè)計與場分布計算[J].生物醫(yī)學(xué)工程學(xué)雜志,2005,22(2):303 -306.
[8] 鄭建斌,霍小林,吳昌哲,于陽.經(jīng)顱磁刺激系統(tǒng)各主要參數(shù)實現(xiàn)的量化設(shè)計[J].中國臨床康復(fù),2004,8(34):7655 -7657.
[9] 崔姍姍.遺傳算法的一些改進及其應(yīng)用[D].北京中國科技大學(xué),2010.
[10] Holland John.Genetic algorithms[J].Scientific American,1992,4:44-50.
[11]趙云珍.遺傳算法及其改進[D].昆明:昆明理工大學(xué),2005.
[12]呂晉君.遺傳算法的改進及其在優(yōu)化上的應(yīng)用研究[D].太原理工大學(xué),2010.
[13] Altenberg Lee,Thomas Back E D,Kaufmann Morgan.Fitness Distance Correlation analysis:an Instruction counterexamole[C].US:IN Proceedings of the 7th International Conference on Genetic Algorithm,1997.
[14]蔣曉鳴,雷霖,王厚軍.一種改進慣性權(quán)重的變異微粒群優(yōu)化算法[J].計算機技術(shù)與發(fā)展,2008,18(6):79 -82.
[15]劉晶晶.粒子群優(yōu)化算法的改進與應(yīng)用[D].武漢:武漢理工大學(xué),2007.
[16]羅德相.粒子群算法改進方法研究[D].南京:廣西民族大學(xué),2009.
[17]陸克中,張秋華,孫蘭娟.一種改進的粒子群優(yōu)化算法及其仿真[J].計算機技術(shù)與發(fā)展,2007,17(11):88 -91.
[18]劉晶晶,吳傳生.一種帶交叉算子的改進的粒子群優(yōu)化算法[J].青島科技大學(xué)學(xué)報,2008,29(1)∶77 -79.
[19]劉冀成.基于改進遺傳算法的生物電磁成像與磁場聚焦應(yīng)用[D].成都:四川大學(xué),2005.
[20] Thomas Back.Evolutionary Algorithms in Theory and Practice[M].New York:Oxford Unuiversity Press,1996.