DOi:10.13590/j.cjfh.2018.01.016
基于空间统计方法的蔬菜中农药残留风险分析
杨蕙1,2,李宁2,计融2,何平1,肖革新2

(1.贵州省疾病预防控制中心卫生监测检验所,贵州 贵阳 550004; 2.国家食品安全风险评估中心,北京 100022)

收稿日期:2017-12-01

作者简介: 杨蕙 女 主管医师 研究方向为食品安全及统计流行病学 E-mail:179393792@qq.com
通信作者:肖革新 男 副研究员 研究方向为空间统计及大数据信息化 E-mail:biocomputer@126.com

基金项目:

摘要:目的 分析蔬菜中农药残留空间格局分布变化及规律,找到聚集热点区域,为蔬菜中农药残留监测、监管提供科学依据。方法 用空间统计方法对某省2014—2016年蔬菜中农药残留监测信息结合地理信息进行空间格局分布及空间聚集性分析。结果 2014—2016年区域内监测的蔬菜中农药残留检出率和超标率两者分布格局较一致,但三年对比变化较大,检出较高区域集中分布在西南部和东北部;超标率涉及区县较少,散在分布在西南部和东北部。2014年监测区域的蔬菜中农药检出率全局存在正相关性(P<0.05,Moran's I=0.410),局部相关主要分布在西南部、北部和南部地区;2016年监测区域的蔬菜中农药超标率全局存在正相关性(P<0.05,Moran's I=0.111),局部相关主要分布在西南部地区;2016年监测区域的蔬菜中农药检出率全局存在负相关性(P<0.05,Moran's I=-0.087),局部相关主要分布在中部地区。时空扫描统计量分析结果显示监测区域内超标率和检出率均有聚集性,超标率有一个聚类区域(LLR=27.11,P<0.05,RR=20.04),分布在东北部区县;检出率有两个聚类区域,一类聚类区域(LLR=43.24,P<0.05,RR=6.9)分布在东北部和南部区县,二类聚类区域(LLR=19.13,P<0.05,RR=4.13)分布在西南部区县。结论 2014—2016年监测区域的蔬菜中农药残留检出率和超标率呈现较为一致的下降趋势。空间自相关性和时空扫描统计量结果显示存在聚集性,聚集区域主要为东北部和西南部邻近区县,且聚集区域的相对风险度RR值均>4,提示这些聚集区县较其他区县农药残留风险高,是监管部门应该关注的热点区域。
关键词: 蔬菜; 空间统计; 农药残留; 风险分析
文章编号:1004-8456(2018)01-0073-06     中图分类号:R155     文献标志码:A    
Risk analysis of pesticide residue in vegetables based on spatial statistical method
    YANG Hui1,2, LI Ning2, JI Rong2, HE Ping1, XIAO Ge-xin2
(1.Guizhou Provincial Center for Disease Control and Prevention, Guizhou Guiyang 550004,China; 2.China National Center for Food Safety Risk Assessment,Beijing 100022,China)
Abstract:Objective To preliminary assess the risk of Norovirus infections by establishing oyster-Norovirus exposure model of Beijing residents from retail to table.Methods Health risk and the influencing factors were estimated using exposure assessment model with 356 quantitative detection data point of Norovirus in oysters collected in Beijing under simulated consumption scenarios and published dose-response model derived from outbreak data.Results The contamination level of GI and GII Norovirus positive samples were 2.62×104 copies/g (95%CI:3.73×103-1.54×105) and 5.02×104 copies/g (95%CI:8.13×103-2.52×105). For positive histo-blood group secretor groups, the risk of consuming one raw oyster potentially contaminated with GI and GII Norovirus were 0.93 (95%CI:0.73-0.98) and 0.95 (95%CI:0.80-0.99) respectively; for negative histo-blood group secretor groups, the risk 0.37 (95%CI:0.04-0.64) and 0.57 (95%CI:0.07-0.99). The estimated risk was mainly related to the assigned contamination level for those non-detected samples with the regression coefficient of 0.49. Uncertainty analysis showed that the detection level of the method used was higher than the ID50/LD50, which was the main source of uncertainty.Conclusion The risk of Norovirus infections after consuming raw oyster was high for Beijing residents; increasing the sensitivity and decreasing the limit of detection could decrease the uncertainty of the risk estimation.
Key words: Wheat flour; geographical; origin; principal component analysis; partial least squares to latent structures-discrimination analysis; orthogonal partial least squares to latent structures-discrimination analysis; multivariate data analysis 
 我国是一个农业大国,据文献[1]报道2004年农药使用量就已居世界首位。中国又是一个传统的蔬菜生产大国,早在1988年我国蔬菜播种面积和总产量已成世界蔬菜生产第一大国。蔬菜作为人们日常饮食中不可或缺的一种食物,我国居民每天对蔬菜的需求量很大[2],所以在蔬菜种植过程中,人们为了防治病虫,大量地使用农药,使得蔬菜农药残留问题日益严重[3]
        随着世界经济的发展,消费者对自身健康的考虑以及环保意识的增强,对蔬菜的安全优质要求越来越大,无农药、无化肥或减农药、减化肥的有机蔬菜己成为一种新的消费时尚,对“绿色蔬菜”的需求将越来越迫切。无污染和保健性蔬菜在欧洲、美国、日本等经济发达国家和地区已经成为一种趋势[4]。目前国内外对于蔬菜中农药残留的研究仅限于监测结果的描述性分析、超标风险、建立监测方法等方面,没有结合数据时空上的属性进行研究,无法深入揭示蔬菜中农药残留在空间上的格局分布变化和时空上的变化规律及相对风险,本研究利用空间统计学方法对区域内蔬菜农药残留风险进行分析,研究其空间格局分布、空间关系、时空变化规律,以探寻聚集和关注的热点区域,为监测、监管部门提供依据。    
1材料与方法     
1.1数据来源
        某省2014—2016年蔬菜中农药残留监测的数据。     
1.2统计学分析
        用R3.2.4和Excel对监测数据进行处理, 采用空间统计方法中的ArcMap对整理出的检出率和超标率进行可视化及空间格局分析,用OpenGeoDa和SaTScan对其进行空间自相关性和时空扫描统计量分析,找到热点聚集区域及相对风险度。    
2结果     
2.1数据整理
        将某省2014—2016年蔬菜中农药监测的19 835条数据进行整理。总体检出率和超标率情况见表1。监测的农药种类主要为杀虫杀螨剂和植物调节剂,共涉及41种,其中检出率较高的是毒死蜱、氯氰菊酯、乙酰甲胺磷、氯氟氰菊酯、甲氰菊酯5类;超标率较高的是毒死蜱、克百威、氯氟氰菊酯3类。监测的蔬菜种类包括非葫芦科茄果类、甘蓝和芸薹类、瓜菜类(葫芦科)、茎类、块根和块茎类、鳞茎类、 鲜豆类、叶菜类、豆芽菜等9类蔬菜。
表12014—2016年蔬菜中农药残留监测情况
Table 1Monitoring of pesticide residues in vegetables 
    from 2014 to 2016    
2.2监测区域检出率和超标率空间可视化及格局分布分析
        监测区域(某省)以地理空间分东、西、南、北、 中五区来进行描述,以下所有空间可视化图上为北,下为南,左为西,右为东,与地图描述一致。2014—2106年监测区域采样覆盖率(80%以上)一致且分布在全区域。     
2.2.1检出率空间可视化及格局分布分析
            将整理好的检出率和相应区域地理信息整合起来连接到ArcMap,结果见图1。由此可见,检出率整体呈下降趋势,检出分布格局三年变化较大,2014年检出率较高区域主要分布在中、南部,2015年检出率较高区域主要分布在东北部和西南部,2016年检出率较高区域主要分布在中部。
 图1某省2014—2016年蔬菜中农药残留检出率空间分布图
     Figure 1Spatial distribution of detection rate of pesticide residues in vegetables in a province from 2014 to 2016        
2.2.2超标率空间可视化及格局分布分析
                将整理好的超标率和相应区域地理信息整合起来连接到ArcMap,结果见图2。由此可见,超标率分布格局三年变化较大,2014年超标率较高区域主要分布在中部和散在分布在南、北部,2015年超标率较高区域主要分布在东北部和西南部,涉及区县比2014年多,而2016年超标率涉及区县变少,其中较高区域散在分布在西南部和东南部的5个区县。     
2.3空间自相关分析
        空间自相关分为全局空间自相关和局部空间自相关。其中,全局空间自相关主要用于描述某区域内某种现象的整体空间分布情况,以判断该现象在空间上是否存在聚集性空间相关性,局部空间自相关分析侧重于研究局部范围内空间对象属性值在的空间相关性[5]。空间相关性用Morans I系数和P值来反映,当显著性检验结果P<0.05时,Morans I系数取值>0时,表明所研究区域存在空间正相关,且取值越接近1,表明空间正自相关性越强,研究对象的值呈聚集分布;当Morans I系数取值<0时,表明所研究区域存在空间负相关,取值越接近-1,表明空间负相关性越强;当Morans I系数接近0时,则表明不存在空间自相关性,空间变量呈随机分布格局[6]。若相邻区域间同一属性值表现出相同或相似的相关程度,即属性值在空间区域上呈现高的地方邻近区域也高,呈现低的地方临近区域也低,则称为空间正相关;若相邻区域间同一属性值表现出不同的相关程度,即属性值在空间区域上呈现高的地方邻近区域低,呈现低的地方临近区域高,则称为空间负相关;若相邻区域间同一属性值不表现任何依赖关系,即呈随机分布,则称为空间不相关。对其蔬菜中农药监测数据进行空间相关性分析,具有空间相关性且差异有统计学意义(P<0.05,Morans I∈[-1,1]且≠0)。 
图2某省2014—2016年蔬菜中农药残留超标率空间分布图
     Figure 2Spatial distribution of pesticide residues in vegetables from 2014 to 2016 in a province    
         P<0.05,Morans I=0.410,提示2014年监测区域的蔬菜中农药检出率全局存在正相关性;局部有7个区县存在高-高正相关性,5个区县存在低-低正相关性,4个区县存在低-高负相关性和1个区县存在高-低负相关性,局部相关主要分布在西南部、北部和南部地区,见图3~5。
图3监测区域蔬菜中农药检出率和超标率Morans散点图
     Figure 3Morans scatter plot of the detection rate of pesticides in vegetables in the monitoring area    
图4监测区域蔬菜中农药检出率和超标率Morans I随机化零假设检验结果
     Figure 4Morans I randomised zero hypothesis test results for detection rate of pesticides in vegetables from 2014 to 2016    
         P<0.05,Morans I=0.111,提示2016年监测区域的蔬菜中农药超标率全局存在正相关性;局部有1个区县存在高-高正相关性,1个区县存在低-低正相关性,2个区县存在低-高负相关性和4个区县存在高-低负相关性,局部相关主要分布在西南部地区,见图3、4和6。
        P<0.05,Morans I=-0.087,提示2016年监测区域的蔬菜中农药检出率全局存在负相关性;局部有1个区县存在低-低正相关性,3个区县存在低-高负相关性和2个区县存在高-低负相关性,局部相关主要分布在中部地区,见图3、4和6。    
2.4时空扫描统计量分析
        时空聚集性分析同时考虑了时间和空间两个因素,主要采用移动窗口法在地理空间上创建扫描窗口(圆柱体),圆柱体的底面对应研究的地理区域,圆柱体的高对应扫描时间间隔,圆柱体的半径对应扫描的风险。SaTScan扫描统计量是窗口所含事件的最大个数, 通过使用局部对数似然函数来估计统计量对数似然比(LLR)的极大值, 后用蒙特卡罗模拟扫描统计量, 进行显著性检验, 最终确定窗口内是否存在异常聚类[7],它可以揭示地理空间上聚集区域随时间变化的规律,同时还能得出聚集区域的相对风险度(RR),更加精确的对空间聚集区域位置进行了定位。当LLR取最大值且P<0.05时,可以认为该区域存在聚集特征。注:A为超标率;B为检出率以监测数据区县为对应研究区域,时间间隔为年,扫描2014—2016年的蔬菜中农药检出率和超标率,得到相应RR值。将扫描结果显示的聚类区域标注到地图上,扫描结果见图7。结果显示区域内超标率和检出率均有聚集性,超标率有一个聚类区域(LLR=27.11,P<0.05, RR=20.04),分布在东北部区县,聚集时间为2015年;检出率有两个聚类区域,一类聚类区域(LLR=43.24,P<0.05, RR=6.9)分布在东北部和南部区县,二类聚类区域(LLR=19.13,P<0.05, RR=4.13)分布在西南部区县,RR值均>4,且差异有统计学意义,其中一类聚集区风险概率大于二类聚集区[8],聚集时间为2016年。
图52014年监测区域蔬菜中农药检出率LISA分析结果
     Figure 5Results of Lisa analysis on the detection rate of 
    pesticides in vegetables in the monitoring area in 2014    
  图62016年监测区域蔬菜中农药超标率和检出率LISA分析结果
     Figure 6Results of LISA analysis on the over-standard rate and detection rate of pesticides in vegetables in the monitoring area in 2016  
 图7某省蔬菜中农药残留超标率和检出率聚类区域图
     Figure 7Cluster regional map of pesticide residue over-standard rate and detection rate in vegetables of a province       
3讨论
        从空间可视化图片可以很直观的看出,三年的检出率和超标率呈现较为一致的下降趋势,检出率高的区域超标率也高,但各年检出率和超标率高低区域分布不一致,2014年检出率较高的区域主要分布在中部和南部,2015年分布变化为东北部和西南部,2016年分布变化为中部较高;2014年超标率较高的区域主要分布在中部和散在分布在南、北部,2015年分布变化为集中在东北部和西南部,2016年散在分布在西南和东南部几个区县。
        空间自相关性和时空扫描统计量结果显示,监测区域内2014—2016年蔬菜中农药残留检出和超标情况存在聚集性,聚集区域主要为东北部和西南部邻近区县,且聚集区域的相对风险度RR值均>4,提示这些聚集区县较其他区县农药残留风险高,是监管部门应该关注的热点区域。
        空间统计学理论在20世纪60年代就被Matheron提出[9],于1990年被引入国内[10],目前空间统计学被应用到多个领域,而还未见引入食品领域,据研究表明[11],90%以上的卫生数据有空间或地理成分,而监测的食品安全数据同样具有空间属性,本研究将空间统计学方法运用到食品安全风险分析中,可以打破以往经典统计分析,从空间、时间、专题属性多维的角度直观、深入的分析及预判食品安全风险,为相关部门提出针对性的预警对策建议具有一定前瞻性意义。
        不足之处在于,监测为一年一次,且数据仅利用了3年数据进行时空扫描,扫描结果比较粗略,没有真正体现出聚集和相对风险在时空上的变化,建议针对监测数据多做此方面的研究,将年限延长为5~10年来进行时空分析,以便更好的发现食品监测数据在时空上的规律。另外,虽然各年监测的蔬 菜种类大类一致,但因区域和季节等不确定原因,难以确保监测的蔬菜细类一致性,每种蔬菜的农药含量和检出限不尽相同,会对检出率和超标率结果有影响,因为蔬菜细类较多且有些蔬菜细分到小类的样品量较少,因此本研究没有进行细类比较分析,不排除某省区域内的蔬菜农药残留空间格局分布因蔬菜种类的不同呈现差异性,建议增加蔬菜细类样品量,进行进一步深入分析。
参考文献
[1]张俊,王定勇.蔬菜的农药污染现状及农药残留危害[J].河南预防医学杂志,2004,15(3):182-186.
[2]安然.蔬菜的农药污染现状及潜在风险[J]. 食品安全导刊,2017,1(3):1.
[3]海力帕木·吾麦尔. 蔬菜农药残留问题[J]. 河南农业,2016,1(2):62.
[4]杨锦秀. 中国蔬菜产业发展的经济学分析[D]. 成都:西南财经大学,2005:211.
[5]吴珣,杨婕,张红.不同空间权重定义下中国人口分布空间自相关特征分析[J]. 地理信息世界, 2017,24(2): 32-38.
[6]熊昌盛,栾乔林,韦仕川.基于空间自相关的耕地质量分布格局[J]. 地域研究与开发, 2016,35(5): 128-132,148.
[7]吴建玮,刘凯,田鹏飞,等. 时空扫描对院前出诊负荷的探索[J].中国保健营养,2017,27(13): 15-16.
[8]李婷,何金戈,杨长虹,等. 基于空间聚集性与时空扫描的肺结核流行特征分析[J]. 中国防痨杂志, 2016,38(12): 1032-1040.
[9]邹亚娟.空间统计学研究进展综述[J].企业技术开发,2014,33(12): 79-80.
[10]戴平生,陈建宝.空间统计学研究应用综述[C].国际应用统计学术研讨会,烟台,2008:1-9.
[11]张喜旺,秦奋,刘剑锋. 地统计学及其在公共卫生领域的应用[J].统计与决策,2006,10(5):129-131.