CN114676636A - 一种综合植被和生境特征的草原区土壤水分快速反演方法 - Google Patents

一种综合植被和生境特征的草原区土壤水分快速反演方法 Download PDF

Info

Publication number
CN114676636A
CN114676636A CN202210338884.6A CN202210338884A CN114676636A CN 114676636 A CN114676636 A CN 114676636A CN 202210338884 A CN202210338884 A CN 202210338884A CN 114676636 A CN114676636 A CN 114676636A
Authority
CN
China
Prior art keywords
data
index
soil moisture
vegetation
reflectance
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.)
Granted
Application number
CN202210338884.6A
Other languages
English (en)
Other versions
CN114676636B (zh
Inventor
魏丹丹
冯春涌
黄熙枝
肖晨超
刘凯
Current Assignee (The listed assignees may be inaccurate. Google has not performed a legal analysis and makes no representation or warranty as to the accuracy of the list.)
Ministry Of Natural Resources Land Satellite Remote Sensing Application Center
Original Assignee
Ministry Of Natural Resources Land Satellite Remote Sensing Application Center
Priority date (The priority date 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 date listed.)
Filing date
Publication date
Application filed by Ministry Of Natural Resources Land Satellite Remote Sensing Application Center filed Critical Ministry Of Natural Resources Land Satellite Remote Sensing Application Center
Priority to CN202210338884.6A priority Critical patent/CN114676636B/zh
Publication of CN114676636A publication Critical patent/CN114676636A/zh
Application granted granted Critical
Publication of CN114676636B publication Critical patent/CN114676636B/zh
Active legal-status Critical Current
Anticipated expiration legal-status Critical

Links

Images

Classifications

    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06F—ELECTRIC DIGITAL DATA PROCESSING
    • G06F30/00—Computer-aided design [CAD]
    • G06F30/20—Design optimisation, verification or simulation
    • G06F30/27—Design optimisation, verification or simulation using machine learning, e.g. artificial intelligence, neural networks, support vector machines [SVM] or training a model
    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06F—ELECTRIC DIGITAL DATA PROCESSING
    • G06F17/00—Digital computing or data processing equipment or methods, specially adapted for specific functions
    • G06F17/10—Complex mathematical operations
    • G06F17/18—Complex mathematical operations for evaluating statistical data, e.g. average values, frequency distributions, probability functions, regression analysis
    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06N—COMPUTING ARRANGEMENTS BASED ON SPECIFIC COMPUTATIONAL MODELS
    • G06N3/00—Computing arrangements based on biological models
    • G06N3/02—Neural networks
    • G06N3/04—Architecture, e.g. interconnection topology
    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06N—COMPUTING ARRANGEMENTS BASED ON SPECIFIC COMPUTATIONAL MODELS
    • G06N3/00—Computing arrangements based on biological models
    • G06N3/02—Neural networks
    • G06N3/08—Learning methods
    • G06N3/084—Backpropagation, e.g. using gradient descent
    • Y—GENERAL TAGGING OF NEW TECHNOLOGICAL DEVELOPMENTS; GENERAL TAGGING OF CROSS-SECTIONAL TECHNOLOGIES SPANNING OVER SEVERAL SECTIONS OF THE IPC; TECHNICAL SUBJECTS COVERED BY FORMER USPC CROSS-REFERENCE ART COLLECTIONS [XRACs] AND DIGESTS
    • Y02—TECHNOLOGIES OR APPLICATIONS FOR MITIGATION OR ADAPTATION AGAINST CLIMATE CHANGE
    • Y02A—TECHNOLOGIES FOR ADAPTATION TO CLIMATE CHANGE
    • Y02A90/00—Technologies having an indirect contribution to adaptation to climate change
    • Y02A90/10—Information and communication technologies [ICT] supporting adaptation to climate change, e.g. for weather forecasting or climate simulation

Landscapes

  • Engineering & Computer Science (AREA)
  • Physics & Mathematics (AREA)
  • Theoretical Computer Science (AREA)
  • General Physics & Mathematics (AREA)
  • Data Mining & Analysis (AREA)
  • Evolutionary Computation (AREA)
  • Mathematical Physics (AREA)
  • Software Systems (AREA)
  • General Engineering & Computer Science (AREA)
  • Artificial Intelligence (AREA)
  • Life Sciences & Earth Sciences (AREA)
  • Computational Mathematics (AREA)
  • Mathematical Analysis (AREA)
  • Pure & Applied Mathematics (AREA)
  • Biomedical Technology (AREA)
  • Biophysics (AREA)
  • Computational Linguistics (AREA)
  • Health & Medical Sciences (AREA)
  • General Health & Medical Sciences (AREA)
  • Molecular Biology (AREA)
  • Computing Systems (AREA)
  • Mathematical Optimization (AREA)
  • Computer Hardware Design (AREA)
  • Bioinformatics & Computational Biology (AREA)
  • Computer Vision & Pattern Recognition (AREA)
  • Evolutionary Biology (AREA)
  • Bioinformatics & Cheminformatics (AREA)
  • Medical Informatics (AREA)
  • Operations Research (AREA)
  • Probability & Statistics with Applications (AREA)
  • Geometry (AREA)
  • Algebra (AREA)
  • Databases & Information Systems (AREA)
  • Image Processing (AREA)
  • Investigating Or Analysing Materials By Optical Means (AREA)

Abstract

本发明提出一种综合植被和生境特征的草原区土壤水分快速反演方法，包括：1、研究区野外数据采集；2、光谱数据预处理、SG滤波；3、利用滤波后的反射率数据计算指示植被生长、冠层水分含量、水分胁迫等的11种指数；4、利用CARS算法进行特征波段选择；5、利用相关性分析的方法，分析植被生物理化参量、冠层含水量指标以及生境特征与土壤水分之间的相关性，筛选出相关性好的几种指数；6、构建数据集，确定BP神经网络结构，构建BP神经网络模型，实现表层土壤水分含量反演。本发明能够挖掘出高光谱数据中表层土壤水分相关的特征，结合相关参量，能够高效地进行干旱半干旱草原区表层土壤水分含量的反演。

Description

一种综合植被和生境特征的草原区土壤水分快速反演方法
技术领域
本发明涉及表层土壤水分遥感应用技术领域，具体的说是一种基于植被冠层参量和生境特征的干旱半干旱草原区表层土壤水分含量快速反演的方法。
背景技术
土壤水分通常指存储在非饱和土壤中的水量，又称土壤含水量或土壤湿度，一般采用土壤重量含水量和土壤体积含水量表示。土壤水分是地表和大气能量交换以及水分交换过程中的重要环节，可以调节空气的温度和湿度。土壤水分是影响地表降雨与径流的重要参量，也是理解地气间以水和能量为媒介的交互反馈作用的关键参量，对该参量的监测对于生态系统健康状况、农业发展等具有重要意义。
土壤水分监测传统的方式一般是通过便于携带且易于安装的仪器进行土壤水分数据的采集，但是此方法在获取土壤水分数据的过程中需要投入大量的人力和物力，同时测量的数据代表范围小，难以满足大范围土壤水分的持续监测。遥感技术具有信息量大、快速易行、无污染、无破坏等优点，尤其是高光谱遥感技术能快速获取土壤的连续光谱信息，探测土壤水分变化的细微差异，为土壤水分的定量反演提供了可能。在分析过程中，通过分析不同土壤类型区供试样方中的光谱形状特征，利用相关系数比较，发现1400nm及1900nm附近的吸收峰是随着土壤水分变化显著的波段。由此说明，通过挖掘高光谱曲线的特征可以找出其与土壤水分之间的规律并实现土壤水分的反演。
发明内容
本发明针对传统土壤水分测量方法测定周期长、成本高、不能实现大范围快速监测等问题，提出一种基于植被冠层参量和生境特征的干旱半干旱草原区表层土壤水分含量快速反演的方法，该方法可进一步应用到卫星影像，为实现区域快速重复监测土壤水分含量，提供理论和技术支持。
这种基于植被冠层参量和生境特征的干旱半干旱草原区表层土壤水分含量快速反演的方法，包括以下步骤：
步骤S1、对研究区进行野外实地调查，采集土壤水分数据、植被冠层反射光谱数据、叶面积指数数据、叶片叶绿素含量数据、土壤温度数据；
步骤S2、对野外实测冠层参量数据以及原始反射率数据进行预处理，随后对反射率数据进行Savitzky-Golay滤波以消除光谱噪声；
步骤S3、利用滤波后的反射率数据进行相关指数(归一化植被指数NDVI、增强型植被指数EVI、比值植被指数SR、红边归一化植被指数NDRE、水分胁迫指数MSI、归一化水分指数NDWI、简单比值水分指数SRWI)的计算；
步骤S4、采用竞争性自适应重加权算法(Competitive Adaptive ReweightedSampling，CARS)对经过滤波后的反射率数据进行波段选择；
步骤S5、利用相关性分析的方法，分析植被生物理化参量(叶面积指数、冠层叶绿素含量、NDVI、EVI、SR、NDRE)、冠层含水量指标(MSI、NDWI、SRWI)以及生境特征(土壤温度)与土壤水分之间的相关性，筛选出相关性好的几种参量；
步骤S6、将筛选出的数据(含特征波段和优选指数)进行训练集和验证集的划分；确定BP神经网络结构，构建BP神经网络模型，从而实现土壤水分的反演。
其中，本发明的技术方案还进一步包括：
执行步骤S1时，主要通过野外实地测量得到研究区的土壤水分数据、植被冠层反射光谱数据、叶面积指数数据、叶片叶绿素含量数据、土壤温度数据，具体步骤如下；
步骤S1.1、土壤水分数据采用TDR-350进行测量，该仪器可以快速简便地测量表层土壤水分，通过新功能的改进，能够对土壤电导率进行测量，修正土壤水分读数。测量样方土壤水分的过程中，同时采集样方四角以及中心点位置的0-10cm深度土壤水分数据及土壤温度数据，随后分别将采集的样方四角以及中心点位置的土壤水分数据以及土壤温度数据进行算数平均作为该样方土壤水分数据以及土壤温度数据并加以记录，其中土壤水分数据单位为％VWC，土壤温度数据单位为℃；
步骤S1.2、植被冠层反射光谱数据采用SVC光谱仪进行采集，其光谱范围为350nm-2500nm。在实测过程中，选择在天气晴朗、无风或风速很小的情况下进行。在每次观测前先进行白板校正，测量过程中，探头对准样方内植被冠层，并保持1m的距离，每个样本点记录5条光谱数据并取平均值；
步骤S1.3、叶面积指数数据采用LAI-2200进行采集，在各样方点观测过程中，一般每个采集点测量一次天空光和4次植被下数据，尽量选择天空光变化较小的时段进行。同时在观测LAI时要变化角度，以确保观测值的精确性，最终对每个像元内的观测值取平均值；
步骤S1.4、叶片叶绿素含量数据使用SPAD-502进行采集，其测量面积为2mm×3mm(厚度不超过1.2mm)，在测量前，先对仪器进行校准，校准后开始对样方内叶片的叶绿素含量进行测量。在测量过程中，尽量避免选择过厚的样品，如果测定有较多叶脉的叶片，可进行多次测量并取其平均值，此时记录的叶片叶绿素含量为SPAD值。
执行步骤S2时，具体步骤如下；
步骤S2.1、将采集的植被冠层反射光谱数据“*.sig”文件转换为操作性强的“*.xls”文件，每个xls文件中包含波段信息及其对应的反射率值；
步骤S2.2、对同一样方内的多条光谱数据进行平均计算，用平均后的光谱数值代表该样方的光谱属性。
步骤S2.3、采用Savitzky-Golay滤波的方法对光谱数据进行多点平滑处理，在较好地保持原始反射率光谱信息的同时，有效地降低和消除数据中的噪声。
步骤S2.4、通过换算，将叶片叶绿素含量换算为冠层叶绿素含量。其计算公式为：
CCC＝LCC×LAI (1)
其中，CCC表示冠层叶绿素含量，单位为μg/cm2；LCC表示叶片叶绿素含量，单位为μg/cm2；LAI表示叶面积指数，单位为m2/m2；
执行步骤S3时，主要包括几种指数的计算。在指数计算过程前，针对每一种指数对波段进行组合以计算该指数对应的所有波段组合，选出与土壤水分相关性最好的波段组合作为后续指数计算的依据。
指数计算的具体过程为：
步骤S3.1、对归一化植被指数的计算，其计算公式为：
Figure BDA0003577979140000031
其中，NDVI表示归一化植被指数；ρ741表示在波长741nm处的反射率；ρ676表示在波长676nm处的反射率；
步骤S3.2、对归一化水分指数的计算，其计算公式为：
Figure BDA0003577979140000032
Figure BDA0003577979140000033
Figure BDA0003577979140000034
其中，NDWI表示归一化水分指数；ρ800表示在波长800nm处的反射率；ρ1137表示在波长1137nm处的反射率；ρ500表示在波长500nm处的反射率；ρ741表示在波长741nm处的反射率；ρ1000表示在波长1000nm处的反射率；
步骤S3.3、对简单比值水分指数的计算，其计算公式为：
Figure BDA0003577979140000041
Figure BDA0003577979140000042
Figure BDA0003577979140000043
其中，SRWI表示简单比值水分指数；ρ800表示在波长800nm处的反射率；ρ1137表示在波长1137nm处的反射率；ρ500表示在波长500nm处的反射率；ρ741表示在波长741nm处的反射率；ρ1000表示在波长1000nm处的反射率；
步骤S3.4、对增强型植被指数的计算，其计算公式为：
Figure BDA0003577979140000044
其中，EVI表示增强型植被指数；ρ741表示在波长741nm出的反射率；ρ676表示在波长676nm处的反射率；ρ408表示在波长408nm处的反射率；
步骤S3.5、对比值植被指数SR的计算，其计算公式为：
Figure BDA0003577979140000045
其中，SR表示比值植被指数；ρ741表示在波长741nm处的反射率；ρ721表示在波长721nm处的反射率；
步骤S3.6、对红边归一化植被指数的计算，其计算公式为：
Figure BDA0003577979140000046
其中，NDRE表示红边归一化植被指数；ρ741表示在波长741nm处的反射率；ρ690表示在波长690nm处的反射率；
步骤S3.7、对水分胁迫指数的计算，其计算公式为：
Figure BDA0003577979140000047
其中，MSI表示水分胁迫指数；ρ1464表示在波长1464nm处的反射率；ρ800表示在波长800nm处的反射率；
执行步骤S4时，主要通过CARS算法进行波段的选择，具体过程如下；
步骤S4.1、利用MATLAB软件，将测得的m×1的土壤水分数据以及m×n的光谱数据作为输入数据，进行CARS计算；
CARS算法是一种结合蒙特卡洛与PLS模型回归系数的特征变量选择方法。CARS算法中，每次通过自适应加权采样保留PLS模型中回归系数绝对值权重较大的点作为新的子集，去掉权值较小的点，然后基于新的子集建立PLS模型，经过多次运算，选择PLS模型交互验证均方根误差(RMSECV)最小的子集中的波长作为特征波长。
步骤S4.2、将CARS算法筛选出的波段进行记录并保存。执行步骤S5时，其具体步骤如下；
步骤S5.1、将测得的土壤水分数据、叶面积指数数据、土壤温度数据，以及计算得到的冠层叶绿素含量数据、NDVI、NDWI、SRWI、EVI、SR、NDRE、MSI数据以样方为单位进行整合；
步骤S5.2、对土壤水分和植被生物理化参量(叶面积指数、冠层叶绿素含量、NDVI、EVI、SR、NDRE、冠层含水量指标(MSI、NDWI、SRWI)以及土壤温度进行相关性分析。
步骤S5.3、筛选出相关系数绝对值较高的指数数据进行记录并保存。
Figure BDA0003577979140000051
其中，σx，y表示变量X，Y的相关系数；μx和μy分别表示变量X，Y的均值(期望)，σx、σy分别表示变量X，Y的标准差。σx，y取值范围为[-1，1]，绝对值越大，相关性越强。执行步骤S6时，主要进行数据集的构建并划分，随后构建BP神经网络结构并进行土壤水分的反演，具体步骤如下：
步骤S6.1、构建含有上述筛选出的光谱波段数据以及各种指数数据的数据集，并按照7：3的比例分别划分测试集和验证集。
步骤S6.2、确定BP神经网络结构，以筛选的波段数据和指数数据作为输入值，以土壤水分作为输出值，构建BP神经网络模型。
步骤S6.3、利用构建的BP神经网络模型进行土壤水分的反演，将反演结果同实测值进行RMSE的计算评价精度，RMSE越小，反演效果越好。
Figure BDA0003577979140000052
本发明的优点及有益效果在于：本发明利用CARS算法对Savitzky-Golay滤波后的光谱数据进行特征波段选择，以选出能够反映表层土壤水分的波段数据；同时利用滤波后的光谱数据计算NDVI、NDWI、SRWI、EVI、SR、NDRE、MSI指数，随后采用相关性分析的方法，计算植被生物理化参量(叶面积指数、冠层叶绿素含量、NDVI、EVI、SR、NDRE)、冠层含水量指标(MSI、NDWI、SRWI)以及土壤温度与表层土壤水分的相关系数，选出相关系数绝对值较高的几种指数。最后通过特征波段与优选指数结合表层土壤水分数据构建BP神经网络模型对表层土壤水分进行反演。本发明快速高效的进行干旱半干旱草原区表层土壤水分的反演，能够大幅提高干旱半干旱草原区表层土壤水分的监测效率，降低成本；同时，本发明涉及的植被冠层参量均可通过波段插值等方法在星载高光谱数据中得以实现，结合热红外遥感技术，可以实现大尺度干旱半干旱草原区表层土壤水分的快速反演。
附图说明
图1为本发明的实施流程图。
图2为Savitzky-Golay滤波前后的对比示例。
图3为土壤水分预测情况。
具体实施方式
为使本发明的技术方案、解决的技术问题更加清楚明白，以下结合具体实施例图1，对本发明做进一步的描述。
本发明提供了一种基于植被冠层参量和生境特征的干旱半干旱草原区表层土壤水分含量快速反演的方法，高光谱技术具有光谱分辨率高、无损等特点，能够快速获取土壤的连续光谱信息，探测土壤水分变化的细微差异，为综合反演土壤水分提供了可能。
步骤S1、主要通过野外实地测量得到研究区的采集土壤水分数据、植被冠层反射光谱数据、叶面积指数数据、叶绿素含量数据、土壤温度数据，具体步骤如下；
步骤S1.1、土壤水分数据采用TDR-350进行测量，该仪器可以快速简便地测量土壤水分，通过新功能的改进，能够对土壤电导率进行测量，修正土壤水分读数。测量样方土壤水分的过程中，采集样方四角以及中心点位置的土壤水分数据及土壤温度数据，随后将采集的5个数据进行平均作为该样方土壤水分数据以及土壤温度数据并加以记录，示例如下：
测量得到的样方土壤水分数据为(％VWC)：12.40、9.60、11.80、9.30、9.70；
样方土壤水分算数平均值计算(％VWC)：(12.40+9.60+11.80+9.30+9.70)÷5＝10.56；
测量得到的样方土壤温度数据为(℃)：34.5、33.9、34.2、34.6、34.8；
样方土壤温度算术平均值计算(℃)：(34.5+33.9+34.2+34.6+34.8)÷5＝34.4；
步骤S1.2、利用SVC进行测量前，先对仪器进行预热，将手簿和光谱仪预备好后，手持光谱仪探头面向太阳，保证样方内没有阴影，测量时光谱仪探头距植被冠层1m。每次测量前后用白板进行校正，每个样方测量5次重复数据，每条光谱不低于10次采样。
步骤S1.3、叶面积指数数据采用LAI-2200进行采集，在各样方点观测过程中，一般每个采集点测量一次天空光和4次植被下数据，尽量选择天空光变化较小的时段进行。同时在观测LAI时要变化角度，以确保观测值的精确性，最终对每个像元内的观测值取平均值；
在测量过程中，尽量避免直射阳光，如果天空有明显的明亮和阴暗区域，测量不均匀冠层时易出现错误；天空最亮方向下的冠层会使测量数据增大。叶面积指数计算示例如下：
叶面积指数算数平均计算：(1.34+1.38+1.32+1.30+1.41)÷5＝1.35；
步骤S1.4、叶绿素使用SPAD-502进行采集，其测量面积为2mm*3mm(厚度不超过1.2mm)，在测量前，先对仪器进行校准，校准后开始对样方内叶片的叶绿素含量进行测量。在测量过程中，尽量避免选择过厚的样品，如果测定有较多叶脉的叶片，可进行多次测量并取其平均值，此时记录的数据为叶片叶绿素含量。叶片叶绿素含量算数平均计算(SPAD)：(35.20+33.40+39.30+34.80+39.90)÷5＝36.5；
步骤S2、野外实测地物光谱由于受环境、仪器和目标自身光谱特性等各种因素的影响，获得的光谱曲线中往往会含有噪声，通过光谱数据的预处理以及Savitzky-Golay滤波可以在一定程度上消除噪声，并突出地物光谱的某些细微差别。
步骤S2.1、野外植被冠层反射光谱数据由SVC采集，采集后的植被冠层反射光谱数据为“*.sig”格式，该格式数据可以通过SVC配套的PC端软件进行查看以及编辑，但是其操作性不够强。为方便后续光谱数据的应用，将SVC光谱仪采集的光谱数据“*.sig”文件转换为操作性强的“*.xls”文件，每个xls文件中包含波段信息及其对应的反射率值；
步骤S2.2、高光谱数据的测量过程存在一定的随机误差，在数据处理的过程中利用平均的方法可以在一定程度上削弱随机误差的影响。平均的方法是指在每个样本的不同位置进行多次重复测量，每次重复采集多条光谱，随后对其数据进行平均计算，用其平均后的光谱值代表该目标的光谱属性。
步骤S2.3、采用Savitzky-Golay滤波的方法对光谱数据进行多点平滑处理，也就是对窗口中奇数个点根据最小二乘方法进行拟合，再将求得的中心点的拟合值作为平滑之后该点的值。使用这种平滑方法能够较好地保持原始反射率光谱信息的同时，有效地降低和消除数据中的噪声。
步骤S2.4、通过换算，将叶片叶绿素含量SPAD值转换为单位为μg/cm2，之后通过下面公式换算为冠层叶绿素含量，单位为μg/cm2。其计算公式为：
CCC＝LCC×LAI (1)
其中，CCC表示冠层叶绿素含量，单位为μg/cm2；LCC表示叶片叶绿素含量，单位为μg/cm2；LAI表示叶面积指数，单位为m2/m2；如：
CCC＝LCC×LAI＝31.4×1.04＝32.66
步骤S3、利用滤波后的反射率数据进行相关指数(归一化植被指数NDVI、增强型植被指数EVI、比值植被指数SR、红边归一化植被指数RedEdgNDVI、水分胁迫指数MSI、归一化水分指数NDWI、简单比值水分指数SRWI)的计算，在指数计算过程前，针对每一种指数对波段进行组合以计算该指数对应的所有波段组合，选出与土壤水分相关性最高的波段组合作为后续指数计算的依据。
步骤S3.1、对归一化植被指数的计算，其计算公式为：
Figure BDA0003577979140000081
其中，NDVI表示归一化植被指数；ρ741表示在波长741nm处的反射率；ρ676表示在波长676nm处的反射率；如：
Figure BDA0003577979140000082
步骤S3.2、对归一化水分指数的计算，其计算公式为：
Figure BDA0003577979140000083
Figure BDA0003577979140000084
Figure BDA0003577979140000085
其中，NDWI表示归一化水分指数；ρ800表示在波长800nm处的反射率；ρ1137表示在波长1137nm处的反射率；ρ500表示在波长500nm处的反射率；ρ741表示在波长741nm处的反射率；ρ1000表示在波长1000nm处的反射率；如：
Figure BDA0003577979140000086
Figure BDA0003577979140000087
Figure BDA0003577979140000088
步骤S3.3、对简单比值水分指数的计算，其计算公式为：
Figure BDA0003577979140000089
Figure BDA0003577979140000091
Figure BDA0003577979140000092
其中，SRWI表示简单比值水分指数；ρ800表示在波长800nm处的反射率；ρ1137表示在波长1137nm处的反射率；ρ500表示在波长500nm处的反射率；ρ741表示在波长741nm处的反射率；ρ1000表示在波长1000nm处的反射率；如：
Figure BDA0003577979140000093
Figure BDA0003577979140000094
Figure BDA0003577979140000095
步骤S3.4、对增强型植被指数的计算，其计算公式为：
Figure BDA0003577979140000096
其中，EVI表示增强型植被指数；ρ741表示在波长741nm出的反射率；ρ676表示在波长676nm处的反射率；ρ408表示在波长408nm处的反射率；如：
Figure BDA0003577979140000097
步骤S3.5、对比值植被指数SR的计算，其计算公式为：
Figure BDA0003577979140000098
其中，SR表示比值植被指数；ρ741表示在波长741nm处的反射率；ρ721表示在波长721nm处的反射率；如：
Figure BDA0003577979140000099
步骤S3.6、对红边归一化植被指数的计算，其计算公式为：
Figure BDA00035779791400000910
其中，NDRE表示红边归一化植被指数；ρ741表示在波长741nm处的反射率；ρ690表示在波长690nm处的反射率；如：
Figure BDA00035779791400000911
步骤S3.7、对水分胁迫指数的计算，其计算公式为：
Figure BDA0003577979140000101
其中，MSI表示水分胁迫指数；ρ1464表示在波长1464nm处的反射率；ρ800表示在波长800nm处的反射率；如：
Figure BDA0003577979140000102
步骤S4、采用竞争性自适应重加权算法(Competitive Adaptive ReweightedSampling，CARS)对经过滤波后的反射率数据进行波段选择。
步骤S4.1、利用MATLAB软件，将测得的m×1的土壤水分数据以及m×n的光谱数据作为输入数据，进行CARS计算；
CARS算法是一种结合蒙特卡洛与PLS模型回归系数的特征变量选择方法。CARS算法中，每次通过自适应加权采样保留PLS模型中回归系数绝对值权重较大的点作为新的子集，去掉权值较小的点，然后基于新的子集建立PLS模型，经过多次运算，选择PLS模型交互验证均方根误差(RMSECV)最小的子集中的波长作为特征波长。
首先采用蒙特卡洛采样法，每次随机从校正集中选择一定数量(一般为80％)的样本进入建模集，剩余的20％作为预测集建立PLS模型。蒙特卡洛的采样次数(N)需要提前设定。记录每一次采样过程PLS模型中的回归系数的绝对值权重：
Figure BDA0003577979140000103
其中，|bi|表示第i个变量的回归系数绝对值；m表示每次采样中剩余的变量数；ωi表示第i个变量的回归系数绝对值权重。
随后利用指数衰减函数强行去除回归系数绝对值权重相对较小的波长。在第i次基于MC采样建立PLS模型时，根据EDF得到保留的波长点的比例Ri为：
Ri＝μe-ki (14)
其中，μ和k实常数。
在每次采样时，都从上次采样时的变量数中采用自适应加权采样选择数量为Ri*n个的波长变量，进行PLS建模，计算RMSECV。在N次采样完成之后，CARS算法得到N组候选的特征波长子集，选择RMSECV最小值所对应的波长变量自己为特征波长。
步骤S4.2、将CARS算法筛选出的波段进行记录并保存。
步骤S5、利用相关性分析的方法，分析植被生物理化参量(叶面积指数、冠层叶绿素含量、NDVI、EVI、SR、NDRE)、冠层含水量指标(MSI、NDWI、SRWI)以及土壤温度与土壤水分之间的相关性，筛选出相关性好的几种指数。步骤S5.1、将测得的土壤水分数据、叶面积指数数据、土壤温度数据，以及计算得到的冠层叶绿素含量、NDVI、NDWI、SRWI、EVI、SR、NDRE、MSI数据以样方为单位进行整合，示例如下：
表1参数整理示例
Figure BDA0003577979140000111
步骤S5.2、对土壤水分和植被生物理化参量(叶面积指数、冠层叶绿素含量、NDVI、EVI、SR、NDRE)、冠层含水量指标(MSI、NDWI、SRWI)以及土壤温度进行相关性分析。
Figure BDA0003577979140000112
其中，σx，y表示变量X，Y的相关系数；μx和μy分别表示变量X，Y的均值(期望)，σx、σy分别表示变量X，Y的标准差。σx，y取值范围为[-1，1]，绝对值越大，相关性越强。
步骤S5.3、筛选出相关系数绝对值较大的指数数据进行记录并保存。
步骤S6、将筛选出的数据(含特征波段和优选指数)进行训练集和验证集的划分；确定BP神经网络结构，构建BP神经网络模型，从而实现土壤水分的反演。
步骤S6.1、构建含有上述筛选出的光谱波段数据以及各种指数数据的数据集，并按照7:3的比例分别划分测试集和验证集。
步骤S6.2、确定BP神经网络结构，以筛选的波段数据和指数数据作为输入值，以土壤水分作为输出值，构建BP神经网络模型。
步骤S6.3、利用构建的BP神经网络模型进行土壤水分的反演，将反演结果同实测值进行RMSE的计算评价精度，RMSE越小，反演效果越好。
Figure BDA0003577979140000113
综上可知，采用本发明一种基于植被冠层参量和生境特征的干旱半干旱草原区表层土壤水分含量快速反演的方法，能够充分挖掘高光谱数据中土壤水分相关的特征，结合热红外遥感技术，能够高效地进行干旱半干旱草原区表层土壤水分含量的反演。
本发明着重强调了一种基于植被冠层参量和生境特征的干旱半干旱草原区表层土壤水分含量快速反演的方法，也可以被认为是针对特定场景下一套完整的操作流程，包括野外数据采集、光谱数据预处理、CARS算法波段选择、相关指数计算及优选、构建土壤水分BP神经网络反演模型。基于本发明的上述具体实施例，本技术领域的技术人员在不脱离本发明原理的前提下，对本发明所作出的任何改进和修饰(包括更换应用场景)，皆应纳入本发明的专利保护范围。

Claims (7)

1.一种综合植被和生境特征的草原区土壤水分快速反演方法，其特征在于，包括如下步骤：
步骤S1、对研究区进行野外实地调查，采集土壤水分数据、植被冠层反射光谱数据、叶面积指数数据、叶片叶绿素含量数据及土壤温度数据；
步骤S2、对野外实测冠层参量数据以及原始反射率数据进行预处理，随后对反射率数据进行Savitzky-Golay滤波以消除光谱噪声；
步骤S3、利用滤波后的反射率数据进行相关指数的计算，包括：归一化植被指数NDVI、增强型植被指数EVI、比值植被指数SR、红边归一化植被指数NDRE、水分胁迫指数MSI、归一化水分指数NDWI及简单比值水分指数SRWI；
步骤S4、采用竞争性自适应重加权算法CARS对经过滤波后的反射率数据进行波段选择；
步骤S5、利用相关性分析的方法，分析植被生物理化参量、冠层含水量指标以及生境特征与土壤水分之间的相关性，筛选出相关性好的几种参量；
步骤S6、将筛选出的数据进行训练集和验证集的划分；确定BP神经网络结构，构建BP神经网络模型，从而实现土壤水分的反演。
2.根据权利要求1所述的一种综合植被和生境特征的草原区土壤水分快速反演方法，其特征在于：执行步骤S1时，通过野外实地测量得到研究区的土壤水分数据、植被冠层反射光谱数据、叶面积指数数据、叶片叶绿素含量数据、土壤温度数据，具体步骤如下：
步骤S1.1、土壤水分数据采用TDR-350进行测量；测量样方土壤水分的过程中，同时采集样方四角以及中心点位置的0-10cm深度土壤水分数据及土壤温度数据，随后分别将采集的样方四角以及中心点位置的土壤水分数据以及土壤温度数据进行算数平均作为该样方土壤水分数据以及土壤温度数据并加以记录；其中，土壤水分数据单位为％VWC，土壤温度数据单位为℃；
步骤S1.2、植被冠层反射光谱数据采用SVC光谱仪进行采集，其光谱范围为350nm-2500nm；在实测过程中，选择在天气晴朗、无风或风速小的情况下进行；在每次观测前先进行白板校正，测量过程中，探头对准样方内植被冠层，并保持1m的距离，每个样本点记录5条光谱数据并取平均值；
步骤S1.3、叶面积指数数据采用LAI-2200进行采集，在各样方点观测过程中，每个采集点测量一次天空光和4次植被下数据，选择天空光变化小的时段进行；同时在观测LAI时要变化角度，以确保观测值的精确性，最终对每个像元内的观测值取平均值；
步骤S1.4、叶片叶绿素含量数据使用SPAD-502进行采集，测量面积为2mm×3mm，厚度不超过1.2mm，在测量前，先对仪器进行校准，校准后开始对样方内叶片的叶绿素含量进行测量；在测量过程中，避免选择过厚的样品，如果测定有多叶脉的叶片，需要进行多次测量并取其平均值，此时记录的叶片叶绿素含量为SPAD值。
3.根据权利要求1所述的一种综合植被和生境特征的草原区土壤水分快速反演方法，其特征在于：执行步骤S2时，具体步骤如下：
步骤S2.1、将采集的植被冠层反射光谱数据“*.sig”文件转换为操作性强的“*.xls”文件，每个xls文件中包含波段信息及其对应的反射率值；
步骤S2.2、对同一样方内的多条光谱数据进行平均计算，用平均后的光谱数值代表该样方的光谱属性；
步骤S2.3、采用Savitzky-Golay滤波的方法对光谱数据进行多点平滑处理，能保持原始反射率光谱信息，降低和消除数据中的噪声；
步骤S2.4、通过换算，将叶片叶绿素含量换算为冠层叶绿素含量；计算公式为：
CCC＝LCC×LAI (1)
其中，LCC表示冠层叶绿素含量，单位为μg/cm2；LCC表示叶片叶绿素含量，单位为μg/cm2；LAI表示叶面积指数，单位为m2/m2。
4.根据权利要求1所述的一种综合植被和生境特征的草原区土壤水分快速反演方法，其特征在于：执行步骤S3时，包括几种指数的计算；在指数计算过程前，针对每一种指数对波段进行组合以计算该指数对应的所有波段组合，选出与土壤水分相关性最好的波段组合作为后续指数计算的依据；
指数计算的具体过程为：
步骤S3.1、对归一化植被指数的计算，其计算公式为：
Figure FDA0003577979130000021
其中，NDVI表示归一化植被指数；ρ741表示在波长741nm处的反射率；ρ676表示在波长676nm处的反射率；
步骤S3.2、对归一化水分指数的计算，其计算公式为：
Figure FDA0003577979130000022
Figure FDA0003577979130000023
Figure FDA0003577979130000024
其中，NDWI表示归一化水分指数；ρ800表示在波长800nm处的反射率；ρ1137表示在波长1137nm处的反射率；ρ500表示在波长500nm处的反射率；ρ741表示在波长741nm处的反射率；ρ1000表示在波长1000nm处的反射率；
步骤S3.3、对简单比值水分指数的计算，其计算公式为：
Figure FDA0003577979130000031
Figure FDA0003577979130000032
Figure FDA0003577979130000033
其中，SRWI表示简单比值水分指数；ρ800表示在波长800nm处的反射率；ρ1137表示在波长1137nm处的反射率；ρ500表示在波长500nm处的反射率；ρ741表示在波长741nm处的反射率；ρ1000表示在波长1000nm处的反射率；
步骤S3.4、对增强型植被指数的计算，其计算公式为：
Figure FDA0003577979130000034
其中，EVI表示增强型植被指数；ρ741表示在波长741nm出的反射率；ρ676表示在波长676nm处的反射率；ρ408表示在波长408nm处的反射率；
步骤S3.5、对比值植被指数SR的计算，其计算公式为：
Figure FDA0003577979130000035
其中，SR表示比值植被指数；ρ741表示在波长741nm处的反射率；ρ721表示在波长721nm处的反射率；
步骤S3.6、对红边归一化植被指数的计算，其计算公式为：
Figure FDA0003577979130000036
其中，NDRE表示红边归一化植被指数；ρ741表示在波长741nm处的反射率；ρ690表示在波长690nm处的反射率；
步骤S3.7、对水分胁迫指数的计算，其计算公式为：
Figure FDA0003577979130000037
其中，MSI表示水分胁迫指数；ρ1464表示在波长1464nm处的反射率；ρ800表示在波长800nm处的反射率。
5.根据权利要求1所述的一种综合植被和生境特征的草原区土壤水分快速反演方法，其特征在于：执行步骤S4时，通过CARS算法进行波段的选择，具体过程如下：
步骤S4.1、利用MATLAB软件，将测得的m×1的土壤水分数据以及m×n的光谱数据作为输入数据，进行CARS计算；
在CARS算法中，每次通过自适应加权采样保留PLS模型中回归系数绝对值权重大的点作为新的子集，去掉权值小的点，然后基于新的子集建立PLS模型，经过多次运算，选择PLS模型交互验证均方根误差RMSECV最小的子集中的波长作为特征波长；
步骤S4.2、将CARS算法筛选出的波段进行记录并保存。
6.根据权利要求1所述的一种综合植被和生境特征的草原区土壤水分快速反演方法，其特征在于：执行步骤S5时，具体步骤如下：
步骤S5.1、将测得的土壤水分数据、叶面积指数数据、土壤温度数据，以及计算得到的冠层叶绿素含量数据、NDVI、NDWI、SRWI、EVI、SR、NDRE、MSI数据以样方为单位进行整合；
步骤S5.2、对土壤水分和植被生物理化参量，即叶面积指数、冠层叶绿素含量、NDVI、EVI、SR、NDRE、冠层含水量指标即MSI、NDWI、SRWI以及土壤温度进行相关性分析；
步骤S5.3、筛选出相关系数绝对值较高的指数数据进行记录并保存；
Figure FDA0003577979130000041
其中，σx，y表示变量X，Y的相关系数；μx和μy分别表示变量X，Y的均值，σx、σy分别表示变量X，Y的标准差；σx，y取值范围为[-1，1]，绝对值越大，相关性越强。
7.根据权利要求1所述的一种综合植被和生境特征的草原区土壤水分快速反演方法，其特征在于：执行步骤S6时，进行数据集的构建并划分，随后构建BP神经网络结构并进行土壤水分的反演，具体步骤如下：
步骤S6.1、构建含有上述筛选出的光谱波段数据以及各种指数数据的数据集，并按照7：3的比例分别划分测试集和验证集；
步骤S6.2、确定BP神经网络结构，以筛选的波段数据和指数数据作为输入值，以土壤水分作为输出值，构建BP神经网络模型；
步骤S6.3、利用构建的BP神经网络模型进行土壤水分的反演，将反演结果同实测值进行RMSE的计算评价精度，RMSE越小，反演效果越好；
Figure FDA0003577979130000042
CN202210338884.6A 2022-04-01 2022-04-01 一种综合植被和生境特征的草原区土壤水分快速反演方法 Active CN114676636B (zh)

Priority Applications (1)

Application Number Priority Date Filing Date Title
CN202210338884.6A CN114676636B (zh) 2022-04-01 2022-04-01 一种综合植被和生境特征的草原区土壤水分快速反演方法

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
CN202210338884.6A CN114676636B (zh) 2022-04-01 2022-04-01 一种综合植被和生境特征的草原区土壤水分快速反演方法

Publications (2)

Publication Number Publication Date
CN114676636A true CN114676636A (zh) 2022-06-28
CN114676636B CN114676636B (zh) 2023-08-29

Family

ID=82076498

Family Applications (1)

Application Number Title Priority Date Filing Date
CN202210338884.6A Active CN114676636B (zh) 2022-04-01 2022-04-01 一种综合植被和生境特征的草原区土壤水分快速反演方法

Country Status (1)

Country Link
CN (1) CN114676636B (zh)

Cited By (3)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN115015258A (zh) * 2022-07-21 2022-09-06 中国科学院西北生态环境资源研究院 一种农作物长势与土壤水分关联确定方法及相关装置
CN115758855A (zh) * 2022-08-18 2023-03-07 中国人民解放军61540部队 一种基于lstm及注意力机制的地表反射率反演方法
WO2024098587A1 (zh) * 2022-11-10 2024-05-16 北华航天工业学院 基于ceemdan的农田土壤水分胁迫甄别与监测方法

Citations (4)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN108663330A (zh) * 2018-04-19 2018-10-16 中国国土资源航空物探遥感中心 一种基于叶片实测光谱的植被覆盖区土壤铜元素反演方法
CN110118742A (zh) * 2019-05-27 2019-08-13 成都理工大学 基于光谱分析的陆表植被冠层生态水含量遥感反演方法
US20200225075A1 (en) * 2019-01-14 2020-07-16 Wuhan University Method and system for optical and microwave synergistic retrieval of aboveground biomass
CN113255874A (zh) * 2021-06-08 2021-08-13 中国农业科学院农业资源与农业区划研究所 一种基于优化的bp神经网络微波遥感反演土壤水分的方法

Patent Citations (4)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN108663330A (zh) * 2018-04-19 2018-10-16 中国国土资源航空物探遥感中心 一种基于叶片实测光谱的植被覆盖区土壤铜元素反演方法
US20200225075A1 (en) * 2019-01-14 2020-07-16 Wuhan University Method and system for optical and microwave synergistic retrieval of aboveground biomass
CN110118742A (zh) * 2019-05-27 2019-08-13 成都理工大学 基于光谱分析的陆表植被冠层生态水含量遥感反演方法
CN113255874A (zh) * 2021-06-08 2021-08-13 中国农业科学院农业资源与农业区划研究所 一种基于优化的bp神经网络微波遥感反演土壤水分的方法

Cited By (4)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN115015258A (zh) * 2022-07-21 2022-09-06 中国科学院西北生态环境资源研究院 一种农作物长势与土壤水分关联确定方法及相关装置
CN115758855A (zh) * 2022-08-18 2023-03-07 中国人民解放军61540部队 一种基于lstm及注意力机制的地表反射率反演方法
CN115758855B (zh) * 2022-08-18 2024-04-12 中国人民解放军61540部队 一种基于lstm及注意力机制的地表反射率反演方法
WO2024098587A1 (zh) * 2022-11-10 2024-05-16 北华航天工业学院 基于ceemdan的农田土壤水分胁迫甄别与监测方法

Also Published As

Publication number Publication date
CN114676636B (zh) 2023-08-29

Similar Documents

Publication Publication Date Title
CN113192117B (zh) 一种基于Sentinel-2卫星影像红边波段改进小麦生长早期叶面积指数估算的方法
CN114676636B (zh) 一种综合植被和生境特征的草原区土壤水分快速反演方法
CN107796764B (zh) 一种基于三波段植被指数的小麦叶面积指数估算模型的构建方法
Li et al. Multi-LUTs method for canopy nitrogen density estimation in winter wheat by field and UAV hyperspectral
Atzberger et al. Comparative analysis of three chemometric techniques for the spectroradiometric assessment of canopy chlorophyll content in winter wheat
CN111044516B (zh) 一种水稻叶绿素含量遥感估测方法
CN111488926B (zh) 一种基于优化模型的土壤有机质测定方法
CN110118742A (zh) 基于光谱分析的陆表植被冠层生态水含量遥感反演方法
CN113252583B (zh) 一种基于枯草植被指数计算高寒枯草覆盖度的方法
CN113466143B (zh) 土壤养分反演方法、装置、设备及介质
CN114372707A (zh) 一种基于遥感数据的高寒湿地退化程度监测方法
CN103940748A (zh) 基于高光谱技术的柑橘冠层含氮量预测与可视化的方法
CN113111799B (zh) 一种基于集合经验模态分解的耕地土壤肥力水平监测方法
CN113065230B (zh) 基于优化光谱指数建立水稻叶片spad的高光谱反演模型
CN111896497B (zh) 一种基于预测值的光谱数据修正方法
CN112924401B (zh) 一种植被冠层叶绿素含量半经验反演方法
CN117152636A (zh) 一种基于双波段关系的浅海底质反射率遥感监测方法
CN111595806A (zh) 一种利用中红外漫反射光谱监测土壤碳组分的方法
CN114755189B (zh) 一种特征优选的自注意力机制高光谱卫星lai反演方法
CN114878508A (zh) 一种基于高光谱成像土遗址文物表面含水量检测方法
CN112229817A (zh) 一种苏打盐碱地重金属定量反演模型建立方法
Zhang et al. Estimation of grassland height using optical and SAR remote sensing data
CN117309780B (zh) 土壤锗元素含量的确定方法
CN117269078B (zh) 土壤阳离子交换量的确定方法
CN117630000B (zh) 土壤含盐量评估方法、装置、电子设备及存储介质

Legal Events

Date Code Title Description
PB01 Publication
PB01 Publication
SE01 Entry into force of request for substantive examination
SE01 Entry into force of request for substantive examination
GR01 Patent grant
GR01 Patent grant