剖面似然比統(tǒng)計及CLs方法在PandaX-II 暗物質實驗數據分析中的應用

阿布都沙拉木·阿布都克力木，諶 勛，吾爾尼沙·依明尼亞孜

(1. 新疆大學 物理科學與技術學院，烏魯木齊 830044； 2. 上海交通大學 物理與天文學院，上海 200240)

近年來，眾多宇宙學和天文觀測證據都強烈暗示著宇宙中存在大量的暗物質[1].理論物理學家對暗物質的本質有若干猜想，其中最為流行的理論之一[2-3]認為暗物質是一種參與某種弱相互作用的大質量粒子(WIMP).在極低本底環(huán)境下，利用高靈敏的粒子探測器對暗物質進行直接探測是研究WIMP的重要手段.

在過去的三十多年時間內，國際上進行了多項暗物質直接探測實驗研究，但是除了個別被廣泛質疑的信號外，大部分實驗都沒有找到暗物質的跡象[4-7].近年來，基于二相型液氙技術的暗物質直接探測實驗將暗物質探測的靈敏度提高了多個數量級，其中由我國主導的PandaX-II 暗物質直接探測實驗更是基于當時世界上最大的5.4×104kg·d的數據集，于2017年給出了對質量大于100 GeV/c2(1 eV/c2=1.782×10-36kg)的暗物質與普通物質自旋無關散射截面的最強限制[6].

自2017年以來，PandaX-II 實驗還分別發(fā)表了關于暗物質與普通物質自旋相關散射、軸子暗物質及非彈性散射暗物質的研究結果[8-10].由于數據分析沒有找到超出本底的信號，因此僅給出了相應的排除曲線.在這些分析中，剖面似然比(PLR)[11]被用于構造檢驗統(tǒng)計，并結合一種被稱為CLs的方法[12]來估計需要檢驗的物理量的可能取值.本文以暗物質與普通物質的非彈性散射過程為例，結合PandaX-II 實驗在2016年獲取的數據，詳細討論如何運用PLR及CLs方法來計算暗物質與核子的等效彈性截面上限.

1 PandaX-II 實驗及其數據分析

PandaX-II 實驗是一個運用二相型(氣/液)氙時間投影室(TPC)技術來進行暗物質直接探測的實驗[13].PandaX-II 使用的TPC大致為一個直徑約60 cm，高約60 cm的圓柱形，內部裝有大約580 kg的液氙.TPC上下各有一個陣列，由55個直徑3 inch(1 inch=2.54 cm)的光電倍增管(PMT)組成，用于探測閃爍光信號.粒子在液氙里面與氙原子核或核外電子發(fā)生散射，這個過程中部分轉移的能量會產生可以被PMT陣列捕獲的閃爍光信號(通常被稱為S1)；還有一部分能量將氙電離，產生自由電子，這些電子在TPC的強電場作用下漂移到氣態(tài)氙中，并產生可以被PMT接收的正比閃爍光信號(通常被稱為S2).采用PMT接收到的光電子(PE)數目來表示光信號的強度.根據一般的暗物質理論，WIMP與原子核散射的截面遠大于其與電子散射的截面，因此它在探測器里面發(fā)生核反沖的概率更大[14].核反沖過程中產生的自由電子更容易發(fā)生重結合，這會增大S1并減少S2；該特征可以用于分辨電子反沖的本底[15].

PandaX-II 實驗的數據分析目標為在理論建議的反沖能量區(qū)間內尋找出可能的超出本底預期的核反沖信號.為了完成這個目標，探測器記錄到的信號需要根據刻度得到的參數進行均勻性修正，并通過一系列數據質量篩選，最終得到待定信號區(qū)間內的候選粒子，并運用統(tǒng)計學的方法對候選粒子進行檢驗.關于數據修正和篩選的細節(jié)，可以參考文獻[4].

在運用統(tǒng)計方法對數據進行檢驗之前，需要理解本底信號的分布.PandaX-II 中的本底主要包括如下來源：① 來自于探測器部件內的放射性同位素雜質釋放出的伽馬光子在液氙內產生的電子反沖事件；② 液氙內部氪、氡等放射性同位素雜質和氙的放射性同位素衰變引起的電子反沖事件；③ 探測器部件內釋放出的中子在液氙內產生的核反沖事件；④ 物理上沒有關聯的S1和S2由于偶然符合形成的“假”事件.使用放射源對探測器進行刻度的分析，結果表明電子反沖和核反沖事例之間存在較大的區(qū)分度(見圖1).此外，本底水平也需要通過其他方法進行初步估計.例如，來自探測器材料的本底貢獻可以基于測得的探測器材料的放射性強度，進行Monte Carlo模擬來進行估計，氙里面的85Kr放射性同位素本底則可以通過測得的85Kr衰變到85Rb過程中的β-γ符合過程的數目來進行估計[4].

圖1 PandaX-II 實驗 2016年核反沖與電子反沖的事例分布Fig.1 The distribution of nuclear recoil events and electron recoil events in PandaX-II 2016 dataset

一般來說，如果在暗物質的直接探測實驗中沒有觀察到明顯超出本底預期的信號，那么數據分析的最終結果是給出暗物質與普通物質發(fā)生相互作用的截面上限.銀河系內的暗物質密度為一個由天文學所給出的固定輸入參數[16]，因此實驗給出的截面上限是一個和暗物質質量相關的值，最終結果表現為在一定的暗物質質量區(qū)間內所給出的一條連續(xù)曲線，即所謂的排除線.目前基于液氙的幾個實驗的最新結果都是基于PLR構造統(tǒng)計量，從而給出了對應的排除線[4-7].

2 剖面似然比方法

粒子物理學家通常運用頻率論的統(tǒng)計檢驗或者貝葉斯檢驗來探索新物理[17].在頻率論的框架下，當處理設置上限問題時，通常將測得數據同時包含本底和信號的假設作為待檢驗的零假設(H0)，并將數據中僅有本底的假設作為備擇假設(H1).一個假設H和觀測到的數據的符合程度是通過計算P值來估計的，即在該假設成立的情況下，得到觀測數據或者與H偏差更大的數據的概率.通常當計算得到的P值低于某一給定閾值時，可以認為假設被排除.

(1)

為計算方便，一般使用的檢驗統(tǒng)計被定義為

tμ=-2lnλ(μ)

(2)

由前文可知，tμ越大意味著數據和假設之間的偏差也越大.用于量化數據和模型之間的差值的P值可以定義為

(3)

式中：f(tμ|μ)是在信號強度μ固定的情況下進行多次實驗得到的tμ的分布，通常需要進行大量的Monte Carlo模擬來得到；tμ,ob則根據實際的觀測數據得到.

在使用剖面似然比方法設置信號強度的上限時，統(tǒng)計檢驗被進一步定義為

(4)

(5)

這樣定義的P也通常被表示為Ps+b，說明這是針對同時存在信號(s)和本底(b)所做的檢驗.

通常，如果要以1-α的置信水平來設置被檢驗的物理量μ的上限，就需要在使得Pμ≥α的μ的可能取值當中找到一個最大的μlimit.

3 CLs方法

一些實驗需要檢驗的信號強度的數值非常小，會得到相近的f(qμ|μ)與f(qμ|0)的分布，從而在接受某個μ值時，其對應的P值與信號不存在(μ=0)時的P值相當接近，導致實驗在該μ值附近已缺乏靈敏度.

為了解決這個問題，CLs方法[17]被用來對P值進行修正.該方法需要計算在僅有本底假設下的P值，即

新修正值被定義為

對于某一個特定的置信水平1-α，若根據某個信號強度μ計算得到CLs≤α,則拒絕這個假設.由于1-Pb<1，CLs值通常大于對應的Ps+b，從而能得到更加保守的結果.由此可見，CLs并非置信水平(Confidence Levels)，它只是一個出于實用主義而給出的名稱，該方法在高能物理實驗里面得到了廣泛運用.PandaX-II 在發(fā)布2016年數據的幾個分析結果都使用了CLs方法來給出對應物理量限制.

4 PandaX-II 對暗物質非彈性散射截面的限制

4.1 WIMP與核子的非彈性散射

為了解釋DAMA/LIBRA實驗中的疑似暗物質信號[18]，一些理論模型考慮暗物質存在多個能級，能級間具有質量劈裂δ，從而暗物質會和原子核發(fā)生非彈性散射，從低能態(tài)躍遷到高能態(tài)[19-21].由于質量劈裂的存在，非彈性散射需要一個最低的暗物質運動速度的存在，相互作用的相空間也被壓縮.在文獻[10]中，PandaX-II 基于2016年獲取的數據來對暗物質與核子的彈性散射等效截面進行了限制.

根據暗物質與核子發(fā)生非彈性散射的微分事例率計算公式，并考慮到部分參數來自天文學的輸入，最終暗物質非彈散射的事例率依賴于3個參數，分別是暗物質粒子的質量mχ, 暗物質粒子的質量劈裂δ以及彈性散射極限下的零動量暗物質與核子的作用截面(以下簡稱等效截面)σn[10]，其中σn為該分析中的目標物理量.選取暗物質質量mχ分別為1和10 TeV/c2作為參考質量，分析結果將得到在不同δ情況下σn的上限.

4.2 似然函數的構造

為了構造對應的檢驗統(tǒng)計，需要寫出對應的似然函數.由于PandaX-II 首批實驗數據包含了不同條件下所得的數據，可根據條件劃分為nset個數據集，因此似然函數可以寫為不同數據集各自似然函數的乘積形式

其中：函數G為高斯函數，表示暗物質的測量誤差δdm和本底的測量誤差δbj的影響，這些誤差也作為待擬合的冗余參數；bj表示不同種類的本底；σdm為信號截面，即我們所關心的等效截面；σbj則表示不同本底的截面，由相關數據分析給出.在PandaX-II 實驗首批數據里面，主要的本底分為5類：① 液氙內溶解的放射性85Kr 產生的電子反沖本底；② 放射性127Xe 產生的電子反沖本底；③ 其他來源的電子反沖本底；④ 無關的S1及S2偶然符合導致的本底；⑤ 來自探測器組件內放射性雜質產生的中子本底.具體來說，每個數據集的似然函數則可寫成

4.3 事例S1和S2的聯合分布

對于暗物質信號來說，可由公式得到其能量分布，而根據能量生成S1和S2的聯合分布需要使用Noble Element Simulation Technique(NEST)模型[22]，該模型描述了液氙內沉積的能量如何轉化為信號S1和S2，在液氙暗物質實驗中得到廣泛應用，模型中的一些自由參數需要使用刻度數據去進行優(yōu)化.對于電子反沖和中子反沖的本底，使用基于Geant4[23]的Monte Carlo模擬得到能譜，同樣運用NEST模型獲得其S1及S2的聯合分布.偶然符合的本底的聯合分布通過將實驗數據分析得到的單S1事例信號和單S2事例信號進行大量隨機匹配而得到.

對于不同的暗物質質量及質量劈裂，信號能譜具有不同的形狀，從而得到的S1和S2分布也不盡相同.不同的等效截面則會影響信號窗口內的總事例數.圖2顯示了在暗物質質量為1 TeV/c2及等效截面為10-40cm2時，不同質量劈裂給出的暗物質事例在信號窗口內的聯合分布；采用歸一化事例密度表示，顏色越深，密度越高，則事例數越多.

mχ=1 TeV/c2, σn=10-40 cm2圖2 不同的質量劈裂對應的暗物質在PandaX-II 探測器中的信號聯合分布Fig.2 The dark matter signal distribution in PandaX-II detector at different mass splitting values

4.4 應用PLR方法擬合候選事例

在首期79.6天內，PandaX-II 實驗記錄到了 24 502 402 個事例.經過挑選后，在3 PE≤S1≤100 PE，100 (raw) PE≤S2≤12 000 PE的信號窗口內保留 716 個事例(其中raw表示原始未修正信號)，其分布如圖3所示.未觀測到明顯超出本底的信號.

圖3 PandaX-II 實驗2016年數據經過篩選后在信號窗口中的分布Fig.3 The distribution of the candidate events from PandaX-II 2016 data in the signal window

假設暗物質質量為1 TeV/c2，在δ=0～300 keV/c2的區(qū)域內平均選取20個點.對于每個被選取的δ，計算預期得到1個暗物質事例所需要的等效截面，在等效截面附近選取一系列掃描點σdm，并針對每個截面掃描點，用所有保留的候選事例來計算其對應的檢驗統(tǒng)計qμ,ob.

為計算Ps+b，需要根據掃描點σdm，擬合得到本底的δbj，從信號及本底的S1和S2聯合分布出發(fā)進行大量Monte Carlo模擬，生成模擬事例，使用這些模擬事例對選取的掃描點進行擬合，以得到信號+本底假設下qμ的分布f(qμ|μ)，從而計算出Ps+b.圖4給出了當暗物質δ=120 keV/c2，假設等效散射截面σdm=7.615×10-42cm2時qμ的概率密度分布；根據數據計算出對應的qμ,ob=3.983，從而得到Ps+b=0.032.

同時需要根據擬合得到的本底δb，使用本底的聯合分布來進行Monte Carlo模擬，生成僅有本底的模擬事例，并使用這些模擬事例對選取的掃描點進行擬合，以得到僅有本底假設下的qμ分布f(qμ|0)，并計算出1-Pb.圖4中同樣給出示例參數在僅有本底假設下的qμ分布，由該分布可以得到對應的1-Pb=0.148.

δ=120 keV/c2, σdm=7.615×10-42 cm2圖4 qμ的概率密度分布Fig.4 The probability density distribution of qμ

這樣，我們可得到每個掃描點σdm，qμ,ob及CLs值之間的對應關系.通過在CLs值接近0.1的附近迭代選取更多掃描點并重復上述過程，使最終計算得到的CLs值與0.1的差值小于某一預設誤差，從而得到在對應δ置信水平為90%的等效截面上限.最終將這些對應不同δ的σ上限繪制在一張圖上，就得到了置信水平為90%的等效截面排除線，如圖5中紅色實線所示.

圖5 PandaX-II 實驗2016年數據對于暗物質與核子非彈散射過程中的等效截面的90%置信水平排除曲線Fig.5 Upper limits at 90% confidence level on the effective cross section of WIMP-nucleon in inelastic scattering set by the PandaX-II 2016 data

使用CLs方法時通常會同時給出誤差帶，計算方式為利用僅有本底假設的Monte Carlo模擬數據，運用上文方法得出置信水平為90%的等效截面.多組模擬得到截面分布的一倍及二倍標準差作為一倍及二倍標準差的誤差帶上下限.圖6給出了對應于暗物質質量劈裂為120 keV/c2時，根據僅有本底假設的Monte Carlo數據所計算得到的90%置信度的截面分布.因為誤差帶是由Monte Carlo模擬給出，而排除線是由數據給出，所以數據的漲落可能導致部分區(qū)域排除線落在一倍標準差之外.

圖6 由僅有本底假設的蒙卡模擬數據所計算得到的90%置信度的截面上限分布(δ=120 keV/c2)Fig.6 The distribution of upper limits of the effective cross section at 90% confidence level calculated from MC simulation with background-only hypothesis (δ=120 keV/c2)

5 結語

PandaX-II 實驗得到的數據包含了S1和S2的值.為了充分利用這些信息，PandaX-II 的數據分析結合應用了PLR和CLs方法，對暗物質的截面性質給出限制.本文基于PLR方法中檢驗統(tǒng)計的構造，具體介紹了PandaX-II 實驗中構造PLR所需的似然函數如何從數據中得到，以及如何計算最終的排除曲線及誤差帶，希望能為其他類似實驗或者需要運用該方法的領域提供參考.

致 謝本文主要工作來源于季向東教授、劉江來教授領導的PandaX-II 暗物質直接探測實驗，該課題得到上海市粒子物理與宇宙學重點實驗室、教育部粒子物理與星系宇宙學重點實驗室的支持.在此表示感謝!