Optimizing agricultural land use: A GIS-based assessment of suitability in the Sana River Basin, Bosnia and Herzegovina

1.1 Background

Land represents one of the most significant natural resources in the world [1]. As a resource, land is defined by the Food and Agriculture Organization (FAO) as “the physical, biotic, environmental, infrastructural, and socio-economic components of a natural land unit” [2]. It is a finite and irreplaceable resource that provides essential ecological, social, and economic benefits, including food security, biodiversity conservation, climate regulation, environmental protection, economic development, and cultural values [3,4,5,6,7,8,9].

From a global perspective, land is a limited resource under increasing pressure [10]. According to Akıncı et al. and Jimoh et al. [11,12], the demand for land following population expansion, urbanization, and industrial development has resulted in land becoming a limited resource at a global level. In the context of using land as a resource for agriculture, urban expansion triggers the expansion of built-up areas, recreational parks, and public infrastructure, gradually occupying the vast agricultural and open areas surrounding a city [13]. According to Khan and Himanchal [14], urban expansion has taken over fertile agricultural land and intensified the issue of food security. Also, according to Beckers et al. [15], urbanization leads to a continuous loss of agricultural land, both directly in the form of land take and indirectly through the use of agricultural land for non-productive rural activities like recreation or hobby farming. Li et al. [16] pointed out that industrial development, along with land use and land cover change, can harm the environment by causing water and land pollution. They highlighted that industrial development leads to increased discharge of industrial waste and strains wastewater capacity, which are primary contributors to water pollution. This development can also adversely impact agriculture, including soil contamination and reduced water availability for irrigation.

Global climate change, inappropriate land use practices, as well as pressure to meet growing food needs, have significantly threatened global agricultural production and food security [1,17,18]. Van Dijk et al. [19] conducted a study based on 57 global food security projection and quantitative scenario studies published in the past two decades, analyzing their methods, underlying drivers, indicators, and projections. Across five representative scenarios encompassing diverse yet plausible socio-economic futures, it was projected that global food demand would increase by 35–56% from 2010 to 2050, with the population at risk of hunger expected to fluctuate between −91 and +8% over the same period. When considering the impact of climate change, these ranges shifted slightly (+30 to +62% for total food demand and −91 to +30% for the population at risk of hunger), but no significant statistical differences were observed overall. Also, according to Godfray et al., Baldos and Hertel, and McKenzie and Wiliams [20,21,22], it is assumed that by 2050, global food demand will double. As mentioned earlier, numerous discussions have emerged in recent times regarding ensuring sustainable agricultural production [23,24]. According to Edrisi et al. and Topuz and Deniz [25,26], prioritizing planned and rational land use as a resource, combined with sustainable agricultural production, should be prioritized as a critical resource in the future, which can lead to greater security in food provision. Proper land use planning ensures that agricultural activities are sustainable, productive, and environmentally friendly. This involves evaluating various factors that impact land suitability for agriculture, such as soil properties, climate change, and topographical conditions. Tashayo et al. [27] emphasize that soil properties, climate change, and variations in topographical conditions are the main reasons for investigating land suitability for agriculture. Xu et al. [28] claim that assessing land suitability for agriculture is of fundamental importance for ensuring food security, maintaining ecological balance, and promoting sustainable agricultural development. In addition, land suitability assessment for agriculture can provide a scientific basis and support to decision-makers such as government agencies and agricultural enterprises, to create an efficient management system for agricultural areas [29,30,31].

Qualitative and quantitative methods, such as the land use planning framework developed by the FAO [32], are employed for these purposes [33]. Moreover, analyzing the assessment of land suitability for agriculture can be used to evaluate scenarios for managing various agricultural practices. These practices contribute to increasing land productivity, reducing chemical inputs in agriculture, preserving natural resources, and mitigating the negative impacts of climate change [34]. Accordingly, land suitability assessment is a process of vital importance for transmitting productive land to future generations to ensure sustainable food production [35]. The multitude of environmental factors involved in land suitability assessment makes it even more complex in the context of sustainable land management without its degradation [11,36]. Hence, agroecological innovations are necessary to develop and determine specific land use methods in line with potentials and limitations, ensuring sustainable food production [37].

In this context, implementing agriculture policies that promote efficient and environmentally friendly practices is essential. In Europe, agriculture and climate change are policy priorities in the European Agenda, which recognizes the need for immediate action [38]. The European Union (EU) has implemented policies such as the European Green Deal and the Common Agricultural Policy (CAP) to promote sustainable, resilient, and profitable agricultural systems. Also, climate-smart agriculture (CSA) has been suggested as an approach to promote sustainable agricultural development for food security under climate change pressures [39]. At a national level, numerous developed as well as developing countries have initiated agricultural policy reforms to preserve and effectively use agricultural land [40]. A good example from a developed country is Switzerland, where progressive reforms in the late 1990s, backed by a substantial majority in a 1996 referendum, shifted agricultural subsidies towards ecological practices. This transformation was facilitated by revisions to the Swiss Federal Agricultural Law of 1992 and subsequent amendments in the Agricultural Act of 2002. These policies now differentiate between three tiers of support based on the sustainability of agricultural practices, focusing on biodiversity conservation, reduced inputs in integrated production, and substantial incentives for organic farming. On the other hand, Bosnia and Herzegovina (B&H), as a developing country, is also working on developing national-level action plans for sustainable agriculture [41,42]. The primary objective of agricultural policies should be to ensure the sustainability of agricultural production [43,44]. The first step in sustainable land use involves adequately planning land use methods, which should also respect and protect sustainable agricultural areas [45,46]. In this regard, it is necessary to ensure land productivity while considering its potential and capacity [47]. Constant investment and promotion of new systems within sustainable agricultural production are essential to ensure land productivity [23,48,49].

Efficient use of land resources in the context of agriculture involves assessing land suitability for agricultural production [50,51], enabling optimization of agricultural land use [52,53], sustainable agricultural development [54], and the implementation of responsible land use planning methods [45,46].

1.2 Problem statement

Over the past few decades, numerous studies have been conducted at both international and national levels on land suitability assessment for agriculture [28]. Contemporary methods for assessing land suitability for agriculture focus on the application of advanced technologies such as geographic information systems (GIS) and remote sensing. According to research results [27,29,55,56,57] GIS and remote sensing “products” in the form of satellite data have great potential in improving the analysis and assessment of land suitability for agriculture.

Studies pertaining to the application of GIS in assessing land suitability for agriculture in B&H are relatively scarce. Witmer and O’Loughlin [58] research analyzed the agricultural land in two study areas (northeast and south) within B&H, characterized by distinct climates, soils, and vegetation. They utilized satellite data methods (30 m Landsat imagery) to produce an initial map of rootable soil depth for B&H, distinguishing three classes of rootable depth within the broader Sana River basin: deep (90–140 cm), medium deep (60–90 cm), shallow (40–60 cm), and very shallow (0–40 cm). Jovanović et al. [59] conducted an analysis of freely available land cover maps for agricultural land use planning at the local level in the Laktaši Municipality within the Republic of Srpska (B&H). They noted that the evaluated Corine land cover (CLC) in the studied area is not a sufficiently precise GIS foundation for local agricultural land use planning. However, it could serve as a valuable starting point for developing sustainable land cover/land use approaches, thereby significantly expediting the planning process. Gekić and Bidžan-Gekić [60] discussed changing trends in agricultural land use in the Uskopaljska Valley within the Mountain-Valley macroregion and the Upper Vrbas-Pliva mesoregion of B&H. According to their analysis, the prevailing issues in agricultural land use include land abandonment and an increase in unused areas. Drašković and Bidžan-Gekić [61] analyzed land use changes over 20 years in East Sarajevo (B&H) using remote sensing methods and the CLC database. Their study revealed a significant expansion of discontinuous urban areas, driven by population migration between municipalities and political entities, as well as ethnic homogenization. In addition, there has been a notable decline in grassland areas and an increase in low-productivity grasslands in agricultural areas. Changes in forest vegetation, which covers a substantial portion of East Sarajevo, include increases in mixed forests and transitional woodland-shrub areas, alongside decreases in coniferous and alpine forests and heaths, partly due to the local economy’s reliance on forest resources. The study by Zurovec et al. [62] investigated the agricultural sector of B&H in the context of climate change, highlighting ongoing challenges and potential opportunities. They emphasized B&H vulnerability to climate change, particularly its impact on agriculture, food security, and socio-economic issues. The study underscores the need for further technological development, including enhancements in weather and climate information systems, agricultural monitoring, crop development, irrigation, and water management.

1.3 Objectives

The main objective of the research is to assess land suitability for agriculture using GIS and remote sensing technologies. This involves a comprehensive analysis of various factors such as soil properties, climate conditions, and topographical variations to determine the most suitable areas for agricultural activities. By leveraging these advanced technologies, the research aims to create detailed maps and models to guide decision-making processes for land use planning and management.

An integral part of this objective is to determine the current level of utilization of suitable land for agricultural purposes. This involves evaluating how much of the land identified as suitable for agriculture is currently being used for agricultural activities and how effectively it is being utilized. This assessment will help identify gaps and opportunities for optimizing land use, improving agricultural productivity, and ensuring sustainable land management practices.

According to the authors’ knowledge, previous studies have scarcely used GIS and remote sensing technologies in B&H for assessing land suitability for agriculture. This research aims to fill this gap by applying these modern technologies to the territory of B&H. The application of GIS and remote sensing in this context is expected to provide a more accurate and detailed understanding of land suitability, which can inform better agricultural policies and practices.

1.4 Significance

The research results have substantial practical significance for various institutions involved in agricultural sector development and sustainable land use planning. By utilizing advanced technologies like GIS and remote sensing, the research provides detailed and accurate assessments of land suitability for agriculture. This information is crucial for government agencies, agricultural enterprises, and other stakeholders to make informed decisions about land management and agricultural practices.

One of the key contributions of this research is the innovative application of modern technologies for assessing land suitability in B&H. Prior studies have rarely employed these technologies within this context, making this research a pioneering effort in the B&H. The use of GIS and remote sensing allows for a comprehensive analysis of land characteristics, enabling a more precise identification of suitable agricultural areas. This not only helps in optimizing land use but also in understanding the inherent potential and limitations of land as a natural resource.

The significance of the research extends to its potential impact on policy-making and strategic planning. By providing a scientific basis for land use decisions, the research supports the development of effective land management strategies that promote sustainable agricultural development. This is particularly important for ensuring food security in B&H, as the efficient use of suitable land can enhance agricultural productivity and sustainability.

Moreover, the research has broader implications for environmental conservation and climate change mitigation. By identifying optimal land use practices and promoting sustainable agriculture, the research contributes to the preservation of natural resources and the reduction of environmental degradation. This aligns with global and regional efforts to address climate change and ensure the long-term viability of agricultural systems.

In summary, the significance of this research lies in its potential to transform land use planning and agricultural practices in B&H. By providing valuable insights and practical tools for decision-makers, the research fosters a more sustainable and productive agricultural sector, ultimately contributing to the well-being of the population and the preservation of the environment.

2.1 Study area

2.1.1 Location and basic characteristics

The study area is the Sana River Basin, situated in the northwestern region of Bosnia and Herzegovina (B&H), between 44.18°N and 45.09°N latitude, and 16.29°E and 17.09°E longitude (Figure 1). The Sana River originates from three karstic springs located on the border of the municipalities of Ribnik and Mrkonjić Grad. It represents the largest tributary of the Una River and belongs to the larger drainage basin of the Sava River. Its spatial contribution to the Sava River Basin is 3.55%. The length of the Sana River is 146 km. The source is situated at an elevation of 414 m a. s. l., while the mouth is at 122 m a. s. l. The total area of the basin according to the HydroSHEDS (https://www.hydrosheds.org/) database [63,64] is 3,470 km². The average elevation of the basin is 505 m, while the average slopes amount to 10.9°.

Figure 1

Study area with locations of meteorological stations used in the validation of satellite data.

2.1.2 Geological and pedological characteristics

The diversity of the lithological composition, duration and frequency of tectonic processes throughout geological history, the geomorphological characteristics of the terrain, and the changeability of hydrometeorological factors and elements created the complex hydrogeological relation within the study area. Accumulations with karst-fracture porosity predominate in the southern (upper) part of the basin, while in the rest of the basin, accumulations with intergranular and/or fracture porosity alternate, along with areas without accumulations [65].

Based on the digitized Basic Soil Map at a scale of 1:50,000 [66] (Figure 2a), the study area exhibits a variety of dominant soil types: cambisols (47.20%), including Calcaric Cambisols (CMc), Chromic Cambisols (CMx), Dystric Cambisols (CMd), Eutric Cambisols (CMe), Ferralic Cambisols (CMo), Humic Cambisols (CMu), and Vertic Cambisols (CMv); leptosols (16.98%), including Eutric Leptosols (LPe), Lithic Leptosols (LPq), Mollic Leptosols (LPm), Rendzic Leptosols (LPk) and Umbric Leptosols (LPu); luvisols (11.10%), encompassing Chromic Luvisols (LVx), Ferric Luvisols (LVf), Haplic Luvisols (LVh), Stagnic Luvisols (LVj), and Vertic Luvisols (LVv); acrisols (9.24%), mainly Ferric Acrisols (ACf); fluvisols (7.75%), such as Calcaric Fluvisols (FLc), Dystric Fluvisols (FLd), Eutric Fluvisols (FLe), and Umbrick Fluvisols (FLu); podzoluvisols (4.96%), including Dystric Podzoluvisols (PDd) and Stagnic Podzoluvisols (PDj); gleysols (2.13%), comprising Calcic Gleysols (GLk), Dystric Gleysols (GLd), Eutric Gleysols (GLe), Mollic Gleysols (GLm), and Umbric Gleysols (GLu); vertisols (0.35%), like Calcic Vertisols (VRk) and Eutric Vertisols (VRe); and podzols (0.25%), characterized by Ferric Podzols (PZf).

Figure 2

Characteristics of the study area: (a) pedology, (b) elevation, (c) slope, (d) land use, (e) average precipitation sum, (f) average air temperature, and (g) hydrographic network.

The diverse soil types present within the study area indicate that each, depending on its characteristics, can more or less influence the agricultural potential of the study area and decisions regarding land management.

2.2 Data preparation for GIS environment

To identify suitable agricultural land, several different geospatial techniques were applied at the study area level. Considering the available data, the following factors were considered: soil characteristics, elevation, slope, land use, hydrological network, and meteorological characteristics (aridity index - satellite data on average annual precipitation and average annual air temperature). The preparation and analysis of these factors were conducted within the GIS environment using QGIS 3.28 ‘Firenze’ software (https://qgis.org/) and the Google Earth Engine (GEE) platform. A systematic methodological presentation of land suitability assessment is provided in Figure 3, which will be further explained.

Figure 3

Algorithm of the methodology for assessing land suitability for agriculture.

The soil map of the Sana River basin was digitized based on the Basic Soil Map (1:50,000) derived according to the FAO classification [66]. The elevation map was created based on the available digital elevation model (EU-DEM – https://www.eea.europa.eu/) with a spatial resolution of 25 m. Using the DEM as input data and through QGIS 3.28 “Firenze” software (https://qgis.org/) and the slope function, a slope map of the study area (expressed in degrees) was created. The land use map was created based on the CLC (https://land.copernicus.eu/) geospatial database from 2018, while the hydrological network map was created based on the HydroSHEDS (https://www.hydrosheds.org/) geospatial database [64].

The aridity map was prepared based on the calculation of the annual values of the [74] aridity index using the formula:

where I _DM is the annual monthly aridity index, P is the mean annual precipitation, and T is the mean annual air temperature. The result of the aridity index is interpreted through the identification of climate types: <10, arid; 10–20, semi-arid; 20–24, mediterranean; 24–28, semi-humid; 28–35, humid; and >35, super-humid [75,76,77,78].

The input data for average annual precipitation and average annual air temperature are satellite based. Satellite precipitation estimation data CHIRPS [69] and satellite air temperature estimation data ERA5 [70] were used. The mentioned data were filtered for a time period of 20 years (2000–2020) and temporally reduced to the vegetation period level. The GEE platform was used for their processing, and prior to their integration into the land suitability assessment, their validity was assessed. According to the methodology of Sabljić et al. [79], the validation of satellite data involved comparing them with data from meteorological stations (MS) based on their average monthly values. In line with the presented methodology, this study conducted a validity assessment of ERA5 satellite data on average air temperature (1981–2023). Meteorological data for the validation process were obtained from the Republic Hydrometeorological Institute of the Republic of Srpska and the Federal Hydrometeorological Institute of the Federation of Bosnia and Herzegovina. Meteorological data from MS located within the basin were considered, as well as data from MS located in the immediate vicinity of the basin (Figure 1 and Table 1). The reason for including MS located outside the basin boundaries is explained by the lack of MS at higher elevations within the basin, as well as the “low” spatial resolution of satellite data (ERA5 ∼ 27.8 km; CHIRPS ∼ 5.5 km) compared to the surface area of the study area.

2.3 Methodology for assessing land suitability for agriculture

The assessment of land suitability for agriculture involves the application of multi-factor analysis. Within this analysis, the analytic hierarchy process (AHP) method was used to generate a synthetic map of land suitability for agriculture. The AHP method was developed by Saaty [80,81]. The goal of AHP is to quantify factors and create a hierarchy of factors according to their priority for assessing suitability [82,83,84]. The method provides systematic and logical decision evaluation [85] and as such enables decision-makers to achieve correct and logical results through consistent assessment of goals, factors, and criteria, modelling in a hierarchy [86,87,88]. The AHP method has broad applications in research aiming to assess land suitability for agriculture [89,90,91,92,93].

According to Liang and Yang [94], the AHP method is susceptible to subjectivity when defining the “weights” of factors, as it relies on the results of previous research in the given field. Accordingly, when defining the hierarchy of priorities based on their importance, user subjectivity can contribute to more objective results. In the process of assessing land suitability for agriculture, not all factors are equally important, and the main reason for using the AHP method is to define the different “weights” of each factor, creating a hierarchy of factors based on their degree of significance for the final result.

The assessment of the importance of factors during pairwise comparisons was conducted using the Saaty Scale [80] of relative importance (1–9), where 1 represents equal importance, 3 represents weak importance, 5 represents strong importance, 7 represents strong, proven importance, and 9 represents extreme importance of one factor over another. Other values (2, 4, 6, and 8) represent intermediate importance levels. Numeric values from 1 to 5 were used in this research. The reason for using a smaller range of values is to mitigate significant differences in weighting coefficients, which would otherwise increase subjectivity in priority assessment [88].

Consistency check of the pairwise comparison matrix, i.e., the assessment of the “weights” of factors, is performed using the consistency index (CI) and the consistency ratio (CR) calculation [95]. CI is calculated using the following formula [96]:

where λ _max is the maximum eigenvalue of the comparison matrix and n is the order of the matrix (size of the comparison matrix) [97].

CR is calculated using the following formula [96]:

where CI is the consistency index and RI is the random index (which depends on the order of the matrix) obtained from the standardized form according to Saaty [98] (Table 2).

If CR ≤ 0.10, the result is considered sufficiently accurate, and the assessments of relative weighting coefficients are deemed acceptable [80,88].

To clearly identify areas of land suitability for agriculture, the result of the land suitability assessment underwent reclassification and k-means cluster analysis. According to Haynes et al. [99], the reclassification process identifies pixels with a specific value, which are then changed to a new value defined by the analyst. Therefore, as mentioned by the authors, this process represents a specific form of map algebra, in which pixel values are evaluated and replaced with a new value. In this regard, since the standardized criterion values within the AHP method range from 1 to 5, the reclassification process of the final land suitability assessment result “moves” within this range with decimal values, thus implying the delineation of land suitability classes as follows: 1 – very unsuitable, 2 – unsuitable, 3 – conditionally suitable, 4 – suitable, and 5 – very suitable. The mentioned process in the research was executed using QGIS 3.28 “Firenze” software (https://qgis.org/), using the Reclassify by the table tool.

The same land suitability assessment result for agriculture was subjected to k-means cluster analysis to identify spatial patterns of land suitability. According to Micić Ponjiger et al. [100], k-means represents a clustering technique based on centroids, where clusters are represented by central vectors. The number of clusters is fixed and denoted as K, and the algorithm finds K cluster centres and assigns objects (O) to the nearest cluster centre, minimizing the squared distances from the cluster [101]. In the research, the land suitability assessment result for agriculture was clustered using the k-means clustering for grids tool within the SAGA GIS software. The Hill-Climbing method [102] was employed. This method operates on the principle of iteratively adjusting cluster centres to maximize similarities within clusters and minimize differences between clusters. Considering that the method performs cluster analysis based on the number of iterations, a total of 15 iterations were defined within the k-means clustering for grids tool.

Finally, the results of the reclassification process and cluster analysis performed on the land suitability assessment for agriculture were compared and analyzed, highlighting the advantages and disadvantages of each approach in interpreting the assessment of land suitability for agriculture at the study area level.

2.4 Methodology for land use identification

The process of land use identification was carried out through the analysis of remote sensing “products” in the form of satellite images. For their processing, the GEE platform based on cloud technology was used [103]. The primary dataset for land use identification was obtained from the Sentinel-2 satellite mission. An overview of the characteristics of the Sentinel-2 satellite mission is provided in Table 3 [104].

Sentinel-2 satellite images, renowned for their high spatial resolution and spectral capabilities, have emerged as indispensable tools in modern remote sensing applications. These images serve as critical input data for the classification process, facilitating the precise delineation of land cover and land use patterns [105]. Leveraging the robust computing capabilities of contemporary platforms like GEE has significantly enhanced the efficiency and accuracy of land cover/use classification using Sentinel-2 data [106,107]. Machine-learning algorithms, such as random forests (RF), k-nearest neighbour (KNN), support vector machine (SVM), and Bayesian methods, have been pivotal in this regard, with the RF classifier notably preferred due to its effectiveness in handling complex datasets [108,109,110].

According to Phiri et al. [105], high accuracies have been reported in the land cover/use classification of Sentinel-2 data, with most classification accuracies exceeding 80%. Additionally, pixel-based methods utilizing the maximum likelihood classifier (MLC) have demonstrated high classification accuracies [111,112,113]. Phiri et al. [105] also compared the accuracies of various machine-learning classifiers and found that RF and SVM outperform other classifiers. By applying the supervised classification process, land use classes are identified at the study area level (Figure 4). The aim is to identify land use types for the year 2023 at the study area level. The image processing consists of two phases: pre-processing and post-processing. During pre-processing, the time period of interest is defined. Since the goal is to identify currently used land potentials for agricultural purposes, supervised classification is conducted for the year 2023. In addition, a crucial step in pre-processing is filtering images based on cloud coverage percentage and spatial extent. In the post-processing phase, land use classes are delineated. For the purpose of classifying land use classes, the year 2023 is divided into three separate temporal sub-periods, considering the spring and autumn sowing and harvesting periods. The goal is to clearly identify and delineate areas used for agricultural purposes (cultivation of crops) from those that are arable but not used for agricultural purposes (e.g., meadows or pastures). A stack (composite) is created based on temporal sub-periods, which is used as a basis for generating training data in the classification process. The mentioned process is executed using the RF algorithm. The following land use classes are extracted: water bodies, forested areas, agricultural land, meadows, and built-up land.

Figure 4

Algorithm for processing Sentinel-2 satellite imagery for land use identification.

The results of the classification process are evaluated for accuracy using confusion matrix. Based on the confusion matrix, global quality metrics such as the overall accuracy (OA) and kappa coefficient (K) were calculated (equations (4) and (5)) to evaluate the impact of composition methods on land cover classification. OA and K were calculated using the following formulas [114]:

where “number of correctly classified samples” refers to the count of samples that were accurately identified or categorized by the model or classifier, and “number of total samples” refers to the total count of all samples in the dataset, including both correctly and incorrectly classified ones.

where “overall accuracy” is the proportion of correctly classified samples out of the total number of samples. It reflects the model’s ability to accurately classify samples. “1” represents the perfect agreement, or the ideal situation where all samples are classified correctly. It serves as the upper limit for the Kappa calculation, and “estimated chance agreement” is the probability of agreement between the observed classifications and the expected classifications by chance alone. It accounts for the possibility that some agreement might occur randomly rather than due to the classifier’s performance.

In addition, based on the confusion matrix, the following parameters were calculated at the land use class level: user accuracy (UA), producer accuracy (PA), and F-1 score (equations (6)–(8)). The F1-score is the harmonic mean between producer’s and user’s accuracies and can be used to evaluate the accuracy at the class level [115]. The formulas for calculating these parameters are presented below [114]:

where “number of correctly classified samples in each class” refers to the total count of samples that were accurately identified as belonging to their true class, and “number of samples classified to that class” refers to the total count of samples that were assigned to that class, regardless of their true class.

where “number of correctly classified samples in each class” refers to the count of samples that were accurately identified as belonging to their true class according to the classification results, and “number of samples from reference data in each class” refers to the total count of samples that actually belong to a specific class according to the reference or ground truth data.

where PA represents the probability that a pixel was correctly classified in a given class. UA represents the probability that a pixel classified in a given class of the map represents that class on the ground [116]. The F1 was found to be the best performance metric and is widely used in the previous research, which gives equal importance to both PA (as a precision) and UA (as a recall) by combining them into a single model performance metric [117,118].

The results of land suitability assessment for agricultural purposes were successively intersected with agricultural areas identified from the supervised classification of land use. The goal of this process is to determine the degree of utilization of suitable land for agricultural purposes.

3.1 Assessment of land suitability for agriculture

The assessment of land suitability for agriculture is based on factors such as pedology, topography (elevation and slope), land use, meteorology (precipitation and temperature), and hydrography. An overview of the basic characteristics of these factors is provided in Section 2.1 (Figure 2).

Before using satellite data on precipitation and air temperatures as meteorological factors within the AHP method framework, it is necessary to assess their validity. Satellite data on precipitation estimation (CHIRPS) have been validated at the level of the study area by Sabljić et al. [119]. According to the results of the mentioned authors, the average amount of precipitation (1981–2023) for the study area matches by 94.95%, while the average annual sum of precipitation matches by 98.19%, making these data valid for research purposes.

According to the earlier described methodology, an assessment of the validity of satellite data on air temperature estimation (ERA5) has been conducted. A comparison of the average air temperature from MS and ERA5 indicates a high degree of validity of the satellite data (Figure 5). Nine of 12 months (March, April, May, June, July, August, September, October, and November) show >90% agreement of satellite data with real data. A lower degree of agreement is observed in the months of January (81.81%), February (83.83%), and December (86.74%). The annual average air temperature (1981–2023) based on MS data is 10.49°C, while based on ERA5 data, it is 10.42°C (with a match of 99.33%).

Figure 5

Comparison of the average annual air temperature MS and ERA5 data (1981–2023).

Taking into account the consequences of climate change [120,121,122], which manifest through changes in precipitation and air temperature quantities during the past decades, for the purpose of calculating the de Marton aridity index, valid CHIRPS and ERA5 data have been reduced to the time period of the last two decades (2000–2020). In addition, the data are temporally limited to the vegetation period. In this regard, for the mentioned time period, the MS data for average annual air temperature amount to 16.61°C, while according to ERA5, they amount to 15.80°C. The agreement for satellite data with real data during the vegetation period is 95.12%. It is concluded that the data are valid for integration within the calculation of the aridity index and can be used as factors in the application of the AHP method.

According to the presented methodology for the previously described factors assessing land suitability for agriculture, a pairwise comparison matrix has been created. The goal of this matrix is to identify the significance of each factor in determining the degree of land suitability for agriculture. During the comparison, values were assigned to the factors based on Saaty’s [98] scale (Table 4). The assignment of values is based on the previous research in this area. The result of the matrix indicates that the highest “weight” or influence on the assessment of land suitability for agriculture is attributed to the pedological characteristics of the study area (0.32), followed by land use (0.19), altitude and slope (equally 0.17), aridity index (0.10), and distance from rivers (0.05).

Based on the pairwise comparison matrix, the calculation of the consistency index (CI), consistency ratio (CR), and random index (RI) was performed:

Using the valuation method, criteria were evaluated within each of the previously mentioned factors for assessing land suitability (Table 5). Within the scale, a criterion with a low “weight” illustrates a small influence, while a criterion with a high “weight” illustrates a significant influence [123] on the degree of land suitability for agriculture. The interval values of the criteria range from 1 to 5, where 1 – very unsuitable, 2 – unsuitable, 3 – conditionally suitable, 4 – suitable, and 5 – very suitable. The valuation process involves standardizing the input criteria to enable mutual comparison both among criteria and among the factors themselves.

In the river basin area of the Sana River, a total of 40 soil types have been identified according to the FAO soil classification. Based on previous research on soils and their suitability [124,125,126], their evaluation was conducted. A slope represents a significant factor in assessing land suitability for agriculture. According to Everest et al. [85], the degree of slope directly affects soil erosion susceptibility, soil tillage methods, types of agricultural machinery used in cultivation, irrigation methods and intensity, plant adaptation, and so on. In addition, steepness, length, and sharpness of the slope directly influence the loss of both soil and water, which can significantly reduce land productivity for agricultural purposes. Gekić and Bidžan-Gekić [60] highlight that steep slopes in the B&H region are not suitable for agriculture. For these reasons, the evaluation of slope was performed in a way that lower degrees of slope were valued as more suitable for agriculture compared to higher degrees of slope. According to the results of Kapluhan [127], elevation influences the duration and characteristics of winter and summer months. With changing climatic conditions due to changes in elevation, physiological and morphological differences may arise in the same plants [128]. In this regard, the vegetation period significantly affects the timing of harvest and profitability in agricultural production [85]. Grujčić et al. [129] emphasize that areas up to 280 m a. s. l. in the northwest of B&H are the most suitable for agricultural production. This area encompasses the Sana River Basin, and accordingly, the evaluation of elevation as a factor was conducted in such a way that lower elevations were valued as more suitable compared to higher elevations. In addition to the aforementioned factors, land use has a significant impact on the assessment of land suitability for agriculture. Within this factor, criteria such as arable land or areas under natural vegetation were valued with a better assessment of suitability compared to criteria such as built-up or water areas. Average annual precipitation and average annual air temperature for the vegetation period (2000–2020) were used as input parameters for calculating the de Marton aridity index. The results of the de Marton aridity index represent different types of climates, and based on their suitability for assessing land for agricultural production, their evaluation was performed. The existing hydrographic network in the river basin area is also an important factor in assessing land suitability for agriculture. Areas in close proximity to rivers may be more financially economical and infrastructurally simpler to irrigate compared to more distant areas. For this reason, the criteria of this factor were evaluated in such a way that areas in close proximity to rivers were valued as more suitable for agriculture compared to those that are more distant.

Based on the evaluation of the criteria for each of the factors for assessing land suitability, a process of reclassification of the factors (1–5) was carried out, and as a result, reclassified maps of the factors were obtained (Figure 6).

Figure 6

Evaluated criteria for assessing favourability: (a) pedology, (b) elevation, (c) slope, (d) land use, (e) aridity index, and (f) distance from rivers.

By applying the AHP method according to the previously described methodology, the factor maps were overlapped together with the predetermined weighting factor of each criterion. As a result, a synthetic map of land suitability assessment for agriculture was obtained (Figure 7). Areas located in close proximity to the main course of the Sana River and its tributaries were identified as suitable for agriculture. These areas are characterized by high-quality soil, gentle slopes and elevations, dominant land use in the form of agricultural surfaces, suitable climate, and proximity to rivers, which is a prerequisite for irrigation. In contrast, the southern part of the river basin is characterized by poor soil for agricultural purposes, as well as steep slopes and elevations. These factors contributed to making the southern part of the river basin mostly unsuitable for agriculture.

Figure 7

Synthesized map of land suitability assessment for agriculture.

To clearly identify land suitable for agriculture, according to the previously described methodology, processes of reclassification (Figure 8a) and k-means cluster analysis (Figure 8b) were conducted based on the results of the suitability assessment. The essence of these processes is to delineate regions/areas with different levels of land suitability for agriculture. The results of the reclassification process (RP) and cluster analysis process (CAP) of land suitability indicate spatial overlaps of areas identified as highly unsuitable for agriculture (1). These overlaps are evident in the southern region of the basin. In addition, according to the CAP results, highly unsuitable land for agriculture (1) is identified in the northern and western regions of the basin. Based on the RP, “smaller” spatial patterns of unsuitable land (2) are identified in the northeastern and western regions of the basin, while clear spatial patterns of this class are identified in the southern region of the basin. The CAP results align with the RP results, further identifying unsuitable land in the northern and central regions of the basin (east and west of the mainstream of the Sana River). These areas are characterized by fewer soil characteristics, steep slopes, higher elevations, and a relatively sparse hydrographic network (limited irrigation potential). Areas conditionally suitable for agriculture (3) according to the RP and CAP results are present in almost the entire basin. More significant spatial patterns of this class are visible in the northeast, west, and central regions of the basin. These areas are characterized by steep slopes and slightly less-suitable soil characteristics, which influenced the identification of these areas as conditionally suitable for agriculture. Significant spatial patterns of suitable land (4) for agriculture are identified in the western, northeastern, and southeastern regions of the basin, thanks to the RP and CAP results. These areas are characterized by gentle slopes and low elevation, very good soil characteristics, suitable climate, proximity to major rivers (Sana and Una), and dominant land use in the form of agricultural surfaces. These factors have led to these spaces having characteristics of suitable land for agricultural production. Highly suitable land (5) for agriculture is predominantly present in the valley of the Sana River, i.e., in the immediate vicinity of its upper, middle, and lower streams. These areas are characterized by high-quality soil, gentle slopes, low elevation, dominant agricultural land use, suitable climate, and proximity to the Sana River.

Figure 8

(a) Reclassified and (b) clustered assessment of land suitability for agriculture.

The areas of the RP and CAP are separated for each class of regions/areas and are compared with each other (Figure 9). The area with the smallest contribution in the river basin is the class of very unsuitable land for agriculture (1). According to the RP results, this class covers an area of 5.32 km² (0.15%), while according to the CAP results, it covers an area of 198.77 km² (5.83%). The class of unsuitable land for agriculture (2) according to the RP results covers an area of 117.47 km² (3.44%), while according to the CAP results, this class of land suitability covers an area of 596.96 km² (17.52%). Conditionally suitable land for agriculture (3) according to the RP results covers an area of 1094.06 km² (32.11%), while according to the CAP results, it covers an area of 969.97 km² (28.47%). The class with the largest spatial participation in the river basin according to the RP and CAP results is the suitable land for agriculture (4). According to the RP results, this class covers an area of 1917.85 km² (56.29%), while according to the CAP results, it covers an area of 973.29 km² (28.57%). The class of very suitable land for agriculture (5) according to the RP covers an area of 271.88 km² (7.98%), while according to the CAP results, it covers an area of 667.59 km² (19.59%).

Figure 9

Comparison of the percentage contribution of the area of each class of land suitability assessment in the total area of the river basin.

Although spatial overlaps between the RP and CAP results are noticeable, there are also certain differences in their outcomes. These differences stem from different methodological approaches, each of which has its own advantages and disadvantages.

RP simplifies data interpolation by grouping pixel values into categories, which facilitates their interpretation and analysis. In addition, this process leads to a clear visual representation of the data, greatly facilitating decision-making in the agricultural sector when planning policies and development programs. However, a drawback of reclassification is the potential loss of detail and data precision since they are grouped into broader categories.

On the other hand, CAP provides insights into spatial patterns and data distribution, greatly facilitating spatial analysis. The advantages of this process lie in identifying classes through “natural” grouping. Therefore, CAP objectively identifies groups of similar data, revealing inherent patterns in the data. In this regard, CAP offers certain advantages over reclassification when issuing classes of land suitability for agriculture. For example, a smaller set of values (1.95–1.99) located in close proximity to a set with higher values (2.00–2.50) would be assigned to class 2 by reclassification, while it would be assigned to class 3 by CAP. Thus, reclassification clearly separates what belongs to each class based on defined criteria, while CAP considers the “broader area.” In the authors’ opinion, in the context of identifying classes of land suitability for agriculture for future planning purposes, the results of CAP are more beneficial, while reclassification results are more suitable for determining the current state on the ground.

3.2 Land use

Within the land use classification process for 2023, a supervised classification process was conducted. The results underwent an accuracy assessment process, which included constructing the confusion matrix. The confusion matrix illustrates the allocation of land use classes for the pixels in the validation set (ground truth) compared to their assignments during the classification process. Diagonal entries in the matrix indicate correctly classified pixels. Entries off the diagonal signify misclassified pixels, suggesting the potential blending of different land use classes. From these data, UA, PA, and F-1 score were calculated. According to the classification results for 2023 (Table 6), PA is 1 for water surfaces, 0.99 for forested areas, 0.97 for meadows, 0.97 for agriculture areas, and 0.98 for built-up areas. UA is 0.99 for water surfaces, 0.99 for forested areas, 0.98 for meadows, 0.98 for agriculture areas, and 0.97 for built-up areas. The F-1 score for all classes exceeds 0.95.

The overall accuracy is 0.99, while the Kappa coefficient is 0.98. This coefficient represents a statistical parameter that measures the agreement between reference and classified data. The Kappa coefficient value is used to assess the validity of the classification. A range of values from 0.81 to 1 indicates almost perfect agreement between the classified and reference data in the classification process [130,131,132]. Considering this, along with the presented accuracy assessment results for the classification at the study area level, it is concluded that the data are valid for further research purposes.

According to the methodology described earlier, a land use map for the year 2023 was created at the study area level (Figure 10). The land use classes identified are as follows: water areas, forested areas, meadows, agriculture areas, and built-up areas. According to the classification results, water areas cover 17.02 km² (0.49%), forested areas cover 2346.73 km² (68.91%), meadows cover 740.13 km² (21.73%), agriculture areas cover 244.49 km² (7.17%), and built-up areas cover 57.04 km² (1.67%).

Figure 10

Land use classes at the study area level (2023).

From the results of supervised land use classification, agricultural areas were identified and successively intersected with reclassified and clustered regions/areas of land suitability for agriculture. It was found that currently, for agricultural purposes, 0.04 km² (RP) or 0.88 km² (CAP) of very unsuitable land (1) is being used. Unsuitable land (2) is used for agriculture on an area of 0.41 km² (RP) or 7.28 km² (CAP). Conditionally suitable land (3) is used for agricultural purposes on an area of 15.75 km² (RP) or 27.52 km² (CAP). Suitable land (4) is used for agricultural purposes on an area of 185.15 km² (RP) or 107.06 km² (CAP), while very suitable land (5) is used for agricultural purposes on an area of 42.99 km² (RP) or 101.65 km² (CAP).

Conditionally suitable, suitable, and very suitable areas represent priority areas for the development of the agricultural sector, and they have been identified through field research (Figure 11a–f).

Figure 11

(a)–(f) Natural characteristics and agricultural areas of the Sana River Basin.