Comparison of stochastic and deterministic methods for mapping groundwater level spatial variability in sparsely monitored basins

In sparsely monitored basins, accurate mapping of the spatial variability of groundwater level requires the interpolation of scattered data. This paper presents a comparison of deterministic interpolation methods, i.e. inverse distance weight (IDW) and minimum curvature (MC), with stochastic methods, i.e. ordinary kriging (OK), universal kriging (UK) and kriging with Delaunay triangulation (DK). The study area is the Mires Basin of Mesara Valley in Crete (Greece). This sparsely sampled basin has limited groundwater resources which are vital for the island’s economy; spatial variations of the groundwater level are important for developing management and monitoring strategies. We evaluate the performance of the interpolation methods with respect to different statistical measures. The Spartan variogram family is applied for the first time to hydrological data and is shown to be optimal with respect to stochastic interpolation of this dataset. The three stochastic methods (OK, DK and UK) perform overall better than the deterministic counterparts (IDW and MC). DK, which is herein for the first time applied to hydrological data, yields the most accurate cross-validation estimate for the lowest value in the dataset. OK and UK lead to smooth isolevel contours, whilst DK and IDW generate more edges. The stochastic methods deliver estimates of prediction uncertainty which becomes highest near the southeastern border of the basin.     

The geostatistic-based spatial distribution variations of soil salts under long-term wastewater irrigation

The purpose of this study was to determine and evaluate the spatial changes in soil salinity by using geostatistical methods. The study focused on the suburb area of Beijing, where urban development led to water shortage and accelerated wastewater reuse to farm irrigation for more than 30 years. The data were then processed by GIS using three different interpolation techniques of ordinary kriging (OK), disjunctive kriging (DK), and universal kriging (UK). The normality test and overall trend analysis were applied for each interpolation technique to select the best fitted model for soil parameters. Results showed that OK was suitable for soil sodium adsorption ratio (SAR) and Na⁺ interpolation; UK was suitable for soil Cl^? and pH; DK was suitable for soil Ca²⁺. The nugget-to-sill ratio was applied to evaluate the effects of structural and stochastic factors. The maps showed that the areas of non-saline soil and slight salinity soil accounted for 6.39 and 93.61 %, respectively. The spatial distribution and accumulation of soil salt were significantly affected by the irrigation probabilities and drainage situation under long-term wastewater irrigation.     

Spatial assessment of soil salinity in the Harran Plain using multiple kriging techniques

The Harran Plain is located in the southeastern part of Turkey and has recently been developed for irrigation agriculture. It already faces soil salinity problems causing major yield losses. Management of the problem is hindered by the lack of information on the extent and geography of the salinization problem. A survey was carried out to delineate the spatial distribution of salt-affected areas by randomly selecting 140 locations that were sampled at two depths (0 to 30 and 30 to 60 cm) and analyzed for soil salinity variables: soil electrical conductivity (EC), soluble cations (Ca^2+, Mg²⁺, Na⁺, and K⁺), soluble anions (SO ₄ ^2? , Cl^?), exchangeable Na⁺ (me 100 g^?1) and exchangeable sodium percentage. Terrain attributes (slope, topographical wetness index) were extracted from the digital elevation model of the study area. Variogram analyses after log transformation and ordinary kriging (OK) were applied to map spatial patterns of soil salinity variables. Multivariate geostatistical methods—regression kriging (RK) and kriging with external drift (KED)—were used using elevation and soil electrical conductivity data as covariates. Performances of the three estimation methods (OK, RK, and KED) were compared using independent validation samples randomly selected from the main dataset. Soils were categorized into salinity classes using disjunctive kriging (DK) and ArcGIS, and classification accuracy was tested using the kappa statistic. Results showed that soil salinity variables all have skewed distribution and are poorly correlated with terrain indices but have strong correlations among each other. Up to 65 % improvement was obtained in the estimations of soil salinity variables using hybrid methods over OK with the best estimations obtained with RK using EC_0–30 as covariate. DK–ArcGIS successfully classified soil samples into different salinity groups with overall accuracy of 75 % and kappa of 0.55 (p?<?0.001).     

Evaluating Kriging as a Tool to Improve Moderate Resolution Maps of Forest Biomass

The USDA Forest Service, Forest Inventory and Analysis program (FIA) recently produced a nationwide map of forest biomass by modeling biomass collected on forest inventory plots as nonparametric functions of moderate resolution satellite data and other environmental variables using Cubist software. Efforts are underway to develop methods to enhance this initial map. We explored the possibility of modeling spatial structure to make such improvements. Spatial structure in the field biomass data as well as in residuals from the map was investigated across 18 ecological zones in the Interior Western U.S. Exploratory tools included directional graphs of summary statistics, three dimensional maps, Moran’s I correlograms, and variograms. Where spatial pattern was present, field and residual biomass were kriged, and predictions made for an independent test set were evaluated for improvement over predictions in the initial biomass map. While kriging has some potential benefit when analyzing the field data and exploring spatial structure, kriging residuals resulted in little or no improvement in the initial biomass map developed using Cubist software. Stationarity assumptions, variogram behavior, and appropriate model fitting strategies are discussed.     

Assessment of groundwater pollution in West Delhi, India using geostatistical approach

The exploration, exploitation, and unscientific management of groundwater resources in the National Capital Territory (NCT) of Delhi, India have posed a serious threat of reduction in quantity and deterioration of quality. The objective of the study is to determine the groundwater quality and to assess the risk of groundwater pollution at Najafgarh, NCT of Delhi. The groundwater quality parameters were analyzed from the existing wells of the Najafgarh and the thematic maps were generated using geostatistical concepts. Ordinary kriging and indicator kriging methods were used as geostatistical approach for preparation of thematic maps of the groundwater quality parameters such as bicarbonate, calcium, chloride, electrical conductivity (EC), magnesium, nitrate, sodium, and sulphate with concentrations equal or greater than their respective groundwater pollution cutoff value. Experimental semivariogram values were fitted well in spherical model for the water quality parameters, such as bicarbonate, chloride, EC, magnesium, sodium, and sulphate and in exponential model for calcium and nitrate. The thematic maps of all the groundwater quality parameters exhibited an increasing trend of pollution from the northern and western part of the study area towards the southern and eastern part. The concentration was highest at the southernmost part of the study area but it could not reflect correctly the groundwater pollution status. The indicator kriging method is useful to assess the risk of groundwater pollution by giving the conditional probability of concentrations of different chemical parameters exceeding their cutoff values. Thus, risk assessment of groundwater pollution is useful for proper management of groundwater resources and minimizing the pollution threat.     

Forest ecological classification and mapping: Their application for ecosystem management in Newfoundland

A prerequisite to sustaining ecosystems is the inventory and classification of landscape structure and composition. Ecological classification and mapping involves the delineation of landscapes into easily recognizable units. Topography, soils, vegetation, physical landscape form, and successional pathways are delineation criteria commonly used.Damman (1967) developed a forest type classification system for Newfoundland using vegetation, soil and landforms as the defining criteria. Damman's forest types were used in combination with mensurational data to assign forest types to timber volume productivity classes. Since each of the Damman forest types is associated with characteristic soils, parent materials, moisture regime and topographic position, the mapping units are similar to Canada Land Inventory (CLI) mapping units. Field work to confirm the correlation between Damman forest types and CLI capability classes was initiated in 1993. CLI maps were recoded in 1994 and Damman forest types were determined; resulting ecosystem-based maps provide a common framework to assess forestry/wildlife interactions in an ecosystem planning process.     

Indicator and probability kriging methods for delineating Cu, Fe, and Mn contamination in groundwater of Najafgarh Block, Delhi, India

Two non-parametric kriging methods such as indicator kriging and probability kriging were compared and used to estimate the probability of concentrations of Cu, Fe, and Mn higher than a threshold value in groundwater. In indicator kriging, experimental semivariogram values were fitted well in spherical model for Fe and Mn. Exponential model was found to be best for all the metals in probability kriging and for Cu in indicator kriging. The probability maps of all the metals exhibited an increasing risk of pollution over the entire study area. Probability kriging estimator incorporates the information about order relations which the indicator kriging does not, has improved the accuracy of estimating the probability of metal concentrations in groundwater being higher than a threshold value. Evaluation of these two spatial interpolation methods through mean error (ME), mean square error (MSE), kriged reduced mean error (KRME), and kriged reduced mean square error (KRMSE) showed 3.52% better performance of probability kriging over indicator kriging. The combined result of these two kriging method indicated that on an average 26.34%, 65.36%, and 99.55% area for Cu, Fe, and Mn, respectively, are coming under the risk zone with probability of exceedance from a cutoff value is 0.6 or more. The groundwater quality map pictorially represents groundwater zones as ??desirable?? or ??undesirable?? for drinking. Thus the geostatistical approach is very much helpful for the planners and decision makers to devise policy guidelines for efficient management of the groundwater resources so as to enhance groundwater recharge and minimize the pollution level.     

Mapping the Spatial Variability of Plant Diversity in a Tropical Forest: Comparison of Spatial Interpolation Methods

Knowledge of the spatial distribution of plant species is essential to conservation and forest managers in order to identify high priority areas such as vulnerable species and habitats, and designate areas for reserves, refuges and other protected areas. A reliable map of the diversity of plant species over the landscape is an invaluable tool for such purposes. In this study, the number of species, the exponent Shannon and the reciprocal Simpson indices, calculated from 141 quadrat sites sampled in a tropical forest were used to compare the performance of several spatial interpolation techniques used to prepare a map of plant diversity, starting from sample (point) data over the landscape. Means of mapped classes, inverse distance functions, kriging and co-kriging, both, applied over the entire studied landscape and also applied within vegetation classes, were the procedures compared. Significant differences in plant diversity indices between classes demonstrated the usefulness of boundaries between vegetation types, mapped through satellite image classification, in stratifying the variability of plant diversity over the landscape. These mapped classes, improved the accuracy of the interpolation methods when they were used as prior information for stratification of the area. Spatial interpolation by co-kriging performed among the poorest interpolators due to the poor correlation between the plant diversity variables and vegetation indices computed by remote sensing and used as covariables. This indicated that the latter are not suitable covariates of plant diversity indices. Finally, a within-class kriging interpolator yielded the most accurate estimates of plant diversity values. This interpolator not only provided the most accurate estimates by accounting for the indices' intra-class variability, but also provided additional useful interpretations of the structure of spatial variability of diversity values through the interpretation of their semi-variograms. This additional role was found very useful in aiding decisions in conservation planning.     

Mapping Soil Gas Radon Concentration: A Comparative Study of Geostatistical Methods

Understanding soil gas radon spatial variations can allow the constructor of a new house to prevent radon gas flowing from the ground. Indoor radon concentration distribution depends on many parameters and it is difficult to use its spatial variation to assess radon potential. Many scientists use to measure outdoor soil gas radon concentrations to assess the radon potential. Geostatistical methods provide us a valuable tool to study spatial structure of radon concentration and mapping. To explore the structure of soil gas radon concentration within an area in south Italy and choice a kriging algorithm, we compared the prediction performances of four different kriging algorithms: ordinary kriging, lognormal kriging, ordinary multi-Gaussian kriging, and ordinary indicator cokriging. Their results were compared using an independent validation data set. The comparison of predictions was based on three measures of accuracy: (1) the mean absolute error, (2) the mean-squared error of prediction; (3) the mean relative error, and a measure of effectiveness: the goodness-of-prediction estimate. The results obtained in this case study showed that the multi-Gaussian kriging was the most accurate approach among those considered. Comparing radon anomalies with lithology and fault locations, no evidence of a strict correlation between type of outcropping terrain and radon anomalies was found, except in the western sector where there were granitic and gneissic terrain. Moreover, there was a clear correlation between radon anomalies and fault systems.     

Development of Ecological Land Classification and mapping in support of forest management in northern Newfoundland,Canada

For the sustainable development of forest land, as recently prescribed by the Canadian Forest Strategy, a land classification project in northern Newfoundland was initiated to support the local forest management activities. The method adopted here is a modification of the Canadian Committee for Ecological Land Classification's (CCELC) system, and it applies various levels of mapping to uniform areas based on geomorphology, soils, vegetation, climate, water, and fauna.In this study, all CCELC levels were mapped; resulting maps were digitized and imported into a Geographic Informations System (GIS). The GIS data base contained the following maps: 1) digital terrain model, 2) bedrock geology, 3) surficial geology, 4) forest inventory, and 5) various levels of the ecological land classification, including Vegetation Types at the lowest level. In addition to the mapping, mensurational data were analyzed to provide stand and stock tables for each of the forest types, including growth curves that could be entered into specific forest growth modelling systems to predict wood supply scenarios based upon different management interventions.     

Multi-scale spatial structure of heavy metals in agricultural soils in Beijing

To effectively investigate the spatial variability of heavy metals in soil, produce a higher quality spatial distribution map, and identify the potential pollution sources of heavy metals, geostatistics was employed to evaluate the effect of scale on spatial variability of heavy metals in Beijing agricultural soils. The results revealed that spatial variability of Cr, Ni, Zn, and Hg was dependent on scale. Validation of the optimality of theoretical semivariance and comparative analysis of the estimation accuracy demonstrated that the multi-scale nested model can reveal the spatial structure of heavy metals effectively and improve the estimation accuracy better than the single-scale method, thereby enabling production a higher quality spatial interpolation map. Thus, the multi-scale kriging nested model is a useful tool for revealing spatial variability of heavy metals in soils, while the spatial distribution maps allow the identification of hot spots with high concentrations of heavy metals.     

Factors influencing acquisition of ecological and exposure information about hazards and risks from contaminated sites

Considerable research indicates that a wide range of socio-economic factors influence attitudes and perceptions about environmental hazards and risks, and that social trust in those who manage a hazard is strongly correlated to judgements about risks and benefits. We suggest that there are three steps that lead to environmental risk perceptions: acquisition of information, interpretation and synthesis of different pieces of information, and understanding of that information in light of previous knowledge, perceptions, or attitudes. In this study we presented 211 college students in the sciences and non-sciences with ecological and exposure information using text, tables and maps to examine the factors that affect information acquisition and interpretation concerning ecological issues at a fictitious hazardous waste site. Students were allowed about an hour to read the materials and answer questions. The percent of students answering each question correctly varied from 4 to 82%, indicating that some questions were extremely difficult to answer. We attributed these differences to variations in the number of places information was presented (in text, tables, maps, or a combination) and the complexity of the information (how many pieces of information were required to answer the question correctly). The correlation between the number of students answering each question correctly and these combined measures (frequency, complexity) was −0.72. Thus, although there were differences in accuracy concerning ecological information as presented in this study, the major differences were accounted for by how the information was presented, and how much information was required. This suggests that risk communicators should carefully determine which ecological information is critical for the target audience, and ensure that it is presented several times (in different forms). That a lower percentage of people correctly answered questions that required synthesizing several pieces of information suggests that this complexity should be reduced where possible, or that the pieces of information should be tied clearly to the conclusion. Self-declaration of effort and perceptions of safety of the site did not markedly influence the relationship between accuracy, difficulty of finding information, and complexity of information. Other possible confounding variables (i.e., science vs non-science major) only accounted for about 27% of the variation in student’s overall score on ecological questions; age, ethnicity, and gender did not enter as significant variables. We cannot manage environmental hazards with appropriate stakeholder input unless we understand how to present environmental information to achieve a full understanding. Protection of human health and the environment requires that people understand ecological and exposure information, particularly on buffer lands.     

Changes in soil microbial biomass and aggregate stability under different land uses in the northeastern Turkey

The characteristics of three neighboring soils from the NE of Turkey were evaluated in order to elucidate the effect of different land-use management on the soil aggregate stability and microbial biomass in Galyan-Atasu dam watershed. Three experimental sites corresponding to three land uses were selected. The first site is a hazelnut orchard (agriculture), the second site is a forest dominated by mature coniferous trees, and the third site is grassland. Soil aggregate stability values for the 1–2-mm aggregates increased from forest (lowest) to agriculture (highest) in the current study. The percentage of clay was highest in agriculture soils with 33.57 %, and overall stability values increased according to soil clay content. The lower aggregate stability in the forest soils probably reflects the highly silty texture soils with 11.95 % compared to agriculture and grassland. However, in our study, there were no significant correlations between aggregate stability and organic C concentrations either in cultivated or forested soils. Aggregate stability depended more on the organic matter content when the organic matter content was greater than 50 or 60 mg g^?1. Below that threshold, aggregate stability may be mainly related to clay content. Furthermore, the results confirmed that higher percentages of Cmic/Corg in agricultural soils are the result of more labile organic substrates maintained in the soil, allowing a higher microbial biomass C per unit of soil organic C. This work gives a better understanding of the relationships between land-use type and soil aggregation and allows to know the soil response to different types of management in humid environments.     

Using of high-resolution topsoil magnetic screening for assessment of dust deposition: comparison of forest and arable soil datasets

Magnetic susceptibility (κ) is an easily detectable geophysical parameter that can be used as a proxy or semi-quantitative tracer of atmospheric industrial and urban dusts deposited in topsoil. An enhanced κ value of topsoil is in many cases also associated with high concentrations of soil pollutants (mostly heavy metals). High-resolution magnetic screening of topsoil in areas of high pollution influx is a useful tool for detection of pollution “hot spots”. General and regional screening maps with a grid density of 10 or 5 km have been performed on the basis of forest topsoil measurement only. The purpose of this study was to perform high-resolution magnetic screening with different grid densities in both forested and agricultural areas (arable land). Our large study area (ca. 200 km²) was located in a relatively more polluted region of the central part of Upper Silesia, and a second (small) one (ca. 100 m²) was located in the western part of Upper Silesia, with considerably lower influx of pollution. In the framework of this study, we applied a statistical comparison of data obtained in forested areas and on arable land. The arable soil showed statistically significantly lower κ values, the result of “physical dilution” of the arable layer caused by annual ploughing. Thus arable soils must be avoided during high-resolution field measurement. From semivariograms, it was clear that the spatial correlations in forest topsoil are much stronger than in arable soil, which suggests that a denser measurement grid is required in forested areas.     

Determination of the compositional change (1999–2006) in the pine forests of British Columbia due to mountain pine beetle infestation

The current mountain pine beetle (Dendroctonus ponderosae Hopkins) outbreak in British Columbia and Alberta is the largest recorded forest pest infestation in Canadian history. We integrate a spatial hierarchy of mountain pine beetle and forest health monitoring data, collected between 1999 and 2006, with provincial forest inventory data, and generate three information products representing 2006 forest conditions in British Columbia: cumulative percentage of pine infested by mountain pine beetle, percentage of pine uninfested, and the change in the percentage of pine on the landscape. All input data were formatted to a standardized spatial representation (1 ha minimum mapping unit), with preference given to the most detailed monitoring data available at a given location for characterizing mountain pine beetle infestation conditions. The presence or absence of mountain pine beetle attack was validated using field data (n?=?2054). The true positive rate for locations of red attack damage over all years was 92%. Classification of attack severity was validated using the Kruskal gamma statistic (γ?=?0.49). Error between the survey data and field data was explored using spatial autoregressive (SAR) models, which indicated that percentage pine and year of infestation were significant predictors of survey error at α?=?0.05. Through the integration of forest inventory and infestation survey data, the total area of pine infested is estimated to be between 2.89 and 4.14 million hectares. The generated outputs add value to existing monitoring data and provide information to support management and modeling applications.     

Nutrient dynamics between sediment and overlying water in the inner part of Izmir Bay, Eastern Aegean

In this article a concept is described in order to predict and map the occurrence of benthic communities within and near the German Exclusive Economic Zone (EEZ) of the North Sea. The approach consists of two work steps: (1) geostatistical analysis of abiotic measurement data and (2) calculation of benthic provinces by means of Classification and Regression Trees (CART) and GIS-techniques. From bottom water measurements on salinity, temperature, silicate and nutrients as well as from punctual data on grain size ranges (0–20, 20–63, 63–2,000 μ) raster maps were calculated by use of geostatistical methods. At first the autocorrelation structure was examined and modelled with help of variogram analysis. The resulting variogram models were then used to calculate raster maps by applying ordinary kriging procedures. After intersecting these raster maps with punctual data on eight benthic communities a decision tree was derived to predict the occurrence of these communities within the study area. Since such a CART tree corresponds to a hierarchically ordered set of decision rules it was applied to the geostatistically estimated raster data to predict benthic habitats within and near the EEZ.     

Predicting saltwater intrusion into aquifers in vicinity of deserts using spatio-temporal kriging

The primary objective of this study was to provide a detailed framework to use the spatio-temporal kriging to model the spatio-temporal variations of salinity data and predict saltwater intrusion into freshwater aquifers in the vicinity of deserts. EC data, measured in extraction wells in the Mahvelat plain located in the Northeastern part of Iran, available from 2007 to 2013, were used to demonstrate the developed framework. The source of data was not a well-designed measurement network. Therefore, to homogenize the data, spatial analysis was used to find EC distribution in the area in each year of study. To conduct the spatial analysis, a guideline and a systematic process were developed to select an appropriate kriging method and optimize its parameters. This process can be applied to different variables. After spatial analysis of EC data for all the years of the analysis period using empirical Bayesian kriging (EBK) method with manually optimized parameters, spatio-temporal and corresponding variogram analysis was conducted using R software. This process was based on a separable product-sum model applied to the data from 2007 to 2012. The data of 2013 and the data available for the years 1999 and 2006 were used for evaluating the performance of the spatio-temporal model. The EC distribution maps, developed for different years until 2021, show a high level of EC in the north, south, and west of the study area and growing saltwater intrusion into the central freshwater aquifer. This result can be attributed to the over-exploitation of the aquifer and hydraulic head and gradient distribution in the area. The framework provided in this study for spatio-temporal analysis of unstructured EC data is useful for groundwater managers in making proper decisions.     

The application of satellite data for the quantification of mangrove loss and coastal management in the Godavari estuary, East Coast of India

The mangrove formations of Godavari estuary are due to silting over many centuries. The estuary covers an area of 62,000 ha of which dense Coringa mangrove forest spread in 6,600 ha. Satellite sensor data was used to detect change in the mangrove cover for a period of 12 years (1992-2004). It was found that an area of about 1,250 ha of mangroves was destroyed by anthropogenic interference like aquaculture, and tree felling etc. It was found that mangrove's spectral response/digital number (DN) value is much lower than non-mangrove vegetation such as plantation and paddy fields in SWIR band. By taking this as an advantage, spectral data was utilized for clear demarcation of mangroves from nearby paddy fields and other vegetation. Simpson's diversity index, which is a measure of biodiversity, was found to be 0.09, showing mangroves dominance. Ecological parameters like mud-flats/swamps, mangrove cover alterations, and biodiversity status are studied in detail for a period of 12 years. The increase in mangrove front towards coast was delineated using remote sensing data. The major advantages of remote sensing data is monitoring of change periodically. The combination of moderate and high-resolution data provided detailed coastal land use maps for implementing coastal regulation measures. The classification accuracy has been achieved is 90%. Overall, simple and viable measures are suggested based on multi-spectral data to sustain this sensitive coastal ecology.     

Land use and land cover (LULC) of the Republic of the Maldives: first national map and LULC change analysis using remote-sensing data

The Maldives islands in recent decades have experienced dramatic land-use change. Uninhabited islands were turned into new resort islands; evergreen tropical forests were cut, to be replaced by fields and new built-up areas. All these changes happened without a proper monitoring and urban planning strategy from the Maldivian government due to the lack of national land-use and land-cover (LULC) data. This study aimed to realize the first land-use map of the entire Maldives archipelago and to detect land-use and land-cover change (LULCC) using high-resolution satellite images and socioeconomic data. Due to the peculiar geographic and environmental features of the archipelago, the land-use map was obtained by visual interpretation and manual digitization of land-use patches. The images used, dated 2011, were obtained from Digital Globe’s WorldView 1 and WorldView 2 satellites. Nine land-use classes and 18 subclasses were identified and mapped. During a field survey, ground control points were collected to test the geographic and thematic accuracy of the land-use map. The final product’s overall accuracy was 85%. Once the accuracy of the map had been checked, LULCC maps were created using images from the early 2000s derived from Google Earth historical imagery. Post-classification comparison of the classified maps showed that growth of built-up and agricultural areas resulted in decreases in forest land and shrubland. The LULCC maps also revealed an increase in land reclamation inside lagoons near inhabited islands, resulting in environmental impacts on fragile reef habitat. The LULC map of the Republic of the Maldives produced in this study can be used by government authorities to make sustainable land-use planning decisions and to provide better management of land use and land cover.     

Landscape Change and Ecosystem Classification in a Municipal District of a Small City (Isernia,Central Italy)

Landscape changes taking place from 1954 to 1992 in the muncipal district of Isernia city (Central Italy) were described in relation to a system of ecosystem classification. Isernia municipal district was selected for study because recent historic changes in this area represent a typical example of landscape transformation similar to many small cities of Italy and other Mediterranean countries. To assess overall changes, three land cover maps (scale 1:25,000) were derived from panchromatic aerial photographs and field surveys. These were then digitalised in a Geographic Information System. A Land Facet (LF) map was derived by combining a phytoclimatic, a lithostatigrafic and a topographic map, and then digitalised as data layers in the same GIS. Results demonstrated two main landscape transformation trends: forest and semi-natural areas increased (8%), whereas agricultural areas decreased (12%). The urban area was relatively small during the entire analysed period, growing from 1% in 1954, to just 5% in 1992. Forest coverage was significant on reliefs, on hillside ecosystems such as limestone and on clay and marl hills LF. Arable land was particularly significant in flat ecosystems with deeper soils, such as on recent alluvial plain LF. These temporal changes were interpreted as being related to the replacement of traditional farming methods (grazing pastures) with more intensive methods (crop fields), especially on alluvial plains.