基于偏差補(bǔ)償遞推最小二乘法的熒光補(bǔ)償方法

王改云

摘要：針對(duì)熒光檢測(cè)技術(shù)中熒光光譜重疊引起熒光值偏差的問(wèn)題，提出了一種基于偏差遞推最小二乘算法辨識(shí)補(bǔ)償矩陣的方法。首先，在多輸入多輸出系統(tǒng)（MIMO）下，利用單染色熒光光譜和多染色熒光光譜實(shí)際測(cè)量的熒光值，通過(guò)遞推最小二乘法進(jìn)行迭代運(yùn)算，推導(dǎo)出參數(shù)估計(jì)值。然后，在其中引入一個(gè)修正項(xiàng)，補(bǔ)償在采集熒光中過(guò)程噪聲引起的誤差。最后，計(jì)算出偏差補(bǔ)償遞推最小二乘法迭代的估計(jì)值。理論分析與仿真表明，該算法在參數(shù)誤差估計(jì)中誤差率小于1%，相對(duì)于遞推最小二乘算法性能提高了50%。所用算法能夠有效提高估計(jì)值精度，同時(shí)也能減小噪聲產(chǎn)生的影響。

Abstract： Focused on the intensity deviation issue caused by spectral overlap in fluorescence detection technology， the bias compensation recursive last squares method for compensation matrix was proposed. Firstly， based on the measured single and multiple staining fluorescence values， the parameter estimates were deduced through recursive least squares method in the multi-input multi-output system （MIMO）. Secondly， by introducing a correction term into the estimated values， the errors caused by noise in the fluorescence acquisition process were compensated. Finally， the estimates was calculated iteratively with bias compensation recursive least squares. The simulation results and theoretical analysis show that with this method， the error rate is less than 1% and the performance is improved by 50% compared with the recursive last squares algorithm. The proposed method can effectively improve the accuracy of estimates， meanwhile reduce the negative effect of noise.

關(guān)鍵詞：熒光檢測(cè)；偏差補(bǔ)償遞推最小二乘；熒光補(bǔ)償；參數(shù)辨識(shí)

Key words： fluorescence detection；bias compensation recursive least squares；fluorescence compensation；parameter identification

中圖分類號(hào)：TP391.9 文獻(xiàn)標(biāo)識(shí)碼：A 文章編號(hào)：1006-4311（2017）06-0114-04

0 引言

熒光檢測(cè)技術(shù)在化學(xué)分析、臨床醫(yī)學(xué)分析、法學(xué)分析等方面具有廣泛應(yīng)用價(jià)值，對(duì)許多學(xué)科領(lǐng)域的發(fā)展產(chǎn)生深遠(yuǎn)的影響[1]。在熒光檢測(cè)過(guò)程中易出現(xiàn)光譜重疊現(xiàn)象，所以需要對(duì)檢測(cè)結(jié)果進(jìn)行補(bǔ)償，熒光補(bǔ)償即是從探測(cè)器除去除匹配熒光以外的熒光信號(hào)的過(guò)程[2]。熒光補(bǔ)償方法的精度對(duì)檢測(cè)結(jié)果起著決定性作用。

常用的熒光補(bǔ)償方法可分為兩類，第一類是通過(guò)可調(diào)電路對(duì)檢測(cè)器的信號(hào)進(jìn)行調(diào)整，以抵消由熒光光譜信號(hào)產(chǎn)生的交叉重疊。通常可調(diào)電路可分為兩種方法：線性放大補(bǔ)償和對(duì)數(shù)放大補(bǔ)償[3]。該類補(bǔ)償方法存在硬件調(diào)節(jié)較為復(fù)雜，靈活度較低，檢測(cè)精度不高等問(wèn)題。第二類是通過(guò)軟件也即運(yùn)用數(shù)學(xué)運(yùn)算進(jìn)行矩陣補(bǔ)償?shù)姆椒?，該類方法更為?jiǎn)單方便且精度更高，所以得到了廣泛的應(yīng)用。該類補(bǔ)償方法主要有全矩陣補(bǔ)償[4]，最小二乘法補(bǔ)償（Least Squares，LS），遞推最小二乘補(bǔ)償?shù)龋≧ecursive Least Squares，RLS）。全矩陣補(bǔ)償在處理信號(hào)噪聲時(shí)存在精度方面不足的問(wèn)題，且補(bǔ)償過(guò)程中由于噪聲信號(hào)的影響，易出現(xiàn)熒光強(qiáng)度補(bǔ)償結(jié)果為負(fù)值的情況。傳統(tǒng)的最小二乘補(bǔ)償算法可以獲得補(bǔ)償矩陣參數(shù)，但是無(wú)法適用于在線辨識(shí)以及實(shí)時(shí)跟蹤[5]。遞推最小二乘補(bǔ)償方法提高了辨識(shí)系統(tǒng)的精度和實(shí)時(shí)性，但是依然存在不能正確補(bǔ)償噪聲所帶來(lái)的誤差的問(wèn)題。

鑒于上述算法的不足，本文提出了一種采用偏差遞推最小二乘算法（Bias Compensation Recursive Least Squares，BCRLS）辨識(shí)補(bǔ)償參數(shù)的方法，通過(guò)在迭代出的估計(jì)參數(shù)上引入一個(gè)修正項(xiàng)來(lái)補(bǔ)償在采集熒光中過(guò)程噪聲引起的誤差。理論分析和仿真表明，本方法可得出補(bǔ)償參數(shù)的無(wú)偏估計(jì)，對(duì)熒光采集過(guò)程中產(chǎn)生的噪聲具有很好的補(bǔ)償效果。

1 熒光矩陣參數(shù)辨識(shí)方法研究

1.1 熒光矩陣參數(shù)

熒光檢測(cè)技術(shù)在使用時(shí)，普遍存在發(fā)射熒光光譜重疊問(wèn)題。以流式細(xì)胞儀為例，為了測(cè)量諸如細(xì)胞的微粒特性，通常使用熒光染料標(biāo)記微粒向熒光染料照射激光束以激發(fā)熒光染料以及測(cè)量由激發(fā)的熒光染料發(fā)出的熒光強(qiáng)度或圖譜[6]。如今這些細(xì)胞大多使用多種熒光染料標(biāo)記微粒，并且使用具有不同接收光波段的多個(gè)光檢測(cè)器如光電倍增管（PMT）或硅光子計(jì)數(shù)器（MPCC）從熒光染料發(fā)出的激光束進(jìn)行多色法測(cè)量。

目前使用的熒光染料在熒光光譜中大多具有重疊的頻段，在使用這些熒光染料組合進(jìn)行多色測(cè)量時(shí)，濾光片能將從各個(gè)熒光染料發(fā)出的熒光分離為不同的頻段。但是，光檢測(cè)器仍然可以接收從別的通道的熒光染料泄露到此通道的熒光，從而產(chǎn)生誤差。

見(jiàn)表1列出了流式細(xì)胞儀中常用的四種熒光染料的激發(fā)和發(fā)射波長(zhǎng)，實(shí)際的各熒光染料的激發(fā)或發(fā)射波長(zhǎng)是正態(tài)或者偏態(tài)曲線，即有很寬的范圍[7]。如圖1所示為FITC、PE、ECD、PE-Cy5的發(fā)射波長(zhǎng)，可以看到四種熒光染料的發(fā)射波長(zhǎng)均為偏態(tài)分布，同時(shí)使用四種熒光染料就會(huì)發(fā)現(xiàn)，四種熒光發(fā)射波長(zhǎng)相互重疊的現(xiàn)象。在流式細(xì)胞儀檢測(cè)光信號(hào)時(shí)，每個(gè)檢測(cè)通道都會(huì)受到其他通道的影響，故此必須要進(jìn)行熒光補(bǔ)償[8]。

在流式細(xì)胞儀中，光檢測(cè)器g1所檢測(cè)的信號(hào)F1（g1）為熒光染料1的熒光強(qiáng)度GN（g1）乘以熒光的泄露矩陣ak加上自身熒光信號(hào)所發(fā)出的自發(fā)熒光的干擾熒光值bk的總和，如式（1）所示，其中ak采用偏差補(bǔ)償遞推最小二乘算法辨識(shí)得到。

以上述四色熒光補(bǔ)償為例，本文分別使用一般最小二乘法、遞推最小二乘法和偏差補(bǔ)償遞推最小二乘法三種方法辨識(shí)補(bǔ)償矩陣的參數(shù)值，把三種算法所辨識(shí)到的補(bǔ)償矩陣參數(shù)用于熒光補(bǔ)償矩陣，最后檢驗(yàn)這三種參數(shù)辨識(shí)算法具體在熒光補(bǔ)償矩陣中的作用。

1.2 遞推最小二乘法原理及設(shè)計(jì)

遞推最小二乘算法的思想可以概括為：

新的參數(shù)估計(jì)值=舊的參數(shù)估計(jì)值+修正值

即新的遞推參數(shù)估計(jì)是在舊的遞推估計(jì)值的基礎(chǔ)上修正而成，這就是遞推的概念[9]。遞推估計(jì)算法無(wú)需存儲(chǔ)全部數(shù)據(jù)，取得一組觀測(cè)數(shù)據(jù)，便可估計(jì)一次參數(shù)，因此所需的計(jì)算量和占用的存儲(chǔ)空間都很小，而且能實(shí)現(xiàn)在線實(shí)時(shí)辨識(shí)。這樣，隨著新觀測(cè)數(shù)據(jù)的逐次引入，一次接一次地進(jìn)行參數(shù)估計(jì)，直到參數(shù)估計(jì)值達(dá)到滿意的精確程度為止。對(duì)于流式細(xì)胞儀來(lái)說(shuō)分析數(shù)據(jù)的速度約為幾萬(wàn)個(gè)細(xì)胞每秒，在進(jìn)行熒光補(bǔ)償時(shí)要求處理數(shù)據(jù)快速、準(zhǔn)確，遞推最小二乘法適合流式細(xì)胞儀熒光補(bǔ)償矩陣的參數(shù)估計(jì)。

2 仿真結(jié)果與分析

要分析每一個(gè)檢測(cè)通道獲得的熒光數(shù)據(jù)必須得到純熒光染料在每個(gè)檢測(cè)通道的校準(zhǔn)數(shù)據(jù)和自動(dòng)熒光控制數(shù)據(jù)，通過(guò)加入空白實(shí)驗(yàn)和校準(zhǔn)實(shí)驗(yàn)可以得出所需熒光數(shù)據(jù)[11]。在采集過(guò)程中，通過(guò)使用從標(biāo)記一種熒光染料的微粒獲得的單染色光譜，然后再分別檢測(cè)其他三種單染色光譜，所采集光譜范圍內(nèi)的10000個(gè)數(shù)據(jù)作為光檢測(cè)器Gn，n=1，2，3，4（10000*4矩陣）所檢測(cè)的熒光值，最后采集標(biāo)記四種熒光染料的微粒獲得多染色的光譜，使用所采集的光譜范圍內(nèi)的10000個(gè)數(shù)據(jù)作為Fn（gn）n=1，2，3，4 （10000*4矩陣）檢測(cè)到的熒光值。

2.1 參數(shù)辨識(shí)結(jié)果分析

根據(jù)所采集的熒光值G（10000*4矩陣）和（10000*4矩陣）的數(shù)據(jù)，利用Matlab使用一般最小二乘算法、遞推最小二乘法和偏差補(bǔ)償遞推最小二乘法三種數(shù)學(xué)模型分別求得補(bǔ)償矩陣KLS，KRLS，KBCRLS。系統(tǒng)中噪聲e（k）零均值方差為σe2， 改變?chǔ)襡2可以影響噪聲信號(hào)，故設(shè)定在不同方差下測(cè)量系統(tǒng)的參數(shù)估計(jì)誤差。當(dāng)σe2=1.002時(shí)，測(cè)量參數(shù)誤差情況如圖2所示，當(dāng)σe2=0.502時(shí)，測(cè)量得出的參數(shù)誤差情況如圖3所示，從三種算法的補(bǔ)償矩陣可以得出誤差參數(shù)δ=，其中θ表示實(shí)際熒光值，表示估計(jì)熒光值，δ為誤差參數(shù)。

對(duì)于一般最小二乘法與遞推最小二乘法可以求得系統(tǒng)的補(bǔ)償矩陣，但是從圖2、3以及表2、3中可以看出一般最小二乘法與遞推最小二乘法辨識(shí)出的矩陣參數(shù)誤差比較大，得出的估計(jì)值是偏離實(shí)際值的。三種方法中，偏差補(bǔ)償遞推最小二乘法的參數(shù)誤差精度最高。

2.2 CV值分析

3 結(jié)語(yǔ)

本文通過(guò)在迭代估計(jì)值上添加一個(gè)補(bǔ)償項(xiàng)，從而得出偏差補(bǔ)償遞推最小二乘法對(duì)熒光矩陣的參數(shù)辨識(shí)。該方法與遞推最小二乘法以及一般最小二乘法比較，其結(jié)果的參數(shù)誤差明顯更小，更加接近準(zhǔn)確值。因此利用偏差補(bǔ)償遞推最小二乘法不僅保留了遞推最小二乘法的特點(diǎn)，還提高了辨識(shí)精度，彌補(bǔ)了遞推最小二乘法中的不足。

參考文獻(xiàn)：

[1]B.Bohn Anja， ，K Bjarne，Moller，S Mikkel，Petersen. Flow cytometry and compensation of highly autofluorescent cells： the example of mesenchymal stem cells[J] .Stem Cell Biology and Research，2015，1： 2-4.

[2]N Aghaeepour，G Finak，TF Consortium，et al.Critial assessment of automated flow cytometry data analysis techniques[J]. Nat.Methods，2013（10）：228-238.

[3]C.Bruce Bagwell， G.Adams Earl. Fluorescence Spectral Overlap Compensation for Any Number of Flow Cytometry Parameters[J].Annals of the New York Academy of Sciences，1993，677：167-184.

[4]S Morishita， H Yokota， H Asama， R Himeno， T Mishima. Compensation method for quantitative observation of multicolor fluorescence with nonlinear mapping [J]. Proceedings of SPIE - The International Society. 2008， 7075：70750J.

[5]C. Stewart Carleton， J. Stewart Sigrid.Four color compensation [J].Cytometry， 1999， 38：161-175.

[6]Yin Zhongbin， Jessica Severin， Michael C. Giddings， Huang Wei-an， Michael S. Westphall，Lloyd M. Smith. Automatic matrix determination in four dye fluorescence-based DNA sequencing[J]. Electrophoresis， 1996， 17： 1143-1150.

[7]Pospichalova Vendula， Dave Zankruti， et al. Simplfied protocol for flow cytometry analysis of fluorescently labeled exosomes and microvesicles using dedicated flow cytometer[J]. Journal of Extracellular Vesicles， 2015， 4： 25530.

[8]JW Tung， DR Parks， WA Moore， LA Herzenberg， LA Herzenberg. New approaches to fluorescence compensation and visualization of FACS data[J]. Clinical Immunology. 2004， 110（3）：277-83.

[9]Shang Xuemei， Dong Mingli， Pan Zhiwen， et al. Identification of the parameters of spectral compensation for flow cytometer based on the recursive least squares algorithm[J]. Computer Simulation， 2015， 32（5）：252-255（in Chinese）.

[10]Zhang Yong， Cui Guimei. Bisa compensation methods for stochastic systems with colored noise[J]. Applied Mathematical Modelling，2011，35：1709-1716.

[11]Zan Liuqin， Zhang Zhenxi， Miao Baogan， et al. Multicolor fluorescence detection in the multiplex quantitative PCR system and spectra crosstalk correction method[J]. Acta Optica Sinica， 2014-01：187-193（in Chinese）.