CN102495948A - 基于mcmc的小麦品种特征参数估算方法 - Google Patents

Abstract

本发明公开了一种基于MCMC的小麦品种特征参数估算方法，其方法主要是通过田间实测小麦品种的生育期及产量与小麦生长模型模拟的生育期及产量构建概率密度函数，利用MCMC方法反演得到小麦品种特征参数。具体过程首先通过测试运算得到小麦品种特征参数的先验概率分布，根据初始参数和先验分布提出候选参数，计算产量和生育期的概率密度函数和释然比，根据M-H准则判断是否接受新的参数，最终得到了品种各特征参数的后验概率密度分布。本发明所述方法估计的结果既准确又高效，在同类模型小麦品种特征参数估算当中具有普遍适用性。

Description

基于MCMC的小麦品种特征参数估算方法

技术领域

本发明属于精准农业中小麦品种特征定量评价领域，涉及小麦品种特征参数估算方法，具体涉及融合了马尔可夫链蒙特卡洛方法(Markov Chain Monte Carlo，MCMC)的小麦品种特征参数估计方法。

背景技术

随着农业技术的发展进步，新的小麦品种不断出现，科学家通过对不同品种的生理生态特性研究，提出了不少可以量化作物品种特性的参数，如株高、千粒重、生理春化天数、灌浆时间、发生分蘖的能力等，这些参数被称为作物品种特征参数，它主要反映作物基因型特征。近半个多世纪以来，国际上出现了不少作物生长模型，如美国的CERES模型、菲律宾国际水稻所的ORYZA模型、澳大利亚的APSIM模型、中国的WheatGrow模型等，它们通过解析“气象-土壤-技术措施”与作物生理生态过程的机理关系，对作物的生长发育过程进行定量的表达，将反映基因型差异的作物品种特征参数融入到具体作物生长发育过程当中，再结合环境因子作为输入变量，达到对作物生长发育过程的动态模拟。科学家在广泛利用模型进行作物生长预测的同时，也在考虑如何逆向地利用作物生长模型进行品种特征参数的反演，以定量评估作物的品种特性，从而指导作物育种。He等利用GLUE方法和CERES-Maize来反演玉米的品种特征参数，金之庆等利用“试错法”和CERES-wheat来反演小麦品种特征参数。但这些方法都存在反演效率低且无法搜寻到全局最优变量等缺点。

发明内容

本发明的目的在于克服上述缺陷，将MCMC方法和小麦生长模型(WheatGrow)相结合，通过构建生育期和产量的目标函数，提出一种小麦品种特征参数估算方法，为定量研究小麦品种特性提供了一种新的途径。

本发明选用了马尔可夫链蒙特卡洛(Markov Chain Monte Carlo，MCMC)方法，其最早被广泛应用于图像处理、统计物理学领域。近年来，随着计算机技术的迅猛发展，MCMC方法在根系水分吸收模型、农业生态模型等诸多复杂模型中得到了广泛应用。该方法是基于贝叶斯理论框架，首先建立平衡分布为π(x)的马尔可夫链，并对其平衡分布进行采样，然后通过不断更新样本信息而使马尔可夫链能充分搜索模型参数空间，最终收敛于高概率密度区；该方法能将一些复杂的高维问题转化为一系列简单的低维问题，因此适用于复杂模型的贝叶斯计算。同时本发明选用了国内建立并具有较好广适性的WheatGrow生长模型，该模型涉及到9个小麦品种特征参数，分别是春化天数(VD)、光周期敏感因子(PS)、基本早熟因子(IE)、温度敏感因子(TS)、灌浆期因子(FD)、叶热间距(LT)、比叶面积(SLA)、收获指数(HI)、小麦分蘖能力(TA)。其中，春化天数(VD)反映了小麦品种必须经历一定天数的春化，才能进行花芽分化，极冬性的小麦品种为60，极春性的小麦品种为0。光周期敏感因子(PS)，该参数反应了小麦能否开花受到的日照长短的影响，用光周期敏感因子反映不同品种对日照长度的敏感性不同。基本早熟因子(IE)，小麦出苗后，即使满足了温光条件，也不会进行花芽分化，这段时间进行营养生长，不同品种表现不同，用基本早熟因子反映该生育阶段的差异。温度敏感因子(TS)，反映了小麦发育受温度高低影响的程度。灌浆期因子(FD)，反映了小麦开花以后不同品种灌浆时间的长短不同。叶热间距(LT)，反映了连续两个叶原基出现的时间间隔。比叶面积(SLA)，单位重量的小麦叶片所包含的叶面积。收获指数(HI)，反映了小麦干物质向穗部分配的比例。小麦分蘖能力(TA)，不同小麦品种产生分蘖的能力不同，受品种基因控制，用分蘖能力反映小麦产生分蘖多少的特征。

具体方案为：

一种基于MCMC的小麦品种特征参数估算方法，包括如下步骤：

1)小麦生长模型即WheatGrow模型的数据获取：通过自动气象站记录小麦生长季田间逐日气象数据，包括日最高气温、日最低气温、日照时数、日降水量；利用五点取样法，将耕层土壤从上到下分4层进行田间土壤取样，然后通过室内物理化学分析，测定土壤容重、饱和含水量、田间持水量、萎焉含水量、实际含水量、土壤全氮、氨态氮和硝态氮；记录小麦生长过程中的田间栽培管理数据；实测小麦主要生育期(包括播种期、拔节期、开花期、成熟期)和收获时的籽粒产量；

2)在小麦品种特征参数范围内选择一套初始品种特征参数组合θ⁽⁰⁾，在未知参数分布的情况下，先假定各个参数的先验分布为均匀分布，分别代表九个品种特征参数(i＝1，2，3，4，5，6，7，8，9)；

3)在

基础上提出一个候选样本

${\mathrm{\θ}}_i^{\mathrm{new}}={\mathrm{\θ}}_i^{k-1}+r\×(\max \left({\mathrm{\θ}}_i\right)-\min \left({\mathrm{\θ}}_i\right))/D$

r是分布在0-1之间的一个随机数；max(θ_i)和min(θ_i)分别是参数θ_i的上下限；D等于5，控制着推荐步长的大小；

4)将田间观测的n年(n＝1，2，3...N)气象数据、土壤数据、田间栽培管理数据以及两套小麦品种特征参数θ^new、θ^(k-1)输入到WheatGrow模型中，运行模型n遍，得到模拟的n年小麦主要生育期和产量结果；

5)通过田间实测的n年小麦生育期及产量数据和模拟的n年小麦生育期及产量结果，计算生育期和产量的释然函数π(θ)，计算公式如下：

$\mathrm{\π}\left(\mathrm{\θ}\right)\∝\exp \{-\frac{1}{2}{[Z\left(t\right)-\mathrm{\θX}\left(t\right)]}^T\mathrm{cov}{\left(e_t\right)}^{-1}[Z\left(t\right)-X\left(t\right)\left]\right\};$

其中，Z(t)为观察值，X(t)是模拟值，cov(et)是e_t的协方差矩阵，θ为品种参数；

6)按照如下公式计算释然比a：

7)对步骤6)计算得到的值，根据M-H准则判断是否接受新的参数，并根据判断的结果对参数进行赋值，同时保存该套参数；

8)重复步骤3)～步骤7)，直到k＝N，N为总抽样次数；

9)判断是否是首次运行，如果不是，则进入步骤11)，如果是，则计算以上步骤8)中所有参数的后验概率分布的平均数E(θ_i)和协方差cov(θ_i，θ_j)，构建一个多元正态分布N(0，cov⁰(θ))，其中cov⁰(θ)是参数θ的协方差矩阵；利用Matlab中的多元正态分布随机抽取样本mvnrnd(0，cov⁰(θ))，按照下面公式更新步骤3)中

用E(θ_i)更新步骤2)中的θ⁽⁰⁾。

$E\left({\mathrm{\θ}}_i\right)=\frac{1}{N-M+1}\underset{k=M}{\overset{N}{\mathrm{\Σ}}}{\mathrm{\θ}}_i^k$

$\mathrm{cov}({\mathrm{\θ}}_i,{\mathrm{\θ}}_j)=\frac{1}{k}{\mathrm{\Σ}}_{k=M}^N({\mathrm{\θ}}_i^k-E\left({\mathrm{\θ}}_i\right))({\mathrm{\θ}}_j^k-E\left({\mathrm{\θ}}_j\right))$

${\mathrm{\θ}}_i^{\mathrm{new}}={\mathrm{\θ}}_i^{k-1}+\mathrm{mvnrnd}(0,{\mathrm{cov}}^0\left(\mathrm{\θ}\right));$

其中M＝N/5，N为达到收敛时的抽样次数，E(θ_i)为θ_i参数的平均值，cov(θ_i，θ_j)为参数之间的协方差；

10)根据新的

重复步骤2)～步骤8)，直到k＝N，N为总抽样次数；

11)重复步骤10)完成3条马尔可夫链；

12)检验所获取的马尔可夫链是否达到收敛的标准，计算参数后验概率密度分布的平均数和标准差。

进一步优选，上述步骤2)所述的品种特征参数是：与小麦生育期相关的品种特征参数，包括春化天数VD、光周期反应的敏感性PS、基本营养生长时间IE、对温度高低的敏感性TS、灌浆时间的长短FD；与产量相关的特征参数，包括相邻两个叶原基出现的时间间隔LT，单位重量的叶片的表面积SLA，品种发生分蘖的能力TA，光合干物质向穗部分配的比例HI。

进一步优选，上述步骤7)：按照M-H准则，将计算出来的释然比a_p、a_y和均匀分布在[0，1]之间的随机变量U进行比较，对于生育期参数：如果a_p≥U，设定

(i＝1，2，3，4，5)，否则设定

对于产量参数，如果a_y≥U，设定

(i＝6，7，8，9)，否则设定

进一步优选，上述步骤12)马尔可夫链收敛的条件是：根据中心极限定理，通过M-H准则的马尔可夫链最终收敛于静态分布，用一种定量收敛判断指标GR_i来终止计算，GR_i的计算基于马尔可夫链链内和链间方差，如下式：

${\mathrm{GR}}_i=\sqrt{\frac{W_i(N-1)-B_i}{W_{i}N}};$

i代表不同的参数，B_i是链之间的变异，N是总抽样次数(k＝1，...，N)，C是链的数目(c＝1，...，C)，W_i是链内的变异，GR_i＜＝1.2通常作为链收敛的标准。

本发明的有益效果：

(1)本发明通过MCMC方法和小麦生长模型(WheatGrow)相结合，有效反演出小麦品种的特征参数，与小麦生育期相关的品种特征参数5个，包含了春化天数(VD)、对光周期反应的敏感性(PS)、基本营养生长时间(IE)、对温度高低的敏感性(TS)、灌浆时间的长短(FD)。与产量相关的特征参数4个，包含了相邻两个叶原基出现的时间间隔(LT)，单位重量叶片的表面积(SLA)，品种发生分蘖的能力(TA)，光合干物质向穗部分配的比例(HI)。在每个品种参数范围内选取一个参数，组成一套品种参数组合θ⁽⁰⁾(共9个参数)。

(2)本发明在计算生育期、产量等释然函数时，并不仅仅将各个生育期看作是独立分布的，而是考虑到了各生育期之间的相互作用，从而确保参数估计的结果更加准确。该过程中应用到蒙特卡洛抽样方法，提高了抽样的效率，Metropolis准则避免了参数陷入局部最优，采用了贝叶斯方法，有效地考虑了参数的先验分布，使得参数估计的结果既准确又高效。

附图说明

图1小麦品种特征参数估算方法的流程图；

图2WheatGrow模型的框架结构图；

图3徐麦25的特征参数抽样路径及概率分布图；

图4随机抽取1000套参数徐麦25品种特征参数的时间序列验证(1000套品种特征参数模拟值用白线表示，1000套品种特征参数模拟结果平均值用红线表示，实测值用蓝线表示)；

图5徐麦25在徐州地区99-03年生育期和产量模拟值与实测值的1∶1关系图。

具体实施方式

下面结合附图与具体实施方式对本发明作进一步详细说明。

利用本发明小麦品种特征参数估算方法在徐州地区进行了应用，用到徐州地区从1999-2003年徐麦25小麦品种的生育期和产量数据，对徐麦25的特征参数进行有效地估计。结合图1的流程进行详细说明。

1)数据的准备，WheatGrow模型(图2)的数据获取，首先通过自动气象站记录小麦生长季田间逐日气象数据，包括日最高气温、日最低气温、日照时数、日降雨量。利用五点取样法，每点从上到下将耕作层分4层进行田间土壤取样，然后通过室内行物理化学分析，测定土壤容重、饱和含水量、田间持水量、萎焉含水量、实际含水量，土壤全氮，土壤氨态氮和硝态氮(表1)。记录小麦生长过程中的田间栽培管理数据(表2)，并实测小麦主要生育期和收获时的籽粒产量(表3)；

表1徐州地区土壤数据

表2徐麦25在徐州地区的栽培管理措施

表3徐麦25在徐州地区实测产量和生育期

2)在每个品种特征参数范围内各随机选取一个初始参数，组成一套初始参数组合θ⁽⁰⁾(表4)，在未知品种参数概率分布的情况下，先假设其先验分布

(k＝1，2，3，...N)为均匀分布。

分别代表9个品种特征参数(i＝1，2，3，4，5，6，7，8，9)。其中与小麦生育期相关的小麦品种特征参数，包含生理春化天数(VD)、对光周期反应的敏感性(PS)、基本营养生长时间(IE)、对温度高低的敏感性(TS)、灌浆时间的长短(FD)；与产量相关的特征参数，包含了相邻两个叶原基出现的时间间隔(LT)，单位重量叶片的表面积(SLA)，品种发生分蘖的能力(TA)，光合干物质向穗部分配的比例(HI)。

表4徐麦25品种特征参数范围、初始参数值和最终估算出的参数平均值

3)在

基础上提出一个候选样本

${\mathrm{\θ}}_i^{\mathrm{new}}={\mathrm{\θ}}_i^{k-1}+r\×(\max \left({\mathrm{\θ}}_i\right)-\min \left({\mathrm{\θ}}_i\right))/D$

r是分布在0-1之间的一个随机数；max(θ_i)和min(θ_i)分别是参数θ_i的上下限(表4)；D等于5，控制着推荐步长的大小。

4)将田间观测气象数据、土壤数据、田间栽培管理数据以及两套小麦品种特征参数θ^new、θ^(k-1)输入到WheatGrow模型中，运行模型5年，得到模型模拟的小麦生育期和产量结果。

5)利用步骤1)获得的实测小麦生育期和实测小麦田间产量(表3)，以及模拟得到的生育期和产量结果，分别计算生育期和产量的释然函数π_p和π_y：

$\mathrm{\π}\left(\mathrm{\θ}\right)\∝\exp \{-\frac{1}{2}{[Z\left(t\right)-\mathrm{\θX}\left(t\right)]}^T\mathrm{cov}{\left(e_t\right)}^{-1}[Z\left(t\right)-\mathrm{\θX}\left(t\right)\left]\right\},$

Z(t)为观察值，θX(t)是模拟值，cov(e_t)是e_t的协方差矩阵，θ为品种参数；

生育期考虑到了拔节期、开花期、成熟期三个生育期，需要将π(θ)扩展为一个三元概率密度函数如下：

$\mathrm{\α}=1+{2\mathrm{\ρ}}_{12}{\mathrm{\ρ}}_{13}{\mathrm{\ρ}}_{23}-{\mathrm{\ρ}}_{13}^2-{\mathrm{\ρ}}_{12}^2-{\mathrm{\ρ}}_{23}^2,$

${\mathrm{\σ}}_i=(1/N)\left[{\mathrm{\Σ}}_{i=1}^N{(x_i-\stackrel{\‾}{x})}^2\right],$

$\mathrm{\ρ}=\frac{{\mathrm{\Σ}}_{i=1}^N(x_i-\stackrel{\‾}{x})(y_i-\stackrel{\‾}{y})}{\sqrt{{\mathrm{\Σ}}_{i=1}^N{(x_i-\stackrel{\‾}{x})}^2}\sqrt{{\mathrm{\Σ}}_{i=1}^N{(y_i-\stackrel{\‾}{y})}^2}},$

$x_j={\mathrm{\Σ}}_{t=1}^S{[O_j\left(t\right)-S_j\left(t\right)]}^2,(j=\mathrm{1,2,3}),$

产量设计一个一元概率密度函数如下：

${\mathrm{\π}}_y\left(\mathrm{\θ}\right)\∝\exp \{-\frac{1}{{2\mathrm{\σ}}_4}{\mathrm{\Σ}}_{t=1}^S{[O_j\left(t\right)-S_j\left(t\right)]}^2\},(j=4),$

σ₁，σ₂，σ₃，σ₄分别是5年小麦拔节期、开花期、成熟期、产量的模拟和实测误差方差。ρ₁₂、ρ₁₃、ρ₂₃是5年小麦拔节期、开花期、成熟期之间的模拟和实测误差相关系数。其中，

表示样本的平均数，N表示样本的数量，xi表示个体。O_j(t)(j＝1，2，3)和S_j(t)(j＝1，2，3)分别是指第t年小麦生育期的观察值和模拟值，j＝1，2，3分别代表着拔节期、开花期、成熟期。O_j(t)(j＝4)和S_j(t)(j＝4)表示第t年小麦产量的观察值和模拟值；

6)计算生育期和产量的释然比a_p和a_v：

7)按照M-H准则，将计算出来的值和均匀分布在[0，1]之间的随机变量U进行比较。对于生育期参数：如果a_p≥U，设定(i＝1，2，3，4，5)，否则设定

对于产量参数：如果a_y≥U，设定

(i＝6，7，8，9)，否则设定

并保存该套参数θ^k；

8)重复3)～7)，直到k＝30000；

9)判断是否是首次运行，如果不是，则进入步骤11)，如果是，则计算以上步骤8)中所有参数的后验概率分布的平均数E(θ_i)和协方差cov(θ_i，θ_j)：

$E\left({\mathrm{\θ}}_i\right)=\frac{1}{N-M+1}{\mathrm{\Σ}}_{k=M}^N{\mathrm{\θ}}_i^k$

$\mathrm{cov}({\mathrm{\θ}}_i,{\mathrm{\θ}}_j)=\frac{1}{k}{\mathrm{\Σ}}_{k=M}^N({\mathrm{\θ}}_i^k-E\left({\mathrm{\θ}}_i\right))({\mathrm{\θ}}_j^k-E\left({\mathrm{\θ}}_j\right))$

其中N＝M/5，E(θ_i)为θ_i参数的平均值，cov(θ_i，θ_j)为参数之间的协方差(表5)。构建参数的多元正态分布N(0，cov0(θ))，这里cov⁰(θ)是参数θ的协方差矩阵；利用Matlab中的多元正态分布随机抽取样本mvnrnd(0，cov⁰(θ))，按照下面公式更新步骤3)中的

用E(θ_i)更新步骤二中的θ⁽⁰⁾；

${\mathrm{\θ}}_i^{\mathrm{new}}={\mathrm{\θ}}_i^{k-1}+\mathrm{mvnrnd}(0,\mathrm{cov}0\left(\mathrm{\θ}\right));$

θ(0)＝{26，1.45，0.007，0.72，0.62，0.0025，0.36，65，0.6}；

表5首次运算得到的小麦品种特征参数的协方差矩阵

10)根据新的

重复步骤2)～步骤7)，直到k＝30000；

11)重复步骤10)完成3条马尔可夫链；

12)九个参数的GR_i均小于1.2，说明抽样结果达到收敛。舍弃每条马尔可夫链的前10000次，计算参数后验概率密度分布的平均数和方差。

图3表示所有参数的抽样路径和概率分布情况以及表6中包含所有参数的平均值和标准差。最终得到徐麦25每个特征参数的概率分布情况，在得到的参数的后验概率密度分布中，将模拟徐麦25各个品种参数与试验获取参数进行对比发现，模拟参数的平均值与试验得到的品种参数比较一致，在每个模拟参数范围内随机抽样1000次，然后输入小麦生长模型中进行模拟计算，并将得到小麦生育期和产量模拟结果的平均值与实测值进行比较分析，如图4和5。结果表明模拟值与实测结果比较一致，说明基于MCMC的参数估算方法可以用于小麦品种特征参数的估算。

表6模拟与实测小麦品种参数对照

Claims (4)

1.一种基于MCMC的小麦品种特征参数估算方法，其特征在于，包括如下步骤：

1)小麦生长模型即WheatGrow模型的数据获取：通过自动气象站记录小麦生长季田间逐日气象数据，包括日最高气温、日最低气温、日照时数、日降水量；利用五点取样法，将耕层土壤从上到下分4层进行田间土壤取样，然后通过室内物理化学分析，测定土壤容重、饱和含水量、田间持水量、萎焉含水量、实际含水量、土壤全氮、氨态氮和硝态氮；记录小麦生长过程中的田间栽培管理数据；实测小麦主要生育期(包括播种期、拔节期、开花期、成熟期)和收获时的籽粒产量；

2)在小麦品种特征参数范围内选择一套初始品种特征参数组合θ⁽⁰⁾，在未知参数分布的情况下，先假定各个参数的先验分布

为均匀分布，

分别代表九个品种特征参数(i＝1，2，3，4，5，6，7，8，9)；

3)在

基础上提出一个候选样本

${\mathrm{\θ}}_i^{\mathrm{new}}={\mathrm{\θ}}_i^{k-1}+r\×(\max \left(\mathrm{\θi}\right)-\min \left(\mathrm{\θi}\right))/D$

r是分布在0-1之间的一个随机数；max(θ_i)和min(θ_i)分别是参数θ_i的上下限；D等于5，控制着推荐步长的大小；

4)将田间观测n年(n＝1，2，3...N)气象数据、土壤数据、田间栽培管理数据以及两套小麦品种特征参数θ^new、θ^(k-1)输入到WheatGrow模型中，运行模型n遍，n＝1，2，3...N，得到模拟的n年小麦主要生育期和产量结果；

5)通过田间实测的n年小麦生育期及产量数据和模拟的n年小麦生育期及产量结果计算生育期和产量的释然函数π(θ)，计算公式如下：

$\mathrm{\π}\left(\mathrm{\θ}\right)\∝\exp \{-\frac{1}{2}{[Z\left(t\right)-\mathrm{\θX}\left(t\right)]}^T\mathrm{cov}{\left(e_t\right)}^{-1}[Z\left(t\right)-X\left(t\right)\left]\right\};$

其中，Z(t)为观察值，X(t)是模拟值，cov(et)是et的协方差矩阵，θ为品种参数；

6)按照如下公式计算释然比a：

7)对步骤6)计算得到的值，根据M-H准则判断是否接受新的参数，并根据判断的结果对参数进行赋值，同时并保存该套参数；

8)重复步骤3)～步骤7)，直到k＝N，N为总抽样次数；

9)判断是否是首次运行，如果不是，则进入步骤11)，如果是，则计算以上步骤8)中所有参数的后验概率分布的平均数E(θ_i)和协方差cov(θ_i，θ_j)，构建一个多元正态分布N(0，cov⁰(θ))，其中cov⁰(θ)是参数θ的协方差矩阵；利用Matlab中的多元正态分布随机抽取样本mvnrnd(0，cov⁰(θ))，按照下面公式更新步骤3)中

用E(θ_i)更新步骤2)中的θ⁽⁰⁾。

$E\left({\mathrm{\θ}}_i\right)=\frac{1}{N-M+1}\underset{k=M}{\overset{N}{\mathrm{\Σ}}}{\mathrm{\θ}}_i^k$

$\mathrm{cov}({\mathrm{\θ}}_i,{\mathrm{\θ}}_j)=\frac{1}{k}{\mathrm{\Σ}}_{k=M}^N({\mathrm{\θ}}_i^k-E\left({\mathrm{\θ}}_i\right))({\mathrm{\θ}}_j^k-E\left({\mathrm{\θ}}_j\right))$

${\mathrm{\θ}}_i^{\mathrm{new}}={\mathrm{\θ}}_i^{k-1}+\mathrm{mvnrnd}(0,\mathrm{cov}0\left(\mathrm{\θ}\right));$

其中M＝N/5，N为达到收敛时的抽样次数，E(θ_i)为θ_i参数的平均值，cov(θ_i，θ_j)为参数之间的协方差；

10)根据新的

重复步骤2)～步骤7)，直到k＝N，N为总抽样次数；

11)重复步骤10)完成3条马尔可夫链；

12)检验所获取的马尔可夫链是否达到收敛的标准，计算参数后验概率密度分布的平均数和标准差。

2.根据权利要求1所述的基于MCMC的小麦品种特征参数估算方法，其特征在于，上述步骤2)所述的品种特征参数是：与小麦生育期相关的小麦品种特征参数，包括了春化天数VD、对光周期反应的敏感性PS、基本营养生长时间IE、对温度高低的敏感性TS、灌浆时间的长短FD；与产量相关的特征参数，包括了相邻两个叶原基出现的时间间隔LT，单位重量的叶片的表面积SLA，品种发生分蘖的能力TA，光合干物质向穗部分配的比例HI。

3.根据权利要求1所述的基于MCMC的小麦品种特征参数估算方法，其特征在于，上述步骤7)：按照M-H准则，将计算出来的释然比a_p、a_y和均匀分布在[0，1]之间的随机变量U进行比较，对于生育期参数：如果a_p≥U，设定

(i＝1，2，3，4，5)，否则设定对于产量参数，如果a_y≥U，设定

(i＝6，7，8，9)，否则设定

4.根据权利要求1所述的基于MCMC的小麦品种特征参数估算方法，其特征在于，上述步骤12)马尔可夫链收敛的条件是根据中心极限定理通过M-H准则的马儿可夫链最终收敛于静态分布，用一种定量收敛判断指标GR_i对终止计算进行判断，该GR_i指标的计算基于马尔可夫链链内和链间方差，如下式：

${\mathrm{GR}}_i=\sqrt{\frac{W_i(N-1)-B_i}{W_{i}N}};$

i代表不同的参数，B_i是链之间的变异，N是总抽样次数(k＝1，...，N)，C是链的数目(c＝1，...，C)，W_i是链内的变异，GR_i＜＝1.2通常作为链收敛的标准。

Images