The study area
The study was conducted in the selected watersheds (Elgeyo and Nyambene) in Kenya. The Elgeyo ecosystem covers 108,194 ha including the state forest (24,354 ha) and adjacent farmland within the five-kilometer buffer zone (83,840 ha). The state forest consists of eight forest blocks, namely Kaptagat, Kipkabus, Kessup, Kapchorua IV, Tingwa Hills, Tumeya, Kapchorua I, and Metkei (Fig. 1). The state forest consists of exotic forest (40%), native forest (38%), and open grassland and shrubland (22%) (KWTA, 2020b). Largely, within Elgeyo Marakwet County with a small section in Uasin-Gishu County with population of 0.5 million and 1.2 million respectively (KNBS, 2019). It extends from 35° 20” to 35° 45” East and 0° 10’ to 0° 20’ North. Precipitation is binomial with a mean of 1200 mm, highest between March and May and lowest in August and October. Temperatures range from11.2 0C to a high 33 0C (KAOP, 2020). The highest point in the ecosystem is estimated at 3350 m above sea level while the slope varies between 2° and 60° (CGoEM, 2018; Elgeyo Marakwet District Plan, 1980; KWTA, 2020a).
The Nyambene ecosystem is part of the Tana and Ewaso-Nyiro watersheds and covers 30,313 ha consisting of the state forest (5,427 ha) and the farmland within the five-kilometer buffer zone (24,886 ha). The state forest is predominantly indigenous and divided into four management blocks including Nyambene, Kilimandingiri, Keiga, and Thuuri (KWTA, 2020c). The Nyambene extends from 0º 17’ N to 0º 8’ N, and from 37º 48’ E to 37 º 52’ E within Meru County and is traversed by the subdistricts of Igembe South, Igembe Central, Tigania East, Tigania West, and Tigania Central (Fig. 2). The five subdistricts have a population of 691,298 (173,743 households) (KNBS, 2019). The precipitation regime is binomial with long rains between March and May and short rains in October and November and a mean of 1700 mm. The altitude of the area ranges from 1000 m to 2,528 m above sea level while temperatures range from 13.7°C to 28.7°C. The ecosystem is endowed with floral diversity, over 200 springs, and a significant number of streams and rivers that serve as water sources for populations within the watershed and further downstream (KWTA, 2020c).
Research And Sampling Design
The study adopted a cross-sectional design with the actual assessment based on ecosystem service type, data, and benefit cohorts. Based on the classification of the Millennium Ecosystem Assessment (MEA) and the TEV framework, the ES data collection was regrouped into three perspectives, socio-cultural, ecological, and economic values. The socio-cultural values comes from households, focus group discussions, participatory mapping tools, expert surveys and the Q methodology. The ecological values using GIS and remote sensing were validated with field studies and substantiated with secondary data, while the economic data leading to monetary allocations used traditional valuation techniques such as market prices, cost-based, stated and revealed preference techniques, and benefit transfer (Baral et al., 2017) to estimate the ES currency unit. However, the study focused on regulatory and support services and used a hybrid approach with both biophysical and socioeconomic attributes(Mengist et al., 2020), as shown.
Data Collection
Assessment of Ecological Values
The valuation of RES involved biophysical quantification and attribution of the monetary unit using non-market valuation techniques to assign monetary values. The evaluation technique used was based on study size, data needs, availability, available resources, topics, and available expertise (Baral et al., 2017; Burkhard et al., 2010, 2012; Häyhä et al., 2015; Paudyal et al., 2015). The assessment began with land use land cover, RES profiling and quantifIcation, attribution of shadow prices to products, and estimation of the grand total.
Land Use Classification
The study used Geographic Information System (GIS) and Remote Sensing (RS) techniques with a spatial resolution of 30 m to generate land cover data for the two ecosystems. The assessment began with image generation, image processing, classification with random forest classifier and creation of coresponding classified maps (LC1990, LC2000, LC2010, LC2020). Four Landsat path/array satellite images from three types of sensors were downloaded from the United States Geological Survey (USGS) website https://earthexplorer.usgs.gov/. The images were taken during the dry season of the year i.e. between January and March to ensure cloud-free and improved image display. The data was processed using ArcGIS 10.7 and R Studio 1.4.1106 and ENVI 5.3. The generated images were projected onto the Universal Transverse Mercator (UTM) coordinate system, datum Arc1960, Zone 36 North, and corrected for geometric errors from the sources using ground control points derived from a 1:50,000 scale topographical map. The other three previuos versions (L5 TM, L7 ETM+, L7 ETM+) of Landsat imagery (1985 TM, 2001 ETM+, and 2010 ETM+) were then each referenced by performing the frame-to-frame registration method using the latest version corrected Landsat 8 OLI /TIRS 2022 image. The IPCC Scheme II classification has been adopted, which considers ten (10) classes, namely; dense forest, moderate forest, open forest, wooded grassland, open grassland, perennial cropland, annual cropland, wetland, open waterbodies and barrens.
The process began by delineating the training site with polygons, encoding the land cover, and enhancing the image features using true and false color composite. Validation of the predefined land cover Landsat imagery training site was performed through field observations from the 100 assessment points per ecosystem, Google Earth imagery, and historical land cover data generated through interviews with the adjacent community. Arandom forest classifier was applier with an accuracy of 0.8 based on the class confusion matrix to create a spectral signature and classification of all pixels in the generated image. Finally, an image filter was applied to smooth the classification results by removing ‘salt’ and ‘pepper’ noise from the classified maps. The final land cover maps were used to generate and analyze the LCLU class area size (ha) using the ‘Tabulae’ area algorithm in ArcGIS version 10.7 which intersects the imagery with the respective study area.
Quantification And Economic Valuation Of Res
Water flow regulation and Water Purification
The study opted for the water storage method with replacement costs as indicated (1), as widely accepted (Langat et al., 2020; Langat, 2016; Xi, 2009) based on the avoided cost principle. Landcover size was determined using 2019 Landsat imagery, while precipitation amount was based on average annual precipitation data obtained on request from the Kenya Metrology Database (MoE&F, 2020). Runoff reduction coefficients were obtained from secondary databases of ecosystemS with similar ecological characteristics (Kateb et al., 2013; Okelo, 2009) and relative land cover coefficients (Blume et al., 2007; Goel, 2011; Karamage et al., 2018; Kauffman et al., 2007). The unit costs of the water regulation were determined by the unit costs of the replacement system (artificial dam) (US$3/m3) based on a replacement cost principle (Eytan & Spuhler, 2020; GOK & World Bank Group, 2005; Wu et al., 2010).
\({\text{V}}_{\text{W}\text{P}}\text{=}\sum _{\text{i=1}}^{\text{N}}{\text{A}}_{\text{LC}}\text{×}{\text{P}}_{\text{C}}\text{×}{\text{RR}}_{\text{Coef.}}\text{×}{\text{C}}_{Sur\text{.}}\) 1
VWP represents the economic value for the watershed; ALC represents the area (ha) of land cover; PC represents the average annual rainfall that the ecosystem receives; RR coef. runoff reduction coefficient of the respective landcover (estimated by the precipitation runoff coefficient of the respective landcover/land use subtracted from the runoff coefficient of the bare area); CSur represents the unit cost per cubic meter of replacement water reservoir.
The function of the water purification ecosystem was based on the avoided water treatment costs according to formula (2). The amount of purified water was based on estimated annual precipitation retained by the two ecosystems. The unit cost of the purification function was based on the unit cost of constructing and maintaining a backup facility (municipal water treatment plant) (Jahanifar et al., 2017). This was based on the assumption that the destruction of the forest ecosystem would result in water quality degradation, which would require the construction of a municipal wastewater treatment plant to replace the ecosystem function.
\({V}_{WQ}\text{=}{Q}_{WC}\text{×ρ}\) 2
Where V WQ represents the economic value to regulate the water quality of the ecosystem; QWC is the amount of water stored and purified by the ecosystem, which can also be represented by total household consumption; ρ represents the unit cost of US$0.3/m3 (Fuente et al., 2015) of the replacement water treatment mechanism
Soil Conservation And Erosion Control
In the study, the relative soil loss of land cover was assumed to be the unit cost of impact mitigation given by formula (3) on the avoided cost principle (Bishop, 1999; Nahuelhual et al., 2007). Land cover size was determined from the 2019 land cover Landsat imagery, while the corresponding land cover soil erosion reduction coefficient was from the secondary database (Hurni, 1988; Kateb et al., 2013; Tessema et al., 2020). The unit cost of the ecosystem’s soil erosion control function was based on the replacement cost of dredged water reservoirs, in this case a hydroelectric power generation dam estimated at US$3.34 per tonne of sediment (Adeogun et al., 2016).
\({\text{V}}_{\text{SC}}\text{=}\sum {\text{LC}}_{\text{A}}\text{×}{\text{SE}}_{\text{RC}}\text{×}{\text{C}}_{\text{Proxy}}\) 3
Where V SC represents the economic value for forest soil protection; LCA is the respective land cover area (ha); SERC is the soil erosion reduction coefficient based on land cover soil erosion coefficients (Hurni, 1988; Tessema et al., 2020); C proxy the proxy unit cost estimated at KES 351 (USD 3.34) per tonne of sediment (Adeogun et al., 2016).
Soil Nutrient Conservation (Dup: Abstract ?)
In the study, it was assumed that the loss of in situ soil minerals (nitrogen (N), phosphorus (P), potassium (K)) is attributed to the relative soil loss across the different land covers and the unit replacement costs formula (4), as is often assumed in similar studies (Nahuelhual et al., 2007). Soil mineral content across different land covers was determined by field sampling and laboratory analysis, while the unit value of the ecosystem soil nutrient protection function was equated to a surrogate (artificial fertilizer) relative unit cost based on a replacement cost principle (Gizaw et al., 2021).
Where: EV SNC is the economic value of soil protection; SLC is soil conserved (kg/ha) of the respective land cover; SNCLC is the soil nutrient content (%) (N, P, K) in the forest soil; QCF is the commercial fertilizer estimated at 150 kg/acre annually in Kenya; and ⸹CF is the ratio of commercial fertilizers (51%, NPK-17-17-17); PCF is the unit price of the commercial fertilizers (KES 60/kg)
Tree Carbon Quantification
The study used a generalized improved pantropical mixed species model (5) to estimate above-ground biomass (AGB). Tree biomass assessment targeted two main carbon pools (stem and root biomass) for each tree with a DBH ≥ 5 cm. Field-based sampling with a nested concentric plot design was used to measure tree dimensions (including tree height, diameter at breast height and crown diameter).The outer circle radius 15 m was used to record and measure trees with DBH ≥ 20 cm, while a 10 m radius was used to measure trees with DBH ≥ 10 < 20 cm, and a radius 5 m was used to record and measure parameters for trees with a DBH ≥ 5 cm while a 2 m radius was used to measure trees with a DBH < 5 cm and seedlings.
\(\text{AGB=0.0673×}{\text{(ρ}{\text{D}}^{\text{2}}\text{H)}}^{\text{0.976}}\) 5
Where AGB is the above-ground weight of the tree (kg), ρ is the wood density D is the diameter at breast height in cm, H is the tree height while α and β are the model coefficients.
Total tree biomass was estimated by multiplying aboveground biomass by 1.25 (Chavan & Rasal, 2010). Aggregated carbon accounts for approximatly 47% of total biomass (IPCC, 2006, 2019). The respective wood densities was obtained from the wood density database (Zanne et al., 2009). Wood specific gravity was a particularly important predictor of AGB, especially when considering a wide range of vegetation types (Chave et al., 2014). The market prices were then used to estimate the economic value of the aggregated ecosystem carbon.
Soil Carbon Quantification
Soil carbon stocks were estimated from both soil organic carbon (SOC) and soil organic matter (SOM) levels. Organic matter (OM) content was determined using the loss on ignition method (LOI) while organic carbon (OC) was calculated using a ratio of 1:0.58 (SOM:SOC). Before the actual processing, a soil sample preparation was carried out. Samples were oven dried, crushed in mortar and pestle for homogenization, then sieved with a 2mm sieve to remove debris and stones, which were weighed separately. After sieving, the soil samples underwent a dry burning process required for carbon analysis to remove residual moisture. Two samples, each weighing 10 grams were placed in a pre-weighed crucible and then burned at 550°C for a minimum of 8 hours and then cooled before their weights were recorded. The difference in weights of the soil before and after heating represented the moisture and organic matter content, while the residue represented the ash. Soil organic carbon (SOC) was estimated by multiplying the weight difference by a factor of 0.58 as given in formula (6), while carbon per unit area was calculated by multiplying the SOC by the respective soil coefficients as given in formula (7).
\({\text{S}}_{\text{OC}}\text{=}\frac{\left({\text{ω}}_{\text{IS}}\text{-}{\text{ω}}_{\text{SR}}\right)\text{×0.58}}{{\text{ω}}_{\text{IS}}}\) 6
\({\text{T}}_{\text{OC}}\text{=(}\text{ρ×D}\text{×}{\text{S}}_{\text{OC}}\text{)100 }\) 7
Where SOC is soil organic carbon (%); ꞷIS is the initial weight of the soil sample; ꞷSR weight of soil residue after inceneration; TOC is total organic carbon (Mg of C per ha); ρ is the bulk density (g/cm); D is the soil tread depth (cm).
In mass calculations, soil samples were weighed for wet weight, air-dried at approximately 40°C for 48 hours, with an aliquot of each sample taken after weighing the air dried samples. The samples were further oven dried at 105°C for twenty-four hours and their weights were recorded. A total of three weights were recorded for each sample (i.e. total soil weight, the weight of the aliquot before oven drying at 105°C, and the weight after oven drying at 105°C) allowing calculation of bulk density.
The study used the market pricing function (Pearce, 2001) as indicated (8) to determine the value of forest carbon sequestration in contrast to the climate change damage function (Ferarro et al., 2011) with potential value overestimation.
\({\text{V}}_{\text{FCR}}\text{=}\sum _{\text{n=1}}^{\text{∞}}{\text{A}}_{\text{LC}}\text{×}{\text{Q}}_{\text{C}}\text{×}{\text{ε}}_{\text{C}}\) 8
Whereby V FCR is the economic value for climate regulation, ALC is area (ha) of the respective land cover, QC is the amount of carbon dioxide sequestered by the respective land cover per unit area, while ꜪC represents the average global carbon market price per unit of carbon.
Prices in global compliance markets currently range from less than US$1/tCO2e to US$30/tCO2e (AU$1–29/tCO2e). While considering the voluntary markets, average prices range from US$1/tCO2e to $5/ tCO2e or (AU$1–6/ CO2e)(World Bank Group, 2020). However, the study chose to use $5 per tonne of CO2 as the prevailing price for carbon traded in Kenya in the Voluntary Carbon Standard (VCS) REDD + market.
Data Analysis
Descriptive statistics in SPSS were used to summarize data in measures of central tendency, spread, and variance. The data were examined to establish agreement with normal distribution assumptions, and using either ANOVA or Friedman’s test depending on the examination outcome, to establish the relationship between study area parameters, for example land use as the independent variable and biomass and soil carbon as dependent variables. If necessary, a logarithmic transformation was carried out in order to comply with the normal distribution assumptions where necessary.