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

王改云

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

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.

關鍵詞：熒光檢測；偏差補償遞推最小二乘；熒光補償；參數辨識

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

中圖分類號：TP391.9 文獻標識碼：A 文章編號：1006-4311（2017）06-0114-04

0 引言

熒光檢測技術在化學分析、臨床醫學分析、法學分析等方面具有廣泛應用價值，對許多學科領域的發展產生深遠的影響[1]。在熒光檢測過程中易出現光譜重疊現象，所以需要對檢測結果進行補償，熒光補償即是從探測器除去除匹配熒光以外的熒光信號的過程[2]。熒光補償方法的精度對檢測結果起著決定性作用。

常用的熒光補償方法可分為兩類，第一類是通過可調電路對檢測器的信號進行調整，以抵消由熒光光譜信號產生的交叉重疊。通常可調電路可分為兩種方法：線性放大補償和對數放大補償[3]。該類補償方法存在硬件調節較為復雜，靈活度較低，檢測精度不高等問題。第二類是通過軟件也即運用數學運算進行矩陣補償的方法，該類方法更為簡單方便且精度更高，所以得到了廣泛的應用。該類補償方法主要有全矩陣補償[4]，最小二乘法補償（Least Squares，LS），遞推最小二乘補償等（Recursive Least Squares，RLS）。全矩陣補償在處理信號噪聲時存在精度方面不足的問題，且補償過程中由于噪聲信號的影響，易出現熒光強度補償結果為負值的情況。傳統的最小二乘補償算法可以獲得補償矩陣參數，但是無法適用于在線辨識以及實時跟蹤[5]。遞推最小二乘補償方法提高了辨識系統的精度和實時性，但是依然存在不能正確補償噪聲所帶來的誤差的問題。

鑒于上述算法的不足，本文提出了一種采用偏差遞推最小二乘算法（Bias Compensation Recursive Least Squares，BCRLS）辨識補償參數的方法，通過在迭代出的估計參數上引入一個修正項來補償在采集熒光中過程噪聲引起的誤差。理論分析和仿真表明，本方法可得出補償參數的無偏估計，對熒光采集過程中產生的噪聲具有很好的補償效果。

1 熒光矩陣參數辨識方法研究

1.1 熒光矩陣參數

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

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

見表1列出了流式細胞儀中常用的四種熒光染料的激發和發射波長，實際的各熒光染料的激發或發射波長是正態或者偏態曲線，即有很寬的范圍[7]。如圖1所示為FITC、PE、ECD、PE-Cy5的發射波長，可以看到四種熒光染料的發射波長均為偏態分布，同時使用四種熒光染料就會發現，四種熒光發射波長相互重疊的現象。在流式細胞儀檢測光信號時，每個檢測通道都會受到其他通道的影響，故此必須要進行熒光補償[8]。

在流式細胞儀中，光檢測器g1所檢測的信號F1（g1）為熒光染料1的熒光強度GN（g1）乘以熒光的泄露矩陣ak加上自身熒光信號所發出的自發熒光的干擾熒光值bk的總和，如式（1）所示，其中ak采用偏差補償遞推最小二乘算法辨識得到。

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

1.2 遞推最小二乘法原理及設計

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

新的參數估計值=舊的參數估計值+修正值

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

2 仿真結果與分析

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

2.1 參數辨識結果分析

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

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

2.2 CV值分析

3 結語

本文通過在迭代估計值上添加一個補償項，從而得出偏差補償遞推最小二乘法對熒光矩陣的參數辨識。該方法與遞推最小二乘法以及一般最小二乘法比較，其結果的參數誤差明顯更小，更加接近準確值。因此利用偏差補償遞推最小二乘法不僅保留了遞推最小二乘法的特點，還提高了辨識精度，彌補了遞推最小二乘法中的不足。

參考文獻：

[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）.