基于随机森林算法的DMSP/OLS数据校正方法、系统及介质

基于随机森林算法的dmsp/ols数据校正方法、系统及介质
技术领域
1.本发明属于遥感图像处理的技术领域，具体涉及一种基于随机森林算法的dmsp/ols数据校正方法、系统及介质。


背景技术：

2.defense meteorological satellite program/operational line scan system(dmsp/ols)卫星雷达的全球夜间灯光(nighttime light，ntl)数据集是唯一可提供1992年-2013年间夜间灯光信息的数据集，其广泛应用于经济总量估算、电力消耗估算、城市建成区面积提取、人口密度估算、自然灾害评估等领域。因此，为研究、补充2013年之前的经济状况或其他相关领域，我们需要有效的dmsp/ols夜间灯光遥感图像数据。dmsp/ols的夜间灯光数据集共有四种：1992-2013年间的原始夜间灯光数据集、稳定夜间灯光数据集、无云覆盖夜间灯光数据集和七年的辐射校准数据。其中，稳定夜间灯光数据集是最常用的数据，包含来自城市、城镇和其他具有持续照明的地点的灯光。
3.dmsp/ols的夜间灯光遥感数据存在三个主要问题，即饱和效应、绽放效应以及年际不一致，在较大程度上影响了数据的准确性和可用性。饱和效应是由于ols无法记录大于63的dn值，在照明强度超过检测上限的城市中心区域无法表示实际灯光强度的空间差异性，限制了夜间灯光数据在进行经济、人口等相关参数的定量性分析中的应用。绽放效应是ols检测到的光照区域大于光源分布的实际空间范围，从而导致对城市区域范围的高估。此外，传感器退化和缺乏机载校准导致了年际不一致的问题，不同年份和不同卫星的传感器采集到的夜间灯光数据不能直接比较，使其无法建立长时间序列。因此，在应用夜间灯光数据前需要对其校正。
4.为了克服年际不一致的问题，elvige等人首先提出了使用不变区域法对全球夜间灯光影像进行校正，使用二次方程模型来实施校正。对于饱和问题，yang等人提出了一种基于回归模型和辐射校准ntl数据的饱和校正方法，以校正年度稳定ntl数据中的饱和效应。cao等人开发了一种自调整模型，基于空间响应函数在不使用任何辅助数据的情况下校正dmsp-ols ntl数据中的光晕效应。总的来看，目前很少有方法能同时且简单有效地解决稳定ntl数据存在的三个主要问题，且目前大多数的校正模型采用的是简单的回归方程模型，都存在过度拟合的问题，而使用随机森林回归算法来建立校正模型可以解决上述问题。
5.机器算法中的随机森林(breiman,2001)是一种集成算法，具有抗过拟合能力高且精准的优点。它通过组合多个弱分类器，最终结果以投票或取均值的方式来决定，使得整体模型的结果具有较高的精确度和泛化性能。


技术实现要素：

6.本发明的主要目的在于克服现有技术的缺点与不足，提供一种基于随机森林算法的dmsp/ols数据校正方法、系统及介质。
7.为了达到上述目的，本发明采用以下技术方案：
8.本发明的一个方面，提供了一种基于随机森林算法的dmsp/ols数据校正方法，包括下述步骤：
9.读取待校正的dmsp/ols夜间灯光图像、辐射标定的dmsp/ols夜间灯光图像并进行预处理；
10.计算待校正的dmsp/ols夜间灯光图像的变异系数，并选择变异系数小于预设值的区域作为参考区域；
11.按照所述参考区域的范围裁剪待校正的dmsp/ols夜间灯光图像，将栅格数据转点；
12.计算所述辐射标定的dmsp/ols夜间灯光图像的局部moran指数，并按照所述参考区域的范围对其进行裁剪，将数据导入所述栅格数据转点的数据中，并将该导入合并后的数据导出；
13.基于随机森林回归算法对所述导入合并后的数据进行夜间灯光图像饱和效应和绽放效应的校正；
14.用经饱和效应和绽放效应校正处理后的dmsp/ols夜间灯光图像替换未校正过的dmsp/ols夜间灯光图像中dn值＞0的部分；
15.对经饱和效应和绽放效应校正、替换处理后的dmsp/ols夜间灯光图像进行多传感器同年校正；
16.对经多传感器同年校正后的dmsp/ols夜间灯光图像进行连续性校正，得到连续性校正后的dmsp/ols夜间灯光图像；
17.用连续性校正后的dmsp/ols夜间灯光图像替换未校正过的dmsp/ols夜间灯光图像dn值＞0的部分。
18.作为优选的技术方案，所述预处理具体为：
19.用研究区域的边界矢量数据对读取的待校正的dmsp/ols夜间灯光图像、辐射标定的dmsp/ols夜间灯光图像进行裁剪，将图像的空间分辨率重采样至1000m，并进行坐标变换。
20.作为优选的技术方案，所述计算待校正的dmsp/ols夜间灯光图像的变异系数，具体如下式：
21.cv＝(sd/mn)
×
100％
22.其中sd是所有年份像元值相加后求得的标准偏差，mn是所有年份像元的dn值相加后求得的平均值。
23.作为优选的技术方案，采用自然断点法对所述变异系数的计算结果进行分类，并选择变异系数小于预设值的区域作为参考区域。
24.作为优选的技术方案，计算所述辐射标定的dmsp/ols夜间灯光图像的局部moran指数，得到辐射标定的dmsp/ols夜间灯光图像的空间分布特征；
25.按照所述参考区域的范围对其进行裁剪，将处理后的dmsp/ols夜间灯光图像的值提取到所述栅格数据转点的数据中，并将得到的数据导出为表格，得到构建校正模型所需原始数据。
26.作为优选的技术方案，所述基于随机森林回归算法对所述导入合并后的数据进行夜间灯光数据饱和效应和绽放效应的校正，具体为：构建校正模型，通过组合多个弱分类
器，将最终结果以投票或取均值的方式决定。
27.作为优选的技术方案，所述多传感器同年校正，具体如下式：
[0028][0029]
式中：dn
(n，i)
为校正后第n年影像中i像元的dn值，和表示校正前2个不同传感器获得的影像中i像元的dn值。
[0030]
作为优选的技术方案，所述连续性校正，具体如下式：
[0031][0032]
其中dn
(n-1，i)
、dn
(n，i)
和dn
(n 1，i)
是第n-1年、第n年和第n 1年的稳定ntl数据的第i个像元的dn值。
[0033]
本发明的另一个方面，提供了一种基于随机森林算法的dmsp/ols数据校正系统，应用于上述的一种基于随机森林算法的dmsp/ols数据校正方法，包括预处理模块、裁剪模块、饱和/绽放效应校正模块、多传感器同年校正模块以及连续性校正模块；
[0034]
所述预处理模块用于读取待校正的dmsp/ols夜间灯光图像、辐射标定的dmsp/ols夜间灯光图像并进行预处理；
[0035]
所述裁剪模块用于计算待校正的dmsp/ols夜间灯光图像的变异系数，并选择变异系数小于预设值的区域作为参考区域，并按照所述参考区域的范围裁剪待校正的dmsp/ols夜间灯光图像，将栅格数据转点，最后计算所述辐射标定的dmsp/ols夜间灯光图像的局部moran指数，并按照所述参考区域的范围对其进行裁剪，将数据导入所述栅格数据转点的数据中，并将该导入合并后的数据导出；
[0036]
所述饱和/绽放效应校正模块用于基于随机森林回归算法对所述导入合并后的数据进行夜间灯光图像饱和效应和绽放效应的校正，并用经饱和效应和绽放效应校正处理后的dmsp/ols夜间灯光图像替换未校正过的dmsp/ols夜间灯光图像中dn值＞0的部分；
[0037]
所述多传感器同年校正模块用于对经饱和效应和绽放效应校正、替换处理后的dmsp/ols夜间灯光图像进行多传感器同年校正；
[0038]
所述连续性校正模块用于对经多传感器同年校正后的dmsp/ols夜间灯光图像进行连续性校正，得到连续性校正后的dmsp/ols夜间灯光图像，并用连续性校正后的dmsp/ols夜间灯光图像替换未校正过的dmsp/ols夜间灯光图像dn值＞0的部分。
[0039]
本发明的另一个方面，提供了一种存储介质，存储有程序，所述程序被处理器执行时，实现上述的一种基于随机森林算法的dmsp/ols数据校正方法。
[0040]
本发明与现有技术相比，具有如下优点和有益效果：
[0041]
(1)在参考区域的选择上，本发明通过计算变异系数来提取稳定像素，方法简便的同时还可以减少主观选择带来的误差。
[0042]
(2)在校正方程的选择上，本发明选择具有较高的精确度和泛化性能的随机森林算法，相比于一元二次方程方程和幂函数，随机森林能避免过拟合的问题，并具有更高的精
确度，使校正后的数据更贴合实际。
[0043]
(3)本发明可以同时且简单有效地解决稳定ntl数据存在的三个主要问题：饱和效应、绽放效应以及年际不一致。
附图说明
[0044]
图1是本发明实施例所述的基于随机森林算法的dmsp/ols数据校正方法的技术流程图；
[0045]
图2是本发明实施例所述的2006年饱和效应最大的五个城市的dmsp/ols稳定夜间灯光图像校正前后对比图；
[0046]
图3是本发明实施例所述的2000年和2013年珠三角地区dmsp/ols稳定夜间灯光图像校正前后对比图；其中图3(a)表示未校正的2000年珠三角地区dmsp/ols稳定夜间灯光图像，图3(b)表示未校正的2013年珠三角地区dmsp/ols稳定夜间灯光图像，图3(c)表示校正后的2000年珠三角地区dmsp/ols稳定夜间灯光图像，图3(d)表示校正后的2013年珠三角地区dmsp/ols稳定夜间灯光图像；
[0047]
图4是本发明实施例所述的校正前图像后计算的中国区域的亮值像元的总数量的对比图；其中图4(a)表示未校正图像的中国区域的亮值像元的总数量，图4(b)表示校正后图像的中国区域的亮值像元的总数量；
[0048]
图5是本发明实施例所述的校正前后图像dn值和省级gdp的拟合度示意图；其中图5(a)表示未校正图像dn值和省级gdp的拟合度示意图，图5(b)表示校正后图像dn值和省级gdp的拟合度示意图；
[0049]
图6是本发明实施例基于随机森林算法的dmsp/ols数据校正系统的结构示意图；
[0050]
图7是本发明实施例的存储介质的结构示意图。
具体实施方式
[0051]
为了使本技术领域的人员更好地理解本技术方案，下面将结合本技术实施例中的附图，对本技术实施例中的技术方案进行清楚、完整地描述。显然，所描述的实施例仅仅是本技术一部分实施例，而不是全部的实施例。基于本技术中的实施例，本领域技术人员在没有做出创造性劳动前提下所获得的所有其他实施例，都属于本技术保护的范围。
[0052]
实施例
[0053]
如图1所示，本实施例提供了一种基于随机森林算法的dmsp/ols数据校正方法，包括以下步骤：
[0054]
步骤一：读取待校正的dmsp/ols稳定夜间灯光图像、辐射标定的dmsp/ols稳定夜间灯光图像(本实施例中选择f12卫星1992-2013年的图像作为待校正的图像，2006年的作为辐射标定图像，均进行去掉了突发事件导致光源处理)，并进行预处理。
[0055]
进一步的，所述预处理具体为：
[0056]
1、读取待校正的dmsp/ols稳定夜间灯光的图像(步骤一中选择f12卫星1992年的图像进行处理)，用研究区域的边界矢量数据对图像进行裁剪，将图像的空间分辨率重采样至1000m，将图像进行坐标变换，转换为china_lambert_conformal_conic(以下称为图像a)。
[0057]
2、读取2006年辐射标定dmsp/ols稳定夜间灯光图像，用研究区域的边界矢量数据对图像进行裁剪，将2个图像的空间分辨率重采样至1000m，将图像进行坐标变换，转换为china_lambert_conformal_conic。
[0058]
步骤二：计算1992-2013年所有的dmsp/ols稳定夜间灯光图像数据的变异系数，如下式：
[0059]
cv＝(sd/mn)
×
100％
[0060]
其中sd是所有年份像元值相加后求得的标准偏差，mn是所有年份像元的dn值相加后求得的平均值，dn值(digital number)是遥感影像像元亮度值，这里指的是dmsp/ols图像的亮度值，dn值越大的像元，其所表示的区域夜间灯光亮度值越高。
[0061]
采用自然断点法对计算后的结果进行分类，选择变异系数小于预设值(本实施例中预设值为23.9)的区域作为参考区域。
[0062]
步骤三：按照步骤二所得参考区域的范围裁剪待校正图像，将栅格数据转点。
[0063]
进一步的，步骤三中，使用arcgis中的按掩膜提取工具用步骤二所得参考区域裁剪待校正图像，将图像使用栅格转点工具转点，将属性表的count值列命名为f101992(按照待校正的夜间灯光图像数据对应名称进行命名)
[0064]
步骤四：计算步骤一所得2006年辐射标定夜间灯光图像数据的局部moran指数，并按照步骤二所得参考区域的范围对其进行裁剪，将数据导入步骤三所得点数据中，数据导出到excel表。
[0065]
进一步的，步骤四中，利用r语言的raster工具包中的moranlocal函数计算步骤一所得2006年辐射标定夜间灯光图像数据的局部moran指数(以下称为图像b)；使用arcgis中的按掩膜提取工具，用步骤二所得参考区域裁剪图像b，对处理后的图像使用提取多值到点，提取到步骤三所得点数据中，将属性表的count值列命名为j2006。将所得栅格图像导入到arcgis中，对处理后的图像使用提取多值到点，提取到步骤三所得点数据中，将局部moran指数属性表的count值列命名为lmn2006，所得结果以csv形式储存在excel表格中。
[0066]
步骤五：将步骤四的数据用随机森林回归算法进行夜间灯光图像数据饱和效应和绽放效应的校正，具体为：
[0067]
构建校正模型，通过组合多个弱分类器，将最终结果以投票或取均值的方式决定，从而使得整体模型的结果具有较高的精确度。
[0068]
进一步的，采用随机森林回归算法工具包中的random forest函数进行夜间灯光图像数据的校正，校正代码如下：
[0069]
[1]library(raster)
[0070]
[2]library(randomforest)
[0071]
[3]r1《-raster()括号内为图像b在电脑中的存储路径
[0072]
[4]re《-raster()括号内为图像a在电脑中的存储路径
[0073]
[5]rr《-read.csv()括号内为步骤四的csv表格在电脑中的存储路径
[0074]
[6]rf《-randomforest(f2006～jf101992 lmn2006,data＝rr)校正不同图像时更改“f101992”对应步骤三的命名
[0075]
[7]re1《-stack(re,r1)
[0076]
[8]names(re1)
[0077]
[9]names(re1)《-c("jf2006","lmn2006")校正不同图像时更改“f101992”对应步骤三的命名
[0078]
[10]prf＜-predict(rel，rf，progress＝
″
txt
″
)
[0079]
writeraster(prf，filename＝
″
moranjzf101992.tif
″
，format＝
″
gtiff
″
，overwrite＝true)校正不同卫星时更改“jf101992”[0080]
步骤六：用校正过的图像替换未校正过的图像中dn值＞0的部分。
[0081]
步骤七：利用公式对步骤六处理后的dmsp/ols夜间灯光影像进行多传感器同年校正，公式如下：
[0082][0083]
(n-1＝1994，1997，1998，
…
，2003)
[0084]
式中：dn
(n，i)
为校正后第n年影像中i像元的dn值，和表示校正前2个不同传感器获得的影像中i像元的dn值。
[0085]
步骤八：利用公式对经步骤七处理后的dmsp/ols影像进行连续性校正，得到校正后的dmsp/ols影像，公式为：
[0086][0087]
(n-1＝1992，1993，
…
，2013)
[0088]
其中dn
(n-1，i)
、dn
(n，i)
和dn
(n 1，i)
是第n-1年、第n年和第n 1年的稳定ntl数据的第i个像元的dn值。
[0089]
步骤九：将步骤七处理所得的图像进行替换未校正过的图像dn值＞0的部分，避免因校正产生新的亮度像素。
[0090]
图2是本发明实施例所述的2006年饱和效应最大的五个城市的dmsp/ols稳定夜间灯光图像校正前后对比图；目测评价表明，本方法在城市级水平上应用效果好，使用本方法校正后图像可很好的减轻图像的绽放效应和饱和效应。
[0091]
图3是本发明实施例所述的2000年和2013年珠三角地区dmsp/ols稳定夜间灯光图像校正前后对比图；其中图3(a)表示未校正的2000年珠三角地区dmsp/ols稳定夜间灯光图像，图3(b)表示未校正的2013年珠三角地区dmsp/ols稳定夜间灯光图像，图3(c)表示校正后的2000年珠三角地区dmsp/ols稳定夜间灯光图像，图3(d)表示校正后的2013年珠三角地区dmsp/ols稳定夜间灯光图像；目测评价表明，本方法在饱和效应和绽放效应严重的地区上应用效果好，使用本方法校正后图像可很好的减轻图像的绽放效应和饱和效应。
[0092]
图4是本发明实施例所述的校正前图像后计算的中国区域的亮值像元的总数量(total light index，tli)的对比图；其中图4(a)表示未校正图像的中国区域的亮值像元的总数量，图4(b)表示校正后图像的中国区域的亮值像元的总数量；结果表明，本方法能减小同年不同传感器数据之间的偏差并生成稳定的时间序列，有效改善数据年际不一致的问题。
[0093]
图5是本发明实施例所述的校正前后图像dn值和省级gdp的拟合度示意图；其中图5(a)表示未校正图像dn值和省级gdp的拟合度示意图，图5(b)表示校正后图像dn值和省级gdp的拟合度示意图；结果表明，相对于未矫正图像，本方法所得图像和省级gdp的相关性有所提高，表明校正后的结果更符合实际。
[0094]
如图6所示，在本技术的另一个实施例中，提供了一种基于随机森林算法的dmsp/ols数据校正系统，该系统包括预处理模块、裁剪模块、饱和/绽放效应校正模块、多传感器同年校正模块以及连续性校正模块；
[0095]
所述预处理模块用于读取待校正的dmsp/ols夜间灯光图像、辐射标定的dmsp/ols夜间灯光图像并进行预处理；
[0096]
所述裁剪模块用于计算待校正的dmsp/ols夜间灯图像的变异系数，并选择变异系数小于预设值的区域作为参考区域，并按照所述参考区域的范围裁剪待校正的dmsp/ols夜间灯光图像，将栅格数据转点，最后计算所述辐射标定的dmsp/ols夜间灯光图像的局部moran指数，并按照所述参考区域的范围对其进行裁剪，将数据导入所述栅格数据转点的数据中，并将该导入合并后的数据导出；
[0097]
所述饱和/绽放效应校正模块用于基于随机森林回归算法对所述导入合并后的数据进行夜间灯光图像饱和效应和绽放效应的校正，并用经饱和效应和绽放效应校正处理后的dmsp/ols夜间灯光图像替换未校正过的dmsp/ols夜间灯光图像中dn值＞0的部分；
[0098]
所述多传感器同年校正模块用于对经饱和效应和绽放效应校正、替换处理后的dmsp/ols夜间灯光图像进行多传感器同年校正；
[0099]
所述连续性校正模块用于对经多传感器同年校正后的dmsp/ols夜间灯光图像进行连续性校正，得到连续性校正后的dmsp/ols夜间灯光图像，并用连续性校正后的dmsp/ols夜间灯光图像替换未校正过的dmsp/ols夜间灯光图像dn值＞0的部分。
[0100]
在此需要说明的是，上述实施例提供的系统仅以上述各功能模块的划分进行举例说明，在实际应用中，可以根据需要而将上述功能分配由不同的功能模块完成，即将内部结构划分成不同的功能模块，以完成以上描述的全部或者部分功能，该系统是应用于上述实施例的基于随机森林算法的dmsp/ols数据校正方法。
[0101]
如图7所示，在本技术的另一个实施例中，还提供了一种存储介质，存储有程序，所述程序被处理器执行时，实现上述实施例的基于随机森林算法的dmsp/ols数据校正方法，具体为：
[0102]
读取待校正的dmsp/ols夜间灯光图像、辐射标定的dmsp/ols夜间灯光图像并进行预处理；
[0103]
计算待校正的dmsp/ols夜间灯光图像的变异系数，并选择变异系数小于预设值的区域作为参考区域；
[0104]
按照所述参考区域的范围裁剪待校正的dmsp/ols夜间灯光图像，将栅格数据转点；
[0105]
计算所述辐射标定的dmsp/ols夜间灯光图像的局部moran指数，并按照所述参考区域的范围对其进行裁剪，将数据导入所述栅格数据转点的数据中，并将该导入合并后的数据导出；
[0106]
基于随机森林回归算法对所述导入合并后的数据进行夜间灯光图像饱和效应和
绽放效应的校正；
[0107]
用经饱和效应和绽放效应校正处理后的dmsp/ols夜间灯光图像替换未校正过的dmsp/ols夜间灯光图像中dn值＞0的部分；
[0108]
对经饱和效应和绽放效应校正、替换处理后的dmsp/ols夜间灯光图像进行多传感器同年校正；
[0109]
对经多传感器同年校正后的dmsp/ols夜间灯光图像进行连续性校正，得到连续性校正后的dmsp/ols夜间灯光图像；
[0110]
用连续性校正后的dmsp/ols夜间灯光图像替换未校正过的dmsp/ols夜间灯光图像dn值＞0的部分。
[0111]
应当理解，本技术的各部分可以用硬件、软件、固件或它们的组合来实现。在上述实施方式中，多个步骤或方法可以用存储在存储器中且由合适的指令执行系统执行的软件或固件来实现。例如，如果用硬件来实现，和在另一实施方式中一样，可用本领域公知的下列技术中的任一项或他们的组合来实现：具有用于对数据信号实现逻辑功能的逻辑门电路的离散逻辑电路，具有合适的组合逻辑门电路的专用集成电路，可编程门阵列(pga)，现场可编程门阵列(fpga)等。
[0112]
上述实施例为本发明较佳的实施方式，但本发明的实施方式并不受上述实施例的限制，其他的任何未背离本发明的精神实质与原理下所作的改变、修饰、替代、组合、简化，均应为等效的置换方式，都包含在本发明的保护范围之内。