A SMAP Supervised Classification of Landsat Images for Urban Sprawl Evaluation

Flavia Di Palma 1, Federico Amato 1,*, Gabriele Nolè 2, Federico Martellozzo 3 and Beniamino Murgante 1 1 School of Engineering, University of Basilicata, 10 Viale dell’Ateneo Lucano, 85100 Potenza, Italy; flaviadipa@libero.it (F.D.P.); beniamino.murgante@unibas.it (B.M.) 2 Italian National Research Council, IMAA C.da Santa Loja, Tito Scalo, Potenza 85050, Italy; gabriele.nole@imaa.cnr.it 3 Department of Methods and Models for Economics, Territory and Finance, University of Rome “La Sapienza” Via Del Castro Laurenziano 9, Roma 00161, Italy; federico.martellozzo@uniroma1.it * Correspondence: federico.amato@unibas.it; Tel.: +39-320-337-1275


Introduction
The environment and settlements have always been the two pans of a balance, both fundamental components of human needs and in persistent imbalance in the management by man [1].Today, more than in the past, this imbalance has become irreversible.For many years, a common approach to territorial governance renounced the vision [2] of the environment as a finite resource, as a common good that we have, but that does not belong to us.Urban growth processes represent the majority of human-induced changes on the environment.Many studies focusing on such phenomena [3,4] highlight their negative impacts on public health [5], climate [6] and the agro-economy [7].In the current context, the environment and settlements are mainly characterized by anthropization and land take.In Italy, for instance, artificial areas (as defined in the CORINE Land Cover nomenclature) increased, inducing a land take change from 2.7% in the 1950s up to 7.0% in 2014, with an increase of 4.3% [8].Land take is now a critical issue at the center of political and scientific debates [9][10][11].However, the complexity of the causes generating this phenomenon and the effects that it induces on ecosystems and environmental services make its assessment extremely complex [12].
Moreover, urban growth is not the only way through which human activity generates negative effects on land use.Just as in the economic and demographic processes, societies are characterized by a succession of different regimes of land use.From natural ecosystems, to the development of local agricultural systems, up to intensive agriculture and the development of so-called urban sprawl, people change with their actions the territory where they live.Even an activity once considered as low impact, such as agriculture, is on the one hand creating an advantage in the short term in terms of food production, while on the other hand, it is producing serious long-term problems in terms of loss of ecosystem services, including some important ones for agriculture itself [13][14][15][16].
Usually, land take is confused with the concept of urban sprawl, defined as low-density urban development [17,18], mainly concentrated in peripherals areas, unplanned and characterized by the simultaneous presence of different land uses [4].Unfortunately, urban sprawl represents only a portion of land take phenomena.Moreover, the processes of the intensive use of soils for agriculture and forestry should be considered, because they cause a loss of the bio-geo-chemical features of soils [8,19].
Although commonly not considered as actual causes contributing to the phenomenon of land take, these processes can still affect the proper functioning of a natural environmental system.Therefore, the evaluation of such factors is a fundamental prerequisite in urban planning, because urban dispersion processes characterize land use forms and their measurement through specific indexes [11,[19][20][21].As an example, analyzing the fragmentation of the territory through the urban sprawl index is a common assessment framework to allow better environmental protection.
The EU has defined objectives to contain land take with the aim of reaching a net zero value by 2050 [5], in order to increase environmental sustainability and improve territorial governance.Several strategies to limit land take have been implemented imposing limitations, pursuing the protection of agricultural soils and stimulating the reuse of abandoned sites.The European Union developed many projects to monitor land use changes, such as the CORINE Land Cover (CLC) and the Land Use and Coverage Area frame Survey (LUCAS).Available data since 1990 indicate that residential urban areas have expanded at a rate four-times greater than the population growth rate, while industrial areas grew more than seven-times faster [22].
It is evident, however, that despite the multitude of available data, the different interpretations do not allow a unified coding.In Italy, for instance, data provided by CLC and LUCAS define different scenarios.According to CLC 2006 data, Italy has 2.8% of sealed soils and 5% of artificial areas; instead, LUCAS 2009 data report sealed soils of 7% and artificial areas of 13%.A time span of three years is not enough to justify this great difference in data.The reason for such different assessments is to be found in the methodological differences of the two monitoring systems and in the different nomenclatures adopted by LUCAS and CLC.Both studies, however, show a trend of the gradual growth of artificial areas in all EU countries.
In 2015, the Agenda for Sustainable Development has been defined, which stimulates United Nations countries to achieve 17 sustainable development goals within the next 15 years.Objective 11 of this agenda concerns the development of sustainable human settlements.In particular, it promotes sustainable territorial planning and management by means of positive economic, social and environmental ties between urban, peri-urban and rural areas, by strengthening planning activities at national and regional levels [23].Italy lags behind in the implementation of European measures [24], due both to a huge variety of urban growth geographies that does not allow a simple transposition of such practices and to the large number of local authorities lacking properly implemented planning instruments as well as, in some regions, due to the lack of real sustainability practices [25].These policies are difficult to pursue because of a context that gives the building sector still a central role in the economy [26].
At the national scale, available data on land use categories and changes over the years are still characterized by more problems and far from a desirable uniformity.Analysis of land cover changes is an important element to assess the impact of human activities on the environment at different spatio-temporal scales.The availability of updated and detailed information on land use and land cover is a key prerequisite for urban areas' management and planning.These data could be fundamental in developing simulation models to support the analysis and the consequences of land use changes.In the absence of such data and simulation models, sustainable urban development can hardly be achieved.Urban sprawl is an important indicator of urban development.The recent growth of the awareness that its monitoring represents a first step in ensuring a reduction of land take is crucial to develop more sustainable planning.However, in Italy, and especially in some southern regions, detailed historical data on land cover are not available.This absence complicates the identification of environmental costs associated with the transitions of land use classes occurring in the past.Furthermore, it is difficult to identify future strategies to develop urban areas both able to ensure greater protection of the ecosystem and the landscape characteristics of an area and to maintain the socio-economic benefits associated with the development of urban areas.This paper proposes the use of remote sensing data analysis techniques as a tool to operate a continuous monitoring of land use class changes over time and, consequently, the amount of consumed soil.The study has been implemented using multi-resolution and multi-spectral data (Landsat) and open and freely available software (QGIS, GRASS).The advantages of this approach are the following: 1.
the opportunity to achieve long-term time series through the use of Landsat data, in orbit since 1972, renouncing the higher spatial resolution provided by the newest generation satellites that would not allow historical evaluations of similar time periods [27,28]; 2.
the possibility to provide national agencies and local authorities with a low-cost dataset.
Analyses have been performed in the area of Vulture Alto-Bradano (Basilicata region, Southern Italy), in order to assess the impact of the absence of territorial policies on landscape protection and on the fragmentation of the most pronounced ecological areas, in the period considered.

Study Area
The study area, the mountain district of Vulture Alto-Bradano, is located in the northeast of the Basilicata region (Figure 1).It extends between Latitude: 40 ˝59 1 45"N and 40 ˝49 1 44"N, Longitude: 15 ˝34 1 55"E and 16 ˝10 1 55"E.The name of this area comes from the mountainous area of Vulture Mountain (1326 m), a high ground isolated from the Apennine ridge, which covers about 45,000 hectares.fundamental in developing simulation models to support the analysis and the consequences of land use changes.In the absence of such data and simulation models, sustainable urban development can hardly be achieved.Urban sprawl is an important indicator of urban development.The recent growth of the awareness that its monitoring represents a first step in ensuring a reduction of land take is crucial to develop more sustainable planning.However, in Italy, and especially in some southern regions, detailed historical data on land cover are not available.This absence complicates the identification of environmental costs associated with the transitions of land use classes occurring in the past.Furthermore, it is difficult to identify future strategies to develop urban areas both able to ensure greater protection of the ecosystem and the landscape characteristics of an area and to maintain the socio-economic benefits associated with the development of urban areas.This paper proposes the use of remote sensing data analysis techniques as a tool to operate a continuous monitoring of land use class changes over time and, consequently, the amount of consumed soil.The study has been implemented using multi-resolution and multi-spectral data (Landsat) and open and freely available software (QGIS, GRASS).The advantages of this approach are the following: 1. the opportunity to achieve long-term time series through the use of Landsat data, in orbit since 1972, renouncing the higher spatial resolution provided by the newest generation satellites that would not allow historical evaluations of similar time periods [27,28]; 2. the possibility to provide national agencies and local authorities with a low-cost dataset.
Analyses have been performed in the area of Vulture Alto-Bradano (Basilicata region, Southern Italy), in order to assess the impact of the absence of territorial policies on landscape protection and on the fragmentation of the most pronounced ecological areas, in the period considered.
In Vulture, it is possible to identify agricultural areas with beautiful landscapes, characterized by a particular mosaic of crops (vineyards, olive trees) growing on the hillside, while the valley is used for intensive cereal cultivation.The system of Norman-Swabian castles in the area of Vulture and the system of paths and connections along the valley in the direction of Potenza municipality still emphasize the close relationship between morphological characteristics, land use and types of settlement.Urban centers are mostly on top of hills, in a defensive and dominant position.
The settlement system of Vulture Alto-Bradano has some peculiarities compared to the remaining part of the region.These features developed in different historical phases as the result of sudden and massive events, such as earthquakes and migration phenomena.Close relationships and short distances occur among the largest municipalities (with more than 10,000 inhabitants).These distances are further reduced by the presence of minor centers that generate a sort of conurbation.
The demographic dynamics of Vulture Alto-Bradano are characterized by a flow towards the major centers from smaller surrounding municipalities, which compensates, in these centers, a general migration to the north of Italy, while the percentage of the population coming from municipalities in neighboring regions is not significant.

Atmospheric Corrections
During the pre-processing phase of the images, an atmospheric correction was performed, because different atmospheric measurements are needed to calculate the reflectance to the ground (ρ).In this way, it is possible to take advantage of image-based techniques without having to take measurements at the site of the acquisition of the image.Among image-based atmospheric correction techniques, DOS1 (Dark Object Subtraction) correction was adopted.This technique is based on Chavez's assumption [29] that within an image, some pixels are in complete shade and their radiances received from the satellite are spurious, due to atmospheric dispersion.In practice, assuming the existence of dark objects with reflecting surfaces approximately equal to zero, the minimum value of the Digital Number (DN) is subtracted from all of the pixels of the image, eliminating the effects created by the interaction with the atmosphere.This assumption is based on the hypothesis that very few surfaces on the Earth are completely black; hence, assuming one percent as the minimum reflectance value is better than zero percent.
The trajectory of the radiance is calculated as: where: L min is the radiance that corresponds to a value of digital calculation where the sum of all pixels with calculated values less than or equal to this value is equal to 0.01% of the totality of the pixels of the image [30].LDO 1 is the radiance of the dark object with the reflectance value of 0.01.L min for Landsat images is expressed by the following equation [31]: while the radiance of a dark object is given by: LDO 1% " 0.01 ˚rpESUNλ ˚cosθs*T z q `Edown s ˚Tv π*d2 The path radiance is: where: 1.
ML is a specific multiplication factor of the band (Landsat metadata) 2.
(RADIANCE_MULT_BAND_x) where x is the band number; 3.
DN min is the minimum digital number; 4.
AL is the specific additive factor (Landsat metadata); 5.
(RADIANCE_ADD_BAND_x) where x is the band number; 6.
θs is the solar zenith angle in degrees, which is equivalent to θs = 90 ˝´θe, where θe is the Sun elevation; 8.
T z is the atmospheric transmittance in the direction of illumination; 9.
E down is the descending irradiance; 10.
T v is the atmospheric transmittance in the direction of the view; 11.
D is the Earth-Sun distance in astronomical units.
In the atmospheric correction technique, DOS1 in particular, it is assumed that: The result is a path radiance equal to: The reflectance at the ground, therefore, is: where L λ is the radiance to the satellite, and the values of ESUN, the mean solar exoatmospheric irradiance, expressed in W/(m 2 ¨µm) for the sensors Landsat 4 [32] and Landsat 5 [33] are defined in Table 1.Assuming that the function of the reflectivity of a body ρ pλq depends on the material that constitutes its surface [34], the function of reflectivity relative to each material is also recognized as its spectral signature.Through the pre-processing Band-Calc tool, the centers of wavelengths have been set according to the characteristics of the Thematic Mapper (TM) for each band, in order to view the composition of bands necessary to distinguish the various local components.

Supervised Classification with the SMAP Algorithm
After the pre-processing phases, it is possible to continue with the supervised classification, an image processing technique that allows the identification of components according to their spectral signatures, distinct for each class of pixels and requiring the supervision of an operator [35].A Region of Interest (ROI) was therefore identified.ROI are polygons automatically created with a region-growing algorithm starting from seed pixels around which it is possible to segment the image, circumscribing each group of pixels belonging to the same land cover class within an area of interest and classified according to their spectral signature.
ROI, homogeneously identified across the whole area, has three parameters to calibrate: ‚ Min size ROI identifies the minimum area of ROI in pixels, if created with the region-growing algorithm; ‚ Max size ROI identifies the length, in pixels, of the side of a square that inscribes the ROI; ‚ a spectral distance value between seed pixels and the surrounding pixels, i.e., the range radius expressed in radiometric units.
Eight land use classes were identified in total (Figure 2).In this phase, orthophotos concerning the considered dates have been used, in order to limit the error level on the areas to classify.The sequential maximum estimation algorithm was used to classify images, which is an a posteriori algorithm (SMAP).It is part of the i.smap GRASS module and represents one of its operating modes with the classification that uses the Maximum Likelihood (ML).


Max size ROI identifies the length, in pixels, of the side of a square that inscribes the ROI;  a spectral distance value between seed pixels and the surrounding pixels, i.e., the range radius expressed in radiometric units.
Eight land use classes were identified in total (Figure 2).In this phase, orthophotos concerning the considered dates have been used, in order to limit the error level on the areas to classify.The sequential maximum estimation algorithm was used to classify images, which is an a posteriori algorithm (SMAP).It is part of the i.smap GRASS module and represents one of its operating modes with the classification that uses the Maximum Likelihood (ML).
The SMAP algorithm is based on the first law of geography by Waldo Tobler [36]: the likelihood that close pixels belong to the same class is higher than the one you would observe with distant pixels.The main reason for the choice of SMAP is the vast extension of the study area and the capability of the algorithm to implement segmentation directly on image regions and not on separately-considered pixels.Generally, Bayesian methods attempt to minimize the average value of wrong segmentation, optimizing the solution of the problem, considering: where C(X, x) represents the value to estimate on the actual segmentation X on the value of the approximate segmentation x.
The   function allows one to gradually assign a higher value to segmentations that include a larger number of misclassified pixels. where: behavior is defined at the coarse scale by K, which contains misclassified pixels.Therefore, if K is the only value:  () ≠  () , but  () =  () for all  > , then: The SMAP algorithm is based on the first law of geography by Waldo Tobler [36]: the likelihood that close pixels belong to the same class is higher than the one you would observe with distant pixels.The main reason for the choice of SMAP is the vast extension of the study area and the capability of the algorithm to implement segmentation directly on image regions and not on separately-considered pixels.
It operates at different scales and resolutions, and it uses assessments at wider scales to drive assessments at scales that are more detailed.It enables one to reduce the spatial dimension of errors and to prevent the homologation of neighboring pixels to the same class, by reducing smoothing operations, if their position in the image changes often [37][38][39].
Generally, Bayesian methods attempt to minimize the average value of wrong segmentation, optimizing the solution of the problem, considering: where C(X, x) represents the value to estimate on the actual segmentation X on the value of the approximate segmentation x.
The C SMAP function allows one to gradually assign a higher value to segmentations that include a larger number of misclassified pixels. where: C SMAP behavior is defined at the coarse scale by K, which contains misclassified pixels.Therefore, if K is the only value: K pKq ‰ x pKq , but X piq " x piq for all i ą K, then: In this way, the total value of the C SMAP function will be equal to: Finally, errors defined by K, at a coarse scale, produce errors in classification at a finer scale; thus, the value of SMAP will be intuitively defined as the width of the largest grouping of misclassified pixels.In this way, it is possible to determine the parameter that minimizes the value of a wrong segmentation, through the following relation: " argmin ÿ L n"0 2 n P ´Xpiq " x piq i ě n ˇˇY " y Since random fields X pnq form a Markov chain, it is possible to recursively calculate this estimation in parameter of n scale.Assuming that x ˆpiq was calculated for i ą n, these results will be used to calculate x ˆpnq .

Accuracy Verification
Once the classification process is completed, in order to evaluate the precision of the image and, therefore, the reliability compared to what has been classified, the accuracy has been verified.
Correct pixels are evaluated through the error matrix, omission errors (pixels of the ground truth incorrectly identified) and commission errors (pixels of the classified map that do not match with the ground truth).
The output is a squared matrix where rows and columns are in equal number compared to the assigned classes of land cover.Accuracy assessment is performed comparing information contained in the map and data related to the ground truth.It is possible to calculate the overall accuracy as the ratio between the number of correctly-classified samples, obtained summing the values of the main diagonal, and the total number of sampling units [40] (the perfect matrix contains non-null values only on the main diagonal).In this phase, in addition to the K coefficient, values to be considered are the producer's accuracy calculated along the columns of the matrix and the user's accuracy calculated along the rows of the matrix.
The overall accuracy is defined by the ratio between the sum of values on the main diagonal and the sum of total values defining the matrix, and it is expressed by: where q is the number of land cover classes, while Pkk the calculated parameters.The user's accuracy defines the probability that a pixel belongs to a given class because the classification algorithm has labelled the pixels in the given class, and it is expressed by: P Ui " P ii {P i` ( 14) where P ii defines the probability of a single cell and P i`i s the sum of marginal elements in the rows.
The producer's accuracy defines the probability that the classification algorithm has labelled a given pixel of the image in a given class because the ground truth detects it as belonging to that class, and it is expressed by: P Aj " P jj {P `j (15) where P jj defines the probability of a single cell and P `j is the sum of marginal elements in the columns.The K coefficient [41] measures the ratio between the probability of expected agreement and disagreement.This measure uses all of the elements of the matrix and not only the ones in the diagonal, and it is expressed by:

Data
Remote sensing technologies implemented in a GIS environment allow developing a detailed analysis of land cover changes during time.The classification step was performed with i.smap GRASS GIS, used in the QGIS software.Specifically, QGIS 2.12.3 Lyon was used.
In the first phase, land cover classification of the study area, Landsat 4-5 TM (Thematic Mapper) satellite images have been used.In order to have greater continuity for all four dates, summer periods were considered: Downloaded images had a .tiffile extension for each Landsat band and an MTL.txt file extension containing metadata information; moreover, images were already georeferenced in WGS84.These images were composed of seven bands with a spatial resolution of 30 m ˆ30 m (except for Band 6, with a spatial resolution of 120 m ˆ120 m), with a 185-km swath, a radiometric resolution of 8-bit and a temporal resolution of 16 days between each repeating cycle.
The spectral resolution of the bands and their corresponding wavelengths' centers according to the characteristics of Thematic Mapper (TM) are represented in the table below (Table 2).
During the ROI definition phase, color and black and white digital orthophotos at the corresponding dates have been used as support for the identification of the areas, with related footprint polygons with the date information of the aerial shot.These images were acquired by the Italian national geoportal with a resolution of 1:10,000.

Supervised Classification with the SMAP Algorithm
Supervised classification was carried out after atmospheric correction DOS1, which provides value changes in bands, since it eliminates distortions caused by atmosphere presence.Subsequently, the centers of wavelengths according to the characteristics of Thematic Mapper (TM) for each band were set.
The purpose was to visualize the composition of bands necessary to distinguish territorial components (RGB natural color 3-2-1, synthesis IR false color 4-3-2, false color image 5-4-3) with the exception of thermal infrared Band 6.
The use of supervised classification allowed a better interpretation of images and a more accurate territory classification into classes.Eight land cover classes have been identified for each date, defined as follows: 1.
Urban areas and artificial lands; 5.
The geographical area under consideration is particularly large; therefore, it was necessary to define a large number of homogeneously-distributed ROI on the entire image, in order to obtain a better response in terms of classification.
Maps generated using the SMAP algorithm (Figure 3) allowed a good representation of the above-mentioned eight classes.
By comparing the classifications for four dates (Figure 4), it is possible to appreciate variations in the development of urban areas and land use in a period of 25 years.Despite some difficulties related to the similarity between spectral signatures predominant in urban areas and bare soil, the distinction between the different areas is easily readable from the images.
The problem of the similarity of the spectral signature, mainly identifiable between urban areas and bare soils, rather than between urban areas and other land covers, is indistinctly present in the images for the four dates (Figure 5).define a large number of homogeneously-distributed ROI on the entire image, in order to obtain a better response in terms of classification.
Maps generated using the SMAP algorithm (Figure 3) allowed a good representation of the above-mentioned eight classes.By comparing the classifications for four dates (Figure 4), it is possible to appreciate variations in the development of urban areas and land use in a period of 25 years.Despite some difficulties related to the similarity between spectral signatures predominant in urban areas and bare soil, the distinction between the different areas is easily readable from the images.The problem of the similarity of the spectral signature, mainly identifiable between urban areas and bare soils, rather than between urban areas and other land covers, is indistinctly present in the images for the four dates (Figure 5).The problem of the similarity of the spectral signature, mainly identifiable between urban areas and bare soils, rather than between urban areas and other land covers, is indistinctly present in the images for the four dates (Figure 5).

Accuracy Verification
From the accuracy verification and observation of the error matrices of each classification, it was possible to highlight the overall accuracy of each classification, as well as the producer's accuracy, user's accuracy and K coefficient (Table 3 summarizes the verifications carried out for all four considered dates).

Accuracy Verification
From the accuracy verification and observation of the error matrices of each classification, it was possible to highlight the overall accuracy of each classification, as well as the producer's accuracy, user's accuracy and K coefficient (Table 3 summarizes the verifications carried out for all four considered dates).Compared to such values and to the proportion of correctly-classified objects over the total amount, among all classifications, the average percentage value of incorrectly-calculated objects was around 9%.A single peak was detected in classification for 1993, where about 13% of items were not correctly classified.In general, percentages describing total classification accuracy never dropped below 80%, thus allowing asserting the trustworthiness of statements, considering errors arising from the presence of objects with very similar spectral signatures (urban and bare soil, urban and plowed soils) and from the minimum percentage of relegated objects in classes different from those in which they were placed.

Change Detection
The assessment of changes occurred in the area between 1985 and 2010, compared to the identified land cover, has highlighted how, in 25 years, urban areas have increased at the expense of agricultural land, rather than bare soil.This datum was already visible observing surfaces (km 2 ) and soil coverage changes in terms of percentages measured for the considered date (Table 4): these variations were most appreciable in a graph (Figure 6).Such an increase occurred at the expense of areas occupied by permanent grassland and in lower percentages by bare soil, as measured by the transition matrices, where there was the migration of objects from one class to another, compared for two dates.Considering other classes, the data on water are interesting, because their percentage decreases from 0.22% in 1985 to 0.12% in 1993, increasing again in 2002 and 2010 up to 0.17% and 0.26%, respectively.
The label assigned to the water class includes permanent water bodies, such as two clearly visible lakes, Vulture and Serracorvo, and non-permanent water bodies, such as small dams and small artificial reservoirs, where the percentage of wet surface is subject to change over time, because it depends also on the amount of precipitation that affects the area.
On these surfaces, an increase in vegetated areas corresponds to a decrease in water surfaces.In other words, these variations in terms of land cover changes are correlated and do not represent sampling errors, as also supported by reliability levels in both classes.The migration of objects belonging to the sparse vegetation class to cultivated soil and vice versa highlights a noticeable increase of cultivated soils in 2002, with occupied areas equal to 533.07 km 2 .This could be explained in two different ways: on the one hand, it describes the actual possibility that sparsely-vegetated soils are converted into productive soils; on the other hand, it takes into account the possibility that these variations may be the result of misclassification due to very similar spectral signatures.Observing the comparison between the extent of both classes in 1985 and 2010 and their transition matrix (Table Considering that the study pays particular attention to the development of urbanized areas, it is interesting to note that the urban class for the four considered dates undergoes a steady increase from 20.62 km 2 in 1985 to 37.41 km 2 in 2010.
Such an increase occurred at the expense of areas occupied by permanent grassland and in lower percentages by bare soil, as measured by the transition matrices, where there was the migration of objects from one class to another, compared for two dates.
Considering other classes, the data on water are interesting, because their percentage decreases from 0.22% in 1985 to 0.12% in 1993, increasing again in 2002 and 2010 up to 0.17% and 0.26%, respectively.
The label assigned to the water class includes permanent water bodies, such as two clearly visible lakes, Vulture and Serracorvo, and non-permanent water bodies, such as small dams and small artificial reservoirs, where the percentage of wet surface is subject to change over time, because it depends also on the amount of precipitation that affects the area.
On these surfaces, an increase in vegetated areas corresponds to a decrease in water surfaces.In other words, these variations in terms of land cover changes are correlated and do not represent sampling errors, as also supported by reliability levels in both classes.The migration of objects belonging to the sparse vegetation class to cultivated soil and vice versa highlights a noticeable increase of cultivated soils in 2002, with occupied areas equal to 533.07 km 2 .This could be explained in two different ways: on the one hand, it describes the actual possibility that sparsely-vegetated soils are converted into productive soils; on the other hand, it takes into account the possibility that these variations may be the result of misclassification due to very similar spectral signatures.Observing the comparison between the extent of both classes in 1985 and 2010 and their transition matrix (Table 5), it is possible to confirm that approximately 17 km 2 of land area have been urbanized and that there was an increase of 1% in agricultural and permanent meadow soils (Figure 7).Data on the development of urban areas during the considered period represent a true estimate, but curtailed by the presence of road infrastructures, which are fully included within the urban conversion of the territory, but are understood as available land.Consequently, a significant portion regarding secondary rural, municipal and local roads has to be added to the percentage, apparently insignificant, which considers an increase of 1% in the territory available to the urban environment, Data on the development of urban areas during the considered period represent a true estimate, but curtailed by the presence of road infrastructures, which are fully included within the urban conversion of the territory, but are understood as available land.Consequently, a significant portion regarding secondary rural, municipal and local roads has to be added to the percentage, apparently insignificant, which considers an increase of 1% in the territory available to the urban environment, between 1985 and 2010.

Discussion
The classification of Landsat images shows the interdependence between the growth of urban areas and the presence of mobility infrastructure (Figure 10).

Discussion
The classification of Landsat images shows the interdependence between the growth of urban areas and the presence of mobility infrastructure (Figure 10).The classification did not allow highlighting these elements in the evaluations, because the large representation scale of the images did not allow their identification, and the study area was composed of small towns having very narrow roads.An expected discrepancy occurred when comparing the increase in urban areas between 1985 and 2010 (Figure 8) with data provided by the Italian National Institute of Statistics, related to the censuses of the resident population in the 19 municipalities of the mountain district of Vulture Alto-Bradano.Although the population decreased from 102,110 inhabitants in 1985 to 96,679 inhabitants in 2010, urbanized areas increased, according to the previously-defined amount (Figure 9).

Discussion
The classification of Landsat images shows the interdependence between the growth of urban areas and the presence of mobility infrastructure (Figure 10).Abandonment of agricultural activity, vehicular congestion in denser urban areas and socio-economic changes lead to the building of housing outside town centers.Problems of accessibility due to urban sprawl, however, urge local authorities to invest in infrastructures and in the completion of newly-developed settlements.This entails new problems of congestion and a new sprawl on soils outside town centers; therefore, the phenomenon continues to repeat itself.This cyclical expansion of urban areas, not necessarily driven by a demand for housing, is one of the major causes of the unsustainability of current land use patterns, and in particular of building activity.
The case of Vulture-Alto Bradano also highlights an issue that is often recurrent in Italian scenarios.Urban planning tools, obsolete or outdated, indeed, usually drive the growth of urban areas.In the Basilicata region, the Law on Urban Planning 1999/23 required the renewal of all municipal planning tools of the region.However, in Vulture-Alto Bradano, only two municipalities have taken steps to upgrade their plans (Figure 11), to date, while in the remaining municipalities, urban development still takes place in accordance with the instruments adopted in the 1980s and, therefore, no longer useful to answer to the changed social and citizenship needs.Abandonment of agricultural activity, vehicular congestion in denser urban areas and socio-economic changes lead to the building of housing outside town centers.Problems of accessibility due to urban sprawl, however, urge local authorities to invest in infrastructures and in the completion of newly-developed settlements.This entails new problems of congestion and a new sprawl on soils outside town centers; therefore, the phenomenon continues to repeat itself.This cyclical expansion of urban areas, not necessarily driven by a demand for housing, is one of the major causes of the unsustainability of current land use patterns, and in particular of building activity.
The case of Vulture-Alto Bradano also highlights an issue that is often recurrent in Italian scenarios.Urban planning tools, obsolete or outdated, indeed, usually drive the growth of urban areas.In the Basilicata region, the Law on Urban Planning 1999/23 required the renewal of all municipal planning tools of the region.However, in Vulture-Alto Bradano, only two municipalities have taken steps to upgrade their plans (Figure 11), to date, while in the remaining municipalities, urban development still takes place in accordance with the instruments adopted in the 1980s and, therefore, no longer useful to answer to the changed social and citizenship needs.
The case of Vulture-Alto Bradano also highlights an issue that is often recurrent in Italian scenarios.Urban planning tools, obsolete or outdated, indeed, usually drive the growth of urban areas.In the Basilicata region, the Law on Urban Planning 1999/23 required the renewal of all municipal planning tools of the region.However, in Vulture-Alto Bradano, only two municipalities have taken steps to upgrade their plans (Figure 11), to date, while in the remaining municipalities, urban development still takes place in accordance with the instruments adopted in the 1980s and, therefore, no longer useful to answer to the changed social and citizenship needs.Similarly, the protection of agricultural and natural landscapes in the Basilicata region is entrusted to landscape plans.These, however, do not study the entire region, but only some of its portions (Figure 12).Specifically, only the area of Mount Vulture and Monticchio volcanic lakes has specific landscape protection rules inside the study area.Similarly, the protection of agricultural and natural landscapes in the Basilicata region is entrusted to landscape plans.These, however, do not study the entire region, but only some of its portions (Figure 12).Specifically, only the area of Mount Vulture and Monticchio volcanic lakes has specific landscape protection rules inside the study area.
The lack of planning tools able to adequately protect natural and environmental systems and to preserve landscape exposes the territory to unregulated and often unsustainable transformations and transitions between land use classes, without an organic vision.The monitoring of only land use changes thus fails to ensure the development of more sustainable usage patterns of the territory.It is therefore necessary for the government to establish updated and consistent urban and regional planning tools [42].

Conclusions
This paper has discussed the importance of land use change monitoring, as a tool for the definition of sustainable government policies for the territory.In particular, a methodology for the development of the monitoring of land use changes through open source software and open data was presented.Therefore, the methodology is characterized by being highly accessible and replicable by research centers and public administrations at a reduced cost.This aspect is very important, considering the limited economical resources of local agencies and authorities.
Land use and land cover changes, as well as the phenomenon of land take have been analyzed through a monitoring method based on remote sensing techniques to perform multi-temporal evaluations to provide important support for planning choices [43][44][45].The use of satellite imagery (Landsat TM 4-5) and remote sensing techniques in a GIS environment allows defining various parameters, through supervised classification, in order to discriminate urban areas from other classes.
The methodology requires some preventive steps in order to clearly distinguish elements composing the image: atmospheric correction (in the specific case DOS1), the perimeter of the area of The lack of planning tools able to adequately protect natural and environmental systems and to preserve landscape exposes the territory to unregulated and often unsustainable transformations and transitions between land use classes, without an organic vision.The monitoring of only land use changes thus fails to ensure the development of more sustainable usage patterns of the territory.It is therefore necessary for the government to establish updated and consistent urban and regional planning tools [42].

Conclusions
This paper has discussed the importance of land use change monitoring, as a tool for the definition of sustainable government policies for the territory.In particular, a methodology for the development of the monitoring of land use changes through open source software and open data was presented.Therefore, the methodology is characterized by being highly accessible and replicable by research centers and public administrations at a reduced cost.This aspect is very important, considering the limited economical resources of local agencies and authorities.
Land use and land cover changes, as well as the phenomenon of land take have been analyzed through a monitoring method based on remote sensing techniques to perform multi-temporal evaluations to provide important support for planning choices [43][44][45].The use of satellite imagery (Landsat TM 4-5) and remote sensing techniques in a GIS environment allows defining various parameters, through supervised classification, in order to discriminate urban areas from other classes.
The methodology requires some preventive steps in order to clearly distinguish elements composing the image: atmospheric correction (in the specific case DOS1), the perimeter of the area of interest and the definition of wavelengths' centers related to the sensor used in remote sensing (in this specific case, the Thematic Mapper).After these preventive steps, it is possible to catalogue land cover classes, defining ROI, and to implement the classification by means of algorithms (in this case, SMAP was used).
The reliability of the method was evaluated in the area of Vulture Alto-Bradano (Basilicata region, Southern Italy), a mountain district located in the northeast of the region, including 19 municipalities with strong historical, natural and landscape connotations.
The possibility to implement a historical evaluation for four dates, using Landsat satellite images, allowed estimating that the increase in urban areas from 1985-2010 was approximately 17 km 2 .Compared to the general trend, this result is lower than the average.However, the mountain district of Vulture Alto-Bradano follows a common trend in Southern Italy; in fact, here, it is possible to observe more than anywhere else that the abandonment of farmlands and small urban centers and at the same time coupled to an increase of the urban sprawl phenomenon [46].
The methodology adopted in this paper highlights a possible way to check land take, using free data and software.The great advantage is that such a methodology can create important support for planning choices at all scales.However, future developments will have to evaluate procedures to overcome some problems evident in the phase of supervised classification, such as the difficulties for classification algorithms to distinguish among artificial surfaces and bare soils, due to the high similarity between the spectral signatures of these elements.Despite these problems, the proposed approach is in any case highly replicable, especially for its low costs, since it is entirely developed using open source software and through open data.

Figure 1 .
Figure 1.Study area location: Italy, the Basilicata region and the Vulture Alto-Bradano area.

Figure 1 .
Figure 1.Study area location: Italy, the Basilicata region and the Vulture Alto-Bradano area.

Figure 2 .
Figure 2. Flowchart of the proposed methodology.DOS, Dark Object Subtraction.

Figure 2 .
Figure 2. Flowchart of the proposed methodology.DOS, Dark Object Subtraction.

Figure 3 .
Figure 3. Classification of Vulture Alto-Bradano in 1985 by means of the SMAP algorithm.

Figure 4 .
Figure 4. (a-d) Classification of the Vulture Alto-Bradano area for four dates.

Figure 4 .
Figure 4. (a-d) Classification of the Vulture Alto-Bradano area for four dates.

Figure 4 .
Figure 4. (a-d) Classification of the Vulture Alto-Bradano area for four dates.

Figure 5 .
Figure 5. Trend of the spectral signatures of the samples of each cover soil class for the four investigated dates.

Figure 5 .
Figure 5. Trend of the spectral signatures of the samples of each cover soil class for the four investigated dates.

Figure 9 .
Figure 9.Comparison of the increase in population and the increase in urbanized areas between 1985 and 2010; elaboration on the data from the Italian National Institute of Statistics.

Figure 9 .
Figure 9.Comparison of the increase in population and the increase in urbanized areas between 1985 and 2010; elaboration on the data from the Italian National Institute of Statistics.

Figure 9 .
Figure 9.Comparison of the increase in population and the increase in urbanized areas between 1985 and 2010; elaboration on the data from the Italian National Institute of Statistics.

Figure 10 .
Figure 10.Infrastructures in the study area and urban areas/artificial lands in 2010.

Figure 10 .
Figure 10.Infrastructures in the study area and urban areas/artificial lands in 2010.

Figure 11 .
Figure 11.Current land use zoning regulation.

Figure 11 .
Figure 11.Current land use zoning regulation.

Table 4 .
Soil covers in 1985Soil covers in  , 1993Soil covers in  , 2002Soil covers in   and 2010.   .

Table 5 .
Transition matrix of objects classified in 1985 and 2010.

Table 5 .
Transition matrix of objects classified in 1985 and 2010.