索瑞鑫 仇玉蘭 王 彤△
微核試驗是公認的檢測染色體或有絲分裂器損傷的一種試驗方法,其最大的優(yōu)點是經(jīng)濟、簡單、快速。該實驗結果通常以嚙齒類動物的骨髓嗜多染紅細胞或者人外周淋巴細胞中的微核率來表示〔1〕。
國際上有多位生物統(tǒng)計學家對微核實驗的統(tǒng)計分析方法進行過討論,但并未形成統(tǒng)一的推薦方法〔2〕。同時我們還注意到,EPA、OECD、ECC、UKEMS等國際上關于微核試驗的指南里也沒有明確統(tǒng)一的數(shù)據(jù)分析方法,大多都是類似于“應使用合適的統(tǒng)計方法分析數(shù)據(jù)”這類陳述〔2〕。為了解目前對這類數(shù)據(jù)所采用的統(tǒng)計分析方法現(xiàn)況,我們利用PubMed-NLM文獻數(shù)據(jù)庫進行了簡單的檢索。限定題目中包含“micronucleus”,年限從1995年到2010年,語言為英文,得到文獻總數(shù)959篇,可得到的全文文章152篇,涉及52種雜志。除去其中關于微核檢測技術、關聯(lián)性研究或僅作統(tǒng)計描述而非組間比較問題的文獻,最終103篇文章中使用t檢驗35篇,方差分析23篇,非參數(shù)秩和檢驗36篇,卡方檢驗13篇,F(xiàn)isher確切概率法7篇,Kastenbaum and Bowman檢驗4篇。這一匯總結果雖然不夠全面,但確實反映出在統(tǒng)計分析方法選擇上的不一致。
通常而言參數(shù)方法的統(tǒng)計功效會高于非參數(shù)方法,但需正確指定其分布假設。而關于微核數(shù)據(jù)的分布假設也有一些不同的討論,如二項分布、泊松分布或負二項分布等〔3〕。由于微核數(shù)的發(fā)生較少(尤其是包含對照組時),所以多數(shù)學者傾向于可以用描述稀有事件發(fā)生數(shù)的泊松分布來假定。然而實際數(shù)據(jù)中泊松分布均數(shù)等于方差的條件有時不一定會被滿足,而負二項分布可以用離散參數(shù)來描述這種過離散現(xiàn)象,并且當離散參數(shù)取極端值的時候,負二項分布就退化為泊松分布,所以將微核數(shù)據(jù)假定為負二項分布可以合理地考慮到動物個體間變異的客觀存在〔4〕。雖然利用隨機效應模型也可以處理這種動物個體間變異,然而,較少的微核數(shù)并不滿足大樣本的漸近性條件(包括卡方檢驗等方法和基于廣義線性混合模型的方法),故而考慮一些確切推斷方法是必要的?;诓此煞植肌⒍椃植嫉拇_切推斷方法存在已久,且在不同方面也有所進展,但基于負二項分布的確切推斷方法最近才被實現(xiàn)〔5〕。故以下我們將介紹基于負二項分布的條件確切概率方法并通過對文獻中常見方法進行統(tǒng)計模擬來比較其Ⅰ類錯誤和檢驗功效。
若Y服從均數(shù)為μ,離散參數(shù)為φ的負二項分布,記為Y~NB(μ,φ),則其概率函數(shù)為:
Y的期望E(Y)=μ,方差 var(Y)=μ+φμ2。當 φ→0時,均數(shù)等于方差,負二項分布則退化為泊松分布。在微核數(shù)據(jù)中,這里的離散參數(shù)φ可描述動物個體間的變異。
在微核實驗中假設每個動物的微核數(shù)Y1,…,Yn相互獨立,服從負二項分布NB(μ=miλ,φ),其中mi表示每個動物的細胞數(shù),λ表示微核率,如果每個動物細胞數(shù)相等(如mi都等于2000),則Y獨立同分布,此時λ的最大似然估計為微核總數(shù)除以細胞總數(shù),與φ無關。
如果所有動物的細胞數(shù)都是一樣的(即mi≡m),那么它們的計數(shù)之和也服從負二項分布,即Z=Y1+…+Yn~NB(nmλ,φn-1)。這里 λ 是充分統(tǒng)計量,λ與φ彼此獨立,可以通過讓(2)式的條件似然值取最大來確定φ值,從而通過(1)式得到推斷的確切概率〔6〕:
若ZtA和ZtB分別是A、B兩組的計數(shù)總和,動物數(shù)分別為 nA和 nB。則在無效假設 ZtA~NB(nkmλ1,),k∈A,B下,可以建立Fisher精確概率的列聯(lián)表,計數(shù)為ZtA+ZtB也是服從負二項分布的隨機變量,我們可以運用公式(1)負二項分布的精確概率計算公式得到比觀測值更極端的概率,即雙側的P值。
我們按照毒理學微核實驗有關的指南文獻提供的依據(jù)對涉及到的模擬參數(shù)進行如下設定:動物數(shù)n的模擬參數(shù)設定為5到10,細胞數(shù)的模擬參數(shù)設定為2000,微核陽性數(shù)的總體差值δ為0到18。所比較的方法包括Poisson精確概率法、Fisher精確概率法、負二項精確概率法、t檢驗及其變換、確切概率秩和檢驗。
由于我們研究的是mi≡m,即每組細胞數(shù)相等的情況,離散參數(shù)φ是獨立的變量,在模擬研究中我們采用了固定離散參數(shù)φ的方法,其值分別設置為1和0.001,當φ→1時,模擬數(shù)據(jù)服從負二項分布;φ=0.001時,模擬數(shù)據(jù)近似服從Poisson分布。
對Ⅰ類錯誤的模擬估計結果見表1,2。
表1 Ⅰ類錯誤及排序(φ→0)
表2 Ⅰ類錯誤及排序(φ→1)
從表1與表2可以看到無論φ為1還是接近0,Poisson精確概率法的Ⅰ類錯誤是最大的,且遠大于規(guī)定的檢驗水準α=0.05。φ→0,即數(shù)據(jù)來自于Poisson分布時,平方根變換t檢驗與平方根反正弦變換t檢驗以及t檢驗的Ⅰ類錯誤大致相同,都僅次于Poisson精確概率法;經(jīng)log變換t檢驗Ⅰ類錯誤是最小的;Fisher精確概率法與負二項精確概率法(NB),及秩檢驗的Ⅰ類錯誤大致相同,均略高于log變換的t檢驗,但遠遠小于Poisson精確概率法的Ⅰ類錯誤。φ→1時,F(xiàn)isher精確概率法的Ⅰ類錯誤較高,僅次于Poisson精確概率法;經(jīng)平方根變換的t檢驗與平方根反正弦變換得到t檢驗的Ⅰ類錯誤大致相同,僅次于Fisher精確概率法;經(jīng)log變換的t檢驗Ⅰ類錯誤僅次于平方根變換的t檢驗;而秩檢驗、t檢驗、負二項精確概率法的Ⅰ類錯誤均小于規(guī)定的檢驗水準α=0.05。由此結果可見,在φ→0時,Poisson精確概率法的Ⅰ類錯誤遠遠超過檢驗水準,故不可接受;在φ→0時,Poisson精確概率法、Fisher精確概率法的Ⅰ類錯誤也遠遠超過檢驗水準α=0.05,為不可接受的方法;其他方法尚可。
當離散參數(shù)φ→0時,檢驗功效結果見表3。
當φ→0時,分布近似于Poisson分布,此時Poisson精確概率法的檢驗功效最高;Fisher精確概率法、負二項精確概率法、t檢驗、平方根反正弦變換的t檢驗(arcsin變換t)、平方根變換的t檢驗、秩檢驗、log變換的t檢驗都隨著n和δ的增高呈現(xiàn)一定的波動,但總得來講,在δ=3到18時各個統(tǒng)計方法的檢驗功效排序趨于穩(wěn)定(圖1),我們從表3中可以得到Poisson,F(xiàn)isher,NB方法的檢驗功效較大,其值也較相近;t檢驗、平方根反正弦變換t檢驗、平方根變換t檢驗、以及l(fā)og變換t檢驗的檢驗功效居中,其大小近似;而秩檢驗的檢驗功效較小。
當離散參數(shù)φ→1時,檢驗功效結果見表4。
當φ→1時,分布近似于負二項分布,此時Poisson精確概率法的檢驗功效最高;Fisher精確概率法的檢驗功效僅次于Poisson精確概率法;負二項精確概率法、t檢驗、平方根反正弦變換的t檢驗(arcsin變換t)、平方根變換的t檢驗、秩檢驗、log變換的t檢驗都隨著n和δ的增高有一定的波動性。我們以總體差別δ為橫坐標,功效為縱坐標作圖,發(fā)現(xiàn)除了Poisson和Fisher在δ為5時趨于穩(wěn)定外,其他方法都隨δ的變化較大,直到δ為12時才相對穩(wěn)定(圖1),故我們以δ=12為例,從表4,圖2中可以發(fā)現(xiàn):負二項精確概率法(NB)的檢驗功效較大;平方根反正弦變換t檢驗、平方根變換t檢驗、以及秩檢驗的檢驗功效居中;而log變換t檢驗、t檢驗的檢驗功效較小。
表3 檢驗功效及排序(φ→0,δ=3)
表4 檢驗功效及排序(φ→1,δ=12)
圖1 δ與功效的關系圖(n=10)
圖2 不同δ值時n與功效示意
雖然 David J.Kirkland 等〔2〕已經(jīng)明確指出過未經(jīng)變換的t檢驗方法很多時候不適用于微核數(shù)據(jù)的研究分析,但用PubMed-NLM檢索近3年的關于微核實驗可獲全文的16篇文章中仍能查到5篇文獻使用此方法。我們的模擬研究結果表明功效隨著動物數(shù)n和容許誤差δ的變化而變化,其中δ的影響略大,但當δ達到一定值時,功效值便趨于穩(wěn)定。為了說明功效的增高是由于δ增大引起還是由于不同統(tǒng)計方法導致,我們選取了功效穩(wěn)定時的結果作為討論依據(jù)。在φ→0時,模擬數(shù)據(jù)服從Poisson分布,此時Poisson分布的精確概率法功效最高,但由于其Ⅰ類錯誤非常大,遠遠超出檢驗水準故不推薦使用;負二項精確概率法與Fisher精確概率法的功效很相近,這是由于二者的算法類似,其值也高且Ⅰ類錯誤均小于0.05;t檢驗在樣本量n小于7時功效很小,在大于7時,功效尚能達到0.8左右,但其Ⅰ類錯誤也增大了,故不可接受,而log變換的t檢驗功效很低故不推薦使用。arcsin變換的t檢驗、平方根變換的t檢驗及秩檢驗在樣本量n為10時,功效較大,他們的Ⅰ類錯誤也小于檢驗水準0.05;故在φ→0時,若動物數(shù)很少,建議使用負二項精確概率法和Fisher精確概率法,若動物數(shù)較大(大于10時)使用Arcsin變換的t檢驗、平方根變換的t檢驗及秩檢驗尚可。在φ→1時,Poisson分布的精確概率法與Fisher精確概率法的效能均很高,但他們的Ⅰ型錯誤很大,遠大于檢驗水準α=0.05,故不推薦使用;t檢驗、log變換的t檢驗功效較低,故不推薦使用;arcsin變換的t檢驗、平方根變換的t檢驗、秩檢驗在動物數(shù)大于9時功效尚可接受,但其Ⅰ型錯誤此時也大于檢驗水準0.05,故不建議使用;負二項分布精確概率法的效能較高,僅次于Fisher精確概率法,其Ⅰ型錯誤也小于檢驗水準0.05。故在φ→1時,建議使用負二項分布精確概率法。綜上表明,在φ→0時,動物數(shù)很少時,建議使用負二項精確概率法和Fisher精確概率法;若動物數(shù)較大(大于10時)使用arcsin變換的t檢驗、平方根變換的t檢驗及秩檢驗尚可。在φ→1時,只推薦使用負二項分布精確概率法。這里需要強調的是,很多時候由于動物間的個體變異存在,負二項分布假定比Poisson分布更合理,應用者應謹慎檢查數(shù)據(jù)分布假設再進行恰當?shù)慕y(tǒng)計分析〔3〕。
1.International Conference on Harmonisation Guidance on Genotoxicity Testing and Data Interprerpretation for Pharmaceuticals Intended for Human Use.Draft January 28th Version 5.3,2008.
2.Kirkland DJ.Statistical evaluation of mutagenicity test data:recommendations of the U.K.Environmental Mutagen Society.Environmental Health Perspectives Supplements vol.102 Suppl,1994,1:43-47.
3.Ludwig Hothorn.Biostatistical analysis of the micronucleus mutagenicity assay based on the assumption of a mixing distribution.Environmental Health Perspectives Supplements vol.102 Suppl,1994,1:121-125.
4.Hothorn LA,Gerhard D.Statistical evaluation of the in vivo micronucleus assay.Arch Toxicol,2009 Jun,83(6):625-634.
5.Robinson MD,McCarthy DJ,Smyth GK.edgeR:a Bioconductor package for differential expression analysis of digital gene expression data,Bioinformatics,2010 Jan 1,26(1):139-140.
6.Robinson MD,Smyth GK.Small Sample Estimation of Negative Binomial Dispersion,with applications to SAGE data.Biostatistics,2008,9(2):321-332.