李 金 楊鵬飛 王菊芳 何進(jìn)鵬 馬 偉 周 珩?

1(蘭州大學(xué) 核科學(xué)與技術(shù)學(xué)院，蘭州 730000)

2(中國(guó)科學(xué)院近代物理研究所 甘肅省空間輻射生物學(xué)重點(diǎn)實(shí)驗(yàn)室 蘭州 730000)

3(中國(guó)科學(xué)院大學(xué) 核科學(xué)與技術(shù)學(xué)院 北京 100049)

引言

纖毛是一種來(lái)源于細(xì)胞基體并突出于細(xì)胞表面的特殊細(xì)胞器，有3 種類型:初級(jí)纖毛(primary cilia)、運(yùn)動(dòng)纖毛(motile cilia)和節(jié)點(diǎn)纖毛(nodal cilia)[1]。 初級(jí)纖毛作為細(xì)胞的“天線”，通過(guò)調(diào)控Hedgehog[2]、Wnt[3]、TGF－β[4]、Notch[5]等多種信號(hào)通路感應(yīng)物理、化學(xué)和生物信號(hào)的刺激，并傳導(dǎo)至細(xì)胞內(nèi)部。 初級(jí)纖毛的結(jié)構(gòu)異常會(huì)導(dǎo)致多種疾病[6-7]，甚至影響腫瘤的發(fā)生和發(fā)展[7]。 通過(guò)化學(xué)藥物抑制初級(jí)纖毛生長(zhǎng)能夠有效控制腫瘤細(xì)胞的惡性增殖[8]。

目前，多采用Photoshop 軟件中的“標(biāo)度尺”對(duì)初級(jí)纖毛進(jìn)行測(cè)量，但存在準(zhǔn)確度低、重復(fù)性差、分析效率低等問(wèn)題，因此，建立一種基于數(shù)字圖像處理技術(shù)的高效分析方法非常必要。

數(shù)字圖像處理技術(shù)可實(shí)現(xiàn)包括噪聲去除[9]、分割[10]、配準(zhǔn)[11]、融合[12]等圖像處理。 噪聲去除是指消除、弱化或壓制圖像中的隨機(jī)變化干擾。 這些干擾可能源于成像過(guò)程中的隨機(jī)效應(yīng)與測(cè)量誤差。分割是將圖像細(xì)分為它的組成要素，細(xì)分水平取決于要解決的問(wèn)題，即感興趣的物體被隔離出來(lái)時(shí)，就應(yīng)該停止分割。 常見的分割方法有基于邊緣分割、基于區(qū)域分割、基于特定理論方法和基于深度學(xué)習(xí)的分割算法等。

初級(jí)纖毛的數(shù)字圖像處理的難點(diǎn)在于噪聲去除與邊界分割。 首先對(duì)相關(guān)醫(yī)學(xué)圖像降噪分割方法進(jìn)行了調(diào)研。 陳炳權(quán)等[13]提出了一種基于離散小波變換和修正中值濾波的醫(yī)學(xué)圖像耦合去噪算法，用于CT 與超聲波圖像降噪。 Sarungbam Bonny等[14]提出了一種非局部均值濾波法用于超聲圖像去除散斑噪聲，該方法能夠保持目標(biāo)邊緣細(xì)節(jié)和結(jié)構(gòu)信息。 Erwin 等[15]介紹了一種改進(jìn)后的Prewitt算子用于視網(wǎng)膜血管分割，取得了較好的分割效果。 另外，還有一些對(duì)邊緣細(xì)節(jié)保持和分割較好的算法如:中值濾波[13]、均值濾波[14]、Prewitt[15]、Sobel[16]、Roberts[17]、高斯－拉普拉斯(Log)[18]等。本研究將上述算法用于對(duì)初級(jí)纖毛形態(tài)特征的定量分析，對(duì)相應(yīng)結(jié)果進(jìn)行了對(duì)比分析，然后，選取了較優(yōu)算法進(jìn)行了二次改進(jìn)，構(gòu)建了基于較優(yōu)算法的初級(jí)纖毛形態(tài)特征的自動(dòng)分析方法。 該方法能自動(dòng)對(duì)每個(gè)初級(jí)纖毛進(jìn)行編號(hào)，自動(dòng)計(jì)算其周長(zhǎng)和面積等形態(tài)參數(shù)。 從而達(dá)到了客觀、準(zhǔn)確、快速、多元化的處理目的，為初級(jí)纖毛圖像數(shù)字特征分析提供了有效的技術(shù)方案。

1 材料與方法

1.1 圖像采集

首先進(jìn)行圖像采集，其質(zhì)量將直接影響后續(xù)操作，如分割算法的難易以及數(shù)據(jù)結(jié)果和分析效率等。 本研究采用電離輻射后的M059K 細(xì)胞(人腦神經(jīng)膠質(zhì)瘤細(xì)胞)的初級(jí)纖毛圖像作為案例，如圖1所示。 X 射線由PXI Precision X-RAY225 系統(tǒng)產(chǎn)生，電壓為225 kV，電流為13.3 mA，照射劑量為4 Gy。 圖像分辨率為2 048×1 536像素，由激光掃描共聚焦顯微鏡(FLUOVIEW FV3000， Olympus)拍攝，圖片保存為jpg 格式。 像素長(zhǎng)度與真實(shí)長(zhǎng)度間的比例尺為11.209 7 pixel/μm。

圖1 X 射線照射M059K 細(xì)胞后初級(jí)纖毛免疫熒光圖像Fig.1 Immunofluorescence image of primary cilia of M059K cells irradiated with X-ray

制片操作過(guò)程:首先利用4%多聚甲醛對(duì)X 射線照射后的細(xì)胞于室溫固定10 min，再使用含0.5%Triton X－100 的磷酸鹽緩沖液(PBS)透膜5 min，然后在含5%山羊血清的PBS 室溫封閉2 h，之后利用Arl13b(1 ∶500)和γ－Tubulin(1 ∶500)一抗室溫孵育2 h，熒光二抗(1 ∶500)室溫避光孵育2 h，最后使用含0.5% Tween－20 的PBS 洗片，并通過(guò)DAPI(4′， 6－diamidino－2－phenylindole)染色細(xì)胞核后封片、采集制作數(shù)字圖像。

1.2 圖像分離與分割

圖像分析設(shè)備為Intel core i5，3 GHz CPU，8 GB內(nèi)存計(jì)算機(jī)，分析軟件為Matlab 2017。 初級(jí)纖毛原始數(shù)字圖像如圖1 所示，采用RGB 彩色模式，每個(gè)像素點(diǎn)由R、G、B 等3 個(gè)字節(jié)構(gòu)成。 經(jīng)Arl13b(Proteintech，＃ 66739)標(biāo)記染色的初級(jí)纖毛呈現(xiàn)紅色；由γ-Tubulin(Sigma，＃SAB4503045)標(biāo)記的細(xì)胞骨架呈現(xiàn)綠色；細(xì)胞核呈現(xiàn)藍(lán)色。 根據(jù)染色特點(diǎn)，可以利用不同的色彩信息將細(xì)胞核與細(xì)胞骨架從圖像中分離出來(lái)。 纖毛與背景中雜質(zhì)的亮度不同，可依據(jù)灰度信息差異將纖毛從背景雜質(zhì)中分割出來(lái)。 以此引入最大類間方差法進(jìn)行分割。

1.3 圖像色彩分離

RGB 圖像可以看作由3 個(gè)灰度圖像形成的“堆棧”，圖像是一個(gè)M×N×3 大小的彩色像素?cái)?shù)組，當(dāng)發(fā)送到彩色監(jiān)視器的紅、綠、藍(lán)輸入端時(shí)，就會(huì)在屏幕上產(chǎn)生彩色圖像。 通過(guò)對(duì)纖毛圖像色彩分析，只需將三維RGB 圖像提取出二維的R 分量，即可將細(xì)胞核與細(xì)胞骨架分離出來(lái)。

1.4 不同濾波方法對(duì)比

濾波是為了突出圖像的空間信息，壓抑或者去除其他無(wú)關(guān)的信息。 濾波方法包括空間濾波和頻域?yàn)V波，空間濾波又包括線性空間濾波和非線性空間濾波。 線性空間濾波基于計(jì)算乘積和(即線性操作)，非線性空間濾波基于涉及領(lǐng)域像素內(nèi)的非線性操作，例如通過(guò)領(lǐng)域內(nèi)的像素最大值代替每個(gè)中心點(diǎn)的響應(yīng)值。 本研究中的濾波結(jié)果如圖2 所示，初級(jí)纖毛圖像(a)分別經(jīng)均值濾波(b)、Log 濾波(c)、Prewitt 濾波(d)、Sobel 濾波(e)、中值濾波(f)以及改進(jìn)型濾波(g)處理。

圖2 不同濾波方法處理后的纖毛圖像。 (a)原始圖像；(b)均值濾波；(c)Log 濾波；(d)Prewitt 濾波；(e)Sobel 濾波；(f)中值濾波；(g)改進(jìn)型濾波Fig.2 The primary ciliary images processed by different filtering methods.(a) Original image；(b)Mean filtering； (c)Log filtering； (d)Prewitt filtering； (e)Sobel filtering； (f)Median filtering；(g)Improved filtering.

線性空間濾波也稱“空間卷積”。 卷積本身是數(shù)學(xué)領(lǐng)域內(nèi)的一種積分變化方法，在數(shù)字圖像處理技術(shù)中，濾波方法如式(1)，濾波原理如圖3 所示。

式中:x，y是像素在圖片中的位置坐標(biāo)；k，l是卷積核中的位置；m[x，y] 是濾波結(jié)果，即圖3 中的結(jié)果40；f[k，l] 是濾波函數(shù)也稱權(quán)重系數(shù)，即圖示中的卷積核；I[x ＋k，y ＋l] 是與w[k，l] 相對(duì)應(yīng)的圖片像素值。

圖3 圖像濾波方法流程Fig.3 The filtering method flow

綜合不同濾波結(jié)果得出結(jié)論:均值與中值濾波能夠過(guò)濾部分細(xì)胞骨架，但也會(huì)導(dǎo)致初級(jí)纖毛邊界模糊；Sobel 濾波算子能夠在一個(gè)方向上增加纖毛的邊界對(duì)比度；Prewitt 濾波算子能夠增強(qiáng)整個(gè)纖毛的對(duì)比度，但當(dāng)纖毛呈豎直排列時(shí)反而會(huì)減弱對(duì)比度；Log 濾波對(duì)比度沒(méi)有Prewitt 濾波算子明顯。 針對(duì)以上問(wèn)題，將Prewitt 處理后的圖片與原圖相加，或者將Log 濾波后的圖像取反與原圖相加。 經(jīng)過(guò)改進(jìn)后的濾波方法使得初級(jí)纖毛圖像的邊界更加清晰，整體亮度提升明顯。

1.5 OTSU′s 閾值處理方法

OTSU′s 最大類間方差法[19-20]是在最小二乘法原理為基礎(chǔ)上進(jìn)行的推導(dǎo)演變。 基本實(shí)現(xiàn)原理依據(jù)閾值k將圖像分為C1 與C2 兩個(gè)部分，最優(yōu)閾值k的選擇標(biāo)準(zhǔn)是C1 與C2 兩部分之間的最大類間方差。 閾值k選取類似于最小二乘法，k的范圍是C1中0 到L－1 內(nèi)的整數(shù)。 OTSU′s 閾值處理方法具有良好的自適應(yīng)分割性能，具體實(shí)現(xiàn)方法如下:

令一幅圖像的直方圖成分由下式所示:

式中，n是圖像中像素的總數(shù)，ni是灰度級(jí)為i的數(shù)目，L是圖像中所有可能的灰度級(jí)數(shù)。 假設(shè)閾值k已經(jīng)選定，C1 是灰度值從0 到k的像素，C2 是灰度值從k＋1 到L－1 的像素，此處的k均為正整數(shù)。

閾值k由Otsu′s 公式選定，即

式中，) 是方差，P1(k)，P2(k) 分別是像素小于等于閾值集合和像素大于閾值集合發(fā)生概率；m1k和m2k分別是像素小于等于閾值集合和像素大于閾值集合的平均灰度；mG全局均值。

如果設(shè)置k＝0，那么擁有為k賦值的任何像素C1 集合的概率為0。

Otsu′s 算法的優(yōu)點(diǎn)是:方差(k) 較大，完全分割一幅圖像的閾值將會(huì)更加接近。 另外，因?yàn)閗是范圍0 到L－1 內(nèi)的整數(shù)，所以更加易于尋找(k)的最大值。 逐步通過(guò)k 的可能個(gè)L 值，計(jì)算每一步的(k) ，然后選擇) 最大值為k值，即最佳閾值。 如果最大值不唯一，則選取的最優(yōu)閾值即所有k值的平均值。

1.6 連通域標(biāo)記與初級(jí)纖毛定量分析

一個(gè)坐標(biāo)(x，y)的像素p具有兩個(gè)垂直和兩個(gè)水平的相鄰像素，它們的坐標(biāo)分別為(x＋1，y)，(x－1，y)，(x，y＋1)和(x，y－1)。 這個(gè)p的相鄰像素集合記為N4(p)。p的4 個(gè)對(duì)角線相鄰像素坐標(biāo)分別為(x＋1，y＋1)，(x＋1，y－1)，(x－1，y＋1)，(x－1，y－1)，這樣p的4 個(gè)對(duì)角像素被記為ND(p)。N4(p)和ND(p)的并集是p的8 個(gè)相鄰像素，記為N8(p)。

現(xiàn)有另一像素q，若q∈N4(p)，則像素p和q稱為4 鄰接。 同樣，若q∈N8(p)，則p和q稱為8鄰接。 若在前景像素p和q之間存在一條完全由前景像素組成的4 連接路徑，則這兩個(gè)前景像素稱為4 連接。 若他們直接存在一條8 連接路徑，則稱為8連接。

由于在分辨率為2 048像素×1 536像素的圖像中初級(jí)纖毛面積大于60 像素，利用最大類間方差法獲得的二值圖像，應(yīng)進(jìn)行形態(tài)學(xué)開運(yùn)算，以便獲取面積大于60 像素的區(qū)域。 本研究采用8 連接處理，并將編號(hào)標(biāo)記在每個(gè)連通域的質(zhì)心，結(jié)果如圖4 所示。 初級(jí)纖毛像素周長(zhǎng)等于圖4 中被標(biāo)記的連通域周長(zhǎng)，像素面積等于連通域所含有像素點(diǎn)的個(gè)數(shù)。纖毛實(shí)際周長(zhǎng)＝纖毛像素面積/比例尺；纖毛實(shí)際面積＝纖毛像素面積/比例尺2。

圖4 連通域標(biāo)記Fig.4 The connected domain markers

1.7 數(shù)據(jù)評(píng)價(jià)指標(biāo)

根據(jù)初級(jí)纖毛定量分析結(jié)果，通過(guò)計(jì)算不同方法下的相對(duì)誤差進(jìn)行評(píng)價(jià)，計(jì)算公式如下:

式中，σA是本模型算法獲取周長(zhǎng)的標(biāo)準(zhǔn)差，σB是手動(dòng)畫線量取長(zhǎng)度的標(biāo)準(zhǔn)差。 考慮兩者應(yīng)在同一水平下對(duì)比，纖毛周長(zhǎng)CA應(yīng)換算成長(zhǎng)度A，換算公式如下:

式中，ε取值范圍0.01～0.05。

2 實(shí)驗(yàn)結(jié)果

2.1 算法性能評(píng)價(jià)

在圖4 中隨機(jī)選取10 個(gè)初級(jí)纖毛，分別采用本研究模型和手動(dòng)畫線方法采集數(shù)據(jù)。 如表1 所示，自動(dòng)欄顯示本研究模型采集并計(jì)算出的10 個(gè)初級(jí)纖毛周長(zhǎng)，手動(dòng)欄顯示直線“標(biāo)度尺”手動(dòng)分段測(cè)量方法采集的10 個(gè)對(duì)應(yīng)初級(jí)纖毛的長(zhǎng)度(單位:pixel)。 根據(jù)表1 數(shù)據(jù)分別算出本模型算法與手動(dòng)畫線量取的檢測(cè)結(jié)果的相對(duì)誤差ω，為了使衡量的結(jié)果更具有說(shuō)服力，選擇ε＝0.01，算出ω＝19.08% 。 從結(jié)果可以得出本模型相比手動(dòng)畫線方法算法離散程度更小，從一定程度可以說(shuō)明測(cè)量結(jié)果更加準(zhǔn)確。 從算法性能方面比較，手動(dòng)畫線方法獲得上述一張圖片的纖毛長(zhǎng)度需要50～80 min，且操作過(guò)程復(fù)雜，精確度較低。 而應(yīng)用本研究模型算法僅需2 s，即可獲取較為準(zhǔn)確的初級(jí)纖毛周長(zhǎng)與面積參數(shù)。

表1 本模型與手動(dòng)畫線量取的初級(jí)纖毛長(zhǎng)度Tab.1 The length of primary cilia measured automatically and manually

3 討論和結(jié)論

本研究提出了一種基于數(shù)字圖像處理技術(shù)的初級(jí)纖毛自動(dòng)分析方法，該方法有助于大幅提高對(duì)初級(jí)纖毛的定量分析精確性。 與手動(dòng)直線“標(biāo)度尺”分段測(cè)量方法相比，本研究通過(guò)二次改進(jìn)Prewitt、Log 算子以及結(jié)合圖像色彩分離、連通域標(biāo)記等技術(shù)，實(shí)現(xiàn)了對(duì)初級(jí)纖毛結(jié)構(gòu)指標(biāo)的全自動(dòng)精確定量分析；由于算法性能的提升，分析時(shí)間較手動(dòng)處理所需的50～80 min 減少至2 s 左右。 另外，本研究的分析離散程度也更低。 綜上所述，本研究在技術(shù)上較大改進(jìn)了初級(jí)纖毛檢測(cè)的準(zhǔn)確度、分析效率、重復(fù)性。 但由于本研究樣本獲取的多樣性和樣本量均不充分，因此仍需更充分的研究予以驗(yàn)證。 為此，相關(guān)后續(xù)研究將補(bǔ)充更多數(shù)據(jù)來(lái)源的圖像，進(jìn)一步更多的訓(xùn)練，來(lái)提高本研究算法的性能。

國(guó)內(nèi)外對(duì)初級(jí)纖毛的研究多集中在纖毛介導(dǎo)的相關(guān)信號(hào)通路響應(yīng)和外界刺激對(duì)纖毛增殖與形態(tài)功能改變的生物學(xué)效應(yīng)影響。 其中，纖毛增殖與形態(tài)學(xué)研究因受限于傳統(tǒng)測(cè)量方法，而影響其快速進(jìn)展。 本研究對(duì)纖毛形態(tài)特征的自動(dòng)分析，快速精確地測(cè)量纖毛的周長(zhǎng)和面積參數(shù)，有助于促進(jìn)初級(jí)纖毛的研究發(fā)展。