基于贝叶斯推断的DNDC模型参数校正与不确定性评价研究.pdf
《基于贝叶斯推断的DNDC模型参数校正与不确定性评价研究.pdf》由会员分享,可在线阅读,更多相关《基于贝叶斯推断的DNDC模型参数校正与不确定性评价研究.pdf(8页珍藏版)》请在淘文阁 - 分享文档赚钱的网站上搜索。
1、第 51 卷 第 2 期土壤学报Vol.51,No.22014 年 3 月ACTA PEDOLOGICA SINICAMar.,2014http:/*中国科学院知识创新工程重要方向项目(KZCX2-EW-QN404)和中国科学院战略性先导科技专项(XDA05050509)资助 通讯作者,E-mail:yczhao 作者简介:秦发侣(1986),男,云南玉溪人,博士研究生,主要从事土壤资源与遥感信息研究。E-mail:dielianhuared 收稿日期:2013 05 20;收到修改稿日期:2013 09 19DOI:10.11766/trxb201305200247基于贝叶斯推断的 DNDC
2、 模型参数校正与不确定性评价研究*秦发侣1,2赵永存1,2史学正1,2于东升1,2徐胜祥1,2(1土壤与农业可持续发展国家重点实验室(中国科学院南京土壤研究所),南京210008)(2中国科学院大学,北京100049)摘要实现土壤有机碳(SOC)动态模拟结果的不确定性定量评价是农田管理决策成败的关键之一。采用贝叶斯推断和马尔科夫链蒙特卡洛方法(MCMC),对 DNDC 模型模拟江苏省宜兴市一个具有 22 年稻-麦轮作历史的田块的 SOC 动态变化,进行模型参数校正和模拟结果的不确定性定量评价。结果表明,DNDC 模型适宜长期监测田块的 SOC 动态变化的模拟,但模型的模拟结果存在一定的不确定性
3、;在 DNDC 模型输入参数数据质量不明的情况下,利用贝叶斯推断和 MCMC 方法能够有效地实现模型输入参数的自动校正和 SOC模拟结果不确定性的定量评价,从而为实现区域或国家尺度农田 SOC 动态模拟的不确定性定量评价提供理论和方法依据。关键词贝叶斯推断;马尔科夫链蒙特卡洛方法;参数校正;不确定性评价;DNDC 模型中图分类号S11+4文献标识码A基于过程的 DNDC(Denitrification-Decomposi-tion)模型已经被广泛应用于中国农田土壤有机碳(SOC)的变化模拟及预测1-5,然而,在模拟和预测过程中不可避免地存在不确定性6-8。一般而言,模型的不确定性是由于模型对复
4、杂世界的抽象、概括而形成的结构误差和模型使用时输入参数的不确定性导致的7-8。对于模型的使用者而言,模型结构是难以改变的,而由于模型输入参数的可获取性、数据质量及模型输入参数本身的变异性等因素所导致的输入参数的不确定性也会随着 DNDC 模拟 SOC 动态变化的过程加以传递,从而对 SOC 动态模拟结果产生深刻影响。因此,在模型结构确定的情况下定量 DNDC 模型输入参数的不确定性及其对 SOC 动态变化模拟结果的潜在影响,对于科学认知 SOC 动态变化规律、建立合理的推荐耕作管理措施,进而实现农田土壤肥力提升和农田土壤固碳具有重要的意义。SOC 动态模拟结果中所包含的不确定性对农田及土壤资源
5、管理决策具有重要影响,提供这些不确定性信息也是使用和交流模型结果的关键环节6-9。然而,目前关于 DNDC 模型模拟 SOC 动态的不确定性定量评价研究还相对较少,已有的研究多集中在模拟温室气体(如 CH4,CO2,N2O 等)排放的不确定性评价方面6,9-10,相关的不确定性评价方法主要包括最敏感因子法(MSF)和蒙特卡洛法(MC),其中 MC 又可分为 DNDC 内置及非 DNDC 内置(用户自定义)两类。如 Li 等10 使用 MSF 和DNDC 内置 MC 方法对比了中国、泰国和美国的三个农田长期定位点所属土壤图斑由于初始 SOC、土壤容重、pH 和黏粒含量这四个输入参数的不确定性所造
6、成的 DNDC 模型模拟温室气体排放的置信区间差异;Hastings 等6 采用自定义的 MC 研究了位于瑞士的一个农田长期定位试验点土壤属性、作物管理措施和气候等的不确定性对 DNDC 模型模拟温室气体排放的影响。DNDC 内置的 MSF 方法假定对模型输出结果影响最大的 4 个参数分别为初始SOC 含量、土壤容重、pH 及黏粒含量,同时假定DNDC 输出结果与这 4 个参数之间的关系为线性正相关或负相关11。MSF 方法将这 4 个参数的最小值和最大值组合成两个参数集,并将它们分别输入248土壤学报51 卷http:/DNDC 模型进行模拟,以输出结果的最小值和最大值来定量模拟结果的不确定
7、性区间。如 Li 等10 认为初始 SOC、pH 的最大值和黏粒含量的最小值组合所得的参数集可以模拟得到 CH4排放量的最大值;相反,初始 SOC、pH 的最小值和黏粒含量的最大值组合所得的参数集可以模拟得到 CH4排放量的最小值。MC 方法则较 MSF 方法更灵活,可以自由选择模型输入参数进行评价,其实施步骤包括:(1)选取关注的模型输入参数;(2)设定输入参数的变化区间及概率分布;(3)依据模型输入参数的概率分布抽样,获得 N 个模型输入参数集;(4)将 N 个模型输入参数集输入模型获得 N 个模型输出结果;(5)对 N 个模型输出结果进行统计分析,用置信区间来表征模型输出结果的不确定性9
8、。DNDC 内置的MC 方法默认输入参数均为均匀分布,即每个参数集对模型输出结果贡献的权重是相等的。用户自定义的 MC 方法则可根据模型输入参数实测数据的统计特征来建立相应的概率分布函数(比如正态分布等)。理论上讲,当模型输入参数的实测数据较为充足时,由于根据实测数据统计推断出的模型输入参数的不确定性区间较为接近输入参数的真实不确定性区间及概率分布,因此能得到较 DNDC 内置的 MC 方法更为可靠的模拟结果置信区间。然而,当模型输入参数的实测数据极为有限时,根据实测数据来统计模型输入参数的可能变化范围变得十分困难。针对此种情况,常用的方法是根据专家经验,或者采用与研究区域类似、实测数据较为齐
9、全的其他区域的输入参数变化区间来代替此区域的输入参数的不确定性区间。此时,所得的输入参数的不确定性区间的准确度并不能得到保证,因此,若基于此输入参数变化区间使用 MSF 和 MC 方法定量 DNDC 模型的输出结果的置信区间,则该置信区间的可靠性有待进一步评估。贝叶斯推断使用实测值来校正和改进模型参数,该方法整合了模型参数的先验概率分布和模拟结果的似然性评价,可以定量后验参数和模型输出结果的不确定性。在模型输入参数的实测数据较少的情况下,以马尔科夫链蒙特卡洛(MCMC)模拟实现的贝叶斯推断提供了一种潜在有效的不确定性评价方法。目前,该方法已被成功用于水文、生态模型的不确定性评价12-14,但在
10、 DNDC 模拟 SOC的不确定性评价方面还未见报道。本研究基于贝叶斯推断,采用 MCMC 模拟实现的贝叶斯推断对DNDC 模型模拟江苏省宜兴市水田长期定位观测点SOC 变化的输入参数进行校正,并定量评价了 SOC模拟结果的不确定性置信区间,旨在为区域或国家尺度农田 SOC 动态模拟的不确定性定量评估提供理论和方法依据。1材料与方法1.1DNDC 模型简介DNDC 模型由美国新罕布什尔州州立大学李长生教授领导的团队开发,具有 SOC 动态模拟能力,是目前被广泛接受的 SOC 动态模拟模型之一。该模型将 SOC 库分为植物残遗体、微生物、活性和惰性等 4 个库,每个库又分为两至三个具有不同分解速
11、率的子库。由用户提供的气候、土壤、作物和耕作管理措施 4 个模块的输入参数共同驱动,完成从作物生长,植物残遗体进入土壤,碳库分配及各子碳库分解过程的模拟。在这些过程中,土壤温度、湿度、pH、Eh 和各种基质(如可溶性有机碳(DOC)等)共同决定了各子碳库的大小及其特殊的分解速率,从而决定了 SOC 的动态变化过程11。关于DNDC 模型的详细介绍,可参考文献 15-17。1.2数据来源DNDC 模型的验证数据来源于江苏省宜兴市的一个土壤类型为漂洗型水稻土的农田长期定位观测点(E 1193000,N 311000)。该观测点记录了 1986 2007 年的 22 年间稻-麦轮作试验数据,包括 1
12、986 年的土壤基本属性(表 1)和 22 年间逐年的SOC 含量(图 1)、水稻小麦产量、化肥及有机肥施用量。气象数据包括 1986 2007 年日最高气温、日最低气温和日降水量均来自中国气象局国家气象信息中 心 构 建 的 气 象 数 据 共 享 服 务 网(http:/ 1 中的实测值及相应的气象数据、作物管理数据输入 DNDC 模型进行 1986 2007 年间该观测点的 SOC 动态变化模拟,然后将模拟结果与实测的逐年 SOC 数据进行对比以检验 DNDC 模拟该观测点 SOC 动态变化的可行性。MCMC 方法的先验数据来源于中国土种志中记载的 55 个漂洗型水稻土剖面信息。为了检验
13、MCMC 方法在输入数据先验知识有限、数据质量不明条件下的 DNDC 模型输入参数校正效果,本文将各输入参数的先验概率密度分布统一设定为最简单的均匀分布,分布的上下限则根据 55 个剖面信息统计所获得的土壤属性的最大值和最小值来设定2 期秦发侣等:基于贝叶斯推断的 DNDC 模型参数校正与不确定性评价研究249http:/(表 1)。以 MCMC 方法获得 SOC 的后验不确定性区间是否覆盖实测的 22 年间逐年的 SOC 含量来检验模型输出结果不确定性区间的可靠性;同时,将MCMC 方法所获得的输入参数的不确定性区间与表 1 的 6 个输入参数实测值进行比较,检验 MCMC方法的输入参数推断
14、结果是否准确。表 1DNDC 模型 6 个输入参数的实测值、先验不确定性区间的统计特征Table 1Statistics of on-site-observed values and priori uncertainty intervals of the six DNDC model input parameters初始 SOC InitialSOC(g kg1)黏粒含量Clay content(%)容重Bulk density(g cm3)pH秸秆还田量Crop residue returnedto soil(%)灌溉指数Irrigation index实测值 Site-observed v
15、alue9.3341.217.050.9先验不确定性区间Priori uncertainty interval最小值 Minimum9.130.985.550.7均值 Mean16.7331.336.7230.85最大值 Maximum24.3621.677.84011.3DNDC 模型模拟 SOC 动态变化可行性检验标准采用相关系数(r)和均方根误差(MSE)来评价模拟结果,其计算方法分别如下:r=ni=1(yoi yo)(ymi ym)ni=1(yoi yo)2ni=1(ymi ym)槡2(1)MSE=ni=1ymi y()oi2槡n(2)式(1)、式(2)中 yoi为第 i 年的 SOC
16、 实测值,ymi为SOC 模拟值,n 为 SOC 的实测和模拟年数,yo和ym分别为 SOC 的 n 年实测和模拟平均值。r 越接近1,模型模拟值与实测值之间的拟合度越高;MSE 越小,模拟结果越准确。1.4贝叶斯推断贝叶斯推断过程可表示为12,18:p(|Y)=p(Y|)p()p(Y|)p()d p(Y|)p()(3)式(3)中,=1m,2m,dm为 DNDC 模型的输入参数集,Y=y1m,y2m,ysm为输出结果集,d为模型的输入参数维数,s 为模型的输出结果维数,m 为 MCMC 方法抽取的参数集的个数。本研究仅针对 DNDC 众多输入参数中的6 个关键参数进行校正,并推断由这 6 个参
17、数的不确定性导致的 SOC 模拟的不确定性区间,因此取 d=6,s=1。p()为输入参数的先验概率密度,p(Y|)为数据似然性的概率密度,p(|Y)为输入参数的后验概率密度。式(3)表示模型输入参数的后验概率密度依比例收敛于()模型输入参数的先验概率密度与数据似然性概率密度的乘积。数据似然性 p(Y|)采用似然函数来描述,本研究中假设其服从正态分布12-14,即:p(Y|)=12槡exp(yo ym)222(4)当有 n 对实测值和模拟值时,有12,14:p(Y|)=ni=112槡iexp(yoi ymi)22i(5)式(5)中,i为第 i 个 SOC 实测值的测量误差,根据国标 HJ-615
18、-2011土壤有机碳测定重铬酸钾氧化-分光光度法,本研究假定 SOC 测量误差为实测值的10%,即 i=0.1yoi。式(4)、式(5)中其他符号的含义同式(1)和式(2)。本研究假定 DNDC 模型输入参数的先验分布为均匀分布,采用贝叶斯推断和 MCMC 方法从输入参数的先验不确定性区间中随机抽样,共建立 m=22 000 个 DNDC 模型输入参数集,将其输入 DNDC模型模拟 1986 2007 年间的 SOC 逐年变化。根据SOC 模拟结果与农田长期定位实测点相应年份的SOC 实测值的似然性(相似程度)来逆向选择满足条件的输入参数集及对应的 SOC 模拟值,最终选择了 1 940 个满
19、足条件的样点集。这些样点集被用于重建 DNDC 输入参数的后验概率分布,并统计 SOC模拟结果的不确定性置信区间,以此实现 SOC 模拟结果的不确定性定量评价。1.5Metropolis-Hastings 算法MH(Metropolis-Hastings)算法通过遍历模型输出结果的样本空间,根据输出结果与实测值的似然性来选择满足条件的样点,并通过已选择的满足条件的输出结果来逆向推断输出参数,从而重构输出250土壤学报51 卷http:/结果和输入参数后验分布的一种简单实用的 MC-MC。MH 算法中,马尔科夫链的起点可以是样本空间中的任意一点,之后按照一定的规则生成新点。在本文中,由点 生成新
20、点 的规则为 q(|)Nd(,S),N 为正态分布,d 为输入参数的维度,S为输入参数的协方差矩阵19。新点被接受的概率按下式来确定:=p(Y|)p()p(Y|)p()(6)按式(6)计算得 后,从(0,1)中产生一个随机数,若,则新点 被接受,设 =,进行下一个样点的生成;若 ,则新点 不被接受,需返回 并重新生成,直至 被接受,并设 =,进行下一个样点的生成19。如此不断地产生新样点,直至由被接受的样点组成的马尔科夫链达到收敛。ob-ert 和 osenthal21 建议通过调节S来使得样点被接受的概率介于 10%40%之间,以便达到高效的马尔科夫链收敛。1.6应用软件MH 算法的实现是一
21、个由起始样点集开始,不断生成并筛选新样点集,直至满足条件的样点集数量达到足以构建平稳的马尔科夫链的一个实时动态过程。每生成一个新样点集后需要将此样点集输入 DNDC 模型进行一次模拟,并提取模型模拟结果,计算模型模拟结果和 22 个逐年 SOC 实测值之间的似然性,以式(6)来决定新样点集可否被保留为后验样点集。MH 算法的实现需要用户编写程序。本研究采用 python 编程语言来实现 MH 算法。对于由 MH 算法所筛选出的 1 940 个样点集所构成的马尔科夫链的后续分析则采用具有数据统计分析功能的 软件来完成。软件中提供了专门分析马尔科夫链是否达到收敛的扩展包(boa)20,采用其中的
22、Heidleberger 和 Welch 方法。2结果与讨论2.1DNDC 模型的适用性评估为检验 DNDC 模型模拟农田长期定位点 SOC动态变化的可行性,将表 1 中的点位实测数据及相应的气象数据和作物管理措施数据输入 DNDC 模型进行模拟,并对比 1986 2007 逐年间的 SOC 模拟值与对应实测值(图 1)。可知 DNDC 模型基本能够重现出各年份 SOC 的动态变化趋势。DNDC 模型输出的 SOC 值与实测的 SOC 值之间的 r 为0.92,MSE 值为 0.7 g kg1,表明 DNDC 模型的 SOC 模拟值与实测值之间的拟合度较高,模型模拟精度可以接受。然而,从图 1
23、 可知 DNDC 模型在重现 SOC的剧烈波动方面还略显不足。主要表现为 DNDC模拟的 SOC 值在部分年份仍与实测值有一定的偏差,最为明显的是在 1994 年至 2007 年间模型的模拟值未能够再现实测值所呈现出的一个较为明显的波峰和波谷。对于这些偏差,可能有模型结构和数据采集分析两方面的原因。模型结构方面可能是由于 DNDC 模型现有的模型结构并没有完全地反映复杂自然环境条件下 SOC 对某些极端条件(如干旱、暴雨等近年来常见的气象事件)的变化响应机制;数据采集分析方面,可能如徐胜祥22 所指出的,长期定位实测试验中 SOC 的测定存在着很多难以克服的困难,如采样位置、深度、测定方法的微
24、小变化以及采样人员的变动等均可能导致 SOC 实测结果较为明显的年际变化。IPCC 曾建议将过程模型用作估算由人类活动引起的 SOC 变化的指导工具之一 23。通过图1 可以发现,模型模拟结果虽然能够在很大程度上反映 SOC的总体演变趋势,但在局部年份并未能完全地再现年际间的波动。因此,简单地将模型的单一模拟值用于SOC 变化的定量评价存在一定的风险。如图 1 中1994 2007 年间,SOC 的模拟值与实测值间尚有一定的偏差。因此,当将过程模型用于 SOC 动态变化估算时,稳健的方法必须要考虑模型模拟结果所存在的不确定性。而为了准确地定量并在一定程度上减小模型模拟结果的不确定性,对于模型输
25、入参数的校正及其不确定性评价就显得尤为重要。图 11986 2007 年间 SOC 的 DNDC 模型模拟值与实测值间的对比Fig.1Comparison of the DNDC model simulatedSOC with the on-site-observed SOC in the periodbetween 1986 and 20072 期秦发侣等:基于贝叶斯推断的 DNDC 模型参数校正与不确定性评价研究251http:/2.2后验马尔科夫链的平稳性和收敛性MH 算法从 6 个输入参数的先验分布区间内共抽取了 22 000 个样点集输入至 DNDC 模型进行模拟,并筛选出了其中符合
- 配套讲稿:
如PPT文件的首页显示word图标,表示该PPT已包含配套word讲稿。双击word图标可打开word文档。
- 特殊限制:
部分文档作品中含有的国旗、国徽等图片,仅作为作品整体效果示例展示,禁止商用。设计者仅对作品中独创性部分享有著作权。
- 关 键 词:
- 基于 贝叶斯 推断 DNDC 模型 参数 校正 不确定性 评价 研究
限制150内