周基陽 ，徐聲偉 ，林楠森 ，蔡新霞

（1.中國科學(xué)院電子學(xué)研究所 傳感技術(shù)聯(lián)合國家重點(diǎn)實(shí)驗(yàn)室，北京 100190；2.中國科學(xué)院大學(xué)電子學(xué)研究所，北京 100190）

神經(jīng)系統(tǒng)通過電生理和電化學(xué)信號(hào)進(jìn)行交流[1]，其中，電生理信號(hào)包括場電位和動(dòng)作電位。利用生物電測量技術(shù)的微電極，可以采集到神經(jīng)元發(fā)放的動(dòng)作電位。

研究動(dòng)作電位，需要把動(dòng)作電位從大量的噪聲信號(hào)中提取出來。目前的動(dòng)作電位的提取方法有閾值檢測，峰值檢測，非線性能量閾值檢測，基于數(shù)學(xué)形態(tài)學(xué)預(yù)處理的神經(jīng)元放電檢測方法[2-3]，其次還有改進(jìn)峰值檢測、小波變換等算法[4-5]。

閾值檢測法是目前最廣泛的檢測方法[6-7]，目前閾值的自適應(yīng)比較困難，計(jì)算復(fù)雜，在電路上沒有很好的實(shí)現(xiàn)。將來把自適應(yīng)閾值算法應(yīng)用到多通道神經(jīng)信號(hào)數(shù)據(jù)采集儀器上去，進(jìn)行動(dòng)作電位的實(shí)時(shí)采集[8-10]，可以濾除干擾，更加準(zhǔn)確地提取動(dòng)作電位信號(hào)，易于后端的進(jìn)一步分析。

本文設(shè)計(jì)了一種自適應(yīng)閾值算法，它可以對(duì)采集到的神經(jīng)信息，經(jīng)過不斷的計(jì)算來調(diào)整閾值，并根據(jù)此閾值，實(shí)時(shí)地提取動(dòng)作電位。而且當(dāng)采集的數(shù)據(jù)中有較大的噪聲波動(dòng)時(shí)，算法能夠根據(jù)數(shù)據(jù)的變化快速的調(diào)整閾值，因此動(dòng)作電位提取結(jié)果不會(huì)受到較大的影響。此外，此算法比一般的自適應(yīng)閾值法計(jì)算量小，在硬件上易于實(shí)現(xiàn)。

1 實(shí)驗(yàn)部分

動(dòng)作電位最明顯的形態(tài)特征是它的幅值，因此，檢測動(dòng)作電位最為簡單的方法是閾值法，當(dāng)信號(hào)中有超過閾值的部分，即認(rèn)為有動(dòng)作電位出現(xiàn)。

目前自動(dòng)設(shè)定閾值的方法主要基于以下公式。信號(hào)標(biāo)準(zhǔn)差的倍數(shù)：

式中：x代表動(dòng)作電位信號(hào)各采樣點(diǎn)的幅值；n為采樣點(diǎn)數(shù)。

信號(hào)的平均灰度值：

式中：x代表動(dòng)作電位信號(hào)各采樣點(diǎn)的幅值；median代表取中值。

根據(jù)上述公式自動(dòng)設(shè)定閾值的方法復(fù)雜，本文將此公式進(jìn)行了簡化，設(shè)計(jì)了一種提取動(dòng)作電位的正向自適應(yīng)閾值算法，并對(duì)算法分別進(jìn)行了Matlab仿真和FPGA電路實(shí)現(xiàn)。

1.1 自適應(yīng)閾值算法

首先，將式（2）做如下修改：

式中：x為大于 0的數(shù)據(jù)；median（x）為求 x的中值；a取 3～9；b是噪聲基線。

上述是用中值法來求閾值，中值法在電路上實(shí)現(xiàn)較為復(fù)雜，為了降低硬件電路復(fù)雜度，對(duì)式（3）做進(jìn)一步的改進(jìn)，用均值代替中值。

式中：mean（x）為求 x 的均值；N 取 1～9，根據(jù)放大的倍數(shù)和實(shí)際的電路決定。當(dāng)b為0時(shí)公式為

本文最后用式 （4），（5）來進(jìn)行仿真和電路實(shí)現(xiàn)，利用均值的方法來計(jì)算閾值，比計(jì)算信號(hào)的標(biāo)準(zhǔn)差的倍數(shù)，信號(hào)的平均灰度值更加簡單，計(jì)算量小，容易實(shí)現(xiàn)。算法如下：

（1）依次判斷采樣數(shù)據(jù)是否大于b。如果是，則把采集到的數(shù)據(jù)賦值給x。

（2）計(jì)算閾值Thr。當(dāng)x中的采樣點(diǎn)數(shù)達(dá)到一定數(shù)值 M（M 預(yù)設(shè)為 200，500，1024 等整數(shù)值）時(shí)，計(jì)算一次 x 的均值，然后根據(jù)式（5）計(jì)算 Thr（Thr1）。不斷的采樣，用新的采樣數(shù)據(jù)覆蓋x，下次采樣到M個(gè)采樣點(diǎn)時(shí)，用此新數(shù)據(jù)計(jì)算x的均值和Thr2，依次計(jì)算 Thr3，Thr4……

（3）在Thr1計(jì)算完畢之前，把所有的采集數(shù)據(jù)全部發(fā)送顯示出來，不做處理。

（4）在 Thr1計(jì)算完畢之后，Thr2計(jì)算完畢之前，用剛剛計(jì)算的Thr1作為閾值，來判斷此段采樣數(shù)據(jù)中是否含有Spike信號(hào)。

（5）在Thr2計(jì)算完畢后，Thr3計(jì)算完畢之前，把這個(gè)新的Thr2作為閾值來判斷此段數(shù)據(jù)是否出現(xiàn)動(dòng)作電位。依次下去，實(shí)現(xiàn)動(dòng)作電位的實(shí)時(shí)提取。

計(jì)算機(jī)計(jì)算閾值的時(shí)間很短，可以忽略不計(jì)。自適應(yīng)閾值算法經(jīng)過一段時(shí)間就調(diào)整一次閾值，以適應(yīng)數(shù)據(jù)的變化，找出最佳的閾值。

1.2 Matlab仿真的實(shí)現(xiàn)

本文首先用Matlab進(jìn)行仿真，測試算法是否能夠準(zhǔn)確地提取動(dòng)作電位。仿真數(shù)據(jù)是一段長0.2 s的實(shí)測數(shù)據(jù)。數(shù)據(jù)的維數(shù)為Q，數(shù)據(jù)分為若干個(gè)（假設(shè)為c）小段，每個(gè)小段包含200（M=200）個(gè)大于噪聲基線b的數(shù)，然后每小段數(shù)據(jù)分別利用式（4）計(jì)算一次閾值，共c個(gè)閾值，利用這c個(gè)閾值判定是否為動(dòng)作電位。即用第一段數(shù)據(jù)計(jì)算第一個(gè)閾值Thr1，用Thr1判定第二段數(shù)據(jù)，用第二段數(shù)據(jù)計(jì)算第二個(gè)閾值Thr2，用Thr2判定第三段數(shù)據(jù)，依次下去，直到最后。模擬算法的步驟，測試此種方法是否有效。

1.3 自適應(yīng)閾值算法的FPGA實(shí)現(xiàn)

仿真后，本文又將算法在FPGA上進(jìn)行了電路驗(yàn)證。

1.3.1 自適應(yīng)閾值法的硬件電路

信號(hào)源信號(hào)首先進(jìn)入前置放大電路，再經(jīng)過帶通濾波放大電路，之后進(jìn)行16位采樣，采樣得到數(shù)據(jù)再經(jīng)過自適應(yīng)閾值算法處理后，通過USB接口發(fā)送到上位機(jī)上去，如圖1所示。

圖1 自適應(yīng)閾值法電路Fig.1 Circuit of adaptive threshold method

信號(hào)源產(chǎn)生的神經(jīng)電生理信號(hào)包括場電位（0.1～100）Hz和動(dòng)作電位（300～1000）Hz。 首先進(jìn)行前置放大，再進(jìn)行帶通濾波，濾波范圍是200 Hz～3 kHz，共放大3800多倍。

16位采樣模塊采用AD7656-1芯片。它是16 bit ADC，滿足神經(jīng)元?jiǎng)幼麟娢徊蓸右?。同時(shí)AD7656-1支持6個(gè)通道的模擬輸入，易于將來進(jìn)行多通道實(shí)驗(yàn)。

FPGA模塊是電路的核心部分，用于實(shí)現(xiàn)自適應(yīng)閾值法、采樣控制、計(jì)時(shí)和控制USB數(shù)據(jù)傳輸?shù)墓δ?。FPGA采用ALTERA FPGA NIOS CYCLONE IV開發(fā)板，它集成了CH376，它支持1.5 Mb/s低速和12 Mb/s全速USB通訊，滿足算法電路驗(yàn)證的要求。

1.3.2 自適應(yīng)閾值法軟件的設(shè)計(jì)

在經(jīng)過仿真后，本文又將算法在FPGA上進(jìn)行了驗(yàn)證。算法主要實(shí)現(xiàn)以下功能：

（1）產(chǎn)生采樣控制信號(hào)，控制AD7656-1完成采樣功能。

控制FPGA產(chǎn)生CONVST轉(zhuǎn)換啟動(dòng)信號(hào)和RD信號(hào)，使每60 μs采樣一次，即采樣率為16.67 kHz。

（2）控制USB的發(fā)送。

控制USB芯片CH376的發(fā)送，實(shí)現(xiàn)當(dāng)一旦檢測到大于閾值，便啟動(dòng)發(fā)送，之后，不再發(fā)送，直到下一次檢測到大于閾值時(shí)才再次發(fā)送。

（3）完成計(jì)時(shí)功能。

每過60 μs，計(jì)數(shù)器加1。程序中定義了一個(gè)32位計(jì)數(shù)器，可以計(jì)時(shí)71 h，滿足實(shí)驗(yàn)需求。

（4）自適應(yīng)閾值法的實(shí)現(xiàn)。

用FPGA實(shí)現(xiàn)自適應(yīng)閾值法，主要有兩個(gè)部分，第1，自動(dòng)閾值的計(jì)算；第2，當(dāng)檢測到信號(hào)大于閾值時(shí)，將此部分動(dòng)作電位提取出來。

第1 計(jì)算自適應(yīng)閾值。在實(shí)驗(yàn)中為了消除負(fù)值影響，把所有采樣數(shù)據(jù)加上十進(jìn)制的32768，即把所有的數(shù)據(jù)提高了 5 V。所以式 （4）中，b=5 V（32768），M 為 1024，N 為 1～9。 自適應(yīng)閾值法的流程圖如圖2所示。

圖2 自適應(yīng)閾值法FPGA實(shí)現(xiàn)流程Fig.2 FPGA flow of adaptive threshold method

首先比較AD采樣得到的輸入數(shù)據(jù)Vin和噪聲基線b，平移5 V，b所以為5 V。在不斷采樣的同時(shí)，不斷的統(tǒng)計(jì)超過噪聲基線（5 V）的值，當(dāng)此值達(dá)到1024時(shí)，啟動(dòng)一次閾值計(jì)算（把sort置為1，sort為啟動(dòng)信號(hào)），開始計(jì)算一次閾值，記為Thr1。計(jì)算完畢之后，sort置為0，結(jié)束計(jì)算。

隨著不斷的采樣，再次重新統(tǒng)計(jì)Vin中大于5 V的數(shù)（再從0統(tǒng)計(jì)超過b的值），當(dāng)超過b的值再次達(dá)到1024個(gè)時(shí)，啟動(dòng)第二次閾值計(jì)算（再把sort置為1，計(jì)算完畢后sort置為0），把第二次計(jì)算的閾值記為Thr2，依次下去，計(jì)算Thr3，Thr4……。

第2 提取動(dòng)作電位的方法主要運(yùn)用了移位技術(shù)，具體的實(shí)現(xiàn)方法如下。

（1）首先定義兩個(gè)16位寬度的32個(gè)數(shù)據(jù)的數(shù)組，兩個(gè)數(shù)組維數(shù)相同，寬度相同。如表1所示。

表 1 Write_Data（Shift_Data）Tab.1 Write_Data（Shift_Data）

Shift_Data數(shù)組用于放置AD采樣的數(shù)據(jù)，Write_Data用于存放檢測到的動(dòng)作電位。

（2）當(dāng)每過 60 μs時(shí)，AD 采樣一次，得到一個(gè)數(shù)據(jù)，把這個(gè)數(shù)據(jù)放到Shift_Data的END位，下一次60 μs時(shí)，把上一次的數(shù)據(jù)向HEAD方向移位一次，把新采集的數(shù)據(jù)放到END位，重復(fù)下去，實(shí)現(xiàn)了數(shù)據(jù)的依次移位，溢出的數(shù)據(jù)刪除。

（3）同時(shí)每過 60 μs，就檢測 Shift_Data 的第 8位一次，當(dāng)此位大于自動(dòng)計(jì)算的閾值時(shí)，就認(rèn)為是動(dòng)作電位，即把這32位數(shù)組Shift_Data全部送給Write_Data，這里存放提取的動(dòng)作電位。當(dāng)還沒有進(jìn)行閾值的第一次計(jì)算，此時(shí)閾值設(shè)為0。當(dāng)Thr1計(jì)算完畢之后，在Thr2計(jì)算完畢之前，用Thr1來判定Shift_Data第8位數(shù)據(jù)，看此位數(shù)據(jù)是否大于Thr1。在Thr2計(jì)算完畢之后，再Thr3計(jì)算完畢之前，用Thr2判斷第8位數(shù)據(jù)，看此位數(shù)據(jù)是否大于Thr2，依次進(jìn)行下去。

（4）為了防止重復(fù)給Write_Data賦值，如果檢測到動(dòng)作電位，將Shift_Data所有位清零，然后再次采樣，移位。

（5）檢測到Write_Data有數(shù)據(jù)一次，立刻發(fā)送一次，即把這里的動(dòng)作電位通過usb接口發(fā)送出去，這樣就實(shí)現(xiàn)了自適應(yīng)閾值法。

1.3.3 上位機(jī)軟件

上位機(jī)軟件采用LabVIEW （laboratory virtual instrument engineering workbench）編寫，實(shí)現(xiàn)實(shí)時(shí)顯示和保存動(dòng)作電位的功能。

上位機(jī)中需要記錄有效的動(dòng)作電位數(shù)據(jù)和計(jì)時(shí)數(shù)據(jù)，來觀測算法是否實(shí)時(shí)有效，同時(shí)存儲(chǔ)數(shù)據(jù)。

2 結(jié)果與討論

2.1 算法的Matlab仿真結(jié)果

本文采用大鼠實(shí)驗(yàn)的實(shí)測數(shù)據(jù)作為仿真數(shù)據(jù)，采樣率為25 k/s。圖3（a）中細(xì)實(shí)線為原始數(shù)據(jù)，粗線為不斷調(diào)整的閾值。 對(duì)應(yīng)式（4），取 N=5，b=0，自適應(yīng)閾值法中M=200。圖3（b）是自適應(yīng)閾值法的結(jié)果。

圖3 自適應(yīng)閾值法Matlab仿真Fig.3 Matlab simulation of adaptive threshold method

根據(jù)算法特點(diǎn)，圖中前一段數(shù)據(jù)為原始數(shù)據(jù)，沒有經(jīng)過閾值提取。不難看出，自適應(yīng)閾值法能有效的抑制噪聲信號(hào)，完成了動(dòng)作電位提取的功能。

在式（4）中，有個(gè)關(guān)鍵參數(shù)N，N的取值對(duì)檢測結(jié)果有直接的影響。以一段數(shù)據(jù)data（[1:35000]）為例，這段數(shù)據(jù)包含55個(gè)動(dòng)作電位，N的取值對(duì)檢測結(jié)果的影響如表2所示。

表2 N的大小對(duì)結(jié)果的影響Tab.2 Influence of N

由表2可以看出，N的取值十分重要，直接影響檢測結(jié)果。N的大小和電路的放大倍數(shù)、動(dòng)作電位的幅度、噪聲、系統(tǒng)的條件等有關(guān)，需要根據(jù)實(shí)際情況進(jìn)行確定。

測試結(jié)果證明，根據(jù)式（4）能夠成功地計(jì)算出閾值，并且能夠有效提取動(dòng)作電位。

2.2 算法的FPGA電路驗(yàn)證結(jié)果

本文采用神經(jīng)信號(hào)發(fā)生器產(chǎn)生的動(dòng)作電位來檢測算法是否成功。經(jīng)過濾波，濾除了神經(jīng)信號(hào)發(fā)生器所產(chǎn)生的場電位，經(jīng)過放大，把300 μV的動(dòng)作電位放大到1 V左右，在此基礎(chǔ)上進(jìn)行實(shí)驗(yàn)。

根據(jù)式（4），在FPGA中根據(jù)AD采樣數(shù)據(jù)的不斷更新，每隔一段時(shí)間自動(dòng)計(jì)算閾值一次，并且根據(jù)此閾值，判定動(dòng)作電位是否出現(xiàn)。式（4）中b取5 V，M為1024，N取7時(shí)，自適應(yīng)閾值法在上位機(jī)上的處理結(jié)果如圖4所示。圖4（a）為自適應(yīng)閾值法得到的部分結(jié)果，可以看出，大于閾值的動(dòng)作電位被發(fā)送了出來，其余的小于閾值的噪聲信號(hào)沒有被發(fā)送出來。圖 4（b）是圖 4（a）中單個(gè)動(dòng)作電位的放大圖形，圖4（c）是動(dòng)作電位疊加的結(jié)果。

圖4 自適應(yīng)閾值法得到結(jié)果Fig.4 Result of adaptive threshold method

當(dāng)N取0時(shí)，相當(dāng)于直接進(jìn)行采樣，沒有進(jìn)行閾值的判定，把所有的動(dòng)作電位和噪聲都發(fā)送出來；當(dāng)N取5時(shí)，把100%的動(dòng)作電位和部分噪聲截取出來了；提高閾值，當(dāng)N取7時(shí)，把100%的動(dòng)作電位和極少部分噪聲截取出來；濾除了絕大部分的噪聲信號(hào)。

3 結(jié)語

經(jīng)過Matlab仿真和FPGA電路驗(yàn)證，本文設(shè)計(jì)的自適應(yīng)閾值法能夠根據(jù)神經(jīng)信號(hào)的數(shù)據(jù)不斷自動(dòng)地計(jì)算閾值，并且可以根據(jù)此閾值實(shí)時(shí)地提取動(dòng)作電位信號(hào)。

本文利用均值計(jì)算閾值的方法簡單，計(jì)算量小，容易在電路上進(jìn)行實(shí)現(xiàn)。此算法隨著數(shù)據(jù)段的變化根據(jù)公式自動(dòng)的調(diào)整閾值，具有較大的靈活性。此種動(dòng)作電位的提取方法，簡化了定閾值的設(shè)定過程，實(shí)現(xiàn)了自適應(yīng)的閾值檢測法。同時(shí)，對(duì)大量的冗余信號(hào)自動(dòng)進(jìn)行了濾除，提高了數(shù)據(jù)的傳輸效率，簡化了后續(xù)的處理過程，這為后續(xù)腦神經(jīng)活動(dòng)的研究提供了有效的依據(jù)。

[1] 袁彩霞.人體神經(jīng)系統(tǒng)的研究與啟示[J].中國醫(yī)學(xué)影像技術(shù)，2003，19（204）：72-73.

[2] 王冬雪.神經(jīng)元峰電位的檢測與分類方法研究[D].安徽：中國科學(xué)技術(shù)大學(xué)，2011：11-19.

[3] 王清波.基于信息壓縮的無線植入式腦機(jī)接口中算法及系統(tǒng)研究[D].浙江：浙江大學(xué)，2011：1-94.

[4] 姚舜，劉海龍，陳傳平，等.一種獲取鋒電位的峰值檢測算法的改進(jìn)方案[J].生物醫(yī)學(xué)工程研究，2005，24（1）：14-17.

[5] Zoran Nenadic，JoelW Burdick.Spikedetection usingthe continuous wavelet transform[J].IEEE Transactions on Biomedical Engineering，2005，52（1）：74-87.

[6] Reid R Harrison.A low-power integrated circuit for adaptive detection of action potentials in noisy signals[C]//Proceedings of the 25’Annual International Conference of the IEEE EMBS Cancun.Mexico，2003：3325-3328.

[7] Paul T Watkins，Gopal Santhanam，Krishna V Shenoy，et al.Validation of adaptive threshold spike detector for neural recording[C]//Proceedings of the 26th Annual International Conference of the IEEE EMBS.San Francisco：CA，2004：4079-4082.

[8] J Morizio，J Parmentier，J Bender，et al.15 Channel Wireless Headstage System for Single Unit Rat and Primate Recordings[P].USA：No.TBME-00565-2005.

[9] Li Huang，Xu Zhang，Ning Guan，etc.Real-time multi-channel system for neural spikes acquisition and detection[C]//New Circuits and Systems Conference（NEWCAS），2012 IEEE 10th International.Beijing：Inst.of Semicond，2012：149-152.

[10]Reid R Harrison，Paul T Watkins，Ryan J Kier，et al.A lowpower integrated circuit for a wireless 100-electrode neural recording system[J].IEEE Journal of Solid-State Circuits，2007，42（1）：123-132. ■