CN107423561A - 一种土壤属性插值的估算方法 - Google Patents
一种土壤属性插值的估算方法 Download PDFInfo
- Publication number
- CN107423561A CN107423561A CN201710561915.3A CN201710561915A CN107423561A CN 107423561 A CN107423561 A CN 107423561A CN 201710561915 A CN201710561915 A CN 201710561915A CN 107423561 A CN107423561 A CN 107423561A
- Authority
- CN
- China
- Prior art keywords
- mrow
- msub
- point
- measured
- value
- Prior art date
- Legal status (The legal status is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the status listed.)
- Pending
Links
- 239000002689 soil Substances 0.000 title claims abstract description 62
- 238000011156 evaluation Methods 0.000 title claims abstract description 9
- 238000000034 method Methods 0.000 claims abstract description 75
- 230000007613 environmental effect Effects 0.000 claims abstract description 47
- 238000009826 distribution Methods 0.000 claims abstract description 37
- 238000005070 sampling Methods 0.000 claims description 24
- 238000011160 research Methods 0.000 claims description 11
- 238000005315 distribution function Methods 0.000 claims description 9
- 238000013507 mapping Methods 0.000 claims description 5
- 238000010606 normalization Methods 0.000 claims description 5
- 238000004364 calculation method Methods 0.000 claims description 4
- 238000004519 manufacturing process Methods 0.000 claims description 3
- 230000008569 process Effects 0.000 abstract description 6
- IJGRMHOSHXDMSA-UHFFFAOYSA-N Atomic nitrogen Chemical compound N#N IJGRMHOSHXDMSA-UHFFFAOYSA-N 0.000 description 4
- 230000000875 corresponding effect Effects 0.000 description 4
- 238000012986 modification Methods 0.000 description 3
- 230000004048 modification Effects 0.000 description 3
- 238000005904 alkaline hydrolysis reaction Methods 0.000 description 2
- 230000000694 effects Effects 0.000 description 2
- 230000004720 fertilization Effects 0.000 description 2
- 238000009499 grossing Methods 0.000 description 2
- 229910052757 nitrogen Inorganic materials 0.000 description 2
- 238000005457 optimization Methods 0.000 description 2
- 230000008520 organization Effects 0.000 description 2
- 239000004016 soil organic matter Substances 0.000 description 2
- 238000012360 testing method Methods 0.000 description 2
- 238000012271 agricultural production Methods 0.000 description 1
- 230000009286 beneficial effect Effects 0.000 description 1
- 238000006243 chemical reaction Methods 0.000 description 1
- 230000002596 correlated effect Effects 0.000 description 1
- 238000013480 data collection Methods 0.000 description 1
- 238000013501 data transformation Methods 0.000 description 1
- 230000007547 defect Effects 0.000 description 1
- 238000011161 development Methods 0.000 description 1
- 238000010586 diagram Methods 0.000 description 1
- 229910001385 heavy metal Inorganic materials 0.000 description 1
- 238000011835 investigation Methods 0.000 description 1
- 238000012886 linear function Methods 0.000 description 1
- 235000021049 nutrient content Nutrition 0.000 description 1
- 239000005416 organic matter Substances 0.000 description 1
- 230000009467 reduction Effects 0.000 description 1
- 238000012732 spatial analysis Methods 0.000 description 1
- 238000007619 statistical method Methods 0.000 description 1
Classifications
-
- G—PHYSICS
- G16—INFORMATION AND COMMUNICATION TECHNOLOGY [ICT] SPECIALLY ADAPTED FOR SPECIFIC APPLICATION FIELDS
- G16Z—INFORMATION AND COMMUNICATION TECHNOLOGY [ICT] SPECIALLY ADAPTED FOR SPECIFIC APPLICATION FIELDS, NOT OTHERWISE PROVIDED FOR
- G16Z99/00—Subject matter not provided for in other main groups of this subclass
Landscapes
- Management, Administration, Business Operations System, And Electronic Commerce (AREA)
Abstract
本发明公开了一种土壤属性插值的估算方法，包括以下步骤：S1：利用克里金法生成土壤属性待测变量的先验概率分布；S2：利用环境因素与预测点的相关关系生成软数据；S3：利用BME法更新土壤属性待测变量的先验概率分布；S4：利用土壤属性待测变量的先验概率分布进行土壤空间属性制图。本发明将操作简单、成熟的克里金法与精确、合理的贝叶斯最大熵法进行结合，使得插值过程更加简便，插值结果更加合理。
Description
技术领域
本发明涉及土壤学研究领域，具体涉及一种土壤属性插值的估算方法。
背景技术
土壤连续属性(土壤养分含量、重金属含量)的分布特征和定量分布信息是进行土壤生产力评价及环境综合评估的基础。近年来，随着精准农业的战略的实施，对土壤属性信息的精度提出了更高的要求。因此如何提高土壤属性的预测精度，为农业生产、生态评价提供精确的土壤属性信息，已成为土壤学研究的热点问题。
以克里金方法为代表的经典地统计学方法在许多学科领域都得到了广泛的应用，但就其方法本身来说还存在一定的缺陷，制约了预测结果的精确度，如：
(1)经典的克里金方法仅考虑变量自身的采样数据，而忽略了与变量相关的其他信息(如环境因素、专家经验、历史资料等)的利用，造成了预测结果的降低，在采样点稀疏的情况下还会造成预测结果出现非常不合理的情况，产生了较大的不确定性；
(2)克里金方法具有一定的平滑效应，针对短距离内变化剧烈的区域，克里金的平滑效应会将降低或完全消除这些变化在空间上的反应，导致信息丢失。
(3)克里金方法要求区域变量满足高斯分布，而在实践中通过数据变换单点的高斯分布容易满足，多点的高斯分布却难以满足。另外数据中往往存在一些离群点与异常点，这些特殊的采样点的存在也将极大的影响模型的精度。
(4)克里金方法为局部估计方法，对估计值的整体空间相关性考虑不够，它仅保证了数据的估计局部最优而不能保证数据估计的全体最优。这就要求进行克里金方法预测的采样数据需要在空间上分布均匀且达到一定的采样密度。
新兴的BME方法虽然在空间属性预测上取得了较好的成果，但是其发展时间较短，缺少相应的可视化软件的支持，如主流的空间分析及制图软件Arcgis和QGIS都不支持该方法。BME法的另一大缺点是在软数据的获取上，要搜集组织出与待测变量相关软数据需要对研究区域和待测变量的特性有较深的了解，因而在同一区域的空间变量预测可能由研究者的对该区域的了解不同造成最终结果的不一致。总而言之，BME方法虽然可以综合利用各方面的数据，但是对这些数据的来源、类型、组织方式等都没有具体的要求，该方法的具体实施还需要对具体研究区域和预测变量进行分析，可以说其实施过程比经典的克里金法要复杂的多，从而限制了BME方法在土壤、环境等领域的利用。
发明内容
本发明的目的在于提供一种土壤属性插值的估算方法，解决目前的土壤属性插值的估算采用单纯的克里金法，估算精度较低的问题。
为解决上述的技术问题，本发明采用以下技术方案：
一种土壤属性插值的估算方法，包括以下步骤：
S1：利用克里金法生成土壤属性待测变量的先验概率分布；
S2：利用环境因素与预测点的相关关系生成软数据；
S3：利用BME法更新土壤属性待测变量的先验概率分布；
S4：利用土壤属性待测变量的先验概率分布进行土壤空间属性制图。
更进一步的方案是，上述S1步骤中利用克里金法生产土壤属性待测变量的先验概率分布的具体方法是：
S101：按照S4步骤中的制图精度要求将土壤研究区域划分为row×col的网格，设土壤研究区域有N个采样点，则第k(k<N)个采样点的属性值可表示为χij，其中i和j为采样点在row×col的网格中的位置，i<row,j<col；对于待测点pij，使用实验变差函数确定待测点周围各点的权重值，计算出每个待测点的属性估值
S102：利用每个待测点的属性估值得到每个待测点的先验概率分布函数。
更进一步的方案是，上述S101步骤中计算每个待测点的属性估值的使用以下公式计算得出：其中X(pn)为待测点pij周围各观测点的属性值；pn为待测点周围各观测点的空间位置；N为待测点周围观测点的数量；λn为观测点属性值X(pn)的权重；n的取值为1～N的正整数。
更进一步的方案是，上述S101步骤中通过以下线性方程组求得观测点属性值X(pn)的权重λn：
其中p0表示待观测点的位置；pi和pj为待测点周围各观测点的空间位置，i、j取值为1～n；表示实验变差函数。
更进一步的方案是，上述S101步骤中设定实验变差函数的取值只与空间中的相对位置有关，因此使用以下公式计算实验变差函数：
其中，N(h)为空间中距离为h的观测变量的对数；pi为观测点在空间中的位置；X(pi)为待测点pi周围各观测点的属性值；
计算出实验变差函数的值后使用变差函数理论模型对实验变差函数值进行拟合，得到一条最优的实验变差函数。
更进一步的方案是，上述S102步骤中利用每个待测点的属性估值得到每个待测点的先验概率分布函数的公式是：
其中，G表示由克里金法得到的概率分布；μ为克里金的预测结果；σ2为预测方差。
更进一步的方案是，上述S2步骤中利用环境因素与预测点的相关关系生成软数据的具体方法是：
S201：将实验区域划分为N个属性值的取值区间，计算第k个属性取值区间对该点上的单个环境因素的隶属度Mk,si,j，计算公式为其中，Ck为属性值属于第k个区间的采样点数；C1为属性值属于第k个区间且该点的第j个环境因子属于第1个环境因子区间的采样点数；
S202：计算环境因子对待测点属性值的决定系数Rj 2，其中是由环境因子j与对应的函数关系确定的属性值估计值；是采样点属性值的期望；χi为采样点属性值的实际值并将决定系数Rj 2进行归一化，
S203：在研究区域随机生成分布均匀的一定数量的软数据点，对于某一个软数据点psi，由S201步骤中确定的隶属度函数确定某个属性值取值区间对该软数据点psi上的单个环境因子的隶属度Mk,si,j，并结合环境因子对待测点属性值的决定系数，确定该软数据点psi的取值范围对某个属性值取值区间的相似度从而计算待测点属性值在软数据点上的分布函数fs(psi)＝((S1,si,N1),…,(Sk,si,Nk),…,(SN,si,NN))，其中fs(psi)表示由环境数据确定属性值概率密度分布；Nk表示属性值第k个取值区间。
更进一步的方案是，上述S3步骤中，利用BME法更新土壤属性待测变量的先验概率分布采用下述公式：
其中，fG为先验概率密度函数；fs为软数据概率密度函数。
更进一步的方案是，上述S4步骤中利用土壤属性待测变量的先验概率分布进行土壤空间属性制图的具体方法是：
根据S3步骤中获得的土壤属性待测变量的先验概率分布函数，确定各待测点的预测属性值，预测属性值由下述公式确定：
其中，即为预测值，它满足当时，函数fK取得极大值，也就说明在为值p0处，随机变量χ0最有可能的取值是
按照上述公式计算出各待测点的预测属性值后，得到属性分布图。
更进一步的方案是，上述S2步骤中，对于描述型环境因子数据，对每一个描述型环境因子赋予唯一的数字编号，根据编号的数量将采样点的属性值划分为与编号数量相等的取值区间，计算每个取值区间对各描述型环境因子的隶属度。
与现有技术相比，本发明的有益效果是：
1、本发明将操作简单、成熟的克里金法与精确、合理的贝叶斯最大熵法进行结合，使得插值过程更加简便，插值结果更加合理。
2、本发明的方法能够有效的利用各方面的数据，从而提高插值精度。
3、本发明的插值过程中考虑到了环境因素对土壤属性值的影响，将其作为软数据参与到插值的过程中来，解决了克里金法在采样点稀疏或分布不均时插值结果不稳定、不合理的问题，降低了对数据的需求。
4、本发明以环境因子与属性值的关系定量的产生软数据，降低了软数据搜集、组织的难度。
具体实施方式
为了使本发明的目的、技术方案及优点更加清楚明白，以下结合实施例，对本发明进行进一步详细说明。应当理解，此处所描述的具体实施例仅仅用以解释本发明，并不用于限定本发明。
具体实施例：
1、使用克里金法生成待测变量的先验概率分布
按照制图的精度要求将研究区域划分为row×col(这里以二维平面插值为例)的格网，设该区域有N个采样点，则第k(k<N)个采样点的属性值可表示为χij，其中i<row,j<col；对于待测点pij，其属性的估计值可以使用一下公式计算得出：
式中：X(pn)为待测点pij周围各观测点的属性值；pn为待测点周围各观测点的空间位置；N为待测点周围观测点的数量；λn为观测点属性值X(pn)的权重。
上式中确定待测点估计值的关键在于确定λn，可以通过解以下线性方程组求得：
式中：p0表示待测点的位置；表示变差函数，在地统计中，我们先前已经假定了变差函数的取值只和空间中的相对位置有关，所以可以使用下式计算出实验变差函数：
式中：N(h)为空间中距离为h的观测变量的对数。pi为观测点在空间中的位置。
计算出实验变差函数的值以后，使用理论模型对实验变差函数值进行拟合，得到一条最优的实验变差函数。
使用变差函数确定待测点周围各点的权重值，从而确定待测点的估计值。使用这种方法得到格网中每个点的属性值及其先验概率密度函数：
式中：G表示由克里金法得到的概率分布；μ为克里金的预测结果；σ2为预测方差。
2、利用相关数据生成软数据
可以利用以下3种方法生成软数据：
(1)利用环境因素与预测点的相关关系
假设预测变量受环境因素的影响，那么在空间点pij的属性分布就可以在一定程度上用该点的环境因
子Eij(Eij,1,…,Eij,m)来表示，按照这一推断我们就可以在区域中利用某一点的环境因子来近似的表达待测点区域化变量的分布。
首先根据制图精度的要求，将采样点的属性值划分为N个区间，对于第k个区间，其属性值含量为：
式中：χmin，χmax为采样点属性值的最小取值和最大取值。
划分出属性值的取值区间后，就可以定量的分析属性取值区间与环境因子的关系。方法是将区域中与待测变量取值有关的环境因素按照同样的方法划分为N个区间，那么对于属性值的第k个区间，它对某一环境因子(例如第j个环境因子)取值区间的隶属度可用以下数据序列表示：
式中：Ck为属性值属于第k个区间的采样点数；C1为属性值属于第k个区间且该点的第j个环境因子属于第1个环境因子区间的采样点数。将上式的计算结果以隶属度为纵坐标，环境因子的区间序号为横坐标用折线表示出来，采用一条平滑的曲线拟合该折线，从而得到采样点属性值范围对第j个环境因子范围的隶属度函数。使用该方法最终得到N×m个隶属度函数。
另外，考虑到每个环境因子对待测点属性值的影响程度不同，这里还需要计算这m个环境因子对待测点属性值的决定系数：
首先由采样点属性值与对应的环境因子的取值建立线性函数关系，按照上式计算环境因子j的决定系数。式中：是由环境因子j与对应的函数关系确定的属性值估计值；是采样点属性值的期望；χi为采样点属性值的实际值。
按照上述方法计算出m个决定系数：
(R1 2,…,Rm 2)
随后对决定系数进行归一化：
这样以来，便可以在研究区域随机生成分布均匀的一定数量的软数据点，对于某一个软数据点psi，可由之前确定的隶属度函数确定某个属性值取值区间对该点上单个环境因子的隶属度Mk,si,j，并考虑到环境因子对属性值的影响程度，使用下式确定该点取值范围对属性值某一取值区间的相似度：
对上式计算的Sk,si进行归一化处理，即软数据点上的属性值分布函数可以用下式表达：
fs(psi)＝((S1,si,N1),…,(Sk,si,Nk),…,(SN,si,NN))
式中:fs(psi)表示由环境数据确定属性值概率密度分布；Nk表示属性值第k个取值区间。
上述关于环境因子的计算主要是数值型环境因子，对于描述型(范畴型)环境因子，需要对该方法进行一些改进，如对土壤碱解氮进行空间统计插值时，根据以往的经验我们可以知道碱解氮的含量与土类或是区域的施肥措施存在一定的联系，而土壤类型和施肥措施属于描述型的环境因子数据，无法直接套用第三项的计算方法，这里的做法是对于每一个土类或者施肥措施，赋予其唯一的数字编号，根据编号的数量将采样点的属性值划分为与编号数量相等的取值区间，这样就可以计算每个取值区间对土类编号的隶属度，从而进行后续的计算。在这一步需要注意的是，由于环境因子中描述性数据的存在，使得我们划分属性取值区间时就不能任意选择区间的数量，土壤属性的取值区间与环境因子的取值区间应统一为描述型环境因子数据的最小编号数量，这样就可以将数值型与描述型的环境因子以相同的方式进行计算，从而拓展了环境因子的选取范围。
(2)使用区域历史图件生成软数据
利用研究区域的历史图件也可生成软数据，例如进行土壤某种属性插值时，往往可以获取到土壤属性的历史分布图件，则可以据此生成软数据点psi，其取值范围为该点所在历史土壤属性分布图的分级范围。以此方法生成区间软数据集。
(3)使用其他与预测变量相关的数据生成软数据
如缺少上述两种方法所需要的数据，也可以使用其他与预测属性值相关的数据生成软数据，如对土壤有机质进行空间统计插值时，由以往的经验可以得知土壤的有机质含量与土地利用\覆被存在一定的联系，因而可以利用土地利用\覆被图或直接利用遥感影像，由土壤有机质含量与土地利用\覆被的关系生成均匀的区间软数据点。
3、使用BME法更新待测变量的先验概率分布
根据贝叶斯条件概率的原理，我们可以在考虑软硬数据的情况下更新由克里金法得到的先验概率密度函数。如预测空间变量χ0在p0位置处考虑软硬数据的完全概率分布为：
式中：fG为先验概率密度函数，由克里金法得到；fs为软数据概率密度函数，由第二步确定。
该公式给出了在空间任意位置pij属性值χij的完全概率分布，基于该公式可以对其进行预测。
4、空间属性制图
有了上一步得到的预测变量完全概率分布，就可以在特定的位置上对属性值进行预测，需要制图时，由于第1步已经将研究区域划分为row×col的格网，那么我们可以将每个格网的节点作为一个像素点，像素值即为预测的属性值。预测的属性值可由下式确定：
式中：即为预测值，它满足当时，函数fK取得极大值，也就说明在为值p0处，随机变量χ0最有可能的取值是按照这种方法计算空间格网中的各个节点的属性预测值，得到属性分布图，这之后就可以按照制图的要求，对结果进行分级以制作各种专题地图。
尽管这里参照本发明的解释性实施例对本发明进行了描述，但是，应该理解，本领域技术人员可以设计出很多其他的修改和实施方式，这些修改和实施方式将落在本申请公开的原则范围和精神之内。更具体地说，在本申请公开和权利要求的范围内，可以对主题组合布局的组成部件和/或布局进行多种变型和改进。除了对组成部件和/或布局进行的变形和改进外，对于本领域技术人员来说，其他的用途也将是明显的。
Claims (10)
1.一种土壤属性插值的估算方法，其特征在于：包括以下步骤：
S1：利用克里金法生成土壤属性待测变量的先验概率分布；
S2：利用环境因素与预测点的相关关系生成软数据；
S3：利用BME法更新土壤属性待测变量的先验概率分布；
S4：利用土壤属性待测变量的先验概率分布进行土壤空间属性制图。
2.根据权利要求1所述的土壤属性插值的估算方法，其特征在于：所述S1步骤中利用克里金法生产土壤属性待测变量的先验概率分布的具体方法是：
S101：按照S4步骤中的制图精度要求将土壤研究区域划分为row×col的网格，设土壤研究区域有N个采样点，则第k(k<N)个采样点的属性值可表示为χij，其中i和j为采样点在row×col的网格中的位置，i<row,j<col；对于待测点pij，使用实验变差函数确定待测点周围各点的权重值，计算出每个待测点的属性估值
S102：利用每个待测点的属性估值得到每个待测点的先验概率分布函数。
3.根据权利要求2所述的土壤属性插值的估算方法，其特征在于：所述S101步骤中计算每个待测点的属性估值的使用以下公式计算得出：其中X(pn)为待测点pij周围各观测点的属性值；pn为待测点周围各观测点的空间位置；N为待测点周围观测点的数量；λn为观测点属性值X(pn)的权重；n的取值为1～N的正整数。
4.根据权利要求3所述的土壤属性插值的估算方法，其特征在于：所述S101步骤中通过以下线性方程组求得观测点属性值X(pn)的权重λn：
<mrow>
<mo>{</mo>
<mrow>
<mtable>
<mtr>
<mtd>
<mrow>
<msubsup>
<mo>&Sigma;</mo>
<mrow>
<mi>i</mi>
<mo>=</mo>
<mn>1</mn>
</mrow>
<mi>n</mi>
</msubsup>
<mover>
<mi>&gamma;</mi>
<mo>^</mo>
</mover>
<mrow>
<mo>(</mo>
<msub>
<mi>p</mi>
<mi>i</mi>
</msub>
<mo>-</mo>
<msub>
<mi>p</mi>
<mi>j</mi>
</msub>
<mo>)</mo>
</mrow>
<msub>
<mi>&lambda;</mi>
<mi>i</mi>
</msub>
<mo>+</mo>
<mi>&mu;</mi>
<mo>=</mo>
<mover>
<mi>&gamma;</mi>
<mo>^</mo>
</mover>
<mrow>
<mo>(</mo>
<msub>
<mi>p</mi>
<mn>0</mn>
</msub>
<mo>-</mo>
<msub>
<mi>p</mi>
<mi>j</mi>
</msub>
<mo>)</mo>
</mrow>
<mo>,</mo>
<mrow>
<mo>(</mo>
<mi>j</mi>
<mo>=</mo>
<mn>1</mn>
<mo>,</mo>
<mn>...</mn>
<mi>n</mi>
<mo>)</mo>
</mrow>
</mrow>
</mtd>
</mtr>
<mtr>
<mtd>
<mrow>
<msubsup>
<mo>&Sigma;</mo>
<mi>i</mi>
<mi>N</mi>
</msubsup>
<msub>
<mi>&lambda;</mi>
<mi>i</mi>
</msub>
<mo>=</mo>
<mn>1</mn>
</mrow>
</mtd>
</mtr>
</mtable>
<mo>,</mo>
</mrow>
</mrow>
其中p0表示待观测点的位置；pi和pj为待测点周围各观测点的空间位置，i、j取值为1～n；表示实验变差函数。
5.根据权利要求4所述的土壤属性插值的估算方法，其特征在于：所述S101步骤中设定实验变差函数的取值只与空间中的相对位置有关，因此使用以下公式计算实验变差函数：
<mrow>
<msup>
<mover>
<mi>&gamma;</mi>
<mo>^</mo>
</mover>
<mo>*</mo>
</msup>
<mrow>
<mo>(</mo>
<mi>h</mi>
<mo>)</mo>
</mrow>
<mo>=</mo>
<mfrac>
<mn>1</mn>
<mrow>
<mn>2</mn>
<mi>N</mi>
<mrow>
<mo>(</mo>
<mi>h</mi>
<mo>)</mo>
</mrow>
</mrow>
</mfrac>
<msubsup>
<mo>&Sigma;</mo>
<mrow>
<mi>i</mi>
<mo>=</mo>
<mn>1</mn>
</mrow>
<mrow>
<mi>N</mi>
<mrow>
<mo>(</mo>
<mi>h</mi>
<mo>)</mo>
</mrow>
</mrow>
</msubsup>
<msup>
<mrow>
<mo>&lsqb;</mo>
<mi>X</mi>
<mrow>
<mo>(</mo>
<msub>
<mi>p</mi>
<mi>i</mi>
</msub>
<mo>)</mo>
</mrow>
<mo>-</mo>
<mi>X</mi>
<mrow>
<mo>(</mo>
<msub>
<mi>p</mi>
<mi>i</mi>
</msub>
<mo>+</mo>
<mi>h</mi>
<mo>)</mo>
</mrow>
<mo>&rsqb;</mo>
</mrow>
<mn>2</mn>
</msup>
<mo>,</mo>
</mrow>
其中，N(h)为空间中距离为h的观测变量的对数；pi为观测点在空间中的位置；X(pi)为待测点pi周围各观测点的属性值；
计算出实验变差函数的值后使用变差函数理论模型对实验变差函数值进行拟合，得到一条最优的实验变差函数。
6.根据权利要求2所述的土壤属性插值的估算方法，其特征在于：所述S102步骤中利用每个待测点的属性估值得到每个待测点的先验概率分布函数的公式是：
<mrow>
<msub>
<mi>f</mi>
<mi>G</mi>
</msub>
<mrow>
<mo>(</mo>
<msub>
<mi>&chi;</mi>
<mrow>
<mi>i</mi>
<mi>j</mi>
</mrow>
</msub>
<mo>)</mo>
</mrow>
<mo>=</mo>
<mfrac>
<mn>1</mn>
<mrow>
<msqrt>
<mrow>
<mn>2</mn>
<mi>&pi;</mi>
</mrow>
</msqrt>
<mi>&sigma;</mi>
</mrow>
</mfrac>
<msup>
<mi>e</mi>
<mrow>
<mo>-</mo>
<mfrac>
<msup>
<mrow>
<mo>(</mo>
<msub>
<mi>&chi;</mi>
<mrow>
<mi>i</mi>
<mi>j</mi>
</mrow>
</msub>
<mo>-</mo>
<mi>&mu;</mi>
<mo>)</mo>
</mrow>
<mn>2</mn>
</msup>
<mrow>
<mn>2</mn>
<msup>
<mi>&sigma;</mi>
<mn>2</mn>
</msup>
</mrow>
</mfrac>
</mrow>
</msup>
<mo>,</mo>
</mrow>
其中，G表示由克里金法得到的概率分布；μ为克里金的预测结果；σ2为预测方差。
7.根据权利要求6所述的土壤属性插值的估算方法，其特征在于：所述S2步骤中利用环境因素与预测点的相关关系生成软数据的具体方法是：
S201：将实验区域划分为N个属性值的取值区间，计算第k个属性取值区间对该点上的单个环境因素的隶属度Mk,si,j，计算公式为其中，Ck为属性值属于第k个区间的采样点数；C1为属性值属于第k个区间且该点的第j个环境因子属于第1个环境因子区间的采样点数；
S202：计算环境因子对待测点属性值的决定系数Rj 2，其中是由环境因子j与对应的函数关系确定的属性值估计值；是采样点属性值的期望；χi为采样点属性值的实际值并将决定系数Rj 2进行归一化，
S203：在研究区域随机生成分布均匀的一定数量的软数据点，对于某一个软数据点psi，
由S201步骤中确定的隶属度函数确定某个属性值取值区间对该软数据点psi上的单个环境因子的隶属度Mk,si,j，并结合环境因子对待测点属性值的决定系数，确定该软数据点psi的取值范围对某个属性值取值区间的相似度从而计算待测点属性值在软数据点上的分布函数fs(psi)＝((S1,si,N1),…,(Sk,si,Nk),…,(SN,si,NN))，其中fs(psi)表示由环境数据确定属性值概率密度分布；Nk表示属性值第k个取值区间。
8.根据权利要求7所述的土壤属性插值的估算方法，其特征在于：所述S3步骤中，利用BME法更新土壤属性待测变量的先验概率分布采用下述公式：
<mrow>
<msub>
<mi>f</mi>
<mi>K</mi>
</msub>
<mrow>
<mo>(</mo>
<msub>
<mi>&chi;</mi>
<mn>0</mn>
</msub>
<mo>)</mo>
</mrow>
<mo>=</mo>
<msub>
<mi>f</mi>
<mi>G</mi>
</msub>
<mrow>
<mo>(</mo>
<msub>
<mi>&chi;</mi>
<mn>0</mn>
</msub>
<mo>|</mo>
<msub>
<mi>&chi;</mi>
<mrow>
<mi>h</mi>
<mi>a</mi>
<mi>r</mi>
<mi>d</mi>
</mrow>
</msub>
<mo>,</mo>
<msub>
<mi>&chi;</mi>
<mrow>
<mi>s</mi>
<mi>o</mi>
<mi>f</mi>
<mi>t</mi>
</mrow>
</msub>
<mo>)</mo>
</mrow>
<mo>=</mo>
<mfrac>
<mrow>
<mo>&Integral;</mo>
<msub>
<mi>f</mi>
<mi>G</mi>
</msub>
<mrow>
<mo>(</mo>
<msub>
<mi>&chi;</mi>
<mn>0</mn>
</msub>
<mo>,</mo>
<msub>
<mi>&chi;</mi>
<mrow>
<mi>h</mi>
<mi>a</mi>
<mi>r</mi>
<mi>d</mi>
</mrow>
</msub>
<mo>,</mo>
<msub>
<mi>&chi;</mi>
<mrow>
<mi>s</mi>
<mi>o</mi>
<mi>f</mi>
<mi>t</mi>
</mrow>
</msub>
<mo>)</mo>
</mrow>
<msub>
<mi>f</mi>
<mi>s</mi>
</msub>
<mrow>
<mo>(</mo>
<msub>
<mi>&chi;</mi>
<mrow>
<mi>s</mi>
<mi>o</mi>
<mi>f</mi>
<mi>t</mi>
</mrow>
</msub>
<mo>)</mo>
</mrow>
<msub>
<mi>d&chi;</mi>
<mrow>
<mi>s</mi>
<mi>o</mi>
<mi>f</mi>
<mi>t</mi>
</mrow>
</msub>
</mrow>
<mrow>
<mo>&Integral;</mo>
<msub>
<mi>f</mi>
<mi>G</mi>
</msub>
<mrow>
<mo>(</mo>
<msub>
<mi>&chi;</mi>
<mrow>
<mi>h</mi>
<mi>a</mi>
<mi>r</mi>
<mi>d</mi>
</mrow>
</msub>
<mo>,</mo>
<msub>
<mi>&chi;</mi>
<mrow>
<mi>s</mi>
<mi>o</mi>
<mi>f</mi>
<mi>t</mi>
</mrow>
</msub>
<mo>)</mo>
</mrow>
<msub>
<mi>f</mi>
<mi>s</mi>
</msub>
<mrow>
<mo>(</mo>
<msub>
<mi>&chi;</mi>
<mrow>
<mi>s</mi>
<mi>o</mi>
<mi>f</mi>
<mi>t</mi>
</mrow>
</msub>
<mo>)</mo>
</mrow>
<msub>
<mi>d&chi;</mi>
<mrow>
<mi>s</mi>
<mi>o</mi>
<mi>f</mi>
<mi>t</mi>
</mrow>
</msub>
</mrow>
</mfrac>
<mo>,</mo>
</mrow>
其中，fG为先验概率密度函数；fs为软数据概率密度函数。
9.根据权利要求8所述的土壤属性插值的估算方法，其特征在于：所述S4步骤中利用土壤属性待测变量的先验概率分布进行土壤空间属性制图的具体方法是：
根据S3步骤中获得的土壤属性待测变量的先验概率分布函数，确定各待测点的预测属性值，预测属性值由下述公式确定：
<mrow>
<mfrac>
<mrow>
<msub>
<mi>df</mi>
<mi>K</mi>
</msub>
<mrow>
<mo>(</mo>
<msub>
<mi>&chi;</mi>
<mn>0</mn>
</msub>
<mo>)</mo>
</mrow>
</mrow>
<mrow>
<msub>
<mi>d&chi;</mi>
<mn>0</mn>
</msub>
</mrow>
</mfrac>
<msub>
<mi>&chi;</mi>
<mn>0</mn>
</msub>
<mo>=</mo>
<msub>
<mover>
<mi>&chi;</mi>
<mo>^</mo>
</mover>
<mn>0</mn>
</msub>
<mo>=</mo>
<mn>0</mn>
<mo>,</mo>
</mrow>
其中，即为预测值，它满足当时，函数fK取得极大值，也就说明在为值p0处，随机变量χ0最有可能的取值是
按照上述公式计算出各待测点的预测属性值后，得到属性分布图。
10.根据权利要求7所述的土壤属性插值的估算方法，其特征在于：所述S2步骤中，对于描述型环境因子数据，对每一个描述型环境因子赋予唯一的数字编号，根据编号的数量将采样点的属性值划分为与编号数量相等的取值区间，计算每个取值区间对各描述型环境因子的隶属度。
Priority Applications (1)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| CN201710561915.3A CN107423561A (zh) | 2017-07-11 | 2017-07-11 | 一种土壤属性插值的估算方法 |
Applications Claiming Priority (1)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| CN201710561915.3A CN107423561A (zh) | 2017-07-11 | 2017-07-11 | 一种土壤属性插值的估算方法 |
Publications (1)
| Publication Number | Publication Date |
|---|---|
| CN107423561A true CN107423561A (zh) | 2017-12-01 |
Family
ID=60427046
Family Applications (1)
| Application Number | Title | Priority Date | Filing Date |
|---|---|---|---|
| CN201710561915.3A Pending CN107423561A (zh) | 2017-07-11 | 2017-07-11 | 一种土壤属性插值的估算方法 |
Country Status (1)
| Country | Link |
|---|---|
| CN (1) | CN107423561A (zh) |
Cited By (9)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN108268977A (zh) * | 2018-01-12 | 2018-07-10 | 武汉智博创享科技股份有限公司 | 降低土壤修复中污染物异常特高值影响的方法及系统 |
| CN109002604A (zh) * | 2018-07-12 | 2018-12-14 | 山东省农业科学院科技信息研究所 | 一种基于贝叶斯最大熵的土壤含水量预测方法 |
| CN109540089A (zh) * | 2018-10-16 | 2019-03-29 | 华南理工大学 | 一种基于贝叶斯-克里金模型的桥面高程拟合方法 |
| CN109541172A (zh) * | 2018-10-25 | 2019-03-29 | 北京农业信息技术研究中心 | 土壤属性值的计算方法及装置 |
| CN110084158A (zh) * | 2019-04-15 | 2019-08-02 | 杭州拓深科技有限公司 | 一种基于智能算法的用电设备识别方法 |
| CN110222911A (zh) * | 2019-06-20 | 2019-09-10 | 中国水利水电科学研究院 | 一种卫星遥感与地面数据协同的雨量站网优化布局方法 |
| CN111693006A (zh) * | 2020-06-12 | 2020-09-22 | 中国科学院地理科学与资源研究所 | 珊瑚沙土壤监测区域的传感器数量和位置确定方法及装置 |
| CN112686497A (zh) * | 2020-12-05 | 2021-04-20 | 辽宁大学 | 基于克里金插值的场地土壤环境损害鉴定实物量化方法 |
| CN115099682A (zh) * | 2022-07-18 | 2022-09-23 | 同济大学 | 一种盾构隧道掌子面软硬分类及开挖风险分级方法 |
Citations (2)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| US20100174514A1 (en) * | 2009-01-07 | 2010-07-08 | Aman Melkumyan | Method and system of data modelling |
| CN106290796A (zh) * | 2016-09-13 | 2017-01-04 | 中国水利水电科学研究院 | 一种基于多重分形的土壤盐分空间变异性表征方法 |
-
2017
- 2017-07-11 CN CN201710561915.3A patent/CN107423561A/zh active Pending
Patent Citations (2)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| US20100174514A1 (en) * | 2009-01-07 | 2010-07-08 | Aman Melkumyan | Method and system of data modelling |
| CN106290796A (zh) * | 2016-09-13 | 2017-01-04 | 中国水利水电科学研究院 | 一种基于多重分形的土壤盐分空间变异性表征方法 |
Non-Patent Citations (4)
| Title |
|---|
| 周健: "基于贝叶斯最大熵的土壤呼吸空间分布研究", 《中国优秀硕士学位论文全文数据库农业科技辑》 * |
| 张楚天: "贝叶斯最大熵方法时空预测关键问题研究与应用", 《中国博士学位论文全文数据库基础科学辑》 * |
| 杨勇等: "基于贝叶斯最大熵的多因子空间属性预测新方法", 《浙江大学学报(农业与生命科学版)》 * |
| 马堃: "基于贝叶斯—克里金的土壤碳通量空间分布研究", 《中国优秀硕士学位论文全文数据库农业科技辑》 * |
Cited By (14)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN108268977B (zh) * | 2018-01-12 | 2020-08-18 | 武汉智博创享科技股份有限公司 | 降低土壤修复中污染物异常特高值影响的方法及系统 |
| CN108268977A (zh) * | 2018-01-12 | 2018-07-10 | 武汉智博创享科技股份有限公司 | 降低土壤修复中污染物异常特高值影响的方法及系统 |
| CN109002604A (zh) * | 2018-07-12 | 2018-12-14 | 山东省农业科学院科技信息研究所 | 一种基于贝叶斯最大熵的土壤含水量预测方法 |
| CN109540089A (zh) * | 2018-10-16 | 2019-03-29 | 华南理工大学 | 一种基于贝叶斯-克里金模型的桥面高程拟合方法 |
| CN109540089B (zh) * | 2018-10-16 | 2021-05-14 | 华南理工大学 | 一种基于贝叶斯-克里金模型的桥面高程拟合方法 |
| CN109541172A (zh) * | 2018-10-25 | 2019-03-29 | 北京农业信息技术研究中心 | 土壤属性值的计算方法及装置 |
| CN110084158B (zh) * | 2019-04-15 | 2021-01-29 | 杭州拓深科技有限公司 | 一种基于智能算法的用电设备识别方法 |
| CN110084158A (zh) * | 2019-04-15 | 2019-08-02 | 杭州拓深科技有限公司 | 一种基于智能算法的用电设备识别方法 |
| CN110222911A (zh) * | 2019-06-20 | 2019-09-10 | 中国水利水电科学研究院 | 一种卫星遥感与地面数据协同的雨量站网优化布局方法 |
| CN110222911B (zh) * | 2019-06-20 | 2021-08-03 | 中国水利水电科学研究院 | 一种卫星遥感与地面数据协同的雨量站网优化布局方法 |
| CN111693006A (zh) * | 2020-06-12 | 2020-09-22 | 中国科学院地理科学与资源研究所 | 珊瑚沙土壤监测区域的传感器数量和位置确定方法及装置 |
| CN111693006B (zh) * | 2020-06-12 | 2021-07-02 | 中国科学院地理科学与资源研究所 | 珊瑚沙土壤监测区域的传感器数量和位置确定方法及装置 |
| CN112686497A (zh) * | 2020-12-05 | 2021-04-20 | 辽宁大学 | 基于克里金插值的场地土壤环境损害鉴定实物量化方法 |
| CN115099682A (zh) * | 2022-07-18 | 2022-09-23 | 同济大学 | 一种盾构隧道掌子面软硬分类及开挖风险分级方法 |
Similar Documents
| Publication | Publication Date | Title |
|---|---|---|
| CN107423561A (zh) | 一种土壤属性插值的估算方法 | |
| CN109142679B (zh) | 基于人工神经网络克里金插值的森林土壤养分的空间预测方法 | |
| CN110263111B (zh) | 基于先验知识的土地利用/覆被信息时空监测方法 | |
| CN108536908B (zh) | 基于非点源氮磷流失风险对流域水环境安全评估的方法 | |
| US10439594B2 (en) | Actually-measured marine environment data assimilation method based on sequence recursive filtering three-dimensional variation | |
| CN105699624B (zh) | 一种基于土壤发生层厚度预测的土壤有机碳储量估算方法 | |
| CN110222911B (zh) | 一种卫星遥感与地面数据协同的雨量站网优化布局方法 | |
| CN113902580B (zh) | 一种基于随机森林模型的历史耕地分布重建方法 | |
| Xia et al. | Model-based analysis of spatio-temporal changes in land use in Northeast China | |
| CN107194054B (zh) | 一种顾及socd空间异质性的克里格插值模型构建方法 | |
| CN102867115A (zh) | 一种基于模糊c均值聚类的农田划分方法 | |
| CN107919983B (zh) | 一种基于数据挖掘的天基信息网络效能评估系统及方法 | |
| CN115618720A (zh) | 一种基于海拔高度的土壤盐渍化预测方法及系统 | |
| CN114398951A (zh) | 一种基于随机森林和众源地理信息的土地利用变化驱动因子挖掘方法 | |
| CN111209968A (zh) | 基于深度学习的多气象因子模式预报温度订正方法及系统 | |
| CN116383589A (zh) | 一种土壤养分的空间插值预测方法 | |
| Yang et al. | Cross-system legacy data applied to digital soil mapping: A case study of Second National Soil Survey data in China | |
| CN106227965B (zh) | 一种顾及时空分布非平稳特征的土壤有机碳空间抽样网络设计方法 | |
| CN113743659B (zh) | 一种基于成分法和马尔可夫元胞自动机的城市布局预测方法及应用 | |
| CN116050460B (zh) | 基于注意力神经网络的气温数据空间插值方法 | |
| Dokoohaki et al. | A novel model–data fusion approach to terrestrial carbon cycle reanalysis across the contiguous US using SIPNET and PEcAn state data assimilation system v. 1.7. 2 | |
| CN102830430B (zh) | 一种层位速度建模方法 | |
| Zhao et al. | Wind resource assessment based on numerical simulations and an optimized ensemble system | |
| CN114239274A (zh) | 多源遥感数据驱动指数滤波模型推求根层土壤水的方法 | |
| CN103955953B (zh) | 一种面向数字土壤制图的地形协同变量选取方法 |
Legal Events
| Date | Code | Title | Description |
|---|---|---|---|
| PB01 | Publication | ||
| PB01 | Publication | ||
| SE01 | Entry into force of request for substantive examination | ||
| RJ01 | Rejection of invention patent application after publication |
Application publication date: 20171201 |
|
| RJ01 | Rejection of invention patent application after publication |