CN108154486B - 基于p范数回归模型的光学遥感图像时间序列云检测方法 - Google Patents

Abstract

本发明提供一种基于p范数回归模型的光学遥感图像时间序列云检测方法，该方法包括步骤：(1)预处理：选择感兴趣的区域；(2)基于p范数的云检测回归建模；(3)通过观测值与回归值的偏差检测云。本发明方法简单，检测结果可靠，能将非常薄的云检测出来。

Description

基于p范数回归模型的光学遥感图像时间序列云检测方法

技术领域

本发明涉及云检测领域，尤其是涉及基于p范数回归模型的光学遥感图像时间序列云检测方法。

背景技术

传统的云检测算法是通过设置像素的阈值、使用像素周围的信息、使用空间信息等方法。但是这些方法过于简单且有局限性，难以广泛的使用。

很多传统的去云算法都是针对特定场景的，当脱离这个场景后，就不适用了，而且没办法检测到不同的云，比如说当云为薄云时，就有可能检测不出来。

发明内容

本发明的目的在于：针对现有技术存在的问题，提供一种基于p范数回归模型的光学遥感图像时间序列云检测方法，解决传统云检测方法有局限性，只能针对特定场景检测的问题。

本发明的发明目的通过以下技术方案来实现：

基于p范数回归模型的光学遥感图像时间序列云检测方法，该方法包括步骤：

(1)预处理：选择感兴趣的区域；

(2)基于p范数的云检测回归建模；

(3)通过观测值与回归值的偏差检测云。

步骤(1)包括：获得固定区域一定时间内如一年的图像，本专利中采用22个时刻采集得到的光学卫星Landsat-8OLI图像，用D＝{S₁，S₂...，S₂₂}表示一个区域的子图像数据集，用D_m，n＝{S₁(m，n)，S₂(m，n)...，S₂₂(m，n)}来表示一个处于(m，n)处22个时刻的像素值序列。

步骤(2)包括：

1)假定：y_i＝f(x_i)+ε_i，i＝1，2...，n，其中f(*)是回归函数，y_i是在采样时刻x_i处给出的一个像素的反射测量，n是采样的数量；

2)假设ε_i遵循p范数分布来模拟云的影响，且服从参数为期望μ，方差σ²和正实数p，则ε_i的密度为：

其中Γ(.)是Gamma函数，E(ε_i)＝μ；

3)假设期望μ等于零，在服从p范数分布的噪声ε_i的影响下，从测量值y_i中找到一种的未知的像素点的变化f(*)。

步骤(2)包括：假设f(*)具有3阶局部平滑性；假设x_i是在x附近的采样时刻，则表达式f(x_i)的泰勒展开式写为：

其中，

β₀＝f(x)，β₁＝f′(x)，β₂＝f″(x)，β₃＝f″′(x)，β＝[β₀，β₁，β₂，β₃]^T，a＝[1，(x_i-x)，(x_i-x)²，(x_i-x)³].因此，有：

又由于ε_i服从p范数分布，因此y_i的密度表示成：

并且样本的似然函数：

则β的似然估计的最大值

由下式计算得到：

通过求解这个优化问题，得到回归系数

步骤(3)包括：像素的真实值通过如下的回归方程估计：

步骤(4)包括：给定回归模型，用s表示误差或内点的标度，通过下式来检测云：

其中T是门限值，

T为2.5。

与现有技术相比，本发明方法简单，检测结果可靠，能将非常薄的云检测出来。

附图说明

图1为本发明的方法流程图；

图2在m，n像素点上的观测时间序列D_m，n和估计的回归函数f(*)；

图3为原始遥感图像；

图4为云检测结果图。

具体实施方式

下面结合附图和具体实施例对本发明进行详细说明。

实施例

如图1～图4所示，本发明方法包括：

1.预处理。

2.基于p范数的云检测回归建模。

3.云检测。

第1个步骤预处理：

我们首先在Landsat-8OLI数据集中对所有图像进行配准，以获得固定区域一年的有云图像。在22个时刻中，用D＝{S₁，S₂...，S₂₂}表示一个区域的子图像数据集,用D_m，n＝{S₁(m，n)，S₂(m，n)...，S₂₂(m，n)}来表示处于(m，n)处22个时刻的像素值序列。

第2个步骤基于p范数的云检测回归建模：

基于云及其阴影总是在时间序列上引起像素反射率突然变化的现象，我们可以将它们视为具有p范数的噪声分布。p范数分布的噪声是单峰并且对称的，通过选择适当的p值，理论模型可以产生只有少量大的非零值和大量接近于零的较小值的噪声。因此对于我们的问题，这种分布的拟合程度要比正态分布更优。

对于给定有云图像的时间序列来寻找每个像素对应的变化模型，这实际上是一种在p范数分布噪声下的鲁棒回归问题。

假定：

y_i＝f(x_i)+ε_i，i＝1，2...，n (1)

其中f(*)是回归函数，y_i是在采样时刻x_i处反射率的观测值，n是时间序列长度。在一个确定的采样时刻，当它被云污染时，测量值yi可能偏离像素真实的反射率。假设ε_i服从期望为μ、方差为σ²的p范数分布来模拟云的影响，p为正实数。那么ε_i的密度函数可以表示为：

其中Γ(.)为Gamma函数，E(ε_i)＝μ。为了简单起见，我们可以假设期望μ等于零。从而问题变为：我们要在服从p范数分布的噪声ε_i的影响下，从测量值y_i中找到一种未知的像素值的变化函数f(*)。

由于像素的反射率总是随着时间序列而缓慢而连续地变化，我们可以假设f(*)具有3阶局部平滑性。假设x_i是在x附近的采样时刻，则表达式f(x_i)的泰勒展开式可写为：

其中，

β₀＝f(x)，β₁＝f′(x)，β₂＝f″(x)，β₃＝f″′(x)，β=[β₀，β₁，β₂，β₃]^T，α＝[1，(x_i-x)，(x_i-x)²，(x_i-x)³].

因此，我们有：

又由于ε_i服从p范数分布，因此y_i的密度函数可以表示成：

并且样本的似然函数为

则β的似然估计的最大值，即

可以由下式计算得到：

通过求解这个优化问题，我们可以得到回归系数

于是，像素的真实值可以通过如下的回归方程估计得到

第3个步骤云检测：

通过观测值与回归值的偏差检测云。对每个像素点上的像素值时间序列与回归值相减得到一个误差序列，该误差可以有效地模拟为高斯噪声，而异常值则是与云对应的。给定回归模型，用s表示正常值的尺度，则云可以通过下式来检测：

其中T是门限值。通常情况下，选择T为2.5，可以将98％的高斯分布点识别为正常值，这也是在我们的实验中使用的。为了有效地区分异常值和高斯噪声，我们必须正确地估计正常值的尺度。我们使用鲁棒的中值尺度估计来获得s的值如下：

其中n是样本点的数量。

图2给出了m，n处的D_m，n序列和估计得到的回归函数曲线结果，其中鲁棒回归函数用水平方向的曲线表示，有云观测值用”*”表示。注意，时刻13、14处像素值分别为79、76，这两个点是异常值，它们对应于厚云；而在时刻17、21处像素值分别为39、38，这两个异常值对应于薄云。图3为有云图像，图4为用本专利方法检测得到的云图像。

以上所述仅为本发明的较佳实施例而已，并不用以限制本发明，应当指出的是，凡在本发明的精神和原则之内所作的任何修改、等同替换和改进等，均应包含在本发明的保护范围之内。

Claims (8)

1.基于p范数回归模型的光学遥感图像时间序列云检测方法，其特征在于，该方法包括步骤：

(1)预处理：选择感兴趣的区域；

(2)基于p范数的云检测回归建模；步骤(2)包括：

21)假定：y_i＝f(x_i)+ε_i，i＝1，2…，n，其中f(*)是回归函数，y_i是在采样时刻x_i处给出的一个像素的反射测量，n是采样的数量；

22)假设ε_i遵循p范数分布来模拟云的影响，且服从参数为期望μ，方差σ²和正实数p，则ε_i的密度为：

其中Γ(.)是Gamma函数，E(ε_i)＝μ；

23)假设期望μ等于零，在服从p范数分布的噪声ε_i的影响下，从测量值y_i中找到一种未知的像素点的变化f(*)；

(3)通过观测值与回归值的偏差检测云。

2.根据权利要求1所述的基于p范数回归模型的光学遥感图像时间序列云检测方法，其特征在于，步骤(1)包括：获得固定区域一定时间内的图像，在22个时刻中，用D＝{S₁，S₂...，S₂₂}表示一个区域的子图像数据集，用D_m，n＝{S₁(m，n)，S₂(m，n)…，S₂₂(m，n)}来表示一个处于(m，n)处的像素的22个时刻值。

3.根据权利要求2所述的基于p范数回归模型的光学遥感图像时间序列云检测方法，其特征在于，从Landsat-8 OLI数据中获得图像。

4.根据权利要求2所述的基于p范数回归模型的光学遥感图像时间序列云检测方法，其特征在于，所述一定时间为一年。

5.根据权利要求1所述的基于p范数回归模型的光学遥感图像时间序列云检测方法，其特征在于，步骤(2)包括：假设f(*)具有3阶局部平滑性；假设x_i是在x附近的采样时刻，则表达式f(x_i)的泰勒展开式写为：

其中，β₀＝f(x)，β₁＝f′(x)，β₂＝f″(x)，β₃＝f″′(x)，β＝[β₀，β₁，β₂，β₃]^T，a＝[1，(x_i-x)，(x_i-x)²，(x_i-x)³].因此，有：

又由于ε_i服从p范数分布，因此y_i的密度表示成：

并且样本的似然函数：

则β的似然估计的最大值

由下式计算得到：

通过求解这个优化问题，得到回归系数

6.根据权利要求5所述的基于p范数回归模型的光学遥感图像时间序列云检测方法，其特征在于，步骤(3)包括：像素的真实值通过如下的回归方程估计：

7.根据权利要求5所述的基于p范数回归模型的光学遥感图像时间序列云检测方法，其特征在于，步骤(3)包括：给定回归模型，用s表示正常值的尺度，通过下式来检测云：

其中T是门限值，

8.根据权利要求7所述的基于p范数回归模型的光学遥感图像时间序列云检测方法，其特征在于，T为2.5。

Images