BBN模型在江蘇居民氰戊菊酯長期膳食暴露評估中的應用*

東南大學公共衛生學院流行病與衛生統計教研室(210009) 田 野 閔 捷 劉 沛 王燦楠

在進行污染物的長期膳食暴露評估時，合理地估計每人每天長期攝入某種物質平均量(usual intake，簡稱日常攝入量)往往是風險評估過程中的重點。日常攝入量的計算涉及該物質在各種食物中的濃度和人群消費各種食物的量，前者可用各種食物污染物檢測的平均濃度估計，不考慮其變異;后者需要采集人群長期消費各種食物的信息，實際工作中是難以實現的，如何利用短期2～7天的24小時膳食回顧調查數據來估計人群長期暴露量就顯得尤為重要。許多專家及機構提出了估計及改進的方法〔1-3〕，目前大家普遍接受的是貝塔二項-正態分布模型(Beta binomial-normal，BBN)，該方法是2005年由De Boer和Van der Voe等人提出的〔4〕。其基本思想是考慮調查人群中，個體間變異可近似反映人群的長期飲食情況，個體內變異則被認為僅反映調查個體在調查的n天內短期膳食結構的波動，采用合適的統計模型去除個體內變異，保留個體間變異的暴露量，近似地作為人群的長期暴露量。本文利用BBN模型對江蘇省各年齡、性別組人群氰戊菊酯的長期膳食安全風險水平進行評估，以期為農藥膳食長期暴露風險評估提供科學依據。

BBN模型原理

BBN模型包括攝入頻率和陽性暴露量兩部分。其中，貝塔二項分布被用來構建攝入頻率模型，對數正態分布被用來構建陽性暴露量模型。

1.攝入頻率模型

假設n和npos分別代表每個調查對象的總調查天數和攝入食物中含有該污染物的天數(即陽性攝入天數)。由于每人每天攝入的食物中可能包含污染物，也可能不包含污染物，因此陽性攝入天數npos服從陽性攝入概率為p，重復數為n的二項分布;其中，陽性攝入概率p的取值在0到1之間，服從貝塔分布;那么，陽性攝入天數就服從貝塔二項分布:

貝塔分布的均數和方差分別為

令π=α/(α+β)，φ=1/(α+β+1)，則該分布的均數和方差可表示為:

通過最大似然估計法擬合貝塔二項分布模型。可以得到π和φ的估計值π和φ，從而反推得到α和β的估查對象在某一調查日食用某種事物的陽性攝入頻率分布即可用貝塔分布

2.陽性暴露量模型

長期膳食暴露評估主要反應人群中某污染物日常攝入情況，假設陽性暴露量即調查對象每人每天膳食污染物暴露量不為零，經過對數變換后近似服從正態分布，變換后的陽性暴露量，通過最大似然估計法擬合線性隨機效應模型(linear random effects model):

其中，transf(yposij)表示經對數變換后的陽性暴露量，ci和uij表示人間效應和人內效應，其二者分布分別滿足正態分布 N(0，σ2B)和N(0，σ2w)。通過最大似然估計法可以得到三個參數的估計值m，分別代表個體平均值、人間效應估計值及人內效應估計值。

3.計算長期暴露量分布

長期暴露量是通過對陽性暴露量分布和攝入頻率分布進行Monte Carlo抽樣，然后進行整合獲得。其具體過程為:

首先通過模擬獲得對數轉換后的陽性暴露量m+zi，其中zi從正態分布N(0)中隨機抽取10萬次，然后反對數變換獲得陽性暴露量:

實例分析

1.資料情況主要涉及的數據有膳食調查數據、污染物常規監測數據、研究人群人口學基本情況數據。膳食調查及人口學數據:來源于2002年全國居民營養與健康狀況調查中的江蘇省數據;污染物監測數據:來源于2000～2006年全國14個省/地區食品污染物監測網數據及2005～2006年海關出口農產品監測數據中的江蘇省數據。本研究資料中氰戊菊酯檢出率為3.10%，按照WHO關于污染物未檢出值處理的原則，本文中未檢出值均采用0替換〔5〕。目前，FAO/WHO食品添加劑聯合專家委員會(JECFA)制定的食品中每日允許攝入氰戊菊酯量(ADI)為2μg/kg bw/day。將模型結果與之進行比較，評價人群的氰戊菊酯攝入水平。

2.結果與分析

應用SAS9.1統計軟件進行數據整理和膳食暴露量的計算，結果如下。

(1)BBN模型參數

利用SAS的Proc nlmixed過程的最大似然估計法迭代獲得貝塔分布參數和的估計值。在對陽性暴露量進行方差分解時，首先將陽性暴露量進行對數變換，使其近似正態分布，然后利用SAS的Proc mixed過程將殘差分解為個體間效應和個體內效應，通過最大似然估計法獲得個體平均值m、人間效應估計值及人內效應估計值。計算的氰戊菊酯膳食暴露BBN模型參數見下表。

表1 氰戊菊酯膳食暴露BBN評估模型參數

(2)氰戊菊酯長期膳食暴露評估結果

如表2可見，氰戊菊酯的長期暴露量在平均值及各百分位數均低于其每日允許攝入量ADI(2 μg/kg bw/day)，說明氰戊菊酯殘留在長期膳食暴露評估中是相對安全的。

表2 不同人群氰戊菊酯長期膳食暴露量(單位μg/kg bw/day)

(3)氰戊菊酯膳食暴露量的不確定性分析

為了解人群氰戊菊酯暴露變異的各百分位數的可信區間，采用Bootstrap方法有放回的從原樣本中進行抽樣，并對抽取的每一個Bootstrap樣本進行Monte Carlo模擬以計算暴露量百分位數，經過n次Bootstrap抽樣計算百分位數分布及可信區間，代表百分位數的不確定性。本文取200次Bootstrap抽樣進行不確定性分析。

由表3可見，各暴露水平的暴露量均位于95%可信區間內，說明暴露量百分位數能夠達到穩定。由于本次研究采用的是全比例抽樣，其不確定性主要存在于攝入頻率p～Beta(α，β)和人間效應zi～N(0)的抽樣中，同時Monte Carlo模擬較好的量化了變異性，故在可信區間內變異較小。

討 論

1.長期膳食暴露評估模型

評價污染物長期膳食暴露有諸多不同的數學模型，其共同思想是考慮在大規模人群短期橫斷面膳食調查情況下，數據資料中變異主要包括個體間變異和個體內變異。不同個體間的變異可以近似反映人群的長期膳食狀況，而個體內變異則被看做僅反映調查個體在橫斷面調查中短期的膳食結構波動，因此，通過采用合適的數學模型來去除個體內變異，而保留個體間變異，從而達到用污染物短期的膳食暴露量近似代表長期的膳食暴露量的目的。本研究在計算長期膳食暴露評估時采用歐盟MCRA(Monte Carlo Risk Assessment)推薦的BBN模型。將長期的污染物膳食陽性攝入頻率用短期陽性攝入頻率的貝塔二項分布模型來表示;同時，用正態分布構建長期陽性暴露量模型。從分布中抽樣，進行Monte Carlo模擬獲得污染物的長期暴露量。利用BBN模型進行長期膳食暴露評估的優點是將個體作為研究對象，通過統計方法進行模擬，使結果更符合實際并獲得人群的暴露量分布，其對膳食調查數據及污染物殘留數據的信息利用程度遠大于點評估模型。國外一些研究進行了模型的比較也證明了這一結論〔6〕。當然，此種暴露評估方法也有其局限性，如模擬次數足夠多才能使最終結果穩定及對數據質量要求較高等。另外該模型適用于污染物檢出率不高的暴露評估，如果檢出率高，例如鉛、鎘等的污染物，則BBN中可不考慮Beta分布。

表3 不同人群氰戊菊酯長期膳食暴露量不確定性分析(單位μg/kg bw/day)

2.關于陽性暴露量的正態性轉換

BBN模型中，需要利用正態分布構建陽性攝入量模型。人群的陽性攝入量的原始數據往往不滿足正態性，常常需要用變量變換的方法將偏態資料正態化。在BBN模型中常用的變量變換有兩種:對數變換和冪變換〔7〕。經過變換后的變量，一般需對其分布進行擬合優度檢驗，常用的檢驗有Anderson-Darling檢驗、Kolmogorov-Smirnov檢驗、Chisquare檢驗等。本研究中，陽性暴露量的原始數據表現為正偏態的基本特征，本研究采用對數變換對陽性攝入量的原始數據進行轉換。

3.Bootstrap抽樣次數的確定

關于Bootstrap抽樣次數的確定，一般根據研究目的的不同取不同的值，如果用來計算均值、方差、標準誤等統計量，200次抽樣即可，但如果計算極端百分位數如99.9th的可信區間，Bootstrap的模擬次數應當增加到至少500次〔8〕。本研究Monte Carlo模擬的10萬個人天數作為一個Bootstrap樣本含量，根據需要估計的百分位數情況，將Bootstrap樣本數定為200次、500次、1000次分別進行模擬，比較抽樣結果發現可信區間寬度沒有明顯差異，最終確定本研究中Bootstrap模擬次數為200次比較合適。需要說明的是，對不確定度更多的應從專業上加以評價，如果風險管理者認為數據的不確定度過大，可以拒絕概率估計的結論，或者采取一定措施(如提高數據質量)來減少不確定度〔9〕。

1.Nusser SM，Carriquiry AL，Dodd KW，et al.A semi-parametric transformation approach to estimating usual intake daily intake distributions.J Am Stat Assoc，1996，91:1440-1449.

2.Slob W.Modeling Long-term exposure of the whole population to chemicals in food.RISK ANAL，1993，13(5):525-530.

3.Slob W.Probabilistic dietary exposure assessment taking into account variability in both amount and frequency of consumption.Food and Chemical Toxicology，2006，44(7):933-951.

4.De Boer WJ，van der Voet H.MCRA，Release 5，a web-based program for Monte Carlo Risk Assessment.Biometris and RIKILT，Wageningen University and Research Centre，2006.

5.王緒卿，吳永寧.食品污染監測低水平數據處理問題.中華預防醫學雜志，2002:837-852.

6.Boon PE，Bonthuis M.Comparison of different exposure assessment methods to estimate the long-term dietary exposure to dioxins and ochratoxin A.Food and Chemical Toxicology，2011，49(9):1979-1988.

7.De Boer WJ，van der Voet H.Comparison of two models for the estimation of usual intake addressing zero consumption and non-normality.Food Additives and Contaminants，2009，26(11):1433-1449.

8.Cullen AC，Frey HC.Probabilistic Techniques in Exposure Assessment:A Handbook for Dealing with Variability and Uncertainty in Models and Inputs.Society for Risk Analysis，1999.

9.Thompson KM.Variability and uncertainty meet risk management and risk communication.Risk Anal，2002，22:647-654.