CN115420690A - 近地表痕量气体浓度反演模型及反演方法 - Google Patents

近地表痕量气体浓度反演模型及反演方法 Download PDF

Info

Publication number
CN115420690A
CN115420690A CN202210465007.5A CN202210465007A CN115420690A CN 115420690 A CN115420690 A CN 115420690A CN 202210465007 A CN202210465007 A CN 202210465007A CN 115420690 A CN115420690 A CN 115420690A
Authority
CN
China
Prior art keywords
data
model
trace gas
value
monitoring
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
Application number
CN202210465007.5A
Other languages
English (en)
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.)
Shaanxi Province Environmental Monitoring Center Station
Zhongyao Environment Xi'an Co ltd
Original Assignee
Shaanxi Province Environmental Monitoring Center Station
Zhongyao Environment Xi'an Co ltd
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 Shaanxi Province Environmental Monitoring Center Station, Zhongyao Environment Xi'an Co ltd filed Critical Shaanxi Province Environmental Monitoring Center Station
Priority to CN202210465007.5A priority Critical patent/CN115420690A/zh
Publication of CN115420690A publication Critical patent/CN115420690A/zh
Priority to CN202211545701.4A priority patent/CN116223395A/zh
Pending legal-status Critical Current

Links

Images

Classifications

    • GPHYSICS
    • G01MEASURING; TESTING
    • G01NINVESTIGATING OR ANALYSING MATERIALS BY DETERMINING THEIR CHEMICAL OR PHYSICAL PROPERTIES
    • G01N21/00Investigating or analysing materials by the use of optical means, i.e. using sub-millimetre waves, infrared, visible or ultraviolet light
    • G01N21/17Systems in which incident light is modified in accordance with the properties of the material investigated
    • G01N21/25Colour; Spectral properties, i.e. comparison of effect of material on the light at two or more different wavelengths or wavelength bands
    • GPHYSICS
    • G01MEASURING; TESTING
    • G01WMETEOROLOGY
    • G01W1/00Meteorology
    • G01W1/02Instruments for indicating weather conditions by measuring two or more variables, e.g. humidity, pressure, temperature, cloud cover or wind speed
    • GPHYSICS
    • G06COMPUTING; CALCULATING OR COUNTING
    • G06FELECTRIC DIGITAL DATA PROCESSING
    • G06F17/00Digital computing or data processing equipment or methods, specially adapted for specific functions
    • G06F17/10Complex mathematical operations
    • G06F17/15Correlation function computation including computation of convolution operations
    • GPHYSICS
    • G01MEASURING; TESTING
    • G01NINVESTIGATING OR ANALYSING MATERIALS BY DETERMINING THEIR CHEMICAL OR PHYSICAL PROPERTIES
    • G01N21/00Investigating or analysing materials by the use of optical means, i.e. using sub-millimetre waves, infrared, visible or ultraviolet light
    • G01N21/17Systems in which incident light is modified in accordance with the properties of the material investigated
    • G01N2021/1793Remote sensing
    • G01N2021/1795Atmospheric mapping of gases
    • YGENERAL 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
    • Y02TECHNOLOGIES OR APPLICATIONS FOR MITIGATION OR ADAPTATION AGAINST CLIMATE CHANGE
    • Y02ATECHNOLOGIES FOR ADAPTATION TO CLIMATE CHANGE
    • Y02A90/00Technologies having an indirect contribution to adaptation to climate change
    • Y02A90/10Information 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)
  • General Physics & Mathematics (AREA)
  • Theoretical Computer Science (AREA)
  • Life Sciences & Earth Sciences (AREA)
  • Mathematical Analysis (AREA)
  • Mathematical Optimization (AREA)
  • Pure & Applied Mathematics (AREA)
  • Environmental & Geological Engineering (AREA)
  • Data Mining & Analysis (AREA)
  • Mathematical Physics (AREA)
  • Computational Mathematics (AREA)
  • Environmental Sciences (AREA)
  • Biochemistry (AREA)
  • Software Systems (AREA)
  • Algebra (AREA)
  • Computing Systems (AREA)
  • Atmospheric Sciences (AREA)
  • Biodiversity & Conservation Biology (AREA)
  • Ecology (AREA)
  • Databases & Information Systems (AREA)
  • Spectroscopy & Molecular Physics (AREA)
  • Chemical & Material Sciences (AREA)
  • Health & Medical Sciences (AREA)
  • Analytical Chemistry (AREA)
  • General Engineering & Computer Science (AREA)
  • General Health & Medical Sciences (AREA)
  • Immunology (AREA)
  • Pathology (AREA)
  • Management, Administration, Business Operations System, And Electronic Commerce (AREA)
  • Geophysics And Detection Of Objects (AREA)

Abstract

一种近地表痕量气体浓度反演模型及基于该模型反演地表痕量气体浓度的方法,反演模型基于遥感数据与站点数据建立反演模型。通过如下步骤建立反演模型,S1,收集地表监测数据、遥感数据、气象数据和其他数据;S2,对步骤S1收集的数据进行预处理,得到预处理后的初始建模数据;S3,对初始建模数据进行特征筛选,得到建模数据集;S4,根据建模数据集搭建模型。反演模型建立后,可通过该模型反演地表痕量气体浓度。本发明结合人工智能与大数据挖掘与技术,建立遥感监测的痕量气体总柱浓度与地面站点监测的痕量气体质量浓度的反演模型,能够准确、高效地反演近地表痕量气体的时空分布特征,支撑大气污染物防治工作的精准实施。

Description

近地表痕量气体浓度反演模型及反演方法
技术领域
本发明涉及气象信息技术领域,特别是涉及一种近地表痕量气体浓度反演模型及通过该模型反演近地 表痕量气体浓度分布的方法。
背景技术
痕量气体包括O3、NO2、SO2、CO。目前痕量气体的监测手段主要为地面监测和遥感监测。其中地面监 测主要有地面监测站点的自动监测、手动取样监测和走航监测,优点是监测的时间频率和精度均较高,缺 点是地面监测站点分布不均匀,监测结果覆盖范围小(只能满足监测站点周围区域),建设和维护成本高, 无法满足区域化大尺度痕量气体浓度实时且准确监控的需求;遥感监测能周期性、重复地从空中乃至宇宙 空间对大范围区域进行对地观测,解决了地面站点覆盖范围小的缺点,但目前的遥感技术针对痕量气体只 能提供不同垂直范围内的总柱浓度、对流层柱浓度或痕量气体剖面产品,针对颗粒物只能提供气溶胶光学 厚度产品,不能提供近地表痕量气体的浓度。
因此,针对现有技术不足,提供一种能准确、高效地反演近地表痕量气体的时空分布特征的近地表痕 量气体浓度反演模型及通过该模型反演近地表痕量气体浓度分布的方法以克服现有技术不足甚为必要。
发明内容
本发明的目的在于避免现有技术的不足之处而提供一种近地表痕量气体浓度反演模型及通过该模型 反演近地表痕量气体浓度分布的方法,能够准确、高效地反演近地表痕量气体的时空分布特征特别是浓 度特征。
本发明的目的通过以下技术措施实现。
提供一种近地表痕量气体浓度反演模型,基于遥感数据与站点数据建立反演模型。结合人工智能与大 数据挖掘与技术,建立遥感监测的痕量气体总柱浓度(以下简称“遥感数据”)与地面站点监测的痕量气 体质量浓度(以下简称“站点数据”)的回归模型(反演模型),准确、高效地反演近地表痕量气体的时空 分布特征,支撑大气污染物防治工作的精准实施。
可选的,上述的近地表痕量气体浓度反演模型,通过如下步骤建立:
S1,收集地表监测数据、遥感数据、气象数据和其他数据;
S2,对步骤S1收集的数据进行预处理,得到预处理后的初始建模数据;
S3,对初始建模数据进行特征筛选,得到建模数据集;
S4,根据建模数据集搭建模型。
3.根据权利要求2所述的近地表痕量气体浓度反演模型,其特征在于:S1中,收集地表监测数据具体 是:
从地方省级环境监测中心站获取地表痕量气体浓度数据,包括O3、NO2、CO、SO2的小时值、8小时均 值和24小时均值,痕量气体使用13:00和14:00的算数平均值,数据中包含了监测站点的经度、纬度和 日期;
收集遥感数据具体是:
痕量气体遥感数据使用的是哨兵5P的L2级总柱浓度产品数据,遥感数据从GoogleEarth Engine(GEE) 平台上下载;
收集气象数据具体是:
气象数据从省级气象局获取,气象条件监测站点密集且分布均匀,平均1km2范围内的站点数大于1个, 气象类型包括平均温度(TEM_Avg)、最大温度(TEM_Max)、最小温度(TEM_Min)、平均相对湿度(RH)、 8:00-8:00累积降雨量(PRE_08)、20:00-20:00累积降雨量(PRE_20)、2分钟平均风速(WIN)、蒸发量(EVP)、 日照时数(SSH),全部为日均值数据;
收集其他数据具体是:
其他用到的数据有数字高程模型(DEM),空间分辨率为30m×30m;地表覆盖类型数据(GLC),空间分 辨率为30m×30m,人口密度数据(PopDen),空间分辨率0.09°×0.09°;
痕量气体具有空间和时间上异质性,将日期数据转化为每年的第几天(DOY),匹配时间异质性,用监 测站点的经纬度(Lon/Lat)匹配空间异质性。
可选的,上述的近地表痕量气体浓度反演模型,对步骤S1收集的数据进行预处理的具体过程是:
S21,数据提取
首先,从空气质量监测的原始文件中提取2018年-2020年目标城市群所有空气质量监测站点的痕量气 体13:00-14:00的监测结果,删除缺失值后,将每天13:00和14:00的监测结果求算术平均值,作为 地面监测结果,得到痕量气体监测数据集。提取过程中同时包含监测站点的经度、纬度和监测日值,根据 监测站点的经度和纬度信息转化为WGS-84坐标系下的点矢量文件,简称“痕量气体点矢量”,用于提取对 应监测站点的气象数据、地理数据、人口密度数据;
其次,将GEE下载的遥感数据重投影至WGS84坐标系,使用最近邻采样方法重采样至1km×1km的分 辨率,使用痕量气体点矢量提取对应点位的遥感数据,得到遥感数据集;
再次,从气象数据原始文件中提取2018-2020年所有监测站点的日均值监测结果,包含气象监测站点 的经纬度和日期,剔除所有缺失值,依次按照日期和气象数据字段信息,按照经纬度和监测结果将各个气 象字段的监测数据转化为WGS-84坐标系下的点矢量文件,然后按照反距离权重插值方法插值到1km×1km 分辨率,进一步矢量转栅格后形成栅格文件,栅格文件的像元值就是气象字段的日均值监测结果,最后使 用痕量气体点矢量提取对每天各字段对应位置的气象结果,得到气象数据集;
最后,DEM、GLC、PopDen均为栅格文件,其中DEM和GLC分辨率为30m×30m,PopDen分辨率为900m×900m, 使用最近邻采样方法重采样至1km×1km,重投影至WGS-84坐标系下,使用痕量气体矢量提取,得到辅助 数据集;
S22,数据结合
数据提取阶段共提取了4个数据集,分别为痕量气体的监测数据集、遥感数据集、气象数据集、辅助 数据集,四份数据集中均包含有经度、纬度和日期,按照这三个字段将四份数据结合成为最终的数据集, 简称“原始数据集”,进一步将原始数据集中的日期转化为每年的第几天,Day Of Year(DOY),原始数据 集中痕量气体的监测结果为标签,其余全部为特征;
S23.异常值剔除
对原始数据集进行异常值检查,检查异常值的方法是绘制箱型图,将含有异常值的样本全部删除;此 外,根据先验知识将气象数据的异常值(如WIN>50m/s,TEM_Avg>40℃,SSH>24h)剔除,即删除极端 天气的影响,同样将含有极端天气的样本删除。
可选的,上述的近地表痕量气体浓度反演模型,S3中对初始建模数据进行特征筛选,具体包括:
S31.绘制散点图研究相关关系
通过绘制散点图的方式,研究数据集中各特征变量与痕量气体之间的相关关系,数据中包含长时间序 列、大范围的气体,散点图可以有效反应痕量气体与其他特征之间的相关关系,根据先验知识将基本无相 关的特征删除;
S32.统计person相关系数和p_value
通过统计特征之间的person相关系数,痕量气体与其他变量之间相关性,将person相关系数为0.00 的特征删除,进一步计算p_value值,将p_value<0.01的变量删除;
S33.机器学习筛选特征
将步骤S31、S32筛选后的数据集带入随机森林模型,使用该模型的内建函数“feature.importance” 计算特征重要性,排序后按照特征重要性得分由小到大排序,依次删除得分较小的特征,查看建模精度变 化,直到模型精度发生较大变化时停止删除特征,得到建模数据集,简称“数据集”。
可选的,上述的近地表痕量气体浓度反演模型,S4搭建模型建模使用的机器学习算法为eXtreme Gradient Boosting(XGBoost),是以决策树为基学习器构建的集成算法,基学习器为Classification And Regression Tree(CART)决策树,使用XGBoost反演近地表痕量气体浓度旨在建立回归模型,模型搭建过 程使用的标签为痕量气体的浓度,特征为数据集中除痕量气体浓度之外的其他特征;
具体过程是:
S41,CART树
对于给定的数据集,选择最优切分特征j(特征信息)与切分点s(特征中所有可能取值;
遍历特征j,对固定的切分特征j扫描切分点s,选择使式1达到最小值的(j,s)对;
Figure BDA0003623562240000021
式1中,R1和R2表示由(j,s)对划分的两个数据集,yi是样本i的标签(锂的检测结果),xi表示数 据集中的(j,s)对,c1是R1上的平均值,c2是R2上的平均值,
用选定的(j,s)对划分区域并决定相应的输出值:
R1(j,s)={x|x(j)≤s},R2(j,s)={x|x(j)>s}
(式2);
Figure BDA0003623562240000031
继续对两个子区域R1和R2调用步骤①,②,直至满足停止条件;
将输入空间(数据集)划分为M个区域(叶子结点)R1,R2,……,Rm,生成决策树:
Figure BDA0003623562240000032
式4中,f(x)表示预测结果,M表示叶子结点的个数,m表示第m个叶子结点,Rm表示叶子结点的集 合,
Figure BDA0003623562240000033
表示叶子结点的平均值(该结点上的预测值),I表示叶子结点的权重(叶子结点的样本个数占总 样本个数的权重);
S42,Boosting
Boosting是一族可将弱学习器提升为强学习器的算法;这族算法的工作机制类似:先从初始训练集训 练出一个基学习器,再根据基学习器的表现对训练样本分布进行调整,使得先前基学习器做错的训练样本 在后续受到更多关注,然后基于调整后的样本分布来训练下一个基学习器;如此重复进行,直至基学习器 数目达到事先指定的值T,最终将这T个基学习器进行加权结合;
S43,GBDT
GBDT专注于回归的树的提升集成模型,按照Boosting算法迭代构建CART决策树,最终的预测结果为:
Figure BDA0003623562240000034
式5中,
Figure BDA0003623562240000035
是最终预测结果,K是树的总数量,k代表第k棵决策树,γk是第k颗树的权重,hk表示第k棵树上的预测结果。
2.XGBoost
XGBoost是GBDT的改进,在
Figure BDA0003623562240000039
上有所不同,是对梯度提升树的进一步优化,目的是为了提高模型运行 效率,防止过拟合,提高模型的泛化能力,对于XGBoost来说,整个模型在这个样本i上给出的预测结果 为:
Figure BDA0003623562240000036
目标函数为:
Figure BDA0003623562240000037
Figure BDA0003623562240000038
式6中,fk表示第k棵决策树的函数,xi表示样本i对应的特征向量,K表示决策树的数量,F表示 所有决策树的集合。
式7中,L(Φ)t表示迭代过程的目标函数,
Figure BDA00036235622400000310
表示前t-1次迭代的预测值,Ω(fx)是防止过拟 合的正则项,γ和λ是正则项系数,防止决策树过于复杂;
模型建立后,进入参数调节步骤;
S42,参数调节具体通过如下过程进行:
首先,把数据集划分两份,一份是训练集,占总数据量的70%,用来调整模型的超参数,一份是测试 集,占总数据量的30%,用来测试模型的泛化能力;
其次,绘制回归模型指标变化的学习曲线,循坏调节每个超参数,选取超参数相对适宜的数值;
最后,对每个超参数选择一个数值范围,使用网格搜索的方法,选择超参数之间的最佳搭配方案。
可选的,上述的近地表痕量气体浓度反演模型,S4还包括模型评估,模型评估具体过程是:
根据搭建的回归模型,采用模型评估指标决定系数R2、平均平方误差MSE、平均绝对误差MAE、均方 根误差RMSE中的至少种参数进行评估;
R2取值范围为0~1,越接近于1,说明模型的预测效果越好,越接近于0,说明模型的预测效果越差, 如果为负值,说明模型的效果非常差,模型基本不可用,计算公式见式8:
Figure BDA0003623562240000041
MAE计算每一个样本的预测值和真实值的差的绝对值,然后求和再取平均值,用于评估预测结果和真 实数据集的接近程度,其值越小说明拟合效果越好,计算公式见式9;
Figure BDA0003623562240000042
MSE计算每一个预测值与真实值差的平方,然后求和再取平均值。该指标计算的是拟合数据和原始数 据对应样本点误差的平方和的均值,其值越小说明拟合效果越好,计算结果见式10;
Figure BDA0003623562240000043
RMSE均方根误差就是在均方误差的基础上再开方,其值越小说明拟合效果越好,计算公式见式11;
Figure BDA0003623562240000044
可选的,上述的近地表痕量气体浓度反演模型,为了防止模型过拟合,使用十折交叉验证的方法做模 型评估,十折是指训练集和验证集随机的分成十份,轮流将其中9份作为训练数据,1份作为测试数据, 总共开展十次,用十次结果的平均值来评估模型精度。
本发明还提供一种基于上述的近地表痕量气体浓度反演模型反演近地表痕量气体浓度的方法,通过如 下步骤进行:
S1,准备数据
除痕量气体浓度之外,需要准备的数据与搭建模型时用到的数据类型一致,需要将数据从文本格式处 理为最终的栅格数据,分辨率为1km×1km,坐标系为WGS-84坐标系,处理方法与建模时的数据预处理方 法一致,得到各个特征的栅格文件;
S2,处理栅格行列数
使用研究区的四至坐标裁剪每一个特征的栅格文件,按照遥感数据的行列号,将其他特征的栅格文件 处理成相同的行列号;
S3,提取数据
按照行列数依次提取各个特征的像元值,每一个像元的所有特征即为一个样本,提取完所有特征之后, 对缺失值进行填补,一般是遥感数据可能有缺失值,根据周围10个像元的平均值进行填补,消除遥感监 测的缺失值;
S4,反演痕量气体浓度
将充填后的所有数据导入已建立的模型中,得到每个样本对应的痕量气体浓度,在将痕量气体浓度处 理成栅格文件,每一个像元即是反演的痕量气体浓度。
本发明的近地表痕量气体浓度反演模型及基于该反演模型反演近地表痕量气体浓度的方法,模型结合 人工智能与大数据挖掘技术,建立了遥感监测的痕量气体总柱浓度(以下简称“遥感数据”)与地面站点 监测的痕量气体质量浓度(以下简称“站点数据”)的(反演模型),能够准确、高效地反演近地表痕量气 体的时空分布特征,支撑大气污染物防治工作的精准实施。
说明书附图
利用附图对本发明作进一步的说明,但附图中的内容不构成对本发明的任何限制。
图1是本发明实施例3中的O3与其他特征之间的相关关系散点图;
图2是O3与其它特征之间的相关系数矩阵;
图3是臭氧模型训练集散点密度图;
图4是臭氧模型训练集空间分布特征;
图5是各城市训练集监测的O3浓度月均值(红线)和预测的O3月均值浓度(蓝线);
图6是臭氧模型验证集散点密度图;
图7是臭氧模型各项评价指标的空间分布特征图,其中图7a为样本量;图7b为均方根误差;图7c为为 决定系数;图7d为平均绝对误差。
图8是各城市测试集监测的O3浓度月均值(红线)和预测的O3月均值。
图9是2020年6月关中城市群模型反演O3的空间分布特征;
图10是2020年6月22日关中城市群监测站点O3浓度空间分布。
具体实施方式
结合以下实施例对本发明作进一步说明。
实施例1。
一种近地表痕量气体浓度反演模型,基于遥感数据与站点数据建立反演模型。结合人工智能与大数据 挖掘与技术,建立遥感监测的痕量气体总柱浓度(以下简称“遥感数据”)与地面站点监测的痕量气体质 量浓度(以下简称“站点数据”)的回归模型(反演模型),准确、高效地反演近地表痕量气体的时空分布 特征,支撑大气污染物防治工作的精准实施。
具体的,该近地表痕量气体浓度反演模型,通过如下步骤建立:
S1,收集地表监测数据、遥感数据、气象数据和其他数据;
S2,对步骤S1收集的数据进行预处理,得到预处理后的初始建模数据;
S3,对初始建模数据进行特征筛选,得到建模数据集;
S4,根据建模数据集搭建模型。
其中,S1收集地表监测数据具体是:
从地方省级环境监测中心站(如陕西省环境监测站、山西省环境监测站等,具体可以根据需要判断的 地理位置选择)获取地表痕量气体浓度数据,包括O3、NO2、CO、SO2的小时值、8小时均值和24小时均值, 痕量气体使用13:00和14:00的算数平均值,数据中包含了监测站点的经度、纬度和日期;
收集遥感数据具体是:
痕量气体遥感数据使用的是哨兵5P的L2级总柱浓度产品数据,遥感数据从GoogleEarth Engine(GEE) 平台上下载;
收集气象数据具体是:
气象数据从省级气象局(如陕西省气象局、河北省气象局,具体跟进需要的地理区域选择)获取,气 象条件监测站点密集且分布均匀,平均1km2范围内的站点数大于1个,气象类型包括平均温度(TEM_Avg)、 最大温度(TEM_Max)、最小温度(TEM_Min)、平均相对湿度(RH)、8:00-8:00累积降雨量(PRE_08)、 20:00-20:00累积降雨量(PRE_20)、2分钟平均风速(WIN)、蒸发量(EVP)、日照时数(SSH),全部为日 均值数据;
收集其他数据具体是:
其他用到的数据有数字高程模型(DEM),空间分辨率为30m×30m;地表覆盖类型数据(GLC),空间分 辨率为30m×30m,人口密度数据(PopDen),空间分辨率0.09°×0.09°;
痕量气体具有空间和时间上异质性,将日期数据转化为每年的第几天(DOY),匹配时间异质性,用监 测站点的经纬度(Lon/Lat)匹配空间异质性。
对步骤S1收集的数据进行预处理的具体过程是:
S21,数据提取
首先,从空气质量监测的原始文件中提取2018年-2020年目标城市群所有空气质量监测站点的痕量气 体13:00-14:00的监测结果,删除缺失值后,将每天13:00和14:00的监测结果求算术平均值,作为 地面监测结果,得到痕量气体监测数据集。提取过程中同时包含监测站点的经度、纬度和监测日值,根据 监测站点的经度和纬度信息转化为WGS-84坐标系下的点矢量文件,简称“痕量气体点矢量”,用于提取对 应监测站点的气象数据、地理数据、人口密度数据;
其次,将GEE下载的遥感数据重投影至WGS84坐标系,使用最近邻采样方法重采样至1km×1km的分 辨率,使用痕量气体点矢量提取对应点位的遥感数据,得到遥感数据集;
再次,从气象数据原始文件中提取2018-2020年所有监测站点的日均值监测结果,包含气象监测站点 的经纬度和日期,剔除所有缺失值,依次按照日期和气象数据字段信息,按照经纬度和监测结果将各个气 象字段的监测数据转化为WGS-84坐标系下的点矢量文件,然后按照反距离权重插值方法插值到1km×1km 分辨率,进一步矢量转栅格后形成栅格文件,栅格文件的像元值就是气象字段的日均值监测结果,最后使 用痕量气体点矢量提取对每天各字段对应位置的气象结果,得到气象数据集;
最后,DEM、GLC、PopDen均为栅格文件,其中DEM和GLC分辨率为30m×30m,PopDen分辨率为900m×900m, 使用最近邻采样方法重采样至1km×1km,重投影至WGS-84坐标系下,使用痕量气体矢量提取,得到辅助 数据集;
S22,数据结合
数据提取阶段共提取了4个数据集,分别为痕量气体的监测数据集、遥感数据集、气象数据集、辅助 数据集,四份数据集中均包含有经度、纬度和日期,按照这三个字段将四份数据结合成为最终的数据集, 简称“原始数据集”,进一步将原始数据集中的日期转化为每年的第几天,Day Of Year(DOY),原始数据 集中痕量气体的监测结果为标签,其余全部为特征;
S23.异常值剔除
对原始数据集进行异常值检查,检查异常值的方法是绘制箱型图,将含有异常值的样本全部删除;此 外,根据先验知识将气象数据的异常值(如WIN>50m/s,TEM_Avg>40℃,SSH>24h)剔除,即删除极端 天气的影响,同样将含有极端天气的样本删除。
S3中对初始建模数据进行特征筛选,具体包括:
S31.绘制散点图研究相关关系
通过绘制散点图的方式,研究数据集中各特征变量与痕量气体之间的相关关系,数据中包含长时间序 列、大范围的气体,散点图可以有效反应痕量气体与其他特征之间的相关关系,根据先验知识将基本无相 关的特征删除;
S32.统计person相关系数和p_value
通过统计特征之间的person相关系数,痕量气体与其他变量之间相关性,将person相关系数为0.00 的特征删除,进一步计算p_value值,将p_value<0.01的变量删除;
S33.机器学习筛选特征
将步骤S31、S32筛选后的数据集带入随机森林模型,使用该模型的内建函数“feature.importance” 计算特征重要性,排序后按照特征重要性得分由小到大排序,依次删除得分较小的特征,查看建模精度变 化,直到模型精度发生较大变化时停止删除特征,得到建模数据集,简称“数据集”。
S4搭建模型建模使用的机器学习算法为eXtreme Gradient Boosting(XGBoost),是以决策树为基学 习器构建的集成算法,基学习器为Classification And RegressionTree(CART)决策树,使用XGBoost反 演近地表痕量气体浓度旨在建立回归模型,模型搭建过程使用的标签为痕量气体的浓度,特征为数据集中 除痕量气体浓度之外的其他特征;
具体过程是:
S41,CART树
对于给定的数据集,选择最优切分特征j(特征信息)与切分点s(特征中所有可能取值;
遍历特征j,对固定的切分特征j扫描切分点s,选择使式1达到最小值的(j,s)对;
Figure BDA0003623562240000061
式1中,R1和R2表示由(j,s)对划分的两个数据集,yi是样本i的标签(锂的检测结果),xi表示数 据集中的(j,s)对,c1是R1上的平均值,c2是R2上的平均值,
用选定的(j,s)对划分区域并决定相应的输出值:
R1(j,s)={x|x(j)≤s},R2(j,s)={x|x(j)>s}
(式2);
Figure BDA0003623562240000062
继续对两个子区域R1和R2调用步骤①,②,直至满足停止条件;
将输入空间(数据集)划分为M个区域(叶子结点)R1,R2,……,Rm,生成决策树:
Figure BDA0003623562240000071
式4中,f(x)表示预测结果,M表示叶子结点的个数,m表示第m个叶子结点,Rm表示叶子结点的集 合,
Figure BDA0003623562240000076
表示叶子结点的平均值(该结点上的预测值),I表示叶子结点的权重(叶子结点的样本个数占总 样本个数的权重);
S42,Boosting
Boosting是一族可将弱学习器提升为强学习器的算法;这族算法的工作机制类似:先从初始训练集训 练出一个基学习器,再根据基学习器的表现对训练样本分布进行调整,使得先前基学习器做错的训练样本 在后续受到更多关注,然后基于调整后的样本分布来训练下一个基学习器;如此重复进行,直至基学习器 数目达到事先指定的值T,最终将这T个基学习器进行加权结合;
S43,GBDT
GBDT专注于回归的树的提升集成模型,按照Boosting算法迭代构建CART决策树,最终的预测结果为:
Figure BDA0003623562240000072
式5中,
Figure BDA0003623562240000077
是最终预测结果,K是树的总数量,k代表第k棵决策树,γk是第k颗树的权重,hk表示第k棵树上的预测结果。
3.XGBoost
XGBoost是GBDT的改进,在
Figure BDA00036235622400000710
上有所不同,是对梯度提升树的进一步优化,目的是为了提高模型运行 效率,防止过拟合,提高模型的泛化能力,对于XGBoost来说,整个模型在这个样本i上给出的预测结果 为:
Figure BDA0003623562240000073
目标函数为:
Figure BDA0003623562240000074
Figure BDA0003623562240000075
式6中,fk表示第k棵决策树的函数,xi表示样本i对应的特征向量,K表示决策树的数量,F表示 所有决策树的集合。
式7中,L(Φ)t表示迭代过程的目标函数,
Figure BDA0003623562240000079
表示前t-1次迭代的预测值,Ω(fx)是防止过拟 合的正则项,γ和λ是正则项系数,防止决策树过于复杂;
模型建立后,进入参数调节步骤;
S42,参数调节具体通过如下过程进行:
首先,把数据集划分两份,一份是训练集,占总数据量的70%,用来调整模型的超参数,一份是测试 集,占总数据量的30%,用来测试模型的泛化能力;
其次,绘制回归模型指标变化的学习曲线,循坏调节每个超参数,选取超参数相对适宜的数值;
最后,对每个超参数选择一个数值范围,使用网格搜索的方法,选择超参数之间的最佳搭配方案。
S4还包括模型评估,模型评估具体过程是:
根据搭建的回归模型,采用模型评估指标决定系数R2、平均平方误差MSE、平均绝对误差MAE、均方 根误差RMSE中的至少种参数进行评估;
R2取值范围为0~1,越接近于1,说明模型的预测效果越好,越接近于0,说明模型的预测效果越差, 如果为负值,说明模型的效果非常差,模型基本不可用,计算公式见式8:
Figure BDA0003623562240000081
MAE计算每一个样本的预测值和真实值的差的绝对值,然后求和再取平均值,用于评估预测结果和真 实数据集的接近程度,其值越小说明拟合效果越好,计算公式见式9;
Figure BDA0003623562240000082
MSE计算每一个预测值与真实值差的平方,然后求和再取平均值。该指标计算的是拟合数据和原始数 据对应样本点误差的平方和的均值,其值越小说明拟合效果越好,计算结果见式10;
Figure BDA0003623562240000083
RMSE均方根误差就是在均方误差的基础上再开方,其值越小说明拟合效果越好,计算公式见式11;
Figure BDA0003623562240000084
可选的,上述的近地表痕量气体浓度反演模型,为了防止模型过拟合,使用十折交叉验证的方法做模 型评估,十折是指训练集和验证集随机的分成十份,轮流将其中9份作为训练数据,1份作为测试数据, 总共开展十次(即十次十折交叉验证),用十次结果的平均值来评估模型精度。这个方法的优势在于,同 时重复运用随机产生的子样本进行训练和验证,每次的结果验证一次。
模型搭建完成之后评价模型在时间上和空间上的表现特征,如果评估效果不理想,需要重新调节参数, 优化模型。
模型调优之后,用测试集验证模型的泛化能力。
该近地表痕量气体浓度反演模型,结合人工智能与大数据挖掘技术,建立了遥感监测的痕量气体总柱 浓度(以下简称“遥感数据”)与地面站点监测的痕量气体质量浓度(以下简称“站点数据”)的(反演模 型),能够准确、高效地反演近地表痕量气体的时空分布特征,支撑大气污染物防治工作的精准实施。
实施例2。
一种基于实施例1的近地表痕量气体浓度反演模型反演近地表痕量气体浓度的方法,通过如下步骤进 行:
S1,准备数据
除痕量气体浓度之外,需要准备的数据与搭建模型时用到的数据类型一致,需要将数据从文本格式处 理为最终的栅格数据,分辨率为1km×1km,坐标系为WGS-84坐标系,处理方法与建模时的数据预处理方 法一致,得到各个特征的栅格文件;
S2,处理栅格行列数
使用研究区的四至坐标裁剪每一个特征的栅格文件,按照遥感数据的行列号,将其他特征的栅格文件 处理成相同的行列号;
S3,提取数据
按照行列数依次提取各个特征的像元值,每一个像元的所有特征即为一个样本,提取完所有特征之后, 对缺失值进行填补,一般是遥感数据可能有缺失值,根据周围10个像元的平均值进行填补,消除遥感监 测的缺失值;
S4,反演痕量气体浓度
将充填后的所有数据导入已建立的模型中,得到每个样本对应的痕量气体浓度,在将痕量气体浓度处 理成栅格文件,每一个像元即是反演的痕量气体浓度。
本实施例的基于近地表痕量气体浓度反演模型反演近地表痕量气体浓度的方法,能够准确、高效地反 演近地表痕量气体的浓度情况。
实施例3。
利用陕西省环境气象监测站、陕西省气象局的数据通过实施例1的方式建立的近地表痕量气体浓度 反演模型,利用此模型对陕西省的近地表痕量气体浓度进行反演。以2020年6月22日近地表臭氧浓度反 演为例,验证本实施例方法的有效性。
1.数据收集
1)地表监测数据
地表O3浓度数据从中国环境监测总站获取,包括了O3的小时值、8小时均值和24小时均值,为了和 遥感数据相匹配,使用13:00和14:00的平均数据搭建臭氧模型,数据中包含了监测站点的经度、纬度和 日期。
2)遥感数据
遥感数据使用的是哨兵5P的L2级O3的总柱浓度产品数据。哨兵5P卫星过境时间为13:00-14:00之间, 搭载的传感器“TROPOMI”是迄今为止技术性能先进、空间分辨率最高的大气监测光谱仪,时间分辨率为1 天,空间分辨率为7km×3.5km。遥感数据从GEE平台上下载。
3)气象数据
气象数据从陕西省气象局获取,全部为日均值文件,气象条件监测站点密集且分布均匀,气象类型包 括平均温度(TEM_Avg)、最大温度(TEM_Max)、最小温度(TEM_Min)、平均相对湿度(RH)、8:00-8:00累 积降雨量(PRE_08)、20:00-20:00累积降雨量(PRE_20)、2分钟平均风速(WIN)、蒸发量(EVP)、日照 时数(SSH),全部为日均值数据。
4)其他数据
其他用到的数据有数字高程模型(DEM),空间分辨率为30m×30m;地表覆盖类型数据(GLC),空间分 辨率为30m×30m,人口密度数据(PopDen),空间分辨率0.09°×0.09°。
O3具有空间和时间上异质性,所以将日期数据转化为每年的第几天(DOY),匹配时间异质性,用监测 站点的经纬度(Lon/Lat)匹配空间异质性。
2.数据整合
提取关中城市群所有地面监测站点(共54个)2018年9月至2019年12月13:00和14:00的O3监测 浓度,按天计算O3的平均浓度作为建模的标签数据。将遥感数据、气象数据、其他数据重投影至WGS-84 坐标系,采用反距离权重法重采样至1km×1km空间分辨率,按照O3监测站点的经纬度提取每日的遥感数 据、气象数据、其他数据的对应值,结合到一起,形成了以O3为标签,以DOY、Lon、Lat、TEM_Avg、TEM_Max、 TEM_Min、RH、PRE08、PRE20、WIN、ENP、SSH、DEM、GLC、PopDen、SP-5O3为特征的数据集,总共包含样 本量25717条样本量。
3.特征选择
根据绘制的散点图,如图1所示,其中蒸发量和日照时数数据明显存在问题,理论上讲关中城市群日 照时数不能超过15小时,蒸发量对O3的影响完全没有规律,将日照时数和蒸发量从数据集中删除。
DOY对臭氧的影响关系明显,每年的5月到9月是臭氧的高发季节;其余时间臭氧的浓度较低。经度 大、纬度高的关中城市群东北部地区(临汾市和运城市)臭氧明显高于西部地区,空间位置对臭氧的浓度 影响较大;温度越高,臭氧浓度越高,温度降低,臭氧浓度相应降低,最高温度、最低温度和平均温度的 影响一致;湿度在30%~60%的范围内,O3浓度最高;降雨在2000mm以内,臭氧的浓度急剧下降,大于2000mm 之后,降雨对臭氧的浓度影响不大;风速大于2m/s,对臭氧的传输效果明显;海拔高度在400m~600m的范 围内,O3浓度最高,这是由于关中城市群海拔低处多为平原地区,人类生产活动强,排放量大,不同的土 地利用类型和地表覆盖类型以及人口密度对O3浓度影响有差别,但影响不明显,且以上变量对臭氧的影响 为非线性影响,所以建模过程中要使用非线性的模型。
所有特征与O3的相关性包括正相关和负相关,没有不相关的特征,如图2所示,但DEM和WIN_S_2mi_Avg 两个特征的p_value均大于0.01,如表1的结果,将这两个特征从数据集中剔除,得到最终建模的数据集。
表1特征的p_value结果表
Figure BDA0003623562240000091
Figure BDA0003623562240000101
4.模型搭建
通过单参数调节、网格搜索等手段结合学习曲线多次尝试,尽可能将模型的验证集精度提高,随后在 保持验证集模型精度基本不变的情况下,通过剪枝操作,将训练集的模型精度下调,最大限度的缩小训练 集和验证集的精度,降低模型的泛化误差,提高模型的鲁棒性,并通过测试集检验模型的泛化能力。
最终确定的模型,训练集R2为0.98,十折交叉验证R2平均值为0.94,测试集R2为0.94。
5模型评价
5.1训练集模型评价
5.1.1训练集模型总体表现
臭氧模型训练集总体表现较好,决定系数R2为0.98,几乎没有偏差,表现出轻微的过拟合现象,表明 XGBoost强大的学习能力。图3展示了建模的结果,训练样本为18001条,均方根误差为6.08μg/m3
5.1.2训练集模型空间表现
为了进一步探索模型在空间上的表现,将训练集中每一个点的模型表现特征做了统计。其中各监测站 点的样本数为300~362个,决定系数R2全部大于0.90,如图4所示,为0.94~1,均方误差MSE为 14.39~59.6μg/m3,平均绝对误差MAE为2.9~5.82μg/m3,均方根误差为3.79~7.72μg/m3,见表2。
从空间上来看,关中城市群临汾、运城、渭南、西安、咸阳、宝鸡这几个城市的站点建模效果最佳, R2为0.98~1,这些城市地形平坦,为关中城市群工业生产发展较多的城市,污染物浓度相对较高,气象条 件相对稳定,所以预测效果最好。其余城市多为山区,气象条件不均匀,污染物浓度相对较低,模型效果 略有下降,但仍在0.94~0.97之间。
以上统计结果表明在训练集数据上,模型在空间上的表现优异。
表2关中城市群各监测站点训练集的模型评估因子统计结果
Figure BDA0003623562240000111
Figure BDA0003623562240000121
5.1.3.训练集模型时间表现
为了研究模型在不同时间的表现,对各个城市的监测站点分别求监测值和模型预测值的月均值,分城 市对比模型在不同时间上的表现。
从时间上来看,每年的5月到9月是臭氧污染物的高发季节,这与高温情况下NO2等臭氧前体物的转 化密切相关,每年的12月次年到2月臭氧浓度最低,关中城市群中各个城市均具有这样的规律,如图5 所示。
从地表监测结果和模型预测结果的差距(图中两线之间的差距)来看,无论是臭氧的高发季节,还是 臭氧浓度较低的季节,监测值和预测值基本一致,尤其是临汾、运城、渭南、西安、咸阳、宝鸡这些城市, 监测值绘制的线与预测值绘制的线基本重合,如图5所示,在其余几个城市上,两条线的重合度也非常高。
在训练集上,模型在时间上的预测精度高,表现优异。
5.2测试集模型表现
5.2.1.测试集模型总体表现
模型在测试集的表现基本上反应了模型在未知数据集上的表现。体现了模型的泛化能力和鲁棒性。
测试集总共7716个样本,来自关中城市群54个监测站点,模型对测试集预测的结果与监测结果之间 的偏差为0.04,决定系数R2为0.94,如图6所示,无论是在低值区域还是在高值区域,点位都比较均匀 的分布在1:1线(黑线)的两侧,收敛效果较好。
模型回归线的斜率为0.92,如图6,截距为6.91μg/m3,表明模型在低值区域,预测的O3浓度可能略 高于监测的O3浓度,相反的,在高值区域,预测的O3浓度可能略低于监测的O3浓度,而且随着浓度的提 升,这一误差可能越明显。模型整体上的预测能力优异,可以用于对面上O3浓度的预测。
5.2.2.测试集模型空间表现
为研究测试集数据上,模型在空间的表现,同样将测试集中每一个点的模型表现特征做了统计。
各监测站点的样本数为112~116个,约有32%的监测站点测试数据超过150个,约有60%的监测站点 测试数据超过140个,约91%的监测站点测试数据超过130个,测试数据的样本量较多,从空间分布来看, 各城市的监测站点分布比较均匀,天水市的样本量略少,如图7a所示。
各监测站点的决定系数R2为0.61~0.98,见表3,平均为0.94,约有70%的监测站点R2≥0.90,约98% 的监测站点R2≥0.80,只有一个监测站点的R2小于0.80。从空间分布来看,临汾、运城、渭南、西安、咸 阳、宝鸡这几个城市所有站点的R2全部大于0.90,铜川市、商洛市、庆阳、平凉这几个城市所有站点的 R2介于0.80只0.90之间,天水市有一个站点的R2仅有0.61,如图7c所示,与其他站点相差较大。
各监测站点的均方误差MSE为80.94~333.43μg/m3,平均约154.17μg/m3,约13%的监测站点 MSE≤100μg/m3,约83%的监测站点MSE≤200μg/m3
各监测站点的平均绝对误差MAE为6.58~12.24μg/m3(表3),平均约9.07μg/m3,约54%的站点 MAE≤8μg/m3,约70%的监测站点MAE≤10μg/m3。从空间分布来看,平均绝对误差最小的监测站点都分 布在临汾、运城、渭南、西安、宝鸡这几个城市,平均误差最大的监测站点分布在商洛(图7d)。
各监测站点的均方根误差为9.00~18.26μg/m3(表3),平均约12.24μg/m3,约31%的监测站点 RMSE≤10μg/m3,约87%的监测站点RMSE≤15μg/m3,从空间分布来看,均方根误差最小的监测站点都分 布在临汾、运城、渭南、西安、宝鸡这几个城市(图7d)。
从空间分布来看,临汾、运城、渭南、西安、咸阳、宝鸡这几个城市的模型表现最好,其次为铜川、 庆阳、平凉,天水和商洛略微差一点。
表3关中城市群各监测站点测试集的模型评估因子统计结果
Figure BDA0003623562240000122
Figure BDA0003623562240000131
5.3测试集模型时间表现
临汾市各时间的模型表现总体较好,预测有偏差的2018年9月、2018年11月,2019年8月,差异 均在5μg/m3范围内,如图8所示。
运城、渭南、西安、咸阳、宝鸡、铜川这几个城市不同时间段的地面监测O3浓度与预测的O3浓度基本 无偏差。
商洛、平凉、庆阳、天水这几个城市不同时间上略有偏差,偏差不大于10μg/m3
通过以上预测,可以看出模型在时维度上表现优异。
6模型应用
收集2020年6月22日的气象等数据(与建模时的数据类型一致),重投影至WGS-84坐标系,反距离 权重法重采样值1km×1km,提取每一个网格的像元值,带入模型进行反演,并将反演的O3浓度(图9)与 站点监测的O3浓度的反距离权重插值结果(图10)进行了对比,总体分布趋势一致,精度大幅度提升。 可见,本实施例的方法能够精确、有效用于近地表痕量气体浓度的反演。
最后应当说明的是,以上实施例仅用以说明本发明的技术方案而非对本发明保护范围的限制,尽管参 照较佳实施例对本发明作了详细说明,本领域的普通技术人员应当理解,可以对本发明的技术方案进行修 改或者等同替换,而不脱离本发明技术方案的实质和范围。

Claims (9)

1.一种近地表痕量气体浓度反演模型,其特征在于:基于遥感数据与站点数据建立反演模型。
2.根据权利要求1所述的近地表痕量气体浓度反演模型,其特征在于:通过如下步骤建立:
S1,收集地表监测数据、遥感数据、气象数据和其他数据;
S2,对步骤S1收集的数据进行预处理,得到预处理后的初始建模数据;
S3,对初始建模数据进行特征筛选,得到建模数据集;
S4,根据建模数据集搭建模型。
3.根据权利要求2所述的近地表痕量气体浓度反演模型,其特征在于:S1中,收集地表监测数据具体是:
从地方省级环境监测中心站获取地表痕量气体浓度数据,包括O3、NO2、CO、SO2的小时值、8小时均值和24小时均值,痕量气体使用13:00和14:00的算数平均值,数据中包含了监测站点的经度、纬度和日期;
收集遥感数据具体是:
痕量气体遥感数据使用的是哨兵5P的L2级总柱浓度产品数据,遥感数据从GoogleEarth Engine(GEE)平台上下载;
收集气象数据具体是:
气象数据从省级气象局获取,气象条件监测站点密集且分布均匀,平均1km2范围内的站点数大于1个,气象类型包括平均温度(TEM_Avg)、最大温度(TEM_Max)、最小温度(TEM_Min)、平均相对湿度(RH)、8:00-8:00累积降雨量(PRE_08)、20:00-20:00累积降雨量(PRE_20)、2分钟平均风速(WIN)、蒸发量(EVP)、日照时数(SSH),全部为日均值数据;
收集其他数据具体是:
其他用到的数据有数字高程模型(DEM),空间分辨率为30m×30m;地表覆盖类型数据(GLC),空间分辨率为30m×30m,人口密度数据(PopDen),空间分辨率0.09°×0.09°;
痕量气体具有空间和时间上异质性,将日期数据转化为每年的第几天(DOY),匹配时间异质性,用监测站点的经纬度(Lon/Lat)匹配空间异质性。
4.根据权利要求3所述的近地表痕量气体浓度反演模型,其特征在于:S1中:对步骤S1收集的数据进行预处理的具体过程是:
S21,数据提取
首先,从空气质量监测的原始文件中提取2018年-2020年目标城市群所有空气质量监测站点的痕量气体13:00-14:00的监测结果,删除缺失值后,将每天13:00和14:00的监测结果求算术平均值,作为地面监测结果,得到痕量气体监测数据集。提取过程中同时包含监测站点的经度、纬度和监测日值,根据监测站点的经度和纬度信息转化为WGS-84坐标系下的点矢量文件,简称“痕量气体点矢量”,用于提取对应监测站点的气象数据、地理数据、人口密度数据;
其次,将GEE下载的遥感数据重投影至WGS84坐标系,使用最近邻采样方法重采样至1km×1km的分辨率,使用痕量气体点矢量提取对应点位的遥感数据,得到遥感数据集;
再次,从气象数据原始文件中提取2018-2020年所有监测站点的日均值监测结果,包含气象监测站点的经纬度和日期,剔除所有缺失值,依次按照日期和气象数据字段信息,按照经纬度和监测结果将各个气象字段的监测数据转化为WGS-84坐标系下的点矢量文件,然后按照反距离权重插值方法插值到1km×1km分辨率,进一步矢量转栅格后形成栅格文件,栅格文件的像元值就是气象字段的日均值监测结果,最后使用痕量气体点矢量提取对每天各字段对应位置的气象结果,得到气象数据集;
最后,DEM、GLC、PopDen均为栅格文件,其中DEM和GLC分辨率为30m×30m,PopDen分辨率为900m×900m,使用最近邻采样方法重采样至1km×1km,重投影至WGS-84坐标系下,使用痕量气体矢量提取,得到辅助数据集;
S22,数据结合
数据提取阶段共提取了4个数据集,分别为痕量气体的监测数据集、遥感数据集、气象数据集、辅助数据集,四份数据集中均包含有经度、纬度和日期,按照这三个字段将四份数据结合成为最终的数据集,简称“原始数据集”,进一步将原始数据集中的日期转化为每年的第几天,Day Of Year(DOY),原始数据集中痕量气体的监测结果为标签,其余全部为特征;
S23.异常值剔除
对原始数据集进行异常值检查,检查异常值的方法是绘制箱型图,将含有异常值的样本全部删除;此外,根据先验知识将气象数据的异常值(如WIN>50m/s,TEM_Avg>40℃,SSH>24h)剔除,即删除极端天气的影响,同样将含有极端天气的样本删除。
5.根据权利要求4所述的近地表痕量气体浓度反演模型,其特征在于:S1中:S3中对初始建模数据进行特征筛选,具体包括:
S31.绘制散点图研究相关关系
通过绘制散点图的方式,研究数据集中各特征变量与痕量气体之间的相关关系,数据中包含长时间序列、大范围的气体,散点图可以有效反应痕量气体与其他特征之间的相关关系,根据先验知识将基本无相关的特征删除;
S32.统计person相关系数和p_value
通过统计特征之间的person相关系数,痕量气体与其他变量之间相关性,将person相关系数为0.00的特征删除,进一步计算p_value值,将p_value<0.01的变量删除;
S33.机器学习筛选特征
将步骤S31、S32筛选后的数据集带入随机森林模型,使用该模型的内建函数“feature.importance”计算特征重要性,排序后按照特征重要性得分由小到大排序,依次删除得分较小的特征,查看建模精度变化,直到模型精度发生较大变化时停止删除特征,得到建模数据集,简称“数据集”。
6.根据权利要求5所述的近地表痕量气体浓度反演模型,其特征在于:S1中:S4搭建模型建模使用的机器学习算法为eXtreme Gradient Boosting(XGBoost),是以决策树为基学习器构建的集成算法,基学习器为Classification And Regression Tree(CART)决策树,使用XGBoost反演近地表痕量气体浓度旨在建立回归模型,模型搭建过程使用的标签为痕量气体的浓度,特征为数据集中除痕量气体浓度之外的其他特征;
具体过程是:
S41,CART树
对于给定的数据集,选择最优切分特征j(特征信息)与切分点s(特征中所有可能取值;
遍历特征j,对固定的切分特征j扫描切分点s,选择使式1达到最小值的(j,s)对;
Figure FDA0003623562230000021
式1中,R1和R2表示由(j,s)对划分的两个数据集,yi是样本i的标签(锂的检测结果),xi表示数据集中的(j,s)对,c1是R1上的平均值,c2是R2上的平均值,
用选定的(j,s)对划分区域并决定相应的输出值:
R1(j,s)={x|x(j)≤s),R2(j,s)=(x|x(j)>S} ((式2);
Figure FDA0003623562230000022
继续对两个子区域R1和R2调用步骤①,②,直至满足停止条件;
将输入空间(数据集)划分为M个区域(叶子结点)R1,R2,……,Rm,生成决策树:
Figure FDA0003623562230000023
式4中,f(x)表示预测结果,M表示叶子结点的个数,m表示第m个叶子结点,Rm表示叶子结点的集合,
Figure FDA0003623562230000024
表示叶子结点的平均值(该结点上的预测值),I表示叶子结点的权重(叶子结点的样本个数占总样本个数的权重);
S42,Boosting
Boosting是一族可将弱学习器提升为强学习器的算法;这族算法的工作机制类似:先从初始训练集训练出一个基学习器,再根据基学习器的表现对训练样本分布进行调整,使得先前基学习器做错的训练样本在后续受到更多关注,然后基于调整后的样本分布来训练下一个基学习器;如此重复进行,直至基学习器数目达到事先指定的值T,最终将这T个基学习器进行加权结合;
S43,GBDT
GBDT专注于回归的树的提升集成模型,按照Boosting算法迭代构建CART决策树,最终的预测结果为:
Figure FDA0003623562230000031
式5中,
Figure FDA0003623562230000032
是最终预测结果,K是树的总数量,k代表第k棵决策树,γk是第k颗树的权重,hk表示第k棵树上的预测结果。
1.XGBoost
XGBoost是GBDT的改进,在
Figure FDA0003623562230000033
上有所不同,是对梯度提升树的进一步优化,目的是为了提高模型运行效率,防止过拟合,提高模型的泛化能力,对于XGBoost来说,整个模型在这个样本i上给出的预测结果为:
Figure FDA0003623562230000034
目标函数为:
Figure FDA0003623562230000035
Figure FDA0003623562230000036
式6中,fk表示第k棵决策树的函数,xi表示样本i对应的特征向量,K表示决策树的数量,F表示所有决策树的集合。
式7中,L(Φ)t表示迭代过程的目标函数,
Figure FDA0003623562230000037
表示前t-1次迭代的预测值,Ω(fx)是防止过拟合的正则项,γ和λ是正则项系数,防止决策树过于复杂;
模型建立后,进入参数调节步骤;
S42,参数调节具体通过如下过程进行:
首先,把数据集划分两份,一份是训练集,占总数据量的70%,用来调整模型的超参数,一份是测试集,占总数据量的30%,用来测试模型的泛化能力;
其次,绘制回归模型指标变化的学习曲线,循坏调节每个超参数,选取超参数相对适宜的数值;
最后,对每个超参数选择一个数值范围,使用网格搜索的方法,选择超参数之间的最佳搭配方案。
7.根据权利要求6所述的近地表痕量气体浓度反演模型,其特征在于:S1中:S4还包括模型评估,模型评估具体过程是:
根据搭建的回归模型,采用模型评估指标决定系数R2、平均平方误差MSE、平均绝对误差MAE、均方根误差RMSE中的至少种参数进行评估;
R2取值范围为0~1,越接近于1,说明模型的预测效果越好,越接近于0,说明模型的预测效果越差,如果为负值,说明模型的效果非常差,模型基本不可用,计算公式见式8:
Figure FDA0003623562230000038
MAE计算每一个样本的预测值和真实值的差的绝对值,然后求和再取平均值,用于评估预测结果和真实数据集的接近程度,其值越小说明拟合效果越好,计算公式见式9;
Figure FDA0003623562230000039
MSE计算每一个预测值与真实值差的平方,然后求和再取平均值。该指标计算的是拟合数据和原始数据对应样本点误差的平方和的均值,其值越小说明拟合效果越好,计算结果见式10;
Figure FDA0003623562230000041
RMSE均方根误差就是在均方误差的基础上再开方,其值越小说明拟合效果越好,计算公式见式11;
Figure FDA0003623562230000042
8.根据权利要求7所述的近地表痕量气体浓度反演模型,其特征在于:S1中:为了防止模型过拟合,使用十折交叉验证的方法做模型评估,十折是指训练集和验证集随机的分成十份,轮流将其中9份作为训练数据,1份作为测试数据,总共开展十次,用十次结果的平均值来评估模型精度。
9.基于如权利要求1至8任意一项所述的近地表痕量气体浓度反演模型反演近地表痕量气体浓度的方法,其特征在于:通过如下步骤进行:
S1,准备数据
除痕量气体浓度之外,需要准备的数据与搭建模型时用到的数据类型一致,需要将数据从文本格式处理为最终的栅格数据,分辨率为1km×1km,坐标系为WGS-84坐标系,处理方法与建模时的数据预处理方法一致,得到各个特征的栅格文件;
S2,处理栅格行列数
使用研究区的四至坐标裁剪每一个特征的栅格文件,按照遥感数据的行列号,将其他特征的栅格文件处理成相同的行列号;
S3,提取数据
按照行列数依次提取各个特征的像元值,每一个像元的所有特征即为一个样本,提取完所有特征之后,对缺失值进行填补,一般是遥感数据可能有缺失值,根据周围10个像元的平均值进行填补,消除遥感监测的缺失值;
S4,反演痕量气体浓度
将充填后的所有数据导入已建立的模型中,得到每个样本对应的痕量气体浓度,在将痕量气体浓度处理成栅格文件,每一个像元即是反演的痕量气体浓度。
CN202210465007.5A 2022-04-29 2022-04-29 近地表痕量气体浓度反演模型及反演方法 Pending CN115420690A (zh)

Priority Applications (2)

Application Number Priority Date Filing Date Title
CN202210465007.5A CN115420690A (zh) 2022-04-29 2022-04-29 近地表痕量气体浓度反演模型及反演方法
CN202211545701.4A CN116223395A (zh) 2022-04-29 2022-12-05 近地表痕量气体浓度反演模型及反演方法

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
CN202210465007.5A CN115420690A (zh) 2022-04-29 2022-04-29 近地表痕量气体浓度反演模型及反演方法

Publications (1)

Publication Number Publication Date
CN115420690A true CN115420690A (zh) 2022-12-02

Family

ID=84196616

Family Applications (2)

Application Number Title Priority Date Filing Date
CN202210465007.5A Pending CN115420690A (zh) 2022-04-29 2022-04-29 近地表痕量气体浓度反演模型及反演方法
CN202211545701.4A Pending CN116223395A (zh) 2022-04-29 2022-12-05 近地表痕量气体浓度反演模型及反演方法

Family Applications After (1)

Application Number Title Priority Date Filing Date
CN202211545701.4A Pending CN116223395A (zh) 2022-04-29 2022-12-05 近地表痕量气体浓度反演模型及反演方法

Country Status (1)

Country Link
CN (2) CN115420690A (zh)

Cited By (2)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN116504330A (zh) * 2023-06-28 2023-07-28 航天宏图信息技术股份有限公司 污染物浓度反演方法、装置、电子设备及可读存储介质
CN117216490A (zh) * 2023-11-08 2023-12-12 中国铁道科学研究院集团有限公司电子计算技术研究所 一种智能大数据采集系统

Cited By (4)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN116504330A (zh) * 2023-06-28 2023-07-28 航天宏图信息技术股份有限公司 污染物浓度反演方法、装置、电子设备及可读存储介质
CN116504330B (zh) * 2023-06-28 2023-09-19 航天宏图信息技术股份有限公司 污染物浓度反演方法、装置、电子设备及可读存储介质
CN117216490A (zh) * 2023-11-08 2023-12-12 中国铁道科学研究院集团有限公司电子计算技术研究所 一种智能大数据采集系统
CN117216490B (zh) * 2023-11-08 2024-01-19 中国铁道科学研究院集团有限公司电子计算技术研究所 一种智能大数据采集系统

Also Published As

Publication number Publication date
CN116223395A (zh) 2023-06-06

Similar Documents

Publication Publication Date Title
CN111859800B (zh) 用于pm2.5浓度分布的时空估算和预测方法
Xu et al. Evaluation of machine learning techniques with multiple remote sensing datasets in estimating monthly concentrations of ground-level PM2. 5
CN112905560B (zh) 一种多源时空大数据深度融合的空气污染预测方法
CN109213964B (zh) 一种融合多源特征地理参数的卫星aod产品校正方法
CN109344865B (zh) 一种多数据源的数据融合方法
CN108227041B (zh) 基于站点实测数据和模式结果的水平能见度预报方法
CN113297527B (zh) 基于多源城市大数据的pm2.5全面域时空计算推断方法
CN114926749B (zh) 基于遥感图像的近地面大气污染物反演方法及系统
CN110751094A (zh) 一种基于gee综合遥感影像和深度学习方法的作物估产技术
Monteil et al. The regional European atmospheric transport inversion comparison, EUROCOM: first results on European-wide terrestrial carbon fluxes for the period 2006–2015
CN115420690A (zh) 近地表痕量气体浓度反演模型及反演方法
CN109782373B (zh) 一种基于改进的Naive Bayesian-CNN多目标分类算法的沙尘暴预测方法
CN113553764B (zh) 一种基于深度学习网络的山火预测方法
Chi et al. Machine learning-based estimation of ground-level NO2 concentrations over China
Yu et al. Deep learning-based downscaling of tropospheric nitrogen dioxide using ground-level and satellite observations
CN114004163A (zh) 一种基于modis和长短时记忆网络模型的pm2.5反演方法
Scheibenreif et al. Toward global estimation of ground-level no 2 pollution with deep learning and remote sensing
He et al. Spatiotemporal high-resolution imputation modeling of aerosol optical depth for investigating its full-coverage variation in China from 2003 to 2020
CN109657988B (zh) 基于hasm和欧氏距离算法的烟叶品质分区方法
CN114882373A (zh) 基于深度神经网络的多特征融合沙尘暴预测方法
Li et al. Generating daily high-resolution and full-coverage XCO2 across China from 2015 to 2020 based on OCO-2 and CAMS data
Liu et al. First satellite-based regional hourly NO2 estimations using a space-time ensemble learning model: A case study for Beijing-Tianjin-Hebei Region, China
CN117219183A (zh) 多云雨地区的高覆盖度近地面no2浓度估算方法及系统
CN111879915B (zh) 一种滨海湿地高分辨率的逐月土壤盐度监测方法及系统
CN116340863B (zh) 空气污染物预测方法、装置、电子设备及可读存储介质

Legal Events

Date Code Title Description
PB01 Publication
PB01 Publication
WD01 Invention patent application deemed withdrawn after publication
WD01 Invention patent application deemed withdrawn after publication

Application publication date: 20221202