Bootstrap非參數(shù)法在γ能譜處理中的應(yīng)用

徐紅鵑 王 鋒 艾憲蕓 李京倫 魏 星 石 磊

1（北京防化研究院 北京 102205）

2（國(guó)民核生化災(zāi)害防護(hù)國(guó)家重點(diǎn)實(shí)驗(yàn)室 北京 102205）

Bootstrap非參數(shù)法在γ能譜處理中的應(yīng)用

徐紅鵑1,2王 鋒1,2艾憲蕓1,2李京倫1,2魏 星1,2石 磊1,2

1（北京防化研究院 北京 102205）

2（國(guó)民核生化災(zāi)害防護(hù)國(guó)家重點(diǎn)實(shí)驗(yàn)室 北京 102205）

為減少統(tǒng)計(jì)漲落對(duì)能譜測(cè)量的影響，提高測(cè)量的精確性，利用數(shù)理統(tǒng)計(jì)理論，研究Bootstrap非參數(shù)法的基本原理及程序流程，并將其應(yīng)用于實(shí)驗(yàn)室γ能譜數(shù)據(jù)處理中。用NaI(Tl)探測(cè)器測(cè)量241Am、137Cs、60Co的γ能譜，對(duì)短時(shí)間測(cè)量的γ能譜采用Bootstrap非參數(shù)法進(jìn)行處理，與同樣條件下長(zhǎng)時(shí)間測(cè)量的γ能譜數(shù)據(jù)進(jìn)行對(duì)比，吻合度較高，得到了滿意的結(jié)果，為實(shí)現(xiàn)γ能譜在低活度水平的快速測(cè)量提供了可行的技術(shù)方案。

統(tǒng)計(jì)自舉方法(Bootstrap)，非參數(shù)法，置信區(qū)間，γ能譜

在γ能譜測(cè)量中，如果被測(cè)樣品活度很低或者探測(cè)器的探測(cè)效率較低時(shí)，為了提高核素分析能力，降低測(cè)量數(shù)據(jù)中的統(tǒng)計(jì)漲落，通常采取增加測(cè)量時(shí)間、加大探測(cè)器尺寸等方法。實(shí)際測(cè)量中，便攜式應(yīng)用中探測(cè)器一般較小，而測(cè)量時(shí)間要求短；在實(shí)驗(yàn)室應(yīng)用中，探測(cè)的樣品放射性比較弱；航空γ能譜測(cè)量中，飛機(jī)飛行速度快測(cè)量時(shí)間短。在這種情況下，通過數(shù)學(xué)的方法對(duì)γ能譜數(shù)據(jù)進(jìn)行處理，降低統(tǒng)計(jì)漲落，有利于快速準(zhǔn)確地分析核素，這是近年來γ能譜數(shù)據(jù)處理的一個(gè)發(fā)展方向。

1 Bootstrap方法

1979年，美國(guó)斯坦福大學(xué)統(tǒng)計(jì)學(xué)教授EfronB[1]在總結(jié)和歸納前人研究成果的基礎(chǔ)上，提出了統(tǒng)計(jì)自舉方法(Bootstrap)。Bootstrap方法能夠很好地處理小樣本及未知分布模型樣本的統(tǒng)計(jì)置信度、中位數(shù)和標(biāo)準(zhǔn)方差等，其優(yōu)點(diǎn)受到了統(tǒng)計(jì)界的重視。隨著計(jì)算機(jī)技術(shù)的快速發(fā)展，該方法已經(jīng)深入到不同的學(xué)科領(lǐng)域，在諸如圖像處理、射線層析成像處理、非線性動(dòng)力學(xué)中相關(guān)維數(shù)的估計(jì)、非平穩(wěn)信號(hào)分析、地震數(shù)據(jù)分析等領(lǐng)域都得到了成功的應(yīng)用。

目前，國(guó)內(nèi)外對(duì)Bootstrap方法用于γ能譜分析的研究較少。20世紀(jì)80年代，美國(guó)新罕布什大學(xué)空間研究中心的Simpson G[2]將Bootstrap非參數(shù)法用于估計(jì)宇宙γ射線源位置及宇宙γ射線特征分析。國(guó)內(nèi)防化研究院的劉穎[3]對(duì)152Eu的121.78 keV短時(shí)間內(nèi)實(shí)測(cè)特征γ峰進(jìn)行了Bootstrap抽樣處理，得到的處理結(jié)果與長(zhǎng)時(shí)間的實(shí)測(cè)能譜相近，但她沒有對(duì)其他能量范圍、其他核素的能譜進(jìn)行研究，缺乏代表性。繼續(xù)深入研究Bootstrap方法在γ能譜分析中的應(yīng)用，對(duì)提高探測(cè)系統(tǒng)的探測(cè)效率、增強(qiáng)放射性核素的識(shí)別能力具有重要意義。

圖1 Bootstrap非參數(shù)法流程圖Fig.1 Flowchart of Bootstrap non-parametric method.

Bootstrap非參數(shù)法的基本原理：假設(shè)總體的分布F未知，但已有一個(gè)容量為n的來自分布F的數(shù)據(jù)樣本x=(x1, x2, x3, …, xn)，從此樣本按放回抽樣的方法抽取一個(gè)容量為n的樣本x=(x1, x2, x3, …, xn)，這種樣本被稱為Bootstrap樣本。相繼地、獨(dú)立地自原始樣本中取得B(B≥1000)個(gè)容量為n的Bootstrap樣本，分別計(jì)算這B個(gè)Bootstrap樣本的均值和標(biāo)準(zhǔn)誤差，并將它們按從小到大的順序排列，按照要求的置信水平，求得均值及標(biāo)準(zhǔn)誤差的置信區(qū)間。對(duì)滿足置信區(qū)間的Bootstrap樣本進(jìn)行整理分析，利用這些樣本對(duì)總體F進(jìn)行統(tǒng)計(jì)推斷，這種方法被稱為Bootstrap非參數(shù)方法[4]，用非參數(shù)方法來求參數(shù)的近似置信區(qū)間的優(yōu)點(diǎn)是不需要對(duì)總體分布的類型作任何的假設(shè)，而且適用于小樣本，且能用于各種統(tǒng)計(jì)量。根據(jù)上述原理，設(shè)計(jì)Bootstrap非參數(shù)法處理γ能譜數(shù)據(jù)的程序流程，如圖1所示。

2 實(shí)驗(yàn)測(cè)量

用NaI(Tl)探測(cè)器分別測(cè)量了核素241Am (3.77μCi)、137Cs (2.87 μCi)和60Co (2.23 μCi)的γ能譜，探測(cè)距離為25 cm，測(cè)量時(shí)間分別為600 s和13000 s，扣除本底后，選取測(cè)量時(shí)間為600 s的59.54keV、661.2 keV、1.173 MeV和1.332 MeV特征峰作為原始小樣本。實(shí)驗(yàn)使用的γ能譜儀由濱松公司生產(chǎn)的7.62 cm×7.62 cm NaI(Tl)探測(cè)器、Inspector 2000和Gennie2000軟件組成，系統(tǒng)的標(biāo)稱能量分辨率為7.5% (0.662 MeV)，探測(cè)射線能量為0-3MeV。

3 用Bootstrap非參數(shù)法處理實(shí)驗(yàn)數(shù)據(jù)

根據(jù)Bootstrap方法原理的實(shí)現(xiàn)要求，采用統(tǒng)計(jì)分析軟件R語(yǔ)言進(jìn)行編程。選用600 s測(cè)量時(shí)間的實(shí)驗(yàn)數(shù)據(jù)作為原始小樣本進(jìn)行Bootstrap非參數(shù)法抽樣分布分析，并與測(cè)量時(shí)間為13000 s的實(shí)測(cè)譜進(jìn)行對(duì)比判斷，探討此方法在實(shí)際測(cè)量中用于復(fù)雜核素識(shí)別的可能性。

3.1 用Bootstrap非參數(shù)法處理241Am和137Cs的γ能譜特征峰

圖2是241Am的實(shí)測(cè)譜，虛線譜的測(cè)量時(shí)間為600 s，實(shí)線譜的測(cè)量時(shí)間為13000 s，后者統(tǒng)計(jì)性明顯好于前者。

圖2 600 s和13000 s241Am實(shí)測(cè)譜對(duì)比Fig.2 600 s and 13000 s measured spectra of241Am.

采用Bootstrap非參數(shù)法處理測(cè)量時(shí)間為600 s的實(shí)測(cè)譜數(shù)據(jù)，求出置信水平為0.9的置信區(qū)間后，進(jìn)行25次放回抽樣，對(duì)滿足置信區(qū)間的樣本進(jìn)行分析后統(tǒng)計(jì)加和，得到計(jì)算譜如圖3所示。由圖3(a)，測(cè)量時(shí)間為600 s的計(jì)算譜與測(cè)量時(shí)間為13000 s的實(shí)測(cè)譜在峰形上很吻合，只是實(shí)測(cè)譜的峰值偏高，適當(dāng)增加抽樣次數(shù)后，得到圖3(b)的計(jì)算譜，此時(shí)的計(jì)算譜與實(shí)測(cè)譜無論是在峰形上還是峰位位置或者峰位計(jì)數(shù)上都非常符合。

用同樣方法處理137Cs的661.2 keV特征峰，得到類似的結(jié)果，即短測(cè)量時(shí)間的實(shí)測(cè)譜的計(jì)算譜與長(zhǎng)測(cè)量時(shí)間的實(shí)測(cè)譜能很好地吻合(圖4)。

圖3 增加抽樣次數(shù)前(a)后(b)241Am實(shí)測(cè)譜和計(jì)算譜對(duì)比Fig.3 Measured and calculated spectrum of241Am before (a) and after (b)increase the number of sampling.□ Measured spectra, · Calculated spectra

圖4 137Cs實(shí)測(cè)譜與計(jì)算譜對(duì)比Fig.4 Measured and calculated spectrum of137Cs.□ Measured spectra, · Calculated spectra

3.2 用Bootstrap非參數(shù)法處理60Co的γ能譜特征峰

對(duì)于測(cè)量時(shí)間為600 s的60Co的小樣本，第一個(gè)峰位于395道，第二個(gè)峰位于447道。由于實(shí)驗(yàn)中用的探測(cè)器是NaI(Tl)探測(cè)器，NaI(Tl)晶體的發(fā)光效率、譜儀系統(tǒng)中的高壓電源、光電倍增管倍增系數(shù)、放大器的放大倍數(shù)以及甄別器的閾值等都會(huì)隨溫度和時(shí)間而變化，造成譜的漂移，從而導(dǎo)致測(cè)量誤差[5]。本實(shí)驗(yàn)中長(zhǎng)時(shí)間測(cè)量13000 s后60Co的γ能譜發(fā)生了漂移，60Co的第一個(gè)峰漂移到了392道，第二個(gè)峰漂移到445道。因此，在使用Bootstrap非參數(shù)法計(jì)算得到小樣本的計(jì)算譜后，要根據(jù)實(shí)際測(cè)量條件和譜儀系統(tǒng)的譜漂情況對(duì)峰位進(jìn)行校正，同時(shí)還要考慮到能量接近的相鄰能譜的干擾。241Am和137Cs的γ能譜在13 000 s時(shí)間內(nèi)的譜漂不顯著,原因是其活度較大，測(cè)量時(shí)光電倍增管的工作電壓相對(duì)高，譜漂的影響小，故不需要峰位校正。圖5是針對(duì)60Co的兩個(gè)特征峰在短測(cè)量時(shí)間實(shí)測(cè)譜的計(jì)算譜與長(zhǎng)測(cè)量時(shí)間實(shí)測(cè)譜的對(duì)比，經(jīng)過峰位校正后，計(jì)算譜和實(shí)測(cè)譜相當(dāng)吻合。

圖5 峰位校正前(a)和峰位校正后(b)60Co實(shí)測(cè)譜與計(jì)算譜對(duì)比Fig.5 Measured and calculated spectrum of60Co before (a) and after (b) peak position correction.□ Measured spectra, · Calculated spectra

表1列出了241Am、137Cs、60Co核素特征峰的實(shí)測(cè)譜與計(jì)算譜的峰位和峰面積，峰面積計(jì)算采用總峰面積法(TPA)。由表1，經(jīng)Bootstrap非參數(shù)法處理的短測(cè)量時(shí)間實(shí)測(cè)譜的計(jì)算譜與長(zhǎng)時(shí)間的實(shí)測(cè)譜相比，峰面積相對(duì)誤差小于5%，其中特征峰241Am (59.54 keV)、137Cs (661.2 keV)及60Co (1.332MeV)的峰面積誤差小于等于1%。對(duì)于60Co的兩個(gè)特征峰，經(jīng)峰位校正后，60Co (1.332 MeV)實(shí)測(cè)譜與計(jì)算譜的峰位一致，而60Co (1.173 MeV)實(shí)測(cè)譜與計(jì)算譜的峰位有1道的誤差，相應(yīng)峰面積的相對(duì)誤差也較大，為4.91%。

表1 實(shí)測(cè)譜與計(jì)算譜峰位與峰面積對(duì)比表Table1 Measured and calculated spectrum peak position and peak area comparison.

4 實(shí)驗(yàn)小結(jié)

通過以上計(jì)算及實(shí)驗(yàn)結(jié)果驗(yàn)證得到以下結(jié)論：

(1) 使用Bootstrap非參數(shù)法可以減少統(tǒng)計(jì)漲落對(duì)能譜測(cè)量的時(shí)間要求。在用Bootstrap非參數(shù)法處理這三個(gè)核素的小樣本數(shù)據(jù)時(shí)，小樣本的測(cè)量時(shí)間都是600 s，為得到與對(duì)應(yīng)的測(cè)量時(shí)間為13000 s的實(shí)測(cè)譜可比擬的計(jì)算譜，三種核素的抽樣次數(shù)稍有不同，得到計(jì)算譜的計(jì)算時(shí)間也不同。影響計(jì)算時(shí)間的主要因素是樣本容量，三個(gè)核素的小樣本數(shù)據(jù)中，樣本容量最大的是60Co，所用抽樣計(jì)算的時(shí)間也最長(zhǎng)(約5 s)，加上小樣本的測(cè)量時(shí)間，實(shí)際用時(shí)約605 s。對(duì)于實(shí)驗(yàn)的實(shí)測(cè)時(shí)間13000 s而言，大大縮短了實(shí)際測(cè)量時(shí)間。

(2) 使用Bootstrap非參數(shù)法可克服γ譜儀系統(tǒng)能譜漂移的不利影響。由于NaI(Tl)閃爍晶體的發(fā)光效率受溫度影響，測(cè)量時(shí)間較長(zhǎng)時(shí)，譜漂的不利影響較大。而利用Bootstap非參數(shù)法對(duì)短時(shí)間測(cè)得的實(shí)驗(yàn)譜處理后，得到的計(jì)算譜的統(tǒng)計(jì)性接近長(zhǎng)測(cè)量時(shí)間實(shí)測(cè)譜的統(tǒng)計(jì)性，即通過數(shù)據(jù)處理的方法縮短測(cè)量時(shí)間，這樣就可能避免長(zhǎng)時(shí)間測(cè)量過程中溫度條件變化所造成的峰位漂移等現(xiàn)象，為快速準(zhǔn)確地進(jìn)行γ能譜分析提供可行的技術(shù)方案。

5 結(jié)語(yǔ)

由于Bootstrap非參數(shù)法的抽樣是放回抽樣，每次的抽樣都是等概率隨機(jī)抽樣，但每次抽取的Bootstrap樣本不一定滿足置信區(qū)間的要求，因此要想得到與實(shí)測(cè)譜在峰值上比較吻合的能譜，每次的抽樣次數(shù)略有不同，應(yīng)進(jìn)一步研究抽樣次數(shù)和置信區(qū)間之間的關(guān)系。

1 Efron B, Tibshirani R J. Bootstrap methods for standard errors, confidence interval, and other measures of statistica1 accuracy[J]. Statistical Science, 1986, 1(1): 54-77

2 Simpson G, Hasselwander H M. Bootstrap sampling: applications in gamma-ray astronomy[J]. Astronomy and Astrophysics, 1986, 162(1-2): 340-348

3 劉穎. 基于γ能譜測(cè)量的軍用放射性核素分析方法研究[D]. 北京: 防化研究院, 2011

LIU Ying. The research in radionuclide analysis method applied to military by gamma energy spectrum measurement[D]. Beijing: Research Institute of Chemical Defence, 2001

4 Mooney C Z, Duval R D. Bootstrapping: A nonpararnetric approach to statistical inference[M]. California: Sage Publication, 1993: 7-95

5 屈國(guó)普, 凌球, 郭蘭英, 等. NaI(Tl)閃爍譜儀譜漂移原因分析[J]. 南華大學(xué)學(xué)報(bào), 2005, 19(1): 47-50

QU Guopu, LING Qiu, GUO Lanying, et al. Analysis of the cause for the spectrum-drifting of NaI(Tl) scintillation spectrometer[J]. Journal of University of South China, 2005, 19(1): 47-50

Application of non-parametric Bootstrap method togamma spectrum analysis

XU Hongjuan1,2WANG Feng1,2AI Xianyun1,2LI Jinglun1,2WEI Xing1,2SHI Lei1,2
1(Research Institute of Chemical Defense, Beijing 102205, China)
2(State Key Laboratory of NBC Protection for Civilian, Beijing 102205, China)

Background: In gamma spectral measurement, if the sample activity or detection efficiency of the detector is low, the most often used method to reduce the statistical fluctuation of the measurement data is to increase the measurement time and the detector dimensions. Purpose: Considering the economic factors as well as the matching problem with other nuclear electronics devices, both the size of the detector and the measurement time are limited. In this case, processing gamma spectrum data by the mathematical method to reduce statistical fluctuations for fast and accurate analysis of radionuclides received widespread attention at home and abroad. Methods: The basic principles of non-parametric Bootstrap method was described and applied to laboratory gamma spectrum data processing. The gamma spectrum of241Am,137Cs and60Co were measured in different time periods by NaI(Tl) detectors. Results: The non-parametric Bootstrap method was used to process the gamma spectra measured in short time and the results were compared with the gamma spectra data measured in long time under the same conditions, and the calculated spectra agreed well with the measured spectra. Conclusion: It provides a feasible technique to quickly measure gamma spectrum at low activity levels.

Bootstrap method, Non-parametric method, Confidence interval, Gamma spectrum

TL816+.2

10.11889/j.0253-3219.2014.hjs.37.010501

徐紅鵑，女，1976年出生，2005年于廈門大學(xué)獲碩士學(xué)位，現(xiàn)為防化研究院在職博士研究生，輻射防護(hù)及環(huán)境保護(hù)

2013-09-22，

2013-11-18

CLCTL816+.2