Land subsidence dynamics and their interplay with spatial and temporal land-use transitions in the Douala coastland, Cameroon

The Douala coastland (DCL), situated within the Douala sedimentary basin along the Gulf of Guinea, is characterised by its low elevation and alluvial geology, making it particularly susceptible to coastal erosion, land subsidence, and relative sea-level rise. The DCL is home to numerous rapidly growing cities, such as Douala, Tiko, and Limbe, which are currently experiencing alarming rates of coastal erosion, frequent flooding, and significant loss of land. Regional and continental investigations have provided evidence of coastal subsidence in this region; however, knowledge of its drivers and impact on the DCL remains limited. To address this knowledge gap, interferometric synthetic aperture radar (InSAR) datasets from the Sentinel-1 C-band satellite were used to quantify vertical land motion (VLM) between 2018 and 2023 with respect to the IGS14 global reference frame and assumed to represents absolute VLM. Digital Elevation Model datasets were used to analyse the elevation of the study area. The results revealed that the rate of VLM ranges from −17.6 to 3.8 mm yr⁻¹ (standard deviation of 0.2 mm yr⁻¹), with a mean and median land subsidence rate of 2.7 and 2.5 mm yr⁻¹. The analysis of land cover datasets from 1992 to 2022 suggests that urbanisation increased fivefold from 1992 to 2022 and that all contemporary urban areas experienced land subsidence, with the highest rates observed in non-residential zones with building heights ranging between 3 and 6 m. Subsidence rates of the DCL are inversely proportional to the time at which a particular land use and land cover (LULC) class changed into an urban area, highlighting the impact of the timing of LULC changes and urban expansion on present-day subsidence. The land subsidence rates decreased with an increase in building height, suggesting shallow land subsidence, where tall building usually with deep foundations are less prone to subsidence caused by groundwater extraction from the phreatic aquifer.

1 Introduction

Land subsidence is a key environmental phenomenon that exacerbates coastal cities' vulnerability worldwide. Gambolati and Teatini (2021) described this geological phenomenon as a sudden or gradual settling of the land surface caused by changes in the stress regime of the subsurface structure driven by a combination of natural processes and human activities (Shirzaei et al., 2021; Candela and Koster, 2022). Recent research has revealed a global distribution of land subsidence processes (Herrera-Garcia et al., 2021; Wu et al., 2022; Davydzenka et al., 2024), ranging from arid to coastal zones, urbanised cities, and emerging cities. In urbanised coastal cities, land subsidence is often driven by a combination of natural processes and human activities, with a predominance of extensive groundwater extraction (Wu et al., 2022; Davydzenka et al., 2024; Huning et al., 2024). For instance, Tokyo, Japan (Cao et al., 2021), Beijing, China (Zhu et al., 2015; Zhou et al., 2019, 2020), California and New York, USA (Parsons et al., 2023; Ohenhen et al., 2024;), Rotterdam and Almere, Netherland (Koster et Al., 2018; Verberne et al., 2023; Maoret et al., 2024), Maceio, Brazil (Mantovani et al., 2024), and Venice, Italy (Teatini et al., 2005; Da Lio and Tosi, 2018), are prominent examples of cities experiencing significant land subsidence. In emerging countries, rapid urbanisation and inadequate infrastructure planning often exacerbate land subsidence, particularly in low-elevation coastal zones. Jakarta, Indonesia (Taftazani et al., 2022), the Mekong Delta, Vietnam (Minderhoud et al., 2017, 2018), and Lagos, Nigeria (Ohenhen and Shirzaei, 2022) are prime examples of such challenges.

From a global perspective, land subsidence can be seen as an environmental hazard that negatively impacts the lives of almost one-fifth of the global population (Herrera-Garcia et al., 2021; Minderhoud et al., 2024). Thus far, land subsidence remains underestimated despite the progressive increase in vulnerable populations and infrastructure, especially in rapidly growing cities (Minderhoud et al., 2024). In low-lying coastal areas, the combination of land subsidence and local sea-level rise often exacerbates coastal cities' exposure to inundation and flooding (Shirzaei et al., 2021). Apart from groundwater extraction, land-use changes and natural compaction of alluvial soils greatly influence the occurrence of land subsidence (Minderhoud et al., 2018; Zoccarato et al., 2018). Urbanisation, agricultural practices, and industrial activities are key factors in land-use change. The conversion of natural landscapes into urban areas typically involves the construction of heavy buildings and infrastructure that can compress underlying soils (De Wit et al., 2021; Parsons et al., 2023). Urbanisation often leads to increased groundwater extraction to meet the demands of a growing population, further contributing to land subsidence (Zhu et al., 2015). In many regions, intensive agriculture requires substantial water for irrigation, often sourced from groundwater (Tzampoglou et al., 2023). Industrial activities, particularly large-scale water use and hydrocarbon extraction, also contribute to the subsidence. For example, the Groningen gas field in the Netherlands has shown significant compaction and subsidence since early 1960, corresponding to the extraction period (Thienen-Visser and Fokker, 2017). Analogous to other deltaic regions globally, the Northern Adriatic has been characterised by land subsidence of both natural and, more recently, anthropogenic origin, with varying rates contingent upon geological events and human activities (Tosi et al., 2010).

Urbanised cities or urbanised coastal landscapes without major recent urbanisation dynamics such as Tokyo, Los Angeles, and Venice have implemented various measures to mitigate subsidence. In contrast, several emerging cities and landscapes that have urbanised recently, such as Douala and Lagos, are grappling with the rapidity of the urbanisation process and inadequate infrastructure. A lack of understanding of subsidence drivers poses additional challenges (Dinar et al., 2021). Coastal cities along the Gulf of Guinea in Africa fall within the second category and, therefore, there is a need for detailed scientific investigations to understand the geomorphological mechanisms going on this zone. The Gulf of Guinea is known to be lowly elevated above sea level (Hauser et al., 2023) and bears major rapidly growing cities such as Douala, Lagos, Cotonou, Accra, Abidjan, and Saint-Louis. Tens of millions of inhabitants live there, and the coastlands bear the countries' respective critical infrastructure (Cian et al., 2019; Dada et al., 2023). The coastline is particularly exposed to erosion and sea-level rise (Ideki and Peace, 2024; Dada et al., 2024), and land subsidence exacerbates these processes, as evidenced by preliminary land subsidence investigations in Nigeria (Ohenhen and Shirzaei, 2022; Ikuemonisan et al., 2023), Ghana (Avornyo et al., 2023; Brempong et al., 2023), and Senegal (Cisse et al., 2024). In contrast, there have been no specific studies on land subsidence mechanisms for the coast of Cameroon, leaving it uncharted.

This study aims to fill the knowledge gaps regarding the land subsidence mechanism in the DCL. It is motivated by the necessity of quantifying land subsidence in the DCL, and understanding the driving mechanisms. Investigating the relations between land subsidence and land-use land cover changes in the DCLtherefore constitutes the first step to unravel the interlinkages between vertical land motion, physical setting and human activities along the coast of Cameroon.

2 Study areas

The coastline of Cameroon, located in the eastern part of the Gulf of Guinea (Fig. 1), covers over 400 km from the border with Equatorial Guinea in the southeast to the border with Nigeria in the northwest (Ondoa et al., 2018). The coast bears several river outlets, including the Sanaga River, Nyong River, Lobe, and Ntem rivers, as well as very important estuaries such as the Ntem estuary (from 3°34^′ N, 9°38^′ E to 2°20^′ N, 9°48^′ E), Cameroon estuary (3°58 N, 9°17^′ E to 3°34^′ N, 9°38^′ E), and Rio Del Rey estuary (from 4°30^′ N, 8°30^′ E to 4°14^′ N, 8°56^′ E) (Ondoa et al., 2018). The Cameroon Estuary, also known as the Wouri Delta, is located at the heart of the Douala sedimentary basin (19 000 km²), and hosts the city of Douala, which is one of the most important cities of Cameroon (Fossi et al., 2023). From a hydrological perspective, the Wouri estuary is drained by three coastal rivers, Mungo, Wouri, and Dibamba, with average annual flows of 8109, 16 109, and 4109 m³ s⁻¹, respectively (Fossi et al., 2023). The Wouri Delta is a critical region in Cameroon because of its rich biodiversity, economic significance, social and cultural value for local communities (Ondoa et al., 2018).

Figure 1 Countries of the Gulf of Guinea region bordering the Atlantic Ocean. The coastline of Cameroon is bounded by the Ntem estuary in the southeast and the Rio Del Rey estuary in the northwest. The inset shows the study area (Douala Coastland), Douala City, and Rio Del Rey. Basemap source: Esri.

The coastland of Cameroon is characterised by low-lying topography and wetlands dominated by mangroves, which are exposed to the constant threat of flooding (Kenneth, 2009; Ellison and Zouh, 2012). The wetlands are geomorphologically fragile zones prone to compaction, and human activities have intensified the alteration of the hydrogeomorphic characteristics of the area over time. In an area where the abundance of surface water in wetlands coexists with a scarcity of land, the challenge of exposure is further exacerbated by increasing population and anthropogenic activities. These changes include the conversion of wetlands for housing, industrial development, and the expansion of highways, whereas natural changes are caused by saltwater intrusion and land subsidence (Kenneth, 2009). In addition, Douala is the largest city and economic capital of Cameroon and bears the most important infrastructure, with a population of approximately 4 million inhabitants in 2022. The city shows a similar geomorphological and urban setting as other African (mega)cities affected by land subsidence; however local evidence of land subsidence does not exist for Douala so far but has only been indicated by studies on global and continental scale (see Cian et al., 2019; Davydzenka et al., 2024).

The low-lying coastal areas of Cameroon, particularly in regions such as the Douala coastal lowland and Rio Del Rey Estuary, are increasingly susceptible to a combination of sea-level rise, coastal erosion, flooding, and land subsidence (Fossi et al., 2019). Relative sea-level rise attributed to climate change and land subsidence exacerbates the impact of coastal flooding, rendering previously infrequent flood events more common and severe (Fossi et al., 2019). This elevated water level contributes to increasing erosion of the coastline, progressively diminishing land and exposing local inhabitants, infrastructure, livelihoods, and ecosystems (Fig. 2). Coastal erosion, intensified by frequent storms, induces shoreline retreat and destabilises coastal communities (Fossi et al., 2019). Traces of the impacts of these hazards along the low-lying coast of Cameroon are shown in Fig. 2.

Figure 2 Environmental hazards along the low-lying coast of Cameroon. (A) Destruction of school infrastructure due to intensified coastal erosion as a result of sea-level rise in Bakassi, Rio Del Rey Estuary (Wantim and Ngeh Ngeh, 2025). (B) Destruction of houses by sea-level rise and coastal erosion in Desbundscha (Wantim et al., 2021). (C) Flooding in coastal city along the Gulf of Guinea (Editor, 2024). (D) Shift of the coastline by 300 m over 30 years on Cap Cameroon Island (Fossi et al., 2019). (E) Shoreline retreat in Douala as a result of relative sea-level rise (Sikem, 2021). Panels (F)–(H) highlight the exposure of Manoka Island to sea-level rise, coastal erosion, and potential land subsidence. For (G), the yellow circle highlights a remaining tree trunk submerged in the waterbody, indicating that this area was previously emerged land. The red line indicates the distance between the current shoreline and the submerged building. (F–H) Photos taken in 2023 by author.

Pluvial flooding, which is a significant concern due to heavy rainfall, can be further aggravated by elevation lowering in areas affected by land subsidence (Miller and Shirzaei, 2019; Huning et al., 2024).

3 Material and Methods

3.1 Delineation of the Douala Coastland

Elevation was used to delineate the study area from the entire littoral region. Applying an elevation threshold allows to focus on the lowly elevated coastal areas which are particularly exposed to relative sea-level rise and related coastal hazards in the context of climate change. We used the digital elevation model (DEM) NASADEM Merged DEM Global 1 arcsec V001 and considered areas of the states bordering the Gulf of Guinea with elevations below 100 m (WGS84/EGM96 datum) as the area of interest (NASA, 2021).

3.2 SAR Data and Interferometric Analysis

We used synthetic aperture radar (SAR) data acquired from the Sentinel-1 C-band satellite to generate a land displacement map for DCL. The analysis processed 313 SAR images acquired with the Interferometric Wide (IW) swath mode in ascending (frame 132, path 1192) and descending (frame 51, paths 575, and 580) orbit geometries between 10 March 2018, and 17 June 2023. Figure 3 and Table S1 in the Supplement show the footprints of the SAR frames and inventory of images, respectively. We utilised GAMMA software (Werner et al., 2000) to generate a set of 1197 interferograms, employing a SBAS-type algorithm with a pair selection strategy optimized through dyadic temporal downsampling and Delaunay triangulation to minimize phase closure errors and temporal decorrelation (Lee and Shirzaei, 2023). The interferometric pairs were constrained to a maximum temporal and perpendicular baselines of 300 d and 80 m, respectively. To improve the signal-to-noise ratio of the interferometric phase, a multi-looking factor of 32 (range) by 6 (azimuth) was used to create a pixel dimension of ∼ 75 m, representing a trade-off between noise suppression over natural terrain and the preservation of spatial detail in urban environments. To generate high-spatial-resolution maps of land displacement over the study area, advanced multi-temporal wavelet-based interferometric synthetic aperture radar (InSAR) analysis was applied to the interferograms (Shirzaei and Bürgmann, 2012; Shirzaei, 2013; Lee and Shirzaei, 2023). The wavelet-based analysis involves identifying and removing noisy pixels, reducing the effects of topographically correlated atmospheric phase delay and spatially uncorrelated DEM error. To this end, we identified and removed noisy pixels following Lee and Shirzaei (2023), by applying a statistical framework that selects elite pixels accounting for distributed and permanent scatterers. Permanent scatterers (PS) were identified using an amplitude dispersion index < 0.3 (Ferretti et al., 2001), while distributed scatterers (DS) were assessed using the mean and variance of coherence values across the interferogram stack. DS pixels with an average coherence > 0.65 and low coherence dispersion (Fig. S4) were retained and assigned to Voronoi cells constructed around adjacent PS. Within each Voronoi cell, a Fisher's F-test was used to evaluate variance homogeneity between DS and PS amplitude histories. DS pixels passing this test were reclassified as permanent-distributed scatterers (DSp). The final set of elite pixels therefore comprised PS and DSp (see Lee and Shirzaei, 2023 for further details). Next, the 1 arcsec (∼ 30 m) Shuttle Radar Topography Mission digital elevation model (Farr et al., 2007) was used to calculate and remove the effects of the topographic phase and flat Earth (Franceschetti and Lanari, 1999). A 2D minimum cost flow phase unwrapping algorithm modified to be applied to a sparse set of elite pixels was to estimate absolute phase changes in each interferogram (Costantini, 1998; Costantini and Rosen, 1999). We corrected for the effects of orbital error (Shirzaei and Walter, 2011), spatially uncorrelated topography error, and the topography-correlated component of atmospheric delay by applying a suite of wavelet-based filters (Shirzaei and Bürgmann, 2012). Finally, the line-of-sight (LOS) time series and velocities for each elite pixel were computed using a local reference point (longitude: 9.61, latitude: 4.22) located in a stable area within the broader Douala coastland region but outside the urban center. The LOS displacement velocity was estimated using the reweighted least-squares method as the slope of the best-fitting line to the associated time series, while the standard deviation of the LOS velocity corresponds the uncertainty of the regression slope. The LOS displacement velocities and associated standard deviation values for the ascending and descending orbit geometries are presented in Figs. S1 and S2, respectively.

Figure 3 Footprint of InSAR images (ascending and descending Sentinel-1A/B) in the Douala coastland, including the position of the local reference point. Basemap source: Esri.

In areas with overlapping spatiotemporal SAR satellite coverage and different viewing geometries (i.e., ascending and descending), surface deformation can be resolved. However, because of the near-polar orbits of the Sentinel-1 satellites, which may cause low sensitivity in the north-south component, only the east-west and vertical components can be reliably retrieved (Wright et al., 2004). Assuming that the deformation in the north-south direction is negligible, we estimate the horizontal (east-west) and vertical components of deformation by jointly inverting the LOS time series of the ascending and descending tracks (Miller and Shirzaei, 2015; Ohenhen and Shirzaei, 2022). Thus, we first identified the co-located pixels of the LOS time series by resampling the pixels from the descending track onto the ascending track to obtain two co-located LOS displacement velocities. Given {LOS_A,LOS_D} and $\mathit{\{}{\mathit{\sigma }}_{\mathrm{A}}^{\mathrm{2}},{\mathit{\sigma }}_{\mathrm{D}}^{\mathrm{2}}\mathit{\}}$ are the LOS displacement and variances for a given pixel, where the subscripts A and D indicate the ascending and descending tracks, respectively. The model computed the horizontal (d_x) and vertical (d_z) components from the LOS displacement as shown in Eq. (1):

$$\begin{array}{}\text{(1)}& \left[\begin{array}{c}{\text{LOS}}_{\mathrm{A}}\\ {\text{LOS}}_{\mathrm{D}}\end{array}\right]=\left[\begin{array}{cc}{C}_{\mathrm{A}}^x& C_{\mathrm{A}}^z\\ C_{\mathrm{D}}^x& C_{\mathrm{D}}^Z\end{array}\right]\left[\begin{array}{c}{d}_x\\ d_z\end{array}\right],\end{array}$$

where C^x and C^z are the unit vectors projecting the 2D displacement (horizontal and vertical) field onto the LOS direction, which are a function of the heading and incidence angles for each pixel (Hanssen, 2001). The solution to Eq. (1) is given by Eq. (2):

$$\begin{array}{}\text{(2)}& X=(G^{\mathrm{T}}\mathbf{P}G{)}^{-\mathrm{1}}{G}^{\mathrm{T}}\mathbf{P}L,\end{array}$$

where X represents the unknowns (d_x and d_z), G is the Green's function given by the unit vectors (C^x and C^z), P is the weight matrix, which is inversely proportional to the variances of the observations (${\mathit{\sigma }}_{\mathrm{A}}^{\mathrm{2}}$ and ${\mathit{\sigma }}_{\mathrm{D}}^{\mathrm{2}}$), and L is the observation (LOS_A and LOS_D). Next, we employ the concept of error propagation (Mikhail and Ackermann, 1982) to obtain the associated parameter uncertainties given the observation errors. Thus, the variance-covariance matrix of the horizontal and vertical components was calculated using Eq. (3):

$$\begin{array}{}\text{(3)}& D_X=(G^{\mathrm{T}}\mathbf{P}G{)}^{-\mathrm{1}}.\end{array}$$

The horizontal velocity, VLM rate, and standard deviation are shown in Fig. S3A and B, respectively (see the Supplement). To transform the VLM from a local frame to a global reference frame, we employed the global VLM model by Hammond et al. (2021), which mostly contains long-wavelength deformation signals (i.e., glacial isostatic adjustment (GIA), hydrological loading, and tectonics), and applied an affine transformation following Blackwell et al. (2020) and Ohenhen et al. (2023) to align the VLM rates with the IGS14 global reference frame. The IGS14 global reference frame was used because the study area has no GNSS station to serve as a local reference point and helps to improve the accuracy, calibration, and interpretation of the InSAR measurements.

3.3 Land use, land cover (LULC) and human settlement data analysis

We investigate trends in and spatiotemporal correlations between land subsidence and past land-use changes, following an approach similar to that of Minderhoud et al. (2018). The analysis included the InSAR data described in the previous section and a collection of Land-Use Land Cover (LULC) data that was prepared (reclassified) for spatial and temporal analysis. To further investigate the relations between land subsidence and urban built-up areas, we integrated additional information characterising built-up settlements.

Information on LULC was obtained from the European Space Agency (ESA) global land cover data for the period of 1992 to 2015 with a resolution of 300 m (ESA Land Cover CCI project team and Defourny, 2019), and Esri Land Cover products for the period of 2017 to 2022 with a resolution of 10 m. We use the Esri's Land Cover (Esri) products for the most recent years as it shows a high accuracy than the ESA land cover dataset (Venter et al., 2022).

Information about the setting of built-up areas was obtained by using the urban built-up surface data (GHS-BUILT-S; Pesaresi and Politis, 2023c), built-up volume data (GHS-BUILT-V; Pesaresi and Politis, 2023d), spatial built-up height data (GHS-BUILT-H; Pesaresi and Politis, 2023b), spatial population raster density (GHS-POP; Schiavina et al., 2023), and morphological settlement zone characteristics (GHS-BUILT-C; Pesaresi and Panagiotis, 2023a). All these datasets were acquired from the Global Human Settlement Layer (GHSL) website. The concept and methodology used to obtain GHSL data were described in European Commission (EC) (2023) for GHS-BUILT-S, GHS-BUILT-C, and GHS-BUILT-V data; Pesaresi and Politis (2023b) for GHS-BUILT-H data; and Freire et al. (2016) for GHS-POP. Recent similar applications of GHSL datasets for land subsidence studies are illustrated by the work of Verberne et al. (2025) and Cigna et al. (2025).

The data were pre-processed to ensure quality and consistency between the LULC classes of different years. This involves aligning spatial data layers to a common coordinate system and resolution, extracting data for the study area, and reclassifying the LULC data to classes appropriate for the goal of this study. The reclassification of land cover data was performed only for the ESA global land cover data from 1992 to 2015 to reduce the number of classes by unifying them and making them similar to Esri's land cover classes. This was performed using ArcGIS Pro 3.2 Spatial Analyst tools. The original 17 land cover classes were reclassified into 7 new classes as shown in Table 1. The LULC classes from 2017 to 2022 were not reclassified as they already met the scope of our analysis.

Table 1 Land-use and land-cover reclassification table for the European Space Agency (ESA) dataset (land-cover data 1992–2015).

We produced LULC maps of the DCL for the year 1992, 2007, 2015, 2017, and 2022 using ArcGIS Pro 3.2. The statistics for each land cover class in each year were retrieved using Excel and Power BI software.

In addition to the large study area, we also performed a local analysis specifically for the urban city of Douala, where ArcGIS Pro 3.2 Spatial Analyst tools were used to extract, preprocess, analyse, and describe spatial patterns of specific LULC type, occupation, and morphological characteristics (urban landscape) from the GHSL data in the year 2020. The outcomes included a map of spatial built-up volume, spatial population density, spatial built-up height, spatial built-up surface, and morphological settlement zone.

LULC change analysis was performed to explore spatio-temporal patterns and shifts from 1992 to 2022. We used the LULC maps processed for 1992, 2007, 2015, 2017, and 2022 and focussed three different periods of 1992–2007, 2007–2015, and 2017–2022. Differences between LULC classes and indices across different periods were quantified using the ArcGIS Pro Raster Calculator tool. The output maps obtained highlight the areas where LULC changes between the initial and final years occurred.

3.4 Spatiotemporal correlation between land subsidence and land-use land cover changes

Land subsidence and LULC change are two significant environmental processes that are often interlinked, leading to various socioeconomic and ecological impacts (Minderhoud et al., 2018). Understanding the spatiotemporal correlation between these processes is crucial for sustainable urban planning and resource management. To investigate the relationship between the two variables in the DCL of Cameroon, we consider spatial patterns and temporal dynamics following the investigations of Minderhoud et al. (2018). We hypothesise that areas experiencing rapid LULC changes are often associated with noticeable land subsidence, and the timing of LULC change can influence the rate and extent of land subsidence. To investigate this hypothesis, we consider LULC classes intersecting with VLM pixels.

We quantified land subsidence and characterised the LULC patterns in the study area by applying spatial analysis tools in the GIS. We visualised the spatial distribution of land subsidence and land-use patterns from 1992 to 2022, and assessed the relationship between land subsidence, LULC, and urban city morphological characteristics statistically. Regression analysis was used to correlate the relationship between land subsidence rates and the various urban morphological characteristics from the GHSL. Correlation coefficients were explored to identify significant associations. Time-series graphs of land subsidence rates were then plotted for different LULC classes and morphological characteristics within the city of Douala.

4 Results

4.1 Elevation of the Douala Coastland

The Douala Coastland was delineated, covering a surface area of 2700 km², with 4.5 % elevation below 5 m, 7 % elevation between 5 and 10 m, and 13 % elevation between 10 and 20 m above sea level (Fig. 4). Within this area, the city of Douala has key infrastructure, such as the Douala International Airport, the Douala seaport, and industrial plants. In addition, the study area includes the Douala Wetland area, Manoka Island, and Cap Cameroon Island, which appear to have very low elevations.

Figure 4 (A) Elevation map of the Douala coastland with the main rivers. (B) Frequency distribution of elevation counts in %. Data source: NASADEM Merged DEM Global 1 arcsec V001, vertically referenced to the EGM96 geoid (NASA, 2021).

4.2 Vertical Land Motion in the Douala Coastland

InSAR analysis revealed that the DCL is experiencing vertical land motion (VLM) rates ranging from −17.6 to 3.8 mm yr⁻¹, with a median of −2.5 mm yr⁻¹ and a mean of −2.7 mm yr⁻¹ (p value < 0.0001 for 99 % confidence interval). The InSAR estimates are referenced to the IGS14 global reference frame and assumed to represents local VLM (Fig. 5). 97.9 % of the 68 680 scattering points in the study area showed land subsidence (VLM < 0 mm yr⁻¹). Statistical analysis of the VLM further revealed that the remaining 2.1 % with positive vertical land motion rate lies between 0 and 4 mm yr⁻¹ (p value < 0.0001 for 99 % confidence interval). A statistical description of the VLM rate range within the study area is shown in Fig. 5.

Figure 5 VLM referenced to IGS1G and associated uncertainty derived from Sentinel-1 InSAR data over the Douala coastland. (A) Mean VLM rates (in mm yr⁻¹), where negative values (yellow to red) indicate land subsidence and positive values (blue) indicate shallow uplift. (B) Standard deviation of the VLM estimates, illustrating the spatial variability of measurement uncertainty.

Across all seven selected sites in Douala in Fig. 6, land displacement between 2018 and 2023 reveals a clear spatial contrast between stable or uplifting natural/coastal zones and strongly subsiding reclaimed or heavily urbanized areas, as shown in Fig. 6. Cap Cameroon is the only location exhibiting sustained uplift, while Manoka Island and the early city zone show relatively minor subsidence, reflecting limited anthropogenic loading and more consolidated ground conditions. In contrast, the airport area, seaport zone, and especially the reclaimed urban expansion zones display progressive and substantial land subsidence in the most recently reclaimed area consistent with ongoing consolidation of fill material, heavy infrastructure loading, and groundwater variability.

Figure 6 Time series of VLM from 2018 to 2023 at the selected locations shown in Fig. 5A.

4.3 Land-use and land cover changes in the Douala Coastland

For the investigation period of 1992 to 2022, we document a constant decrease in the tree cover class, from 50.6 % of the total study area in 1992 to 38.7 % in 2015 (Fig. 7). High-resolution land-use data show a decrease in tree cover from 86.6 % in 2017 to 80.9 % in 2022. The reduction of forested areas occurs in favour of built-up, cropland, and rangeland areas as the respective LULC classes record substantial gains within the same period (Table 2). The cropland class increased from 14.4 % in 1992 to 20.6 % in 2015. Between 1992 and 2022, the built area class increased from 1 % to 4.7 %, reflecting the rapid urbanisation rate of Douala City, which has already been documented by previous studies (UN-Habitat, 2022).

Figure 7 Land use and land cover (LULC) of Douala coastland in (A) 1992, (B) 2007, (C) 2015, (D) 2017, and (E) 2022. Panels (A)–(C) were extracted from the European Space Agency (ESA) Global Land Cover Maps, Version 2.0.7 datasets with a resolution of 300 m × 300 m. Panels (D)–(E) were extracted from the Sentinel-2 10-Meter LULC and Esri Land Cover (Esri) datasets. Tree cover definitions differ between ESA (A–C) and ESRI (D, E): it represents significant clusters of tall (∼ 15 m or higher) dense vegetation, typically with a closed or dense canopy, for the former, and areas dominated by trees with > 10 % tree cover for the latter.

Table 2 Distribution (in %) of Land use and land cover (LULC) classes for the years 1992, 2007, 2015, 2017, and 2022 in the Douala coastland. The percentage of classes cannot be summed to 100 because cloud-covered areas were excluded.

The outcomes of the land-use change analysis revealed that the main land cover classes undergoing significant dynamics included tree cover, cropland, rangeland, and built-up areas. Figure 8 shows the different LULC change maps for 1992–2007, 2007–2015, and 2017–2022. To further explore the different LULC changes and significance between 1992 and 2022, a LULC change matrix for each change interval was developed. The transition matrices for the periods 1992–2007 (Table 3), 2007–2015 (Table 4), and 2017–2022 (Table 5) provide quantitative changes between different classes.

Figure 8 Land use and land cover (LULC) changes in Douala coastland over the following time periods: (A) 1992–2007, (B) 2007–2015, and (C) 2017–2022.

Table 3 Land use and land cover (LULC) changes from 1992 to 2007 expressed as a percentage of the total study area, excluding cloud-covered regions. A color-coded scale from green to red is used to highlight null to large values, respectively. LULC remained unchanged in 89 % of the study area, see the values in italics located along the diagonal.

Table 4 Land use and land cover (LULC) changes from 2007 to 2015 expressed as a percentage of the total study area, excluding cloud-covered regions. A color-coded scale from green to red is used to highlight null to large values, respectively. LULC remained unchanged in 94 % of the study area, see the values in italics located along the diagonal.

Table 5 Land use and land cover (LULC) changes from 2017 to 2022 expressed as a percentage of the total study area, excluding cloud-covered regions. A color-coded scale from green to red is used to highlight null to large values, respectively. LULC remained unchanged in 90 % of the study area, see the values in italics located along the diagonal.

In Douala City, analyses of urban characteristics reveal a heterogeneous and dynamic city with varying population density, built-up height, and built-up surface (Fig. 9). In 2023, Douala City had an average spatial population density of 200 inhabitants per 10 000 m² (ranging from 1 to 1072 inhabitants). At the same spatial scale, the average and median built-up surface percentages were 32 % and 35 % of the total built-up area, respectively. The mean and median building heights were 3 and 4 m, respectively. This shows that the level of building construction in terms of building height in Douala remained low. Combining these urban characteristics provides preliminary information on factors that could determine urban-differential soil consolidation.

Figure 9 Landscape characteristics of Douala urban zone in 2020: (A) built-up surface intensity, (B) building height, (C) built-up volume, and (D) population density map. The values refer to a raster resolution of 100 m × 100 m. Basemap source: Esri.

4.4 Relation between land subsidence rate and LULC change

We investigated whether the documented LULC change patterns from 1992 to 2022 correlate with the measured land subsidence rates, considering four main categories of LULC change trajectories that potentially affect land subsidence. Each trajectory characterises the LULC change scenario occurring in a class from 1992 to 2022 (Table 6).

Table 6 Statistical description of VLM with respect to various categories of LULC changes and urbanisation between 1992 and 2022. A color-coded intensity was used, independently in each column and within each path, to highlight land subsidence (negative VLM) amount. Darker colours indicate more severe land subsidence.

The different paths are as follow:

• Path 1 (C_1 to C_3) represents cropland areas in 1992 were progressively converted into built-up areas;

• Path 2 (R_1 to R_3) represents rangeland areas in 1992 were converted progressively into built-up areas;

• Path 3 (W_1 to W_3) represents wetland areas in 1992 were progressively converted to built-up areas;

• Path 4 (T_1 to T_4) represents the tree cover areas in 1992 were progressively converted into built-up areas.

We used as a reference the areas that were urbanised in 1992 and remained so throughout the investigated period. Inspection of Table 6 reveals that tree cover areas that recently (2017) changed to built-up areas are characterised by the highest mean and median rate of land subsidence. We observed similar results in all categories, except for the one in which the initial class in 1992 was cropland (Path 2), where no correlation between LULC change trajectories and subsidence rate was observed. In contrast, the areas that remained in the built-up area throughout 1992–2022 had the lowest mean and median subsidence rates (Fig. 10).

Figure 10 Relationship between vertical land motion (based on InSAR for 2018–2023) and LULC change trajectories (1992–2022) as defined in Table 6. Built-up zones resulting from recent (i.e., over the period from 2015 to 2017) land-use conversion (trajectories C_3, W_3, and T3) show substantially higher land subsidence than zones urbanized earlier, indicating a correlation between built-up period and land subsidence.

Further analysis was carried out to investigate the possible correlation between land subsidence and landscape within the urban city of Douala. The following landscape characteristics were considered: vegetation density, wetland surfaces, road surfaces, building functional use, and building height (Fig. 11). It appears that the mean subsidence rate in vegetated areas is relatively low at approximately 2.5 mm yr⁻¹, irrespective of low-, medium-, and high-vegetation surfaces. The highest subsidence rate was observed in non-residential (industrial) zones, with buildings heights ranging from 3 to 6 m. This range corresponds to the mean and medium ranges of building height. In contrast, comparatively tall buildings in both residential and non-residential zones show very low subsidence, and the subsidence rate decreased with increasing building height.

Figure 11 Douala city: VLM rates derived by InSAR between 2018 and 2023 of various landscapes (wetlands, vegetation, roads, residential and non-residential zones) for various built-up landscapes. Notice that land subsidence shows a negative correlation to building heights (BH).

5 Discussion

5.1 Observed rates of vertical land motion

Our study reveals that the Douala Coastland experienced InSAR-estimated VLM ranging from −17.6 to 3.8 mm yr⁻¹ between 2018 and 2023, resulting in a mean land subsidence rate of approximately 3.0 mm yr⁻¹ (Fig. 5). A comparable value was associated with Douala based on the results obtained from a global study by Davydzenka et al. (2024), in which they implemented deep learning approaches combined with worldwide subsidence measurements to predict global land subsidence-susceptible areas based on 23 environmental parameters. Our results also align with the findings of Cian et al. (2019), who used InSAR to map land subsidence in 18 African coastal cities and determined that 17 out of 18 experienced subsidence, including Douala. This significant rate poses substantial threats to the region, exacerbating the impacts of sea-level rise, a critical concern for coastal communities worldwide (Nicholls et al., 2021). Furthermore, we find that 97.8 % of the coastal lowlands of Douala exhibit subsidence. Evidence of land subsidence in the Douala coastland is not an isolated case in the Gulf of Guinea, and similar investigations have highlighted the exposure of other coastal cities in Ghana (Avornyo et al., 2023; Brempong et al., 2023), Nigeria (Ohenhen and Shirzaei, 2022; Ikuemonisan et al., 2023), and Senegal (Cisse et al., 2024). This suggests that most urban cities along the Gulf of Guinea are undergoing subsidence, indicating an uncertain future for the Douala coastland and Gulf of Guinea in general. The occurrence of inundation and building collapse that has been observed in Douala (Tchamba and Bikoko, 2016) may be attributed to land subsidence, as reported in Nigeria by Ohenhen and Shirzaei (2022), who revealed that approximately 4000 buildings are exposed to a high to very high risk of long-term collapse in Lagos.

The combined effects of land subsidence and sea-level rise increase the vulnerability of coastal infrastructure, ecosystems, and populations, necessitating urgent attention and intervention strategies. Therefore, it is essential to provide a comprehensive understanding of land subsidence dynamics in the Douala coastland, examining the values in the context of other coastal cities worldwide that face similar challenges. Jakarta in Indonesia is one of the most severely affected cities by land subsidence, with contemporary rates reaching 10 cm yr⁻¹ (Taftazani et al., 2022). This rapid subsidence can be primarily attributed to excessive groundwater extraction. In Italy, Venice faces subsidence rates of approximately 1–2 mm yr⁻¹, compounded by rising sea levels (Teatini et al., 2005). Comparative analysis shows that the current coastland subsidence in Douala, averaging 3 mm yr⁻¹, is significantly smaller than those experienced by other cities, such as Jakarta. However, the rate is comparable to the sea-level rise. This is particularly significant given the low elevation of Douala above sea level (Fig. 4; Kemgang et al., 2024).

Moreover, this study offers a snapshot of current VLM rates, which could change rapidly particularly under large-scale groundwater exploitation (Gambolati and Teatini, 2015; Minderhoud et al., 2020). The consequences of land subsidence for Douala may be far-reaching, affecting not only infrastructure and human settlements but also the health and sustainability of mangrove ecosystems, which are essential for coastal protection and biodiversity. Hence, further detailed investigations are required to elucidate the factors causing land subsidence in the Douala coastland, mitigate its occurrence, predict and prevent larger magnitudes in the future.

In addition to the uncertainty linked to other factors influencing land subsidence, we observed that 2.1 % of the SAR scattering points have positive vertical land motion with values ranging from 0 to 4 mm yr⁻¹. Positive values of VLM can result from the effect of unaccounted horizontal motions in the north-south (NS) direction, which may entirely map onto the vertical signal due to limitations in the processing procedure. These scattering points are in a natural wetland area with little or no human activities. Another potential explanation is that the uplift results from sediment accumulation because of river-ocean interaction. Further investigations are needed to understand why this area shows uplift.

5.2 Landuse land cover (LULC) change

The Douala coastland has experienced significant LULC change dynamics (Fig. 7). An increasing trend was observed in areas allocated to croplands, built-up areas, wetlands, and rangelands, accompanied by a substantial decrease in tree cover (Table 2). The five-fold expansion of built-up areas is attributed to population growth and economic development, a phenomenon common in many rapidly growing and developing urban centers. These changes have far-reaching implications for environmental sustainability, ecosystem services, and the overall resilience of the region. The expansion of the anthropogenic footprint correlates with the study by Nfotabong-Atheull et al. (2013), which revealed that 39.9 % of mangrove forests disappeared in the peri-urban setting of Douala from 1974 to 2003, with a net loss of 22.1 % between 2003 and 2009 alone. Nfotabong-Atheull et al. (2013) highlighted a substantial increase in settlements (60 %), agricultural areas (16 %), non-mangrove areas (193.3 %), and open water (152.9 %) around Douala over a 35 year period from 1974 to 2009.

The transition matrices used to investigate the temporal dynamics between LULC classes revealed that the LULC change mechanism is not uniform over time. Between 1992 and 2007, there was a significant increase in cropland and rangeland, mainly because of a decrease in tree cover. The increase in the built-up area during the same interval was due to changes in three different classes: cropland, rangelands, and tree cover (Table 3). These changes were less significant between 2007 and 2015, which is probably also due to the shorter time range considered (eight years compared to 15 years). Tree cover areas remain the most dynamic class among the LULC classes. The major changes were from tree cover to rangeland and cropland. Between 2017 and 2022, when higher-resolution land cover datasets were used, the dynamics between classes were notably (Table 5). Approximately 6 % of the study area lost tree cover, and transformed the land into rangelands (3.3 %), wetlands (1.7 %), built areas (0.9 %), and croplands (0.8 %).

5.3 Land-use-driven land subsidence

Anthropogenic activities related to land-use practices, such as agricultural activities, groundwater extraction for irrigation, industries, and heavy infrastructure, directly influence subsidence (Galloway and Burbey, 2011; Gambolati and Teatini, 2021). We conducted a statistical analysis of the land-use change trajectories and subsidence rates within specific LULC classes to assess the correlation between the two processes. The results indicate that the mean and median land subsidence rates in each category, except cropland, increased with recent changes in the class to built-up areas. Relatively low subsidence was observed in the trajectories where the change to built-up areas occurred earlier. Indeed, areas that remained built-up from the beginning (1992) to the end (2022) of the investigated period exhibited the lowest average subsidence rate (1.9 mm yr⁻¹), while areas that remained tree cover from 1992 to 2017, and subsequently changed to built-up areas over the last five years demonstrated the highest average subsidence rate (7.1 mm yr⁻¹). This observation also applies to wetland and rangeland categories (Table 6 and Fig. 10). This trend can be attributed to the fact that recently urbanised areas are still undergoing consolidation due to the recent loading increase exerted by structures on the land surface. Moreover, tree root degradation during land-use change potentially contribute to land subsidence (Huning et al., 2024). Croplands, built-up areas, wetlands, and rangelands have notable implications for soil consolidation, which is a critical factor that affects the stability and sustainability of an area. Therefore, an appropriate land-use management mechanism must be considered as an adaptive or mitigating measure for controlling land subsidence.

With reference to the urban city of Douala, we investigated the influence of various urban morphological (landscape) characteristics, including wetlands, vegetation zones, roads, and residential and non-residential zones on land subsidence rates. The results revealed significant variations in land subsidence between the residential and non-residential zones. In residential zones, the highest subsidence rate was observed when building heights ranged between 3 and 6 m, whereas taller buildings exhibited relatively lower land subsidence. A similar trend was observed in the non-residential zones (Fig. 11). Taller buildings experienced less subsidence, potentially because of their deep foundations, which made them more stable and unaffected by shallow soil compaction, such as auto-compaction and lowering of the groundwater level in the phreatic aquifer system (Minderhoud et al., 2025). De Wit et al. (2021) reported similar observations in the Mekong Delta, as did Koster et al., (2018), Verberne et al., (2023), and Maoret et al., (2024) in cities in The Netherlands. Nevertheless, we acknowledge that the ability to distinguish deformation patterns between tall and low-rise buildings could be improved using a Persistent Scatterer InSAR (PS-InSAR) processing chain. The InSAR dataset employed in this study provides first-order insights into spatial deformation trends. However, the processing approach used introduces spatial smoothing, which may obscure localized variations in displacement, particularly in densely built-up areas. Future studies applying PS-InSAR methods with higher temporal and spatial sampling could therefore better resolve these fine-scale deformation differences and strengthen the interpretation of structural–geotechnical interactions across the urban landscape. In this context, this study provides a foundation for future detailed investigations of the vulnerability of the Douala coastland to the compound effects of environmental hazards associated with differential land subsidence and consequent relative sea-level rise.

In addition, groundwater pumping could significantly contribute to land subsidence especially in Douala, a rapidly growing coastal city with increasing groundwater demand for population and industrial activities. However, groundwater extraction data are scarce for the study area, limiting our possility to quantify aquifer depletion and its contribution to land subsidence. Nonetheless, observations from the GRACE (Gravity Recovery and Climate Experiment) satellite mission indicate a progressive decline in terrestrial water storage between 2017 and 2023, suggesting ongoing groundwater loss at the regional scale. However, the GRACE data have a coarse spatial resolution (∼ 300–400 km), which constrains their capacity to resolve localized groundwater variations within the Douala coastal plain. As a result, while GRACE provides useful large-scale trends, it does not capture the fine-scale dynamics necessary to establish a direct spatial correspondence between subsidence patterns and groundwater extraction. This calls for the need to establish in-situ monitoring network for data collection.

6 Conclusion

We quantified vertical land motion in the Duala coastland between 2018 and 2023 using interferometric synthetic aperture radar from Sentinel-1 and investigated its relations to LULC change during the period 1992 to 2022. We found that VLM in the DCL ranges between −17.6 to 3.8 mm yr⁻¹, while the DCL on average subsides 2.7 (median 2.5 mm yr⁻¹). We documented a rapid, fivefold increase in urbanisation and that all contemporary urban areas experienced land subsidence, with the highest rates of 2.3 to 7.1 mm yr⁻¹ occurring in non-residential zones with building heights between 3 and 6 m. In contrast, lower rates (1.9 to 4.4 mm yr⁻¹) were found in older built-up areas, as most consolidation from building loads had occurred prior to the observation period. Subsidence rates of the DCL are inversely proportional to the time at which a particular LULC class changed into an urban area, highlighting the impact of the timing of LULC changes and urban expansion on present-day subsidence. The subsidence rates decreased with an increase in building height, suggesting the potential influence of foundation type on land subsidence as tall buildings are often associated with deeper foundations which likely makes these buildings less susceptible to shallow subsidence processes.

Land subsidence in the DCL presents a growing challenge, potentially driven by both natural and anthropogenic processes, analogous to the issues faced by other coastal cities globally. Over the past three decades, rapid urbanisation and LULC changes have likely contributed to substantial land subsidence rates, particularly in newly developed non-residential areas, specifically industrial zones and reclaimed lands. To address ongoing subsidence and its potential to exacerbate coastal erosion, relative sea-level rise and inundation, effective monitoring and mitigation strategies should be developed and implemented. However, it is imperative to acknowledge that multiple factors can cause land subsidence, and for the Douala coastland and city, more comprehensive studies are needed to thoroughly understand the ongoing and potential future land subsidence mechanisms. Investigating the effects of subsidence drivers such as groundwater extraction, soil compaction, building loading, and surface hydrology variation is essential to safeguard the rapid growth of urban zones on the Douala coastland and ensure their resilience against the adverse impacts of land subsidence, relative sea-level rise, and climate-related hazards such as flooding.