Land Use Recognition by Applying Fuzzy Logic and Object-Based Classification to Very High Resolution Satellite Images

The earth observation (EO) sector has seen notable technological progress in recent years, especially in the satellite sector, where increasingly higher-resolution images have been put on the market. Satellite resolutions rapidly went from values lower than one meter, for example, the SPOT 6 and 7 satellites with a Ground Sample Distance (GSD) of 1.5 m on the ground in the panchromatic image, defined as medium–high resolution. Up to a GSD of 31 cm, like those of the acquisitions of the WorldView-3 and WorldView-4 satellites, or more generally, acquisitions below one meter always refer to the panchromatic image, defined as having very high resolution. In addition to the very high resolution, these products are available with return times, i.e., the possibility of having an image of the same area of the planet in less than a week in some cases. All these advances have attracted more and more interest in the development of methodologies for soil classification, more or less automated, aimed at facilitating many applications, from the creation of maps to various uses in precision agriculture, monitoring the health of crops, to the observation of cities and their development. Thanks to the very short revisit times, the monitoring of an area by applying change detection techniques to the images is increasingly used even in the most extreme situations, in which, following natural disasters, it can provide support in estimating the damage and planning the necessary interventions to restore the territory [1]. Our world is changing more and more rapidly, and satellite images have definitively become a valid tool for promptly observing and monitoring this change. At the same time, increasingly high-performance classification methodologies have been developed, from pixel-based analysis to object-based analysis (OBIA), deep learning algorithms, and, although computationally more expensive, convolutional neural networks (CNNs) have become more popular in recent years. The state of the art sees the use of such algorithms and methods for most applications, not only in the classification of the scene in terms of land use but also in the previous phases, such as image fusion, segmentation, and registration of images. These methods, although very effective, must be some more, some less trained in some way to perform their function for the classification of land use by providing known image elements, elements already classified that are sometimes not available in large quantities. In the work presented below, an object-based approach combined with the use of fuzzy logic is proposed; radiometric and geometric aspects of the objects segmented in the image are used simultaneously in a logical scheme designed to reproduce what could be the line of thought of a human observer looking at the image. In this way, the training phase is eliminated, which, while maintaining the logical descriptive scheme of the class, is currently replaced by the calibration of the parameters, which, applied to the individual membership functions, duly combined with each other, define the class of an object. Mainly, the proposed method has the advantage of exploiting a logical description of the classes, which, once defined, is not strictly linked to the context analyzed, but precisely because it is developed on the basis of what could be the considerations of a human observer, it has greater versatility. In addition to this, you do not need to sacrifice for training, if present, part of the ground truth; this phase is replaced by a necessary calibration of the parameters that describe the various classes. However, it must be taken into consideration that this operation does not require much time, as the geometric parameters can be easily adapted depending on the resolution of the image available. Radiometrically, the description of the classes and the main factors involved are the indices and the standard deviation; the indices already have normalized values and are therefore suitable in contexts like the one analyzed, while the standard deviation values can be adjusted in various attempts or by evaluating the correlation with the histograms of the image bands, with images already classified. The ultimate objective of our work is to obtain a soil permeability map, attributing permeability values to the identified classes, through a classification method accurate enough that we need the ground truth only for what is strictly necessary to validate the results since, if not already available, creating ground truth is a very time-consuming activity. In this work, a step forward has been made, using a portion of the Pavia image, for the development of the logical scheme to identify the class of fields and to test the logical scheme describing the water class developed in previous work to effectively distinguish it from the shadows in an urban context, a problem well known in the literature [2]. In both works, the image used is the same, but different areas were analyzed; specifically, the image was acquired by the WorldView3 satellite [3] in March 2021 over the entire area of the municipality of Pavia, a city located in Italy in the Po Valley. It has the ability to acquire a panchromatic image at a ground resolution of 30 cm with a revisit frequency of 4.5 days at 20° off-nadir or less. The raw data are composed of a panchromatic band with a GSD of 31 cm and eight multispectral bands: Coastal Blue, Blue, Green, Yellow, Red, Red Edge, Nir1, Nir2, with a GSD of 124 cm. From these, by carrying out a pansharpening operation based on the Gram–Schmidt orthogonalization [4], using the ArcGIS Pro software, it was possible to analyze and classify an image with a GSD of 31 cm and an estimation of the spectral information of the eight bands. The developed algorithm is mainly composed of three phases: segmentation, classification, and refinement; the result is then compared with the ground truth created specifically in our laboratory by an operator, who patiently identified and manually classified a good number of objects before the classification had been performed. The area manually classified to constitute the ground truth is equal to 975.092,3 square meters, or 15.64% of the analyzed area. The result obtained satisfied us because, beyond the overall accuracy value achieved, 95.32% is already very good in itself. Given the heterogeneity of the scene, it makes us think we can proceed with the classification of the entire image in our possession, keeping the methodology and parameters used unchanged and obtaining a good result. Tests carried out on other portions of the image do not reveal any particular limitations related to the method and seem to confirm its effectiveness, at least visually, by returning a good classification. However, at the moment, we do not have the possibility of rigorously showing its effectiveness in terms of accuracy as we do not have sufficient truth to the ground in those areas.

The area analyzed is a portion of the image acquired by the WorldView3 satellite in March 2021 of the entire municipal area of the city of Pavia. The product purchased from the Geomatics Laboratory of the University of Pavia includes for the entire area an 8-band multispectral image with a lower ground resolution, equal to approximately 120 cm for each pixel, and a panchromatic image, which, in a gray scale, shows the entire area with a ground resolution of approximately 30 cm, as summarized in Figure 1. In addition to the panchromatic and multispectral components, the WorldView3 satellite is capable of acquiring another 8 short-band infrared (SWIR) bands, which have an even lower resolution of approximately 320 cm on the ground [3]. The SWIR component of the image was not purchased because we are interested in seeing what can be obtained from a high-resolution image, and the data resulting from this component have a high but not sufficient spatial resolution. For this reason, only the panchromatic and multispectral components of the image of the Pavia area were purchased and used for classification.

The two components, panchromatic and multispectral, were combined by applying the pansharpening algorithm based on Gram–Schmidt orthogonalization [5]. Within ESRI’s ArcGIS Pro^TM v.3.0.2 software, the details of the technique used by the software are described in this patent [6], obtaining an image with the spectral information of the 8 bands (Coastal Blue, Blue, Green, Yellow, Red, Red Edge, Nir 1, Nir 2, and the degree of detail of the panchromatic image) at an almost 30 cm resolution on the ground [4]. Clearly, the result obtained, distorts the basic multispectral information; various algorithms could have been explored to enhance the image quality and minimize spectral errors [7]. The achieved level of detail, even without additional algorithms, is noteworthy. This detail proves crucial for distinguishing diverse elements within the urban area. However, the advantage of working on a much more detailed image is that it allows a significantly higher degree of accuracy to be achieved for classification purposes.

The area analyzed is located slightly east of the historic center of the city, as shown in Figure 2, and features a modest-sized industrial area and small residential areas; the rest is made up of agricultural areas and isolated farms. The choice was influenced by the fact that this area involves various elements, in particular large industrial buildings, which in past classifications were sometimes erroneously classified as uncultivated fields, small lakes and waterways, and isolated buildings in the middle of the countryside, all located in an area adjacent to the residential center (see Figure 3). This proximity allows us to see if the algorithm and the methodology developed are suitable for dealing with a mixed context and, therefore, applicable to the entire image without having to modify the parameters specifically. This occurs through a dedicated segmentation process, using the multiresolution segmentation algorithm [8], set in a first phase to identify very large objects, the size of fields, and after having identified them as fields or not, a detailed segmentation is necessary to continue in image classification. Both in the segmentation and classification phases, two indices are exploited: the Normalized Difference Vegetation Index (NDVI) [9] and the Normalized Difference Water Index (NDWI) [10], which are very useful for highlighting and, therefore, more easily identifying vegetation and water, respectively. The result obtained from the classification is ten refined using algorithms that allow the fusion of adjacent objects if they satisfy certain criteria and algorithms that allow the classified objects to grow in neighboring pixels under certain conditions, which are exceptionally effective for improving the edges of objects classified as water. Overall, the area studied is composed of 64,877,220 pixels and covers an area of 5.8389 square kilometers. In the area, there are no large watercourses, only small streams and lakes, which are identifiable during the construction of the ground truth in some cases only because the trees above them do not yet have fully developed foliage in early spring. Approximately half of the fields present have crops in an advanced state, while others have either recently passed the sowing phase or have still been left in an uncultivated state, depending on the type of harvest. Paradoxically, we will see in the next chapters how some small wooded areas will be classified as uncultivated fields precisely because the foliage of the trees is not yet developed and no other vegetation has grown on the ground because it is covered by dead leaves.

The methodology applied is basically the same as already presented in a previous article [2] where the approach used for the classification of urban land was exposed. Focusing mainly on the correct distinction during the classification phase between shadows and water is a well-known problem in the literature [11,12,13,14], during the classification of urban areas. The main steps have remained unchanged, and in the order, the useful indices are first calculated, then a first segmentation is carried out to identify and subsequently classify larger objects such as fields, and a second segmentation to identify and classify smaller objects, water, and vegetation. Finally, the classified objects are refined to better define their edges, and what has not yet been classified is attributed to the impermeable class. Even if some parameters have been slightly adapted during the segmentation phase, the main feature remains unchanged, namely the creation of a logical scheme that exploits both the geometric and spectral characteristics of the segmented objects. A scheme is derived on the basis of what could be the logical processes that a human observer could generally adopt to distinguish the various elements in the image.

Figure 14, shown below, displays the results of the classification, with impermeable objects shown in red, uncultivated fields in yellow, cultivated fields in light green, water in blue, and vegetation in green. For comparison purposes, a shapefile was specifically created within ArcGIS Pro, containing several objects belonging to the impermeable fields, water, and vegetation classes. These objects were manually identified by a user, who observed the image and determined the edges of the objects as accurately as possible. During the comparison phase, cultivated fields were merged into the vegetation class, and individual pixels belonging to the ground truth were compared with the corresponding pixels in the classified image. From this comparison, a confusion matrix was derived, summarizing the achieved accuracy. In the lower part of the image, what immediately catches the eye is the presence in the ground truth of a small stream of a lake mainly used for agricultural purposes, which the developed classification algorithm is not able to identify. This is mainly due to the small size of the objects in that area and the presence of trees, which largely cover the path of the stream. Although the trees are bare enough to allow the operator to see the stream and identify it, radiometrically, they reflect enough light in the near-infrared wavelength to mask it from the algorithm’s eyes. This condition was not present in the area used for the development of the logical scheme of the water class since it was an urban environment; therefore, we are not surprised that this error was made, while the distinction between water and shadows, the objective for which it was developed, is carried out very well.

The result is that these objects, due to this particular interaction, as the water class has been described, are not suitable to be part of it and, at the same time, do not have high enough NDVI index values due to the thin foliage of the trees to be classified as vegetation. Unfortunately, this ultimately leads to the attribution of the impermeable soil class to objects, for which it is obvious that this should not be the case. In cases where the vegetation is sufficiently abundant, the small portions of the river are incorporated into larger objects classified as cultivated fields or into objects of the vegetation class. However, it is a less impactful scenario, given our ultimate goal of obtaining a soil permeability map from the classification result. From the comparison between the results obtained from the classification by applying the developed algorithm and the ground truth identified in the study area, a confusion matrix and an accuracy matrix were derived and reported, respectively, in Table 1 and Table 2 below.

In the first table for each column, we see the corresponding class in which the same pixel was classified. As an example, of the 251,804 ground truth pixels that were cataloged as water, 60,357 pixels were classified by the software as water class, 59,486 pixels as vegetation class, and 108,647 pixels as uncultivated fields, while in the second table, the corresponding accuracies are shown. Ultimately, an overall accuracy of 95.32% and a Cohen’s kappa coefficient of 0.93 were achieved in the analyzed area.

As we have seen, there are difficulties in classifying small watercourses, for which the presence of tree foliage is a significant problem, due to the acquisition period, which we remember took place in March 2021. Fortunately, for the applications we are aiming for, this defect does not have major repercussions. The problem is mainly linked to the classification of water, which, if present in small-sized elements, is difficult to identify and is mainly linked to the fact that the image resolution is not sufficient. If we consider that the “real” data used have 120 cm of resolution on the ground, we realize that those rivers are only 2/3 pixels wide, and if there are trees, the radiometric data are significantly altered, considerably complicating correct identification. Furthermore, as mentioned in the belief, this error in the failure to identify water is due to the fact that the logical scheme developed to distinguish between water and shadows in an urban context. Therefore, this result can be seen from two points of view; on the one hand, it highlights the need to improve the water class description if there is a need to at least partially identify these small elements. On the other hand, the effectiveness of what has been developed is confirmed; in the area studied, there are no shadows, neither in the urban area nor outside it, incorrectly classified as water, which is a non-negligible result.

More impactful is the incorrect classification of the field adjacent to the small watercourse, which involves the attribution of the impermeable class to a fairly large surface, highlighted above in Figure 15. In the validation, that particular field was not taken into consideration; however, we can make some considerations regarding the accuracy obtained, in particular underscoring that the extension of the fields plays a fundamental role both in the classification and in the final validation. In the specific case over the entire area, the greatest classification error occurred precisely in this case, and even if the erroneously classified area is considerable (39,167.46 m²), compared with the overall extension of the uncultivated fields present in the classified (2,072,011.59 m²) image, it is reduced to a very low error percentage of 1.89%. Therefore, it is not considered necessary to modify the description currently adopted for the identification of cultivated and uncultivated fields to try to fix this error, which, as demonstrated by the results, has an excellent success rate in identifying these elements. It must also be considered that there are cases in which the result obtained is very satisfactory, where the result obtained reflects almost perfectly what the operator identified in the ground truth. As shown in Figure 16, the parts of non-permeable soil and those that are permeable in the form of vegetation are almost perfectly identified; in this case, the software also outlines the roads and edges of the square serving the shed very well. Therefore, after taking these cases into account, we can state that we are satisfied with the result obtained from the sole application of the developed algorithm, considering that errors of a gross nature, such as the incorrect classification of the fields shown in Figure 16, can always be easily fixed by the user in a review phase of the result obtained. Another aspect that can certainly be improved is the image creation process through pansharpening; currently, in the ArcGIS^TM software, it is possible to perform pansharpening by applying the Gram–Schmidt orthogonalization only on four bands (R, G, B, Nir).

Based on the sensor that acquired the image, the software proposes a different weight for each band to adjust its contribution to creating a low-resolution panchromatic. These weights depend on various factors and are available by default for the four bands listed above. For the other four bands at our disposal—Coastal Blue, Yellow, Red Edge, and Nir 2—the same weights have been adopted, knowing that this necessarily introduces an error in the final result when these four bands are shown at high resolution. This necessarily presupposes a recalibration of the parameters involved when the pansharpening technique is changed or the weights with which to combine the various multispectral bands to form the low-resolution panchromatic are optimized. The technique used at this stage of development clearly presupposes the presence of an operator ready to recalibrate the parameters involved in the description of the classes as the analyzed image changes. The logical process that links and combines the membership functions between them can be considered valid, and for the water class, it has been proven valid, where there is only the need to improve the identification of small watercourses in the case in which they are partially covered by vegetation. Making the parameter calibration process automatic would constitute a notable step forward towards the application of this methodology for the segmentation phase using segmentation procedures already present in the literature [18]. While for the optimization of the values to be attributed to the right and left limits of the membership functions, for the average values of the bands, it could be interesting to use auto-threshold algorithms or, alternatively, to statistically analyze the values adopted in the classification of the image in this work and contextualize them in the overall image analyzed. Subsequently, test on a new image to see whether, by fixing the parameters a priori and applying the same correlations to the radiometric aspects of the new image, the results obtained are good or not. From what emerged from the application of the water logic scheme, it seems that this operation is not necessary if only the region studied is changed. It would certainly be interesting to test the defined parameters on a more current image acquired from the same satellite, and only through a rigorous comparison would we have confirmation or not of the real degree of dependence on the calibration of the parameters involved in the proposed method. It must also be taken into consideration that the classification obtained is absolutely independent of the ground truth, and this represents an advantage compared to machine and deep learning techniques. Creating ground truth is very time-consuming, whether you create it manually or extract it using different methods [20], even more so if part of this information must be used as examples to train this kind of algorithm. Also, in this case, the algorithms must be trained on the basis of the radiometric information of the image, and therefore, the process must be repeated if the sensor that acquires the image changes.

In conclusion, the work carried out allowed us to evaluate the effectiveness of the description of the water class previously developed in another work on a different area of the city of Pavia, discovering limitations due to the presence of vegetation but confirming its effectiveness in distinguishing between water and shadows. Furthermore, the logical scheme describing the field class has been improved, producing excellent results without misclassifying any industrial building, which had been our main problem in previous tests. The distinction between cultivated and uncultivated fields was added later, given the simplicity with which it is possible to identify them, even if it does not constitute essential information for the creation of soil permeability maps. We are satisfied with the results obtained, considering that we are still working to improve the classification result and to make it less dependent on the region studied, reducing its dependence on the recalibration of the parameters involved from time to time. I can imagine an improvement in the classification, in particular, given the critical issues with small watercourses, to further improve the potential of the methodology adopted in identifying land use. Similarly, although not fundamental for the creation of soil permeability maps, the distinction of impervious elements could be evaluated and addressed as a separate problem to distinguish in detail between the portions of soil representing buildings or roads and further enhance the ability to classify land use. Ultimately, the proposed methodology allowed us to achieve an overall accuracy equal to 95.32%, a result that I consider very satisfactory if we consider that this result takes into consideration the small quantity of water present in the scene, which, with a success rate of 23.97%, significantly lowers the overall accuracy achieved on the classification of the entire area. Future developments can certainly concern an improvement in the phase of the creation of the high-definition multispectral image through pansharpening. In my opinion, this constitutes the most critical step because we know for sure that it introduces the Coastal Blue, Yellow, Red Edge, and Nir2 bands, an error because the weights adopted for the creation of the low-resolution panchromatic, even if they are an approximation, are not correct. Ideally, I would like to finish developing an algorithm in Matlab to create the image in a single step through the pansharpening process that adopts the Gram–Schmidt orthogonalization algorithm to have more control over the individual factors that influence the result and have a final product in which the distortion of the original data of the multispectral bands in the final high-resolution image is minimal. Further progress can be obtained by inserting another step in the refinement phase, as highlighted in this work, and trying to extract, if possible, from objects classified as cultivated fields or vegetation any small and partially hidden waterways to improve further the quality of the result obtained. In the last few years, experiments have been carried out in various test areas. The one shown in this work is only the latest, and the results achieved suggest that the methodology adopted can give good results in very different situations, such as urban, rural, and mixed contexts, and therefore, it is suitable for carrying out the classification of the entire area of the municipality of Pavia.