盧九星 韓德民 張 羅

(首都醫(yī)科大學(xué)附屬北京同仁醫(yī)院耳鼻喉頭頸外科北京市耳鼻咽喉科研究所教育部耳鼻咽喉頭頸科學(xué)重點(diǎn)實(shí)驗(yàn)室，北京100730)

鼻腔位于呼吸道的首端，功能有通氣、過(guò)濾、清潔、加溫、加濕、共鳴、反射及嗅覺(jué)，無(wú)一不依賴于呼吸氣流通過(guò)鼻腔而加以實(shí)現(xiàn)。鼻腔鼻竇結(jié)構(gòu)復(fù)雜，在人類進(jìn)化過(guò)程中逐漸形成，具有物種特異性，并非簡(jiǎn)單的氣流通道。以呼吸氣流運(yùn)動(dòng)為切入點(diǎn)，找出氣流運(yùn)動(dòng)與鼻腔解剖結(jié)構(gòu)、功能狀態(tài)之間的規(guī)律性關(guān)系，進(jìn)而幫助評(píng)估鼻腔生理、病理狀態(tài)，對(duì)上呼吸道的功能判斷、疾病預(yù)防以及矯正手術(shù)的量化分析均有一定價(jià)值。

鼻腔空氣動(dòng)力學(xué)興起于20世紀(jì)80年代［1］，目前的主要研究方法是數(shù)值模擬，利用有限元的方法對(duì)上呼吸道進(jìn)行氣流場(chǎng)分析、計(jì)算和預(yù)測(cè)。Tan等［2］根據(jù)20例鼻竇計(jì)算機(jī)斷層掃描(computed tomography，CT)圖像，三維重建后對(duì)氣流場(chǎng)進(jìn)行網(wǎng)格劃分，求解納維葉-斯托克斯(Navier-Stokes)方程，得到中國(guó)人正常鼻腔在吸氣、呼氣2個(gè)時(shí)相，穩(wěn)態(tài)(氣流速率不隨時(shí)間變化)條件下的氣流運(yùn)動(dòng)參數(shù)。Liu等［3］采用相似的方法，分析鼻中隔偏曲對(duì)氣流的影響，發(fā)現(xiàn)偏曲位置會(huì)影響氣流，若合并下鼻甲代償性肥大，偏曲側(cè)之凹面氣流會(huì)出現(xiàn)明顯變化。鼻腔空氣動(dòng)力學(xué)在不侵入氣道的情況下得到氣流場(chǎng)的全部信息，具有任何器械檢查不可比擬的優(yōu)勢(shì)，但人為控制邊界條件，使得模擬結(jié)果只具有參考性而不能反映真實(shí)情況，例如穩(wěn)態(tài)條件的設(shè)定就簡(jiǎn)化了呼吸過(guò)程。且據(jù)目前文獻(xiàn)［1-3］，往往一個(gè)樣本都設(shè)置相同的邊界條件，雖然如此設(shè)置可以充分顯示不同解剖結(jié)構(gòu)對(duì)氣流的影響，但卻不能個(gè)體化反映當(dāng)時(shí)的氣流情況。為了逐漸將鼻腔空氣動(dòng)力學(xué)應(yīng)用于臨床，將氣流運(yùn)動(dòng)情況加入評(píng)估體系，本研究綜合利用CT、鼻阻力儀和鼻聲反射儀3種檢查手段，進(jìn)行上呼吸道氣流場(chǎng)的數(shù)值模擬，力圖對(duì)受試者進(jìn)行個(gè)體化測(cè)量，真實(shí)反映氣流場(chǎng)情況。

1 資料與方法

1.1 研究對(duì)象

受試者10名，漢族，19～21歲，平均年齡(19.5±10.97)歲。其中男性5名，平均(19.6±1.68)歲;女性5名，平均(19.4±2.32)歲。受試者身體質(zhì)量指數(shù)(body mass index，BMI)在18～25之間，平均21.31±0.97，其中男性平均21.67±1.11，女性20.96±0.87。身體健康，五官端正，無(wú)先天性顱面部畸形，齒列整齊，無(wú)急慢性鼻炎、鼻竇炎、鼻外傷、鼻出血、鼻囊腫、明顯鼻中隔偏曲穿孔、鼻塞、鼻異物、鼻畸形等，近3月未患急性鼻科疾病、無(wú)鼻腔局部用藥史、無(wú)鼻堵。試驗(yàn)符合人體試驗(yàn)倫理學(xué)標(biāo)準(zhǔn)，受試者嚴(yán)格簽署知情同意書(shū)。

1.2 醫(yī)學(xué)檢查

1)鼻功能檢查:鼻阻力儀(NR6，英國(guó)GM公司)、鼻聲反射儀(A1，英國(guó) GM公司)。測(cè)試房間室溫20℃ ～25℃，相對(duì)濕度50% ～60%，受試者靜坐20 min后端坐位進(jìn)行檢查。鼻噴減充血?jiǎng)┖笤僮?次鼻阻力、鼻聲反射檢查，每次2項(xiàng)檢查在5 min內(nèi)完成。

2)鼻竇CT掃描:鼻功能檢查后30 min內(nèi)行CT(Brilliance 64，荷蘭飛利浦公司)檢查。掃描上呼吸道，層厚0.67 mm，軟組織窗，向上包括額竇、向下包括甲狀軟骨，斷層圖像數(shù)據(jù)的保存格式為數(shù)字影像和通信標(biāo)準(zhǔn)(digital imaging and communications in medicine，Dicom)。

1.3 鼻腔空氣動(dòng)力學(xué)分析

將CT斷層圖像讀入Mimics10.0(比利時(shí)Materialise集團(tuán))軟件，建立三維模型，以Stl(Stereo Lithography)文件格式導(dǎo)出(圖1)。

圖1 Mimics軟件界面Fig.1 Mimics software interface

采用Icem-cfd(Ansys12.1，美國(guó)Ansys公司)軟件劃分流場(chǎng)網(wǎng)格。讀入Stl模型文件，在面部前方設(shè)置外流場(chǎng)，與呼吸面罩近似。整個(gè)上呼吸道流場(chǎng)劃分非結(jié)構(gòu)性網(wǎng)格，八叉樹(shù)(Octree/Robust)方法，呼吸道網(wǎng)格最大邊長(zhǎng)1 mm，平滑網(wǎng)格至98%網(wǎng)格質(zhì)量大于0.6 (圖2)。

圖2 Icem-cfd軟件界面Fig.2 Icem-cfd software interface

網(wǎng)格文件讀入Fluent(Ansys12.1，美國(guó)Ansys公司)軟件。視氣流為連續(xù)介質(zhì)。外流場(chǎng)設(shè)為入口，取值標(biāo)準(zhǔn)大氣壓;呼吸道末端為出口，瞬態(tài)條件，代入氣流速率曲線方程(圖3);鼻腔壁無(wú)滑移邊界條件。采用剪切壓力傳輸(shear stress transport，SST)K-ω 模型［4］計(jì)算呼吸氣流的三維流動(dòng)，收斂殘差0.001。

連續(xù)方程

動(dòng)量方程

湍流能量方程

能量方程

標(biāo)準(zhǔn)K-ε方程轉(zhuǎn)換為

SST 由此可得，φ1 為 K-ω，φ2 為 K-ε

圖3 呼吸模式圖Fig.3 Patterning breathing，hypotheses for asymmetric sine curve

運(yùn)算結(jié)果進(jìn)入Cfd-post(Ansys12.1，美國(guó)Ansys公司)軟件行三維可視化處理，得到1個(gè)完整呼吸過(guò)程中，氣流在整個(gè)上呼吸道中的流動(dòng)物理參數(shù)(氣壓、速率、壁面剪切力)，及在整體和局部的詳細(xì)分布圖像。1個(gè)呼吸周期中，近似平均取9個(gè)時(shí)間點(diǎn)，順序說(shuō)明呼吸過(guò)程中氣流的變化(圖3);每例呼吸道，從鼻孔前至鼻咽，依據(jù)解剖特征取9個(gè)截面，說(shuō)明流經(jīng)時(shí)的氣流狀態(tài)(圖4)。

1.4 統(tǒng)計(jì)學(xué)方法

應(yīng)用醫(yī)學(xué)統(tǒng)計(jì)分析軟件SPSS 17.0(美國(guó)SPSS公司)對(duì)氣流運(yùn)動(dòng)參數(shù)進(jìn)行統(tǒng)計(jì)。計(jì)量數(shù)據(jù)用均數(shù)±標(biāo)準(zhǔn)差(±s)表示，配對(duì)數(shù)據(jù)比較用配對(duì)t檢驗(yàn)，數(shù)值模擬參數(shù)和儀器實(shí)際測(cè)量值比較采用t檢驗(yàn)。以P＜0.05為差異具有統(tǒng)計(jì)學(xué)意義。

圖4 上呼吸道截面示意圖Fig.4 The schematic section of upper respiratory tract

2 結(jié)果

2.1 驗(yàn)證有限元模型

利用鼻聲反射儀測(cè)量受試者的鼻腔最小截面積和其距前鼻孔的距離，0～5 cm、2～5 cm、0～7 cm的鼻腔容積［6］，同一受試者兩側(cè)鼻腔形態(tài)差異無(wú)統(tǒng)計(jì)學(xué)意義(P＞0.05)。同時(shí)測(cè)量有限元模型的鼻閾截面積(圖4，C截面)、C截面距前鼻孔的距離，距前鼻孔5 cm、2～5 cm之間、7 cm的模型體積。經(jīng)比較兩者之間的差異無(wú)統(tǒng)計(jì)學(xué)意義(均P＞0.05)(表1)。

表1 受試者的鼻腔形態(tài)數(shù)據(jù)Tab.1 The data of volunteers’nasal cavity

2.2 驗(yàn)證氣流計(jì)算結(jié)果

利用鼻阻力儀，測(cè)量男性、女性受試者在鼻壓差150 Pa、75 Pa時(shí)的鼻腔單側(cè)阻力值。同一受試者兩側(cè)鼻腔阻力之間的差異無(wú)統(tǒng)計(jì)學(xué)意義(P＞0.05)。鼻阻力取決于鼻腔的氣流量和氣流產(chǎn)生的壓力，鼻阻力=壓力差/流量。

在鼻前庭處(圖4，B截面)、后鼻孔處(圖4，H截面)取1個(gè)呼吸周期的壓力值，制成曲線，求出經(jīng)鼻壓差，再與流量相比得到鼻阻力值。計(jì)算得到的鼻阻力值和實(shí)際測(cè)量得到的鼻阻力值相比較，2者間的差異無(wú)統(tǒng)計(jì)學(xué)意義(P＞0.05)(表2)。

2.3 氣流的運(yùn)動(dòng)情況

呼吸過(guò)程中，氣流主要流經(jīng)總鼻道和中、下鼻道，以層流為主，平靜呼吸時(shí)未形成明顯渦流，在呼吸道邊緣，如鼻前庭頂、底、側(cè)部，嗅裂、鼻底、鼻腔后上部、鼻甲表面，咽鼓管圓枕前、后等多處散在有小的反向流動(dòng)或渦流形成，部位、大小不定。

表2 受試者鼻腔阻力Tab.2 Volunteers’nasal resistance (Pa·cm-3·s-1)

將呼吸氣流在各截面處的壓強(qiáng)、速率依時(shí)間順序繪成曲線，詳見(jiàn)圖5。A～Ⅰ分別對(duì)應(yīng)鼻孔前、鼻前庭、鼻閾至鼻咽中部9個(gè)截面(圖5)。A為氣流進(jìn)入、呼出鼻孔的壓強(qiáng)速率變化。B至G為左右鼻腔內(nèi)氣流的壓強(qiáng)速率變化。整個(gè)呼吸周期內(nèi)，鼻前庭(B)的氣流壓力都為負(fù)壓。固有鼻腔(C～G)內(nèi)的氣流壓強(qiáng)在吸入時(shí)相為負(fù);呼出時(shí)相為正，且形成明顯的平臺(tái)期。呼吸道內(nèi)的氣流隨著深入體內(nèi)而逐漸平緩，速率變小。

呼吸周期內(nèi)Ⅰ～Ⅸ時(shí)刻，氣流在呼吸道各截面處的壓強(qiáng)速率變化，詳見(jiàn)圖6。Ⅴ是吸呼轉(zhuǎn)換時(shí)相，之前是吸氣時(shí)相，之后是呼氣時(shí)相。吸氣時(shí)相(Ⅰ～Ⅳ)，氣流速率進(jìn)入鼻腔后逐漸減小;壓力為負(fù)，接近體表而不斷變小。呼氣時(shí)相，固有鼻腔內(nèi)的氣流壓力為正。

圖5 一個(gè)呼吸周期內(nèi)，呼吸道A～I(xiàn)截面處的壓強(qiáng)、速率曲線Fig.5 One respiratory cycle，the pressure，velocity curve in the respiratory tract section A～Ⅰ

3 討論

鼻腔空氣動(dòng)力學(xué)可以排除儀器檢查中的各種不確定干擾因素，捕捉現(xiàn)有檢查手段難以觀測(cè)到的現(xiàn)象;其價(jià)值突出表現(xiàn)在能描述呼吸道內(nèi)任意部位的氣流運(yùn)動(dòng)。因此，驗(yàn)證數(shù)值模擬的真實(shí)性、結(jié)果的可信性極其必要，否則鼻腔空氣動(dòng)力學(xué)無(wú)異于美輪美奐的空中樓閣。

圖6 Ⅰ～Ⅸ時(shí)刻，氣流流經(jīng)呼吸道A～I(xiàn)截面時(shí)的壓強(qiáng)、速率曲線Fig.6 Time pointⅠ ～Ⅸ，the pressure，velocity curve in the respiratory tract section A ～ Ⅰ

本研究以Dicom連續(xù)數(shù)據(jù)集為圖像基礎(chǔ)，采用人-機(jī)交互的方式將上呼吸道和面部從背景圖像中分割出來(lái)。圖像分割技術(shù)是三維重建的基石［7］，目前的技術(shù)只能逼近而無(wú)法完全反映真實(shí)的解剖結(jié)構(gòu)，對(duì)特定的解剖結(jié)構(gòu)，計(jì)算機(jī)自動(dòng)分割也很難達(dá)到視覺(jué)思維的準(zhǔn)確性。因此，全自動(dòng)圖像分割很難實(shí)現(xiàn)，人-機(jī)合作十分必要，主要體現(xiàn)在自動(dòng)分割算法的選擇和結(jié)果的修正上。圖像分割的準(zhǔn)確性與操作者密切相關(guān)，結(jié)果具有不可重復(fù)性。本研究為保持一致性，盡量采用計(jì)算機(jī)自動(dòng)分割。Mimics軟件生成的三維模型是用三角面片(Stl格式)來(lái)描述的，一例呼吸道三維模型的三角面片數(shù)量平均2×105左右，局部細(xì)節(jié)的修整可以通過(guò)刪改三角面片來(lái)完成，改動(dòng)誤差在10-6m2。

在Icem-cfd軟件中對(duì)三維模型進(jìn)行網(wǎng)格劃分，四面體網(wǎng)格最大邊長(zhǎng)設(shè)為1 mm，能夠捕捉呼吸道大于1 mm的解剖特征。通過(guò)鼻聲反射儀的檢測(cè)數(shù)據(jù)，驗(yàn)證模型的形態(tài)準(zhǔn)確性。本研究得到的有限元模型在形態(tài)上能較真實(shí)地反映上呼吸道的解剖結(jié)構(gòu)特征，與實(shí)際的結(jié)構(gòu)輪廓有很好的相關(guān)性，應(yīng)用此模型對(duì)呼吸氣流進(jìn)行數(shù)值模擬，其結(jié)果能夠較好的符合實(shí)際。

鼻孔入口是一個(gè)曲面，很難進(jìn)行準(zhǔn)確的圖像分割，外鼻及周圍面部形態(tài)不僅影響呼吸氣流進(jìn)入鼻孔的方式，還會(huì)影響氣流在呼吸道內(nèi)的運(yùn)動(dòng)情況［8］。因此，本研究在鼻孔前加設(shè)外流場(chǎng)，形似呼吸面罩，直徑9 cm左右，邊界條件取一個(gè)大氣壓。如此設(shè)置，呼吸氣流呈自然狀態(tài)進(jìn)入鼻孔;因受試者鼻周形態(tài)不同，而對(duì)氣流的不同影響也都會(huì)被保留。

本研究比較由鼻阻力儀實(shí)際測(cè)得的鼻阻力值和由數(shù)據(jù)模擬結(jié)果計(jì)算得來(lái)的鼻阻力值，來(lái)驗(yàn)證數(shù)值模擬結(jié)果是否合理。男女受試者測(cè)定鼻阻力值時(shí)，佩戴面罩，呼吸力度較弱，潮氣量較小，得到的氣流平均流速也較小，由此計(jì)算出的氣流速率曲線方程偏于正常值范圍下游。據(jù)文獻(xiàn)［4］，實(shí)際測(cè)量的呼吸速率方程曲線近似矩形，呼氣時(shí)相長(zhǎng)于吸氣時(shí)相。類似研究中，呼吸模式多設(shè)定為兩時(shí)相相等或不相等的三角形［9］、正弦曲線［10］。本研究將呼吸曲線設(shè)置為呼吸時(shí)相/吸氣時(shí)相=1.2/1，曲線下面積相等的半正弦曲線(圖3)，較為符合真實(shí)情況。

模型和數(shù)值模擬得到驗(yàn)證后，就可系統(tǒng)分析呼吸氣流的運(yùn)動(dòng)情況。本研究選擇的是正常鼻腔，對(duì)稱的結(jié)構(gòu)產(chǎn)生相似的氣流，兩側(cè)的壓力、速率曲線基本伴行。氣流進(jìn)入鼻腔，至鼻閾處壓力陡然變大，占?jí)毫ψ兓档?5% ～63%;速率在鼻閾處也達(dá)到峰值;呼氣時(shí)，氣流壓力下降至鼻閾后略有抬升，在鼻孔前接近大氣壓;氣流速率平穩(wěn)增高，在鼻閾處有一坡值，隨后下降出鼻腔。呼吸周期內(nèi)，壓力、速率曲線明顯成雙峰狀，吸氣時(shí)相壓力、速率曲線的峰值和坡度大于呼氣時(shí)相，說(shuō)明呼氣氣流變化小于吸氣時(shí)相，吸氣時(shí)氣流變化更劇烈。

本研究得到的呼吸道模型在形態(tài)上能真實(shí)的反映上呼吸道的解剖特征，應(yīng)用此模型對(duì)上呼吸道流場(chǎng)進(jìn)行數(shù)值模擬，結(jié)果較好的符合實(shí)際。

鼻腔空氣動(dòng)力學(xué)是生物力學(xué)的一個(gè)分支，主要研究呼吸過(guò)程中，鼻腔在同氣體作相對(duì)運(yùn)動(dòng)情況下的受力特性，氣體流動(dòng)規(guī)律和伴隨發(fā)生的物理化學(xué)變化;是在醫(yī)學(xué)和流體力學(xué)的基礎(chǔ)上，隨著計(jì)算流體動(dòng)力學(xué)的發(fā)展而成長(zhǎng)起來(lái)的一個(gè)學(xué)科。其指導(dǎo)思想貫穿著醫(yī)學(xué)唯物觀:“沒(méi)有功能的結(jié)構(gòu)是沒(méi)有意義的，沒(méi)有結(jié)構(gòu)的功能是不存在的”。鼻腔結(jié)構(gòu)是鼻腔功能的物質(zhì)基礎(chǔ)，鼻腔功能是鼻腔結(jié)構(gòu)的外在表現(xiàn)。結(jié)構(gòu)變異首先導(dǎo)致鼻腔氣流變化，功能隨之改變;反之，氣流變化會(huì)影響鼻腔功能的實(shí)現(xiàn)，同時(shí)也作用于鼻腔結(jié)構(gòu)，并在時(shí)間維度上逐漸積累，最終導(dǎo)致結(jié)構(gòu)的改變［11］。結(jié)構(gòu)和功能兩者矛盾運(yùn)動(dòng)，由量變達(dá)到質(zhì)變，最終產(chǎn)生病理癥狀。“結(jié)構(gòu)-功能-癥狀”三者的辯證關(guān)系體現(xiàn)了鼻腔疾病從發(fā)生到發(fā)展、量變到質(zhì)變、微觀到宏觀的因果關(guān)系。呼吸氣流是切入點(diǎn)，充分研究氣流的運(yùn)動(dòng)規(guī)律，就可掌握鼻腔“結(jié)構(gòu)-功能-癥狀”的辯證關(guān)系。

［1］Leong S C，Chen X B，Lee H P，et al.A review of the implications of computational fluid dynamic studies on nasal airflow and physiology［J］.Rhinology，2010，48(2):139-145.

［2］Tan J，Han D，Wang J，et al.Numerical simulation of normal nasal cavity airflow in Chinese adult:a computational flow dynamics model［J］.Eur Arch Otorhinolaryngol，2012，269(3):881-889.

［3］Liu T，Han D，Wang J，et al.Effects of septal deviation on the airflow characteristics:using computational fluid dynamics models［J］.Acta Otolaryngol，2012，132(3):290-298.

［4］Liu Y，Edgar A M，Matida，et al.Numerical simulation of aerosol deposition in a 3-D human nasal cavity using RANS，RANS/EIM，and LES［J］.J Aerosol Sci，2007，38:683-700.

［5］嚴(yán)碧涯，張有為，郁秀芳，等.我國(guó)健康成人呼吸氣流速率曲線圖的研究［J］.北京醫(yī)學(xué)，1980，2(5):257-260.

［6］Hilberg O，Pedersen O F.Acoustic rhinometry:recommendations for technical specifications and standard operating procedures［J］.Rhinol Suppl，2000，16:3-17.

［7］蘇秀云，劉蜀彬.Mimics 軟件臨床應(yīng)用 計(jì)算機(jī)輔助外科入門技術(shù)［M］.北京:人民軍醫(yī)出版社，2011:73.

［8］Taylor D J，Doorly D J，Schroter R C.Inflow boundary profile prescription for numerical simulation of nasal airflow［J］.J R Soc Interface，2010，7(44):515-527.

［9］Zhang J，Liu Y X，Sun X Z，et al.Computational fluid dynamics simulations of respiratory airflow in human nasal cavity and its characteristic dimension study［J］.Acta Mech Sin，2008，24(2):223-228.

［10］Ishikawa S，Nakayama T，Watanabe M，et al.Visualization of flow resistance in physiological nasal respiration:analysis of velocity and vorticities using numerical simulation［J］.Arch Otolaryngol Head Neck Surg，2006，132(11):1203-1209.

［11］里曉紅，楊克林.慢性鼻竇炎60例變應(yīng)性因素分析［J］.臨床誤診誤治，2010，23(11):1028-1029.