{'fraction_non_alphanumeric': 0.028198494142271577, 'fraction_numerical': 0.01107590336928399, 'mean_word_length': 4.552865213882163, 'pattern_counts': {'":': 0, '<': 9, '<?xml version=': 0, '>': 5, 'https://': 0, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 4, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
Fast automatic deforestation detectors and their extensions for other spatial objects

This paper is devoted to the problem of detection of forest and non-forest areas on Earth images. We propose two statistical methods to tackle this problem: one based on multiple hypothesis testing with parametric distribution families, another oneon non-parametric tests. The parametric approach is novel in the literature and relevant to a larger class of problems– detection of natural objects, as well as anomaly detection. We develop mathematical background for each of the two methods, build self-sufficient detection algorithms using them and discuss numerical aspects of their implementation. We also compare our algorithms with those from standard machine learning using satellite data.


Introduction
Within the last decade, there have been some attempts to create remote deforestation detectors. Contemporary literature on the topic seems to be geared towards machine learning approaches. Just to name a few examples we refer to [10] that employed deep learning for deforestation detection in the Brazilian Amazon and [15] that presented k-nearest neighbour classification to track deforestation in Peru. Bayesian approach has also seen use, e.g. in [12] on time series data from Fiji. The present work aims to create methods that not only show excellent performance, but are also, unlike machine learning methods, fast, simple and interpretable. To achieve these characteristics we employ sample means and covariances and empirical characteristic functions (ECFs) of colour intensities of image samples. The idea behind it is that forests seen from above, like many other natural objects, show rather simple homogeneous structure. Moreover, in many cases natural objects are distinguishable due to their colour features rather than their spatial patterns. Thus, operating on distributions and their ECFs and completely erasing all the spatial dependencies seems to be reasonable.
The general statistical setting of the problem is as follows. In visible spectrum images of geographical objects are represented as 3 matrices. Each of them consists of intensities of one of the red, green or blue colours. We employ supervised learning where tags (forest / non-forest) are assigned to entries of the training data and learn features of the distributions of the matrices in those entries. During classification, test images are split into much smaller square sub images for the purpose of obtaining better resolution, and these sub-images are tested against reference images that are known to contain forest. We intend to investigate two different approaches to the problem: non-parametric multiple hypothesis testing and parametric one. In the case of parametric multiple hypothesis testing, we assume parametric distributions for the colour intensities of pixels and use training data to estimate the parameters of said distributions. New data is then tested against this distribution. In the case of the non-parametric setting, we do not assume any known distribution for the colour intensities but instead test if the multivariate distributions of the colour intensities of the training data and new data have equal means and covariance matrices.
We build both the parametric and the non-parametric algorithms in such a way that they work on images in visible spectrum, in particular, RGB composites. Thus, the algorithms do not require any certain pre-processing technique for the data and can be applied to a broad range of sources of images of Earth. There are at least two ways of obtaining large datasets of Earth images: satellite data obtained through open access Earth observation missions, e.g. Landsat [13] and Sentinel-2 [4] and aerial photo imaging. In this work, Sentinel-2 images are used because they provide good spatial and temporal resolution. Sentinel-2 mission provides data from visible, near infrared and short-wave infrared sensors in 13 corresponding spectral bands. Although we use RGB composites obtained via a post-processing procedure of Sentinel-2 original tiles, our methods could work on all of the data in those 13 bands.

A practical example
In order to give more intuition of what task we tackle, we present a real-life example of deforestation detection. Figure 1a shows a Sentinel-2 image of a piece of land near the city of Arkhangelsk in Russia suffering from logging. Figure 1b is a black-and-white image representing the detection result. White represents forest and black is used for all types of terrains not covered by it.  One can see that the detector highlighted the areas with logging traces (black rectangles) as well as natural objects: lakes, bare ground and pieces of rivers.
The rest of the paper is organized as follows. Section 2 briefly mentions the key steps in data preparation and discusses what features forest data possess. Sections 3 and 4 describe correspondingly novel non-parametric and parametric methods. Section 5 discusses theoretical properties of the two methods, while the last section compares the two with well-known methods from machine learning using Satellite data.
Our main contribution contains the application of the known theory of Hotelling T 2 statistic to the problem of deforestation detection, which occupies Section 3. There, we prove Theorem 1 and using it, build classification and training algorithms in 3.2 and 3.3. Another major contribution is Section 4 (excluding 4.1), where we develop a new setting of hypothesis testing for stable models based on empirical characteristic functions in 4.2 and use this setting to build algorithms for detection of forest in 4.3 and 4.4.

Data and exploratory analysis
To obtain data in visible spectrum we use a simple preprocessing procedure, collecting red, green and blue bands from Sentinel-2 Level-1C products and combining them into one GeoTIFF file. It is then handled by GDAL [5], which creates three numeric matrices (one per each of the red, green and blue channels) with values between 0 and 1. The procedure is described in appendix A in more detail.
Initial inspection of the empirical density of the forest data gave us some additional ideas for potentially suitable distributions. We considered distributions including, but not limited to, stable, Normal, Log-Normal, Gamma, Cauchy, Weibull and Levy. We estimated parameters for each of the colour intensity distributions and calculated the root-mean-square-error (RMSE) wheref (·) is the probability density function using the estimated parameters, f * (·) is the empirical density and x i is colour intensity of data point i = 1, . . . , n. To choose the best fitting distributions we then plotted the densities of the three distributions with the lowest RMSE's together with the empirical densities and compared. These distributions appeared to be Normal, Gamma and Stable consistently throughout the set of images. An example of a forest image with such comparison can be seen in Figure 2, where we can observe the three best fitting distributions for a forest image with all parameters estimated. Further from the plots we see that the Normal distributions do not seem to fit very well while the Gamma and Stable distributions are quite close, with Stable edging out Gamma especially for the intensity of the colour red.
Since the focus in this work is to be able to detect if an area has forest or not, we describe detection of presence or lack of this type of natural objects. However, our exploratory data analysis shows us that the framework developed could be used to detect other natural objects from satellite images as well. In Figures 11 and 12 in appendix B, we can see promising results of using the stable distribution to describe the colour intensities of images of mountains and sea. When it comes to man-made objects such as, for example a city, our approach does not work, and we have to at least assume a multimodal mixture of some kind to be appropriate since satellite images contain distinctly different parts, i.e., an orange roof, green grass and grey pavement.

Detection with a non-parametric method
In this section, the non-parametric method for forest detection is developed. Later, in Section 4, the ideas presented here are used to create the theory and implementation for the parametric one.

Comparing images with Hotelling T-square statistic
The Mahalanobis distance D was originally introduced as a distance between a given multivariate point p and a distribution. For a point p and a distribution with mean vector µ and covariance matrix Σ; the distance is given by Where, in practice, the mean and covariance matrix of the distribution would be estimated by the sample mean and covariance matrix. Further extensions can be made to calculate the distance between two multivariate populations with WhereX 1 andX 2 are the sample mean vectors of the two populations. The variableΣ −1 is the inverse of the pooled sample covariance matrix, given byΣ = (n 1Σ1 + n 2Σ2 )/(n − 2) for samples with respective sample sizes n 1 , n 2 , total sample size n = n 1 + n 2 and sample covariance matriceŝ Σ 1 ,Σ 2 . An interested reader can find more on the topic in [3].
The Mahalanobis distance is often used in cases where there is correlation present in the considered multivariate distribution. As the reader might have noticed, the Mahalanobis distance expression is closely related to the statistical concept of standardization, i.e, subtracting the mean and dividing by the standard deviation. The idea is to decorrelate and scale each dimensions variance to be 1 and then calculate the Euclidean distance. Theorem 1. Let X 1 and X 2 be independent samples of sizes n 1 × p and n 2 × p composed from i.i.d p-dimensional random vectors from distributions with mean µ and covariance matrix Σ. Then, for D from (1) we have that n1n2 n1+n2 D 2 d − → χ 2 (p).
Proof. By the multidimensional central limit theorem After introducing A := (γΣ) 1/2 andÂ := (γΣ) 1/2 we can write Taking into account (3) and thatÂ where the later is a matrix of constants, by a property of weak convergence Then, we introduce a new quantity The multivariate continuous mapping theorem applied to (4) and using the following chain of transformations we obtain Remark 1. Note, that X 1 and X 2 do not have to be drawn from the same distribution. If they are drawn from different ones, but possessing the same mean and covariance matrix, the result of the theorem still holds.

Hypothesis testing based on Mahalanobis distances between images
Suppose we have an RGB image of forest ( Figure 3). It consists of three matrices corresponding to each of the red, green and blue channels. The matrices contain values of colour intensities of pixels in the range [0, 1]. Then, each matrix is split into column-vectors that are subsequently stack together to form one column-vector of the length p 2 . The three obtained large vectors are bound to form a (p 2 ) × 3 matrix X. Since the same operation is applied to every matrix, the triplets of pixel intensities are kept unchanged, so it is possible to analyse dependencies between colours. On the contrary, matrix X does not contain any information about spatial locations of the pixels. From X, a three-dimensional mean vector and a 3 × 3 covariance matrix are calculated.  The assumptions of the theorem hold for data of type of X: since all entries are bounded between 0 and 1, the means and covariances will always be finite. This naturally leads us to a hypothesis testing framework. Having two images-one picture that is known to depict forest and one test image, we can use the theorem to accept or reject the null hypothesis that the test image is forest. This is done by first converting RGB composites into matrices X 1 and X 2 , obtaining their sizes, means and covariance matrices, and second, computing the Mahalanobis distance between X 1 and X 2 , and checking whether or not T := n1n2 n1+n2 D 2 is lower than a certain threshold.

T min T min < Thres
Non-forest

Forest
Training forest images Test image False True Figure 4: Classification in the non-parametric model A scheme for multiple hypothesis testing framework is depicted on Figure 4. When it comes to classification of a certain image, rescaled Mahalanobis distances T 1 , T 2 , . . . between the image and every forest picture in the training set are computed. Out of all T 's, the minimal one is chosen, and is then compared to a threshold obtained during model training. If the minimal T is smaller that the threshold, the test image is tagged forest, if greater, then-non-forest.
Remark 2. This setting is equivalent to multiple hypothesis testing with the thresholds set to be equal. This restriction can be dropped, in which case the precision of predictions will generally increase. But the downside of such a generalization would be that the additional computational complexity during both training and classification would increase linearly with the number of forest images.
Remark 3. This algorithm describes dependences between colours by means of the pooled sample covariance matrix. On the contrary, all the information about spatial features are erased. Also, distributions of colour intensities in each channel are described only by their location parameters.

Training and cross-validation of the non-parametric model
Training of the non-parametric model is schematically depicted on Figure 5. One can observe two sets of images there: reference forest images in the left column (set 1) and images of both types in the upper row (set 2). For every pair i-th image from set 1 and j-th image from set 2, statistiĉ T ij = n1n2 n1+n2 D 2 ij is computed, where D ij is the Mahalanobis distance between the two pictures.
Then, for every image from set 2, the minimal statistic min iTij is found. Lastly, the threshold T is chosen to minimize misclassification error when classifying images from set 2 according to a result of evaluation of Boolean expression min iTij < T . As previously, True corresponds to forest, while False-to non-forest. Aside from simple training on a set of images, we adopt the setting of k-fold cross-validation for this model. Initially, we shuffle pictures and group them into k folds. In this procedure, in order to preserve internal features, we don't split or merge pictures. We follow a similar procedure to what was described above. Once a fold is chosen, set 1 consists of all forest pictures of a training set, while set 2 is a holdout set. Then statistics min iTij are computed for each of the k holdout sets and T is chosen to minimize the misclassification error on all the combinations. Cross-validation is summarized in Algorithm 1. for each image in fold i do 6: ComputeT = n1n2 n1+n2 D 2 from Theorem 1, between the image and all images labelled forest in the remaining k − 1 folds.

7:
Keep only the smallestT for the image. For all images and allT 's, tag images withT < t as forest andT > t as non forest.

12:
Compute the accuracy for the given threshold t, Acc(t). 13: end for 14: return arg max t Acc(t) and n's, means and covariances of forest images.

Detection with a parametric method
The non-parametric approach described in Section 3 is simple and universal in terms of underlying data, it doesn't make any assumptions on the data coming from RGB composites. On the other hand, it possesses a theoretical drawback-it doesn't take advantage of the distribution of the data, only its mean and covariance matrix. The parametric approach to the detection was designed basing on the opposite ideas-to capture the actual distributions of the data and parametrize them in order to make use of parametric theoretical machinery.

Stable distributions for signal processing
The second classification method we develop in this study is based on univariate stable distributions as they proved to fit well to the empirical data (Section 2). In this short section we give some context for use of stable distributions in the field of signal processing, describe how they are fit to satellite data in the present paper and how we estimate sampled data.
The stable distribution family is a class generalising the Gaussian one. In signal processing it appears due to two factors. First, stable distribution is a result of the Generalised Central Limit Theorem which says that if there exists any limit of a mean of i.i.d. variables with possibly unbounded variations, then this limit must have a stable distribution. Thus, this class is a natural choice for modelling signals with potentially large deviations. Second, stable distributions might be asymmetric, and it is controlled by one of the parameters. This family of random variables has already seen use in remote sensing areas such as ship detection [16] and video background subtraction [2]. While it is unclear to us what role the Generalised Central Limit Theorem plays in shaping the distributions of forest colour intensities, it is certain that parametric control of their shape provides accurate fitting (Figure 2) when the general stable distribution is assumed. Definition 1. A random variable X is said to have a stable distribution if there are parameters 0 < α ≤ 2, σ ≥ 0, −1 ≤ β ≤ 1 and δ ∈ R such that its characteristic function has the following form: Where sign t = 0 if t = 0. The parameter α is known as the index of stability and controls the decay of the tails. The parameters β is not the classical statistical skewness but it indicates how the distribution is skewed, indicating left skewness if β < 0, right skewness if β > 0 and symmetry if β = 0. The parameter σ is a scale parameter controlling variability and δ a location parameter which shifts the distribution. The location parameter only corresponds to the mean when α > 1, otherwise the mean is undefined. In fact, for α < 2 the p-th absolute moment E[|X| p ] < ∞ iff p < α.
Apart from a few cases, the closed form probability function is unknown for stable laws. There exists a method for computing stable densities numerically [9] but it is computationally demanding. That is why we developed a mathematical foundation of our parametric method based on CF's and ECF's. The most famous special cases where the probability density function can be expressed explicitly is when the stable distribution follows a Gaussian, Cauchy or Lévy distribution. A stable distribution is Gaussian when α = 2.

Parameter estimation
This section will cover estimation of the four parameters of the univariate stable distribution used in this article. There is a number of techniques available for estimation of these parameters, including quantile-based estimators, Maximum-likelihood-based estimators and those based on the method of moments. For this work, the Koutrouvelis regression-type estimator [6] was chosen, as it shows a great balance between precision and computational performance. Estimators which performed well but were very slow were also studied, including ones based on the generalized method of moments. Similarly, estimators which were faster than the Koutrouvelis one were assessed, but these did not perform well enough for our purposes: the McCulloch estimator [8] based on quantiles and the one developed by Press [11]. A more in-depth study on the performance of different parameter estimators may be of interest to individuals interested in improving the algorithm proposed in this thesis.
The Koutrouvelis regression-type estimator is quite fast, but one big downside is that the estimator of the location parameter δ is often unreliable. For our purposes this turns out not to be a big issue though because experimentally we have seen that the parameter α > 1 in all our forest data and when this is the case the location parameter δ equals the mean of the distribution [14, Property 1.2.16]. Which implies that the sample mean is a consistent estimator for the location parameter δ for our data.

A test for empirical characteristic functions of stable distributions
In this part, we develop a statistical test that rejects a null hypothesis that a sample comes from a stable distribution defined by (9) with certain parameters. First, we provide a limit theorem for vectors composed of parts of ECFs of stable random variables. For a sequence {X k } k=1,... of iid stable random variables the ECF is given by and define new random vectors and the elements of the matrix are given by Remark 4. Note that when X j , j = 1, . . . , n are also symmetric and centered around the origin Σ z reduces to the form This case was described in Section 3 of [11].
Proof of proposition 1. First, rewrite u j and v j : Since both of them are bounded, then Z n (t) is a mean of random vectors with components having the second moment. Thus, the multidimensional central limit theorem immediately applies. Therefore, it is only left to find the entries of the matrix Σ z (13).
Cov(u j , v j ) = Cov Var Having substituted covariances in (17) using (18) one has Other elements of the covariance matrix Var(u j ) and Var(v j ) are computed in the same manner: Proposition 1 provides a limiting distribution for a two-dimensional vector Z n − Z 0 , making the result inconvenient for statistical inference. To address that, we define an analogy of Mahalanobis distance-a statistic characterising how far Z n is from Z 0 . The following theorem provides a limiting distribution for this statistic and plays the key role in the parametric inference.
Proof. a) According to Proposition 1 and properties of weak convergence Hence, by Continuous mapping theorem, b) follows from independence of Z (1) n and Z (2) n .

Hypothesis testing based on ECFs of colour intensities
Now, let us return to the data extraction ( Figure 6). Assume we have an image of forest and a test image. Having extracted matrix X, one can estimate parameters of the stable distributions of the forest picture channels and compute vectors Z 0 r (t), Z 0 g (t) and Z 0 b (t). Similarly, after obtaining X from the test image, empirical characteristic functions are computed and subsequently Z n r (t), Z n g (t) and Z n b (t) are obtained. Note that we use CFs and ECFs of each channel separately, so dependencies between colours are ignored so as spatial dependencies. Hypothesis testing in the parametric setting is depicted on Figure 7. A test image is compared with each of the forest images in the training set. When comparing a test image and a forest one, [Z 0 ] are used to test separately that the colour intensities of the test image are those of forest: It implies that in order to be tagged as forest, a test image must have all three colour samples passing their corresponding test. Since the limiting distributions in (25) are the same, thresholds for each channel are set to be equal.
For each pair the test image and a forest image, T tf := min(T r , T g , T b ) is computed and then T min := min f (T tf ) is obtained and compared with the threshold.

Training and cross-validation of the parametric model
The training procedure starts with findingΣ z;r 's, Z 0 's and Z n 's for forest images and Z n 's for nonforest ones in formulas (25). After that, training proceeds as in 3.3: T statistics for pairs of images are computed and the threshold is set to minimise the misclassification rate. Throughout classification, Σ z;r 's and Z 0 's will be available in advance, so that only Z n 's will be left to compute to perform hypothesis testing based on (25). For that reason we do not build the parametric framework on Theorem 2 b) withΣ z;r 's obtained from both the test image and a training forest image: estimation of parameters of stable distribution is quite computationally demanding. Although, this approach would provide better accuracy. Cross-validation for the parametric model is summarised in Algorithm 2.

Comparison of the methods
There are features typical for both of our methods. First, if by accident there is a non-forest image in the training forest set, test images of its type start to be tagged as forest because they pass the corresponding test. As a consequence, additional false positive results will occur. Specially designed cross-validation may be exploited to eliminate this problem. Second, the setting is multiple hypothesis testing, so it could naturally be re-set to classify forests images into multiple types such as rain forest, coniferous forest, etc. Another thing is that clouds and sunbeam scattering affect the algorithms because they introduce noise to the distributions of objects. It could be cured by atmospheric correction [7]. Also, both models increase accuracy with growth of the numbers of points in the training data and test data. Improvement of prediction accuracy could be achieved in three ways. Increasing resolution, e.g. by switching from satellite data to aerial imagery, would add more pixels to objects under consideration while adding new tagged images to the training set would provide better training of thresholds and more tagged classes. In practice, merely increasing the size of an image with fixed resolution (adding more area) works well only in the case of training pictures: enlarging test images would imply observing pictures with objects belonging to different Algorithm 2 Parametric cross validation procedure 1: Inputs: Data set of images labelled forest or non-forest, number of folds k, max threshold T max , parameter for CF's (t in (9)). 2: For every forest image compute Z 0 , parameters and the covariance matrix. 3: For every non-forest image compute Z n . 4: Randomly shuffle images and group them into k folds. 5: for i = 1 to k do 6: for each image in fold i do 7: Compute 3T 's from (25), between the image and all images labelled forest in the remaining k − 1 folds.

8:
Keep only the smallestT for the image. For all images and allT 's, tag images withT < t as forest andT > t as non forest.

13:
Compute the accuracy for the given threshold t, Acc(t). 14: end for 15: return arg max t Acc(t) and Z 0 's, parameters and the covariance matrices of forest images.
classes and tagging some areas incorrectly because of that. And lastly, unlike classical machine learning methods, the models developed in the present article are robust to data transformations such as rotation and reflection.
There is also a number of differences between the two models. An advantage of the parametric model is that it can work on incomplete data: e.g. when only 1 channel is available or different channels are available at different times. The non-parametric method needs at least two channels simultaneously. On the other hand, the disadvantage of the parametric model is that the algorithm generally requires more pixels per image in training as it involves parameter estimation and because we chose to use theorem 2 a). The estimation procedure also makes the parametric method slower than the other one. A con of the non-parametric method is that the model does not take into account the actual distributions of colour intensities, which can contribute to additional false positives.
To summarise the comparison, we point out that at the current stage of research it is not possible to claim superiority of one method over the other as they both have pros and cons.

Empirical application
Sentinel-2 satellite images are used to illustrate how the new Mahalanobis distance classifier (MDC) and stable distribution classifier (SDC) perform in comparison to some other commonly used classification methods. The methods used for comparison are support vector machine (SVM), logistic regression (LR) and deep neural network (DNN). Three experiments with real data are shown here: two comparing accuracy of the methods and one giving visual representation of performance. As the source of terrain images, the region Occitanie in France is used, more precisely, vicinity of the city of Albi.
In each of the experiments SVM is implemented with the radial basis function kernel and the DNN has four hidden layers of sizes 350, 200, 100 and 50 together with ReLu as activation function. The feature space for the comparison methods consists of all pixels and their RGB colour values.
In practise, one would typically use some dimensionality reduction technique, e.g. PCA, before implementing the classification methods, but since no such technique is used for our new methods we want to be consistent with using all of the available data. Considering that the images are rather small this should not cause any issues.
In the first two experiments we look at how the methods perform using a training and test split of the available dataset. In the first one the training set consists of 510 images of sizes 20×20 pixels with an equal number forest and non-forest images. The test set is also equally balanced between forest and non-forest images but it consists of 200 images of sizes 20 × 20 pixels. The second experiment is similar, except that each 20 × 20 image is split into four equal pieces, making 2040 images of sizes 10 × 10 pixels in the training set and 1000 10 × 10 images in the test set. It should be noted that all of the images are collected from the same region and under similar conditions. This means that the models are not trained for different seasons or different types of forests. However, it would be straightforward to make the classification models more general by including more diverse data.
The performance measure used to evaluate the methods is the error rate, i.e., the number of incorrectly classified images in the test set divided by the total number of test images. Figure 8 illustrates how the methods perform when the training set varies between 10 % and 100 % of the total training set using images of sizes 10 × 10 pixels and 20 × 20 pixels, respectively. Even though the training set size varies, an equal number of forest and non-forest images are always used for training. The test set is kept unchanged.  It is clear in Figure 8 that the new methods outperform all of the comparison methods. The error rate drops to zero almost immediately in the 10 × 10 pixel case when using SDC whereas the error rate corresponding to MDC gradually decreases from 2.4% and reaches zero when the training set size is about 1800. All of the standard classification methods perform worse and their error rates do not fall below 2%. The best performing methods among the comparison ones are DNN and SVM. They show similar results apart from the fact that SVM is more robust to changes in performance when adding new data to the training set. Analogous results are presented in the 20 × 20 pixel case. Here, MDC generally has a lower error rate, and it drops to zero after the size of the training set reaches 300. On the contrary, the error rate of SDC is higher than in the previous case. Moreover, when the size of the training set is lower than 62, the error rate is higher than those of DNN and SVM. When the size of the training set reaches 62, the error sharply drops and becomes lower or comparable to those two.
In the third experiment we visualise how the new methods and DNN (one of the best performing comparison methods) behave when every sub-image of size 10 × 10 pixels is classified of an image of size 2000 × 2000 pixels. Sub-images that are classified as forest are marked as white and non-forest sub-images are marked as black. This makes it possible to get a better feeling for how the methods differ and identify their strengths and weaknesses. The same (complete) training set as in the first part is used to calibrate the methods. This training set does not contain images that are part of the 2000 × 2000 pixels image but they were collected form a nearby region under similar conditions. The result of the classification is illustrated in Figure 9. Generally, all three models performed well, though there are differences in their outcomes. One can see that unlike MDC and SDC, DNN is not so sensitive to the intensity of the green colour. It is best illustrated by performances of the models on forest covered by a shadow of the cloud in the right part of Figure 9. MDC and SDC incorrectly treat the shadowed vegetation as non-forest, while DNN tags it as forest. Shades from clouds generally could be accounted for by cloud detection, which is a solved problem in Earth image processing (see e.g. [7]). In some other cases DNN tags non-forest vegetation as forest. In particular, in Figure 10 there is forest rounded by non-forest vegetation. Unlike MDC and SDC, DNN fails to detect some non-forest parts of the circle. A Data pre-processing technology We have built a preprocessing procedure of creating RGB composites from Sentinel satellite images using Python programming language. First, we download Sentinel 2 Level-1C products from PEPS [1] with the means of EODAG. The PEPS (Plateforme d'Exploitation des Produits Sentinel) platform, redistributes the products of Sentinel satellites, S1A, S1B, S2A and S2B, S3A and S3B from COPERNICUS, the European system for Earth monitoring. EODAG (Earth Observation Data Access Gateway) is a command line tool and a Python package for searching and downloading remotely sensed images while offering a unified API for data access regardless of the data provider.
Second, for each product we do a preprocessing procedure of creating RGB composites from Sentinel 2 Level-1C products. Using rasterio package 1 we extract red, green and blue bands from channels B04, B03 and B02 respectively and combine them in a single file using GeoTIFF file format 2 .
Third, we compress RGB-composites with lossless compression algorithm to decrease the amount of storage space needed for data using gdat translate utility [5].
Lastly, in order to obtain images with sufficient brightness, we normalize values in the channels to range from 0.0 to 1.0 by dividing them by the constant 2000. We truncate values with normalised intensity more than 1.0 by setting them to 1.0, which primarily affects bright white clouds and conforms with typical "crop and scale" image processing procedure. Deforestation detection algorithms use this normalised values.

B Distributions of colour intensities of other natural objects
In this section we explore two RGB images: one-of a mountain and another one-of sea. Each picture is accompanied by empirical density functions plotted for all channels. For comparison, best fitting stable, gamma and normal distributions are shown.