=Paper=
{{Paper
|id=Vol-2391/paper7
|storemode=property
|title=Applying doubly stochastic filters to evaluate the dynamics of object sizes on satellite image sequences
|pdfUrl=https://ceur-ws.org/Vol-2391/paper7.pdf
|volume=Vol-2391
|authors=Vitalii Dementev,Dmitrii Kondratev
}}
==Applying doubly stochastic filters to evaluate the dynamics of object sizes on satellite image sequences ==
https://ceur-ws.org/Vol-2391/paper7.pdf
Applying doubly stochastic filters to evaluate the dynamics of
object sizes on satellite image sequences
V E Dementyev1, D S Kondratyev1
1
Ulyanovsk State Technical University, ul. Severny Venets, 32, Ulyanovsk, Russia, 432027
e-mail: dve@ulntc.ru
Abstract. One of the important tasks facing the regional authorities is to monitor the condition
of roads and power lines. In the Ulyanovsk region more than 8 thousand km of power lines and
more than 9 thousand km of roads (including rural). A significant part of these facilities is
located outside the settlements in places with medium and low availability. In many such
places there is a problem of uncontrolled forest overgrowth. This work is devoted to solving
the problem of automated satellite monitoring of such areas. For this purpose, it is proposed to
use a modified convolutional neural network that processes time sequences of multispectral
satellite images and allows to allocate territories occupied by forest and undergrowth with high
accuracy. This approach allows us to assess the dynamics of overgrowth of the territory and
perform the appropriate forecast with sufficient accuracy for practice.
1. Introduction
One of the important tasks of the satellite image processing is its thematic mapping, i.e. division of
image into identifiable areas containing the similar visual, correlation or texture characteristics of
pixels. The use of standard segmentation algorithms [1-3] for thematic mapping of satellite images
usually leads to significant errors caused by two reasons. Firstly, these algorithms are largely
incapable of taking into account the multi-zonal nature of remote sensing (RS), so each satellite image
contains the results of the Earth's surface registration in different spectral ranges. Some works [4-6]
suggest the possibility of processing hyperspectral images. Thus, the authors based their theory on
criterion of uniformity for reception of connected areas of such hyperspectral image [4], modification
and generalization of algorithm K-means [5] and use of physical properties of a satellite data [6].
Secondly, the existing approaches unable to use data on the observed territory received at previous
points of time for image segmentation. But using such data can significantly improve the quality of
processing at the expense of a fundamentally larger amount of information, but it is fraught with
difficulties of aggregation.
It is possible to overcome the mentioned disadvantages by using neural network procedures of
segmentation and classification of multidimensional data. In work [9] the variant of the U-NET neural
network with full-connected layers (FCN) modification is presented. Herewith the input layer of the
network consisting from spectral layers of the multizonal image has been extended by three auxiliary
2d halftone images obtained from the original using NDVI, EVI, SAVI transformations and two 2d
arrays , representing the segmentation results at the previous point of time and one year ago. The use
of two such reference markings allows one's to reduce the error of classification in case of rapid
changes of the terrain due to the change of year time and in case of marking array absence at the
V International Conference on "Information Technology and Nanotechnology" (ITNT-2019)
Image Processing and Earth Remote Sensing
V E Dementyev, D S Kondratyev
previous point of time in connection, for example, with cloudiness. The quality analysis of such an
algorithm shows the processing comparable accuracy to the qualified operator results.
It should be noted that the results of the satellite images thematic processing provides information
related to the peculiarities of the Earth's surface in the previous and current moments of time.
However, the necessary element of the RS processing system is the tool that allow to forecast the state
of certain objects and form recommendations for responsible persons. For example, let's consider the
task related to the forecasting about the dangerous convergence possibility of various natural and
technogenic objects. An example of such convergence may be the gradual growth of the forest area
along roads or power lines or landslide processes leading to the destruction of various infrastructures.
2. Satellite Image Processing
Let us formulate the given problem as follows. Let there be an aggregate of segments describing a
certain extended object ๐๐จ (a road, an electric power transmission line etc.). The given aggregate is
usually described by a vector object having geographical reference in absolute coordinates. Let us
separate on this extended object a set of points ๐๐จ๐ข = (๐ฑ๐๐ข , ๐ฒ๐จ๐ข ) having the distance between them
equal to โ๐จ . Let us assume that at a certain time instant ๐ญ next to the object there lies a certain extended
domain ๐๐ ๐ญ defined by a set of points each corresponding to a pixel in the original halftone image. For
each point ๐๐จ๐ข let us construct a perpendicular to the segment [๐๐จ๐ข , ๐๐จ๐ข+๐ ]. Let us find the point
๐๐๐ข ๐ญ = ๏ฟฝ๐ณ๐๐ฑ๐ข
๐ญ ๐ญ
, ๐ณ๐๐ฒ๐ข ๏ฟฝ of intersection of this perpendicular and the domain ๐๐ ๐ญ. Obviously the set of
points {๐๐๐ข ๐ญ , ๐ข = ๐, . . , ๐๐จ } describe a conditional boundary of the domain ๐๐ ๐ญ as viewed from the
extended object. It enables to use estimates of the points coordinates {๐๐๐ข ๐ญ , ๐ข = ๐, . . , ๐๐จ } obtained at
different times ๐ญ = ๐, . . , ๐ as a source of information on the domain ๐๐ ๐ญ dynamics for the purpose of
forming a prediction about its boundaries at future times ๐ญ > ๐ป as well. (Figure 1).
Figure 1. Geometric interpretation of the distance estimate between the stationary object and the
movable area.
In view of errors arising at satellite images recording and processing we get the following relations
(Eq. 1), in which ๐๐๐ฌ๐๐ and ๐๐๐ฌ๐๐ โ white noise samples with zero mean and variance ๐๐๐ .
๐ก ๐ก ๐ก ๐ก ๐ก ๐ก
๐ง๐ธ๐ฅ๐ = ๐ฅ๐ธ๐ + ๐๐ธ๐ฅ๐ , ๐ง๐ธ๐ฆ๐ = ๐ฆ๐ธ๐ + ๐๐ธ๐ฆ๐ , ๐ = 1, . . , ๐๐ , ๐ก = 1, . . , ๐, (1)
๐ ๐ ๐
where ๐๐ฌ๐๐ and ๐๐ฌ๐๐ โ white noise samples with zero mean and variance ๐๐ .
Direct measurements based on the results of satellite material image-type related mapping show
that MSE ๏ฟฝ๐๐๐ง is approximately equal to ๐. ๐๐๐ฑ๐ฒ , where ๐๐ฑ๐ฒ โ resolution of the original images.
Let us suppose that the boundary of the domain ๐ฎ๐น ๐ can move non-uniformly. Thus, for example,
the area of a precipice or a ravine can increase by tens of centimeters per each year and at a certain
moment its rate might increase exponentially. Then we will use doubly stochastic (DS) model [3,10]
to describe unknown coordinates in the following form (Eq. 2). In this case ๐ซ๐๐ฑ , ๐ซ๐๐ฒ โ scalar parameters
determining change potential for the accelerations ๐๐ญ๐๐ฑ๐ข and ๐๐ญ๐๐ฒ๐ข ; ๐๐ญ๐๐ฑ๐ข , ๐๐ญ๐๐ฒ๐ข โ independent normal
random variables with zero mean and variance ๐๐๐ :
๐ก ๐กโ1 ๐กโ2 ๐ก ๐กโ1 ๐กโ2
๐ฅ๐ธ๐ = 2๐ฅ๐ธ๐ โ ๐ฅ๐ธ๐ + ๐๐ธ๐ฅ๐ ๏ฟฝ๐ฅ๐ธ๐ โ ๐ฅ๐ธ๐ ๏ฟฝ,
V International Conference on "Information Technology and Nanotechnology" (ITNT-2019) 55
Image Processing and Earth Remote Sensing
V E Dementyev, D S Kondratyev
๐ก ๐กโ1 ๐กโ2 ๐ก ๐กโ1 ๐กโ2
๐ฆ๐ธ๐ = 2๐ฆ๐ธ๐ โ ๐ฆ๐ธ๐ + ๐๐ธ๐ฆ๐ ๏ฟฝ๐ฆ๐ธ๐ โ ๐ฆ๐ธ๐ ๏ฟฝ, (2)
๐ก ๐กโ1 ๐ก ๐ก ๐กโ1 ๐ก
๐๐ธ๐ฅ๐ = ๐๐๐ฅ ๐๐ธ๐ฅ๐ + ๐๐๐ฅ๐ , ๐๐ธ๐ฆ๐ = ๐๐๐ฆ ๐๐ธ๐ฆ๐ + ๐๐๐ฆ๐ ,
The model (2) can be rewritten in the form:
๐๏ฟฝ๐ธ๐
๐ก ๐ก
= ๐๐ธ๐ฅ๐ (๐๏ฟฝ๐ธ๐
๐กโ1
) + ๐๐ฅ๐ฬ
๐ก , ๐๏ฟฝ๐ธ๐
๐ก ๐ก
= ๐๐ธ๐ฆ๐ (๐๏ฟฝ๐ธ๐
๐กโ1 ฬ
๐ก ,
) + ๐๐ฆ๐
where
๐กโ1 ๐ก
๐ฅ๐ธ๐ ๐ฃ๐ธ๐ฅ๐ 0
๐ก ๏ฟฝ ๐กโ1
๐๐ธ๐ฅ๐ ๏ฟฝ๐๐ธ๐ ๏ฟฝ = ๏ฟฝ 0 ๐กโ1 ๐ก
๐ฃ๐ธ๐ฅ๐ (1 + ๐๐ธ๐ฅ๐ ) 0 ๏ฟฝ,
๐กโ1
0 0 ๐๐ธ๐ฅ๐ ๐๐๐ฅ
๐กโ1 ๐ก
๐ฆ๐ธ๐ ๐ฃ๐ธ๐ฆ๐ 0
๐กโ1 ๐ก
๐ก
๐๐ธ๐ฆ๐ ๏ฟฝ๐๏ฟฝ๐ธ๐
๐กโ1
๏ฟฝ=๏ฟฝ 0 ๐ฃ๐ธ๐ฆ๐ (1 + ๐๐ธ๐ฆ๐ ) 0 ๏ฟฝ,
๐กโ1
0 0 ๐๐ธ๐ฆ๐ ๐๐๐ฆ
0 0
ฬ
๐ก = ๏ฟฝ 0 ๏ฟฝ, ๐๐ฅ๐
๐๐ฅ๐ ฬ
๐ก = ๏ฟฝ 0 ๏ฟฝ.
๐ก ๐ก
๐๐๐ฅ๐ ๐๐๐ฆ๐
The assigned notations enable to apply DS nonlinear filtering to observations and for construction
of predictions of behaviour of the region ๐น๐ in reference to the object๐ป๐ . In doing so let us introduce
๐ฟ ๏ฟฝ ๐โ๐
๏ฟฝ ๐ะญ๐ฌ๐ = ๐๐๐ฌ๐๐ (๐ฟ ๏ฟฝ๐ ๐ ๏ฟฝ ๐โ๐
๐ฌ๐ ) and ๐๐ฌ๐ = ๐๐ฌ๐๐ (๐ฟ๐ฌ๐ )โ extrapolated predictions for the point ๐๐๐ข coordinates at
๐ญโ๐
the time point ๐ญ based on the preceding observations ๐ณ๐๐ฑ๐ข ๐ญโ๐
and ๐ณ๐๐ฒ๐ข . Denote by ๐๐ฑ๐ข ๐ญโ๐ , ๐๐ฒ๐ข ๐ญโ๐ โ filtering
๐ ๐
error covariance matrices at the time point (๐ญ โ ๐), ๐๐ฑ๐๐ข ๐ญ = ๐{๐๏ฟฝ๐ญ๐ฑ๐ข ๐๏ฟฝ๐ญ๐ฑ๐ข }, ๐๐ฒ๐ ๐ข ๐ญ = ๐{๐๏ฟฝ๐ญ๐ฒ๐ข ๏ฟฝ๐๐ญ๐ฒ๐ข } โ diagonal
covariance matrices for random increments๐๏ฟฝ๐๐๐ . Then error covariance matrices for such extrapolation
have the following form:
๐
๐ ๐ก = ๐ ๏ฟฝ๏ฟฝ๐๏ฟฝ ๐ก โ ๐๏ฟฝ ๐ก ๏ฟฝ๏ฟฝ๐๏ฟฝ ๐ก โ ๐๏ฟฝ ๐ก ๏ฟฝ ๏ฟฝ = ๐๐ก โฒ(๐๏ฟฝ ๐กโ1 )๐ ๐กโ1 ๐๐ก โฒ(๐๏ฟฝ ๐กโ1 )๐ + ๐ ๐ก ,
ะญ๐ฅ๐ ะญ๐ธ๐ ะญ๐ธ๐ ะญ๐ธ๐ ะญ๐ธ๐ ๐ธ๐ฅ๐ ๐ธ๐ ๐ฅ ๐ธ๐ฅ๐ ๐ธ๐ ๐ฅ๐ ๐
๐
๐ะญ๐ฆ๐ ๐ก = ๐ ๏ฟฝ๏ฟฝ๐๏ฟฝะญ๐ธ๐
๐ก
โ ๐๏ฟฝะญ๐ธ๐
๐ก
๏ฟฝ๏ฟฝ๐๏ฟฝะญ๐ธ๐
๐ก
โ ๐๏ฟฝะญ๐ธ๐
๐ก ๐ก
๏ฟฝ ๏ฟฝ = ๐๐ธ๐ฆ๐ โฒ(๐๏ฟฝ๐ธ๐
๐กโ1
)๐๐ฆ๐ ๐กโ1 ๐๐ธ๐ฆ๐
๐ก
โฒ(๐๏ฟฝ๐ธ๐ ) + ๐๐ฆ๐ ๐ ๐ก ,
๐กโ1 ๐
If we denote by ๏ฟฝ ๐๐ะญ๐ฌ๐ ,๐
๏ฟฝ๐ะญ๐ฌ๐ - the first elements of the vectors ๐ฟ ๏ฟฝ ๐ะญ๐ฌ๐ and ๐
๏ฟฝ ๐ะญ๐ฌ๐ ; ๐ฉ๐๐ ๐ =
โ๐ โ๐
๐ทะญ๐๐ ๐ ๐ช๐ ๐ป ๐ซ๐๐๐ ; ๐ฉ๐๐ ๐ = ๐ทะญ๐๐ ๐ ๐ช๐ ๐ป ๐ซ๐๐๐ ; ๐ซ๐๐๐ = ๐ช๐ ๐ทะญ๐๐ ๐ ๐ช๐ ๐ป + ๐๐๐ ; ๐ซ๐๐๐ = ๐ช๐ ๐ทะญ๐๐ ๐ ๐ช๐ ๐ป + ๐๐๐ ;
ะก๐ = ะก๐ = (๐ ๐ ๐), then we can write the following relations (Eq. 27) for DS coordinate filters:
๐๏ฟฝ ๐ก = ๐๏ฟฝ ๐ก + ๐ต ๐ก ๏ฟฝ๐ง ๐ก โ ๐ฅ๏ฟฝ ๐ก ๏ฟฝ, ๐๏ฟฝ ๐ก = ๐๏ฟฝ ๐ก + ๐ต ๐ก ๏ฟฝ๐ง ๐ก โ ๐ฆ๏ฟฝ ๐ก ๏ฟฝ.
๐ธ๐ ะญ๐ธ๐ ๐ฅ๐ ๐ธ๐ฅ๐ ะญ๐ธ๐ ๐ธ๐ ะญ๐ธ๐ ๐ฆ๐ ๐ธ๐ฆ๐ (3)
ะญ๐ธ๐
The filtering error variance at each step is determined by the matrices ๐ท๐๐ ๐ = (๐ฌ โ ๐ฉ๐๐ ๐ )๐ทะญ๐๐ ๐ ,
๐ท๐๐ ๐ = (๐ฌ โ ๐ฉ๐๐ ๐ )๐ทะญ๐๐ ๐ .
Special attention must be given to the fact that due to the necessity for the distance from the
extended object ๐ป๐ to the domain ๐น๐ to be surveyed, it is possible to simplify the boundaries
๐ญ
coordinates filtering process ๐ ๐ by processing only one coordinate, namely, the distance ๐ฑ๐ซ๐ข from the
point ๐๐จ๐ข belonging to ๐๐จ up to the boundary ๐ ๐ at the time point ๐ญ. An aggregate of similar
๐ญ
observations ๐ณ๐๐ฑ ๐ซ๐ข
can be processed by a technique identical to the above-described one.
As an illustration of such a technique in figures below series satellite images fragments for the
forest tract in Cherdakly district of the Ulyanovsk region for the period 2001-2017 years (figure 2) and
Milanovsky opencast colliery on riverbank of the Volga in the northern part of the city of Ulyanovsk
for the period 2013-2017 years (figure 3) are presented. Here for convenience of color image
perception and its recovery overlapping of visible spectral bands and superposition of the segmented
image fragment and normals to the object to be monitored is carried out. In the first case the number of
multispectral images to be processed amounted to 42 snapshots, in the second case - 32 snapshots. The
minimal time interval for satellite information production amounts to 14 days.
V International Conference on "Information Technology and Nanotechnology" (ITNT-2019) 56
Image Processing and Earth Remote Sensing
V E Dementyev, D S Kondratyev
July 2003 June 2011 July 2017
Figure 2.Forest tract satellite images fragments and the results of these images processing.
August 2015 April 2016 May 2016
Figure 3. Milanovsky opencast colliery satellite images fragments and the results of these images
processing.
The above-mentioned groups of multispectral images were separated into a training and an
operating samples. The training samples were used to specify the filtering parameters, in particular, to
estimate the parameters ๐ซ๐๐ฑ and ๐๐๐ . The operating part of the sample was processed by three
๐ญ+๐
algorithms enabling to carry out the distance prediction ๐ฑ๐ซ๐ข basing on the preceding observations.
๐ญ+๐ ๐ญ ๐ญ
The first algorithm (I) involves constructing a simple prediction ๐ฑ๏ฟฝ๐ซ๐ข = ๐๐ณ๐ซ๐ข โ ๐ณ๐ซ๐ข wherein only the
variable ๐ฑ๏ฟฝ๐ซ๐ข change speed is taken into account for the time interval (๐ญ โ ๐, ๐ญ). The second algorithm
๐ญ
(II) assumes linear Kalman filtering of the observations ๐ณ๐ซ๐ข and construction of the extrapolated
๐ญ+๐
๏ฟฝ
prediction ๐ฑะญ๐ซ๐ข based on the results of the processing. The third algorithm (III) is in the above-
๐ญ
described doubly stochastic filtering of an aggregate of the observations ๐ณ๐ซ๐ข and construction of the
๏ฟฝ ๐ญ+๐
vector ๐ ะญ๐๐ข . In the table 1 below the values of prediction average errors depending on the object kind
are presented.
On average the DS filter provides prediction accuracy 6% higher than in case of using Kalman
filter and 58% higher than in case of simple linear predictions. It enables to estimate coordinates and
rate change dynamics for boundaries of the domain to be monitored by using the DS filter. It is
V International Conference on "Information Technology and Nanotechnology" (ITNT-2019) 57
Image Processing and Earth Remote Sensing
V E Dementyev, D S Kondratyev
essential that the DS filter enables to quicker respond to abrupt rate change of the processes
determining the object behaviour. As an illustration we provide the estimates behaviour for the
๏ฟฝ ๐๐
distance from opencast colliery to one of the points to be monitored (figure 4a) and the parameter ๐
estimate (figure 4b).
Table 1. Prediction average errors.
Average error Average error Average error
for the for the for the algorithm
algorithm I algorithm II III
The forest tract snapshot. October 2014
6.7 m 2.7 m 2.6 m
The forest tract snapshot. May 2015
10.7 m 3.9 m 3.7 m
The forest tract snapshot. June 2016
6.2 m 3.6 m 3.3 m
Milanovsky opencast colliery snapshot. May
2014 7.1 m 3.8 m 3.6 m
Milanovsky opencast colliery snapshot. May
2015 7.3 m 3.9 m 3.8 m
Milanovsky opencast colliery snapshot. May
2016 6.9 m 3.9 m 3.7 m
Milanovsky opencast colliery snapshot. April
12.4 m 8.9 m 7.8 m
2016. Beginning of the avalanche processes.
Milanovsky opencast colliery snapshot. April 30.7 m 32.8 m 12.6 m
2016. Continuation of the avalanche processes.
Milanovsky opencast colliery snapshot. May 20.3 m 18.1 m 17.3 m
2016. Cessation of the avalanche processes.
Figure 4. Dependence of the filtering results versus survey time.
3. Conclusion
Direct analysis of the given results in comparison with the data of objective monitoring (solid line)
indicates superiority of the DS filter over conventional linear Kalman filter in filtering accuracy. As it
takes place, this superiority makes itself evident in the most distinct manner in case of abrupt change
of the rock collapsing process (and the corresponding reduction of the distance between the opencast
V International Conference on "Information Technology and Nanotechnology" (ITNT-2019) 58
Image Processing and Earth Remote Sensing
V E Dementyev, D S Kondratyev
colliery and the point to be monitored). This change corresponds to a significant change of the
parameter ๐๐๐ estimate which enables to register considerable changes in the opencast colliery domain
state basing only on this estimate dependence nature versus survey time.
4. References
[1] Vasiliev K K 2016 Optimum processing of signals in discrete time (Moscow: Radiotehnika)
[2] Krasheninnikov V R and Vasil'ev K K 2018 Multidimensional image models and processing
Intelligent Systems Reference Library 135 11-64
[3] Vasiliev K K and Dementyev V E 2017 Presentation and processing of satellite multi-zone
images (Ulyanovsk)
[4] Fursov V A, Bibikov S A and Baida O A 2014 Thematic classification of hyperspectral images
by conjugacy index Computer Optics 38(1) 154-158
[5] Zimichev E A, Kazansky N L and Serafimovich P G 2014 Spatial classification of hyperspectral
images using the k-means ++ clustering method Computer Optics 38(2) 281-286
[6] Blokhinov Yu B, Gorbachev V A, Rakutin Yu O and Nikitin A D 2018 Development of an
algorithm for semantic segmentation of real-time aerial photographs Computer Optics 42(1)
141-148 DOI: 10.18287/2412-6179-2018-42-1-141-148
[7] Belim S V and Larionov S V 2016 Image segmentation algorithm based on searching
communities on graphs Computer Optics 40(6) 904-910 DOI: 10.18287/2412-6179-2016-40-6-
904-910
[8] Dementyev V E, Kondratyev D S and Frenkel A G 2018 Application of convolutional neural
networks for thematic mapping of sequences of satellite multispectral images (Moscow:
Radiotehnika)
[9] Dementyev V E 2017 The use of double-stochastic autoregressive models for processing
satellite images (Moscow: Radiotehnika)
[10] Andriyanov N A and Dement'ev V E 2018 Application of mixed models of random fields for
the segmentation of satellite images CEUR Workshop Proceedings 2210 219-226
V International Conference on "Information Technology and Nanotechnology" (ITNT-2019) 59