A multi-class, multi-temporal crop and Land cover mapping framework for Morocco using Sentinel-1/2 monthly composites and advanced machine learning ensembles

Timely and accurate crop-type mapping is fundamental for sustainable agricultural management and food security in semi-arid regions, where climate variability and fragmented landscapes present persistent challenges. This study develops and validates an operational, multi-temporal framework for classifying five key agricultural classes (i.e., soft wheat, durum wheat, barley, trees, and other crops) across diverse Moroccan agroecosystems. By integrating monthly Sentinel-1 Synthetic Aperture Radar and Sentinel-2 optical time series spanning six growing seasons (2018–2025), we extracted 156 features comprising 13 spectral indices across 12 monthly composites. Ground truth data from the national Al Moutmir database, strategically balanced to address natural class imbalances, supported comprehensive training and validation of six machine learning models (i.e., Random Forest, Extra Trees, XGBoost, LightGBM, Voting Ensemble, and Stacking Ensemble). Our findings showed that the LightGBM and Stacking Ensemble achieved the highest performance with 88.04% overall accuracy, followed closely by XGBoost (87.93%). Feature importance analysis revealed that monthly temporal resolution significantly outperformed traditional phenological-stage approaches, with March and April indices (particularly Normalized Difference Vegetation Index and Normalized Difference Red Edge) contributing most to class discrimination. Notably, early-season radar features (Vertical-Vertical polarization in September) provided valuable complementary information when optical data were limited. The framework demonstrated robust generalization through 10-fold cross-validation while explicitly quantifying a 12.56% overfitting gap (train-CV difference), acknowledging a non-negligible overfitting risk. Offering transparent performance assessment. Error analysis identified persistent confusion between spectrally similar cereals, particularly durum and soft wheat, highlighting priority areas for future sensor integration. This scalable, cloud-based pipeline directly supports Morocco’s Green Generation strategy by providing a reproducible, high-accuracy solution for annual crop inventories, with transferable applications across similar Mediterranean and semi-arid agricultural systems.

Over the last decade, the field of remote sensing (RS) has undergone a profound transformation driven by the rapid growth of open-access and high-resolution satellite imagery (Choukri et al., 2024). Initiatives such as the European Space Agency’s Copernicus program, which provides free and systematic coverage through Sentinel-1 and Sentinel-2, and NASA’s Landsat program have dramatically expanded the temporal, spectral, and spatial coverage of earth observation data worldwide (Sudmanns et al., 2020). This unprecedented availability of multi-source imagery often referred to as the era of “open-access data” has universalised access to remote sensing and the development of operational applications across multiple disciplines, including environmental monitoring, hydrology, disaster management, and agricultural production systems (Kussul et al., 2019; Orynbaikyzy et al., 2019). For agriculture in particular, near-real-time monitoring of vegetation dynamics, crop phenology, and land use changes has become an essential tool for improving food security, minimizing climate impacts, and guiding national and regional policy decisions (Wijitdechakul et al., 2016).

The integration of satellite remote sensing and artificial intelligence (AI) has fundamentally changed agricultural monitoring (Sharma et al., 2024). However, the application of these advanced methodologies, especially the combination of high-resolution and machine learning, remains notably limited in the African context (El Bouanani et al., 2025). While the fusion of Sentinel-1 and Sentinel-2 data with machine learning has proven effectiveness for crop-type mapping in different contexts (Orynbaikyzy et al., 2019; Dobrinić et al., 2021), an overfitting gap persists in translating these methodologies into operational, national-scale frameworks adapted to specific agro-ecosystems. This is particularly remarkable in the semi-arid, fragmented landscapes of North Africa, where smallholder dominance, irregular field geometries, and heterogeneous cropping calendars present unique challenges (Choukri et al., 2024). Although previous studies in Morocco have demonstrated the potential of remote sensing for crop and soil monitoring (Khechba et al., 2021; Misbah et al., 2022), most have been limited in spatial extent, temporal duration, or focused on a single sensor. Consequently, there is a critical lack of scalable, robust, and multi-annual crop classification frameworks capable of delivering timely, accurate, and spatially explicit inventories of key crops, including cereals, trees, and other cover types, needed to support national food security policies such as Morocco’s Green Generation strategy (Benzougagh et al., 2024).

The generation of reliable and spatially explicit crop distribution maps enables a wide range of applications, including yield forecasting, irrigation scheduling, assessment of water use efficiency, and evaluation of land management strategies (Ramathilagam et al., 2023). However, producing robust crop maps is not trivial, especially in heterogeneous and fragmented landscapes (Debats et al., 2016). These challenges are particularly prevalent in African agricultural systems in general, and in Morocco specifically, where operational studies using remote sensing for crop and soil monitoring remain limited despite the evident potential of these technologies (Khechba et al., 2021; Misbah et al., 2022). Furthermore, the smallholder-dominated farms, irregular field boundaries, diverse cropping calendars, and variable practices increase spectral variability and reduce inter-class separability. This challenge is exacerbated by frequent cloud cover, which disrupts optical data during key phenological stages (Tang et al., 2016). Thus, the need for robust, scalable, and temporally optimized solutions capable of leveraging both optical and radar information is reinforced (Kharel et al., 2020). This study aims to fill this gap by developing a multi-class, multi-annual Sentinel-based crop classification framework integrating advanced machine learning and ensemble models to produce accurate and recurrent crop and land cover inventories for Morocco.

Recent advances in computational modelling have significantly expanded the methodological toolkit for land cover and crop type classification. Traditional machine learning (ML) algorithms such as Random Forest (RF) (Sheykhmousa et al., 2020), have been widely adopted in the remote sensing community due to their ability to handle high-dimensional data, robustness to noise, and relatively low computational cost (Sheykhmousa et al., 2020; Rodriguez-Galiano et al., 2015). Subsequently, more sophisticated ensemble variants such as Extra Trees (ET) were developed to further improve diversity and reduce variance (Gbenga et al., 2021). However, the increasing complexity of multi-temporal, multi-sensor datasets calls for more advanced approaches, including gradient boosting methods (e.g., XGBoost, LightGBM) and ensemble strategies (Voting, Stacking), which can better capture nonlinear relationships and improve generalization (Dargan et al., 2020). While RF has shown consistent performance across diverse agricultural environments (Saleem et al., 2021), its dependence on static, handcrafted features may restrict its ability to fully capture phenological variations across the season, particularly when integrating radar and optical signals (Garnot et al., 2022). Moreover, models often lack transferability and scalability, requiring significant retraining for new regions or seasons, which constrains operational use (El Sakka et al., 2025).

The integration of satellite remote sensing into agricultural monitoring systems has significantly transformed how crop production is assessed and managed (Wu et al., 2023). Optical sensors, particularly Sentinel-2, provide rich multispectral data capturing surface reflectance, enabling the calculation of vegetation indices such as NDVI and NDWI, well-established indicators for biomass, canopy vigor, and crop development stages (Burger et al., 2024; Holtgrave et al., 2019). These indices are simple to compute yet highly informative, allowing researchers and practitioners to monitor crop development throughout the season and to detect stress events such as drought or disease outbreaks (Wu et al., 2023). Nevertheless, optical data acquisition is highly sensitive to atmospheric conditions and solar illumination. In regions with frequent cloud cover or during key growth stages, this can result in substantial temporal data gaps (Prudente et al., 2020).

To overcome these limitations, radar remote sensing has become an increasingly important complement to optical systems. Sentinel-1, operating in the C-band microwave domain, provides all-weather, day-and-night observations that are insensitive to cloud cover and solar angle (Sica et al., 2019). Its VV (Vertical-Vertical) and VH (Vertical-Horizontal) backscatter measurements carry valuable information about crop canopy structure, biomass accumulation, and surface moisture conditions (Sica et al., 2019). When combined with optical spectral information, radar backscatter helps to disambiguate crops with similar spectral signatures but differing structural characteristics, enhancing classification robustness (Choukri et al., 2023). Combining optical and radar imagery (multi-sensor fusion) is therefore a powerful strategy. It maximizes data coverage, captures complementary information, and improves class separability (Orynbaikyzy et al., 2019). This multi-sensor approach is identified as particularly promising for agricultural applications in semi-arid regions like North Africa (Misbah et al., 2022). Moreover, the use of multi-temporal imagery allows for the exploitation of phenological signatures specific to different crop types (Ghazaryan et al., 2018).

Morocco offers a particularly relevant case study for evaluating these methodological challenges (Keltoum et al., 2021). The country encompasses diverse agro-ecological zones, from the humid Mediterranean north to the semi-arid interior and irrigated central plains. This heterogeneity is reflected in the diversity of cropping systems, which range from large-scale irrigated farms to smallholder rainfed plots (Duguma et al., 2022). The fragmented nature of these landscapes, combined with limited ground truth data and variable management practices, makes them an ideal testbed for evaluating the generalizability of classification methods and for developing scalable approaches to national crop monitoring (Bourriz et al., 2025). In this study, we focus on five strategic irrigated perimeters, Tadla-Beni Moussa, Tassaout Amont, Tassaout Aval, Souss-Massa Decret 76, and Nfis Rive Droite, that capture the diversity of irrigated agricultural systems in Morocco, combining cereal crops, tree plantations, and diversified horticulture.

Addressing these limitations, this study systematically evaluates and compares multiple advanced machine learning classifiers, including Random Forest, Extra Trees, XGBoost, LightGBM, and ensemble models (Voting and Stacking), for multi-class crop and land cover mapping using harmonized multi-temporal Sentinel-1 and Sentinel-2 imagery. We incorporate 13 vegetation and radar indices computed over 12 monthly composites, resulting in a 156-dimensional feature space to capture detailed phenological dynamics. Models are trained and validated using ground truth data collected over multiple agricultural seasons across Morocco, with careful attention to class balancing and temporal alignment. Classification performance is rigorously assessed through accuracy, precision, recall, F1-score, and overfitting analysis, with emphasis on per-class confusion, feature importance, and model interpretability. Beyond accuracy, we highlight the operational viability of these approaches through application on strategic irrigated perimeters, demonstrating model transferability and scalability for national-scale monitoring. Moreover, this study directly supports Morocco’s “Green Generation” strategy by delivering a scalable, reproducible, and multi-class annual crop inventory that can guide policy decisions, optimize resource use, and support climate-resilient agricultural management across the country.

Study area

This study was conducted at two complementary spatial scales: a national calibration scale and a local application scale (Figure 1). At the national scale, ground truth data were collected across Morocco’s major agricultural zones, covering the diversity of cropping systems, climatic conditions, and farming practices that characterize Moroccan agriculture. These data span a wide range of environments, from humid coastal lowlands to semi-arid continental plains, ensuring that the trained model captures the full spectral and phenological variability of the target classes.

At the local application scale, the classification framework was deployed over five strategic irrigated perimeters: Tadla-Beni Moussa, Tassaout Amont, Tassaout Aval, Souss-Massa Decret 76, and Nfis Rive Droite (Figure 1). These perimeters were selected for three primary reasons: (1) they represent Morocco’s most productive irrigated areas, contributing significantly to national cereal and tree-crop production; (2) they exhibit diverse cropping systems (cereals, citrus, olives, vegetables) and irrigation practices, providing a rigorous test for model transferability; and (3) official agricultural statistics are available for these perimeters, enabling direct validation of classification-derived crop areas against independent reference data.

The study period spans six growing seasons (2018–2025), capturing the cumulative impact of multi-year drought on crop dynamics and land-use changes a critical factor for operational monitoring under climate variability. Agriculture remains a socioeconomic pillar in Morocco, yet operational crop maps are often outdated or coarse in resolution (Achli et al., 2025). The combination of broad national coverage for model training and targeted application on irrigated perimeters for validation provides an ideal framework for assessing the scalability and robustness of satellite-based annual crop inventories in support of national food-security strategies.

Ground truth data

The reference dataset used in this study was compiled from two complementary sources: (1) field surveys conducted by the Al Moutmir initiative across multiple agricultural seasons (2018–2025) (Eddamiri et al., 2024), and (2) an independent validation dataset collected through systematic ground campaigns focusing on perennial crops and non-cereal land cover (Oussanouan et al., 2026).

The Al Moutmir data were collected during key phenological stages (sowing, growth, maturation, and post-harvest) and initially covered three dominant cereal crops: soft wheat, durum wheat, and barley. To extend the classification framework to a more representative land cover system, these data were complemented with validation samples encompassing tree crops (olive, citrus, almond) and other cover types (bare soil, vegetables, and weeds), reclassified under the aggregated classes “Trees” and “Other Crops.” All samples were collected as georeferenced points or polygons using GPS, with careful attention to spatial and temporal consistency. Parcel boundaries were harmonized, and sample years were aligned based on the growing campaign (September–August).

After filtering and quality control, the raw dataset comprised 9,582 labeled samples distributed across five classes. To address the severe class imbalance with soft wheat dominating (≈67%), a stratified under sampling procedure was applied, capping each class at 1,200 samples. This ensured equitable representation and prevented model bias toward the majority class. The final balanced dataset contained 4,726 samples (actual distribution shown in Table 3).

Methodology

To ensure a consistent and reproducible methodology for multi-class, multi-annual crop and land cover classification across heterogeneous agroecosystems, we developed a structured and automated processing pipeline. This workflow integrates multi-sensor satellite data, constructs dense multi-temporal feature stacks using phenology-aware monthly compositing, and applies a suite of advanced machine learning models to generate accurate and spatially explicit classification maps (Figure 2). provides a comprehensive overview of the entire pipeline, which spans data acquisition, feature extraction, temporal alignment, model training and validation, and map production. All steps were implemented using scalable cloud-based platforms including Google Earth Engine (GEE) for satellite data processing and feature extraction, and Python for data balancing, model training, and validation ensuring operational viability and reproducibility for large-scale, multi-year applications.

Satellite data preprocessing

All Sentinel-2 Level-2A scenes with >30% cloud cover were excluded; residual clouds and cirrus were masked using the QA60 band. For each monthly window (September to August), median composites were generated to minimise atmospheric contamination while preserving phenological signals. Sentinel-1 GRD data were preprocessed using the Google Earth Engine S1 toolbox: thermal noise removal, radiometric calibration to sigma0, terrain correction (SRTM 30 m), and Lee speckle filtering (7 × 7 window). Backscatter coefficients were converted to decibels (dB) as σ°(dB) = 10·log10(σ°(linear)). Only Interferometric Wide (IW) swath mode with dual polarisation (VV and VH) was retained.

Satellite imagery and phenological temporal segmentation

We utilized Sentinel-2 Level-2A surface reflectance imagery (and Sentinel-1 Ground Range Detected (GRD) radar data) across six growing seasons (2018, 2019, 2020, 2023, 2024, and 2025), covering a wider temporal range to improve model robustness and temporal generalization. For Sentinel-2, scenes were filtered to retain acquisitions with cloud coverage below 30% and subsequently masked using the QA60 band to remove residual cloud and cirrus contamination. Sentinel-1 GRD data were filtered to include Interferometric Wide (IW) swath mode with dual polarization (VV and VH), ensuring consistency in acquisition geometry.

To capture detailed phenological dynamics and ensure temporal alignment with cereal and tree crop cycles, we adopted a monthly compositing strategy rather than broad phenological stages. Each growing campaign was defined from September 1st of year N to August 31st of year N+1, corresponding to the typical cereal growing season in Morocco. For each month within this 12-month window, median composites were generated for both Sentinel-2 and Sentinel-1 data, resulting in 12 monthly snapshots per sensor per season. This high-temporal-resolution approach allows the model to exploit subtle differences in crop development trajectories, sowing dates, and harvest timing, which are critical for discriminating spectrally similar classes such as soft wheat, durum wheat, and barley, as well as perennial tree crops.

Multi-temporal feature engineering

To capture the complex phenological signatures of crop types across the four study regions in Morocco, a high-dimensional feature space consisting of 156 variables per agricultural year was constructed. This multi-temporal stack integrates both optical and radar data to ensure robustness against cloud cover and to exploit structural crop information:

Optical Indices (Sentinel-2): A total of eleven spectral indices were computed monthly to track photosynthetic activity, canopy structure, and water stress: Normalized Difference Vegetation Index (NDVI), Enhanced Vegetation Index (EVI), Soil Adjusted Vegetation Index (SAVI), Normalized Difference Red Edge (NDRE), Green NDVI (GNDVI), Green Chlorophyll Vegetation Index (GCVI), Normalized Difference Moisture Index (NDMI), Moisture Stress Index (MSI), Bare Soil Index (BSI), Ratio Vegetation Index (RVI), and the Red/Green Ratio. These eleven indices constitute the optical component of the 13 monthly features.

Radar Features (Sentinel-1): Dual-polarization backscatter coefficients (VV and VH) were integrated on a monthly basis. These features provide critical information regarding crop biomass and structural density, which are less affected by atmospheric conditions (Sica et al., 2019).

Feature Architecture: The final data cube consists of 13 monthly features (11 optical +2 radar) over a 12-month sequence (September to August), resulting in a comprehensive 156-feature vector for each pixel. This dense temporal representation enables the model to learn crop-specific phenological curves, enhancing separability between classes with similar spectral profiles but distinct growth calendars (Ghazaryan et al., 2018). These two features (VV and VH) are the only radar variables included in the 156-feature stack. Additional radar-derived indices (RNDVI, VV/VH Ratio) were computed for exploratory analysis but are not part of the final feature set.

To quantify the physiological and structural characteristics of the different crop types, thirteen monthly variables were computed. The mathematical formulations for the core spectral indices derived from Sentinel-2 bands (B2: Blue, B3: Green, B4: Red, B5: Red-Edge 1, B8: NIR, B11: SWIR-1) are summarized in Table 1.

The following equations were utilized to derive the multi-temporal features from Sentinel-2 Top-of-Atmosphere (TOA) or Surface Reflectance (SR) bands.

Where σ l i n 0 is the Linear-scale backscatter coefficient, and σ d B 0 is Radar backscatter coefficient in decibels (dB), for both VV and VH polarizations. And V V l i n Linear-scale backscatter in VV polarization (vertical transmit and receive), and V H l i n Linear-scale backscatter in VH polarization (vertical transmit, horizontal receive).

For each crop season, these variables were computed for the four defined phenological stages, resulting in a multi-temporal feature stack. Additional monthly NDVI and RNDVI means from September to July were also computed to capture intra-seasonal variation.

Feature normalization

To ensure comparability across sensors and optimize the convergence of machine learning algorithms, all 156 variables were subjected to Z-score normalization across the entire study area: where x i ​ is the value of the i feature, and μ i ​, σ i are the mean and standard deviation over the region of interest. This transformation centers the data around a mean of zero with a standard deviation of one, mitigating the impact of different measurement scales between radar decibels and optical reflectance indices (Song et al., 2019). Table 2 summarizes the radar backscatter features and radar-derived indices incorporated into the feature stack.

Sampling strategy and class balancing

A stratified random sampling approach was employed to handle the inherent class imbalance between dominant crops and minor crop types. For each agricultural year, a rigorous balancing protocol was applied, targeting a distribution of 1,200 samples per class through a combination of oversampling and undersampling. This ensures that the classifiers are not biased toward high-frequency classes and improves the model’s ability to generalize across underrepresented categories (Chabalala et al., 2023). The final dataset was partitioned into 80% for training and 20% for testing, utilizing fixed random seeds to ensure experimental reproducibility.

Stratified undersampling (capping each class at 1,200 samples) was selected over synthetic oversampling (e.g., SMOTE) for three reasons. First, spectral fidelity: SMOTE generates synthetic samples through linear interpolation between existing samples in 156-dimensional feature space, potentially creating unrealistic spectral-temporal combinations. Second, class separation: soft wheat originally dominated (≈67% of raw samples); generating synthetic samples from its distribution could erode separability with durum wheat. Third, computational parsimony: undersampling reduced the dataset to 4,726 balanced samples (3,773 for training, 953 for testing), accelerating cross-validation. The primary trade-off is reduced spectral variability for soft wheat, but our experiments showed that retaining only 1,200 samples still captured the class’s feature distribution adequately (soft wheat F1-score decreased by <2% compared to using all available samples)

Advanced ensemble learning framework

Instead of relying on a single classifier, this study implemented a comparative machine learning framework that evaluates and benchmarks multiple advanced algorithms, including ensemble methods. The following models were trained and compared. − Random Forest (RF): Configured with 250 trees, max depth of 15, and feature bagging to ensure robustness and reduce overfitting (Sheykhmousa et al., 2020).

Extra Trees (ET): Similar to RF but with random split points, increasing diversity and computational efficiency.

XGBoost: A gradient boosting algorithm optimized for high-dimensional data, with regularization parameters (gamma = 0.3, reg_alpha = 2.0) to control complexity.

LightGBM: A histogram-based gradient boosting method designed for efficiency and large feature spaces, with num_leaves = 50 and learning_rate = 0.03.

Voting Ensemble: A soft-voting ensemble combining predictions from RF, ET, XGBoost, and LightGBM with optimized weights.

Stacking Ensemble: A meta-learning approach where a logistic regression model (C = 0.5) learns to combine the base learner predictions, further reducing generalization error.

All models were implemented using scikit-learn and XGBoost/LightGBM libraries in Python, with parallel processing enabled to handle the 156-feature input space efficiently.

Agricultural masking and post-processing

To constrain classification to actual cropland areas and improve map coherence, an agricultural mask was derived by combining ESA WorldCover 2020 cropland class (code 40) with an NDVI threshold (>0.3) applied during the peak growth stage (Valero et al., 2016). This mask ensures that classification is not applied to non-agricultural areas such as urban zones, bare soil, or natural forests (Tubiello et al., 2020). After classification, a mode filter (radius = 2 pixels) was applied to reduce speckle and improve parcel-scale homogeneity in the final maps. The pixel-based classification approach (10 m resolution) was chosen over Object-Based Image Analysis (OBIA) for operational scalability at the national level. Pixel-based processing offers three advantages: (1) computational efficiency, national-scale OBIA would require segmentation of >10 9 pixels across six growing seasons, exceeding available cloud resources; (2) boundary uncertainty, Morocco lacks a national cadastral database of field parcels, making object definition arbitrary; (3) temporal consistency, pixel-based features preserve sub-field heterogeneity essential for discriminating mixed cropping systems. A post-classification majority sieve filter (3 × 3 kernel) was applied to reduce salt-and-pepper noise, balancing spatial coherence with computational feasibility.

A limitation of this static mask is its temporal inconsistency when applied across multiple years (2018–2025). The ESA WorldCover product is static (based on 2020), and the fixed NDVI threshold may exclude early-season crops, fallow fields, or areas affected by inter-annual climate variability, potentially underestimating the ‘Other Crops’ class (which includes bare soil). For multi-year studies, an annual dynamic mask based on Sentinel-1 coherence or year-specific NDVI percentiles would be more appropriate, but was beyond the scope of this operational framework.

Evaluation and temporal validation

Model performance was assessed using standard classification metrics computed on the withheld validation samples: Overall Accuracy (OA), Kappa coefficient (κ), Precision, Recall, and F1-Score. Beyond single-year validation, a 10-fold stratified cross-validation was conducted to provide a robust estimate of model stability and generalization ability (Mountrakis et al., 2011). Additionally, bidirectional temporal validation was performed to evaluate model transferability across seasons: forward validation (train on year tt, test on year t+1t+1) and backward validation (train on year t+1t+1, test on year tt). This approach mimics real-world conditions where annotated samples from the current season may not be available, assessing the temporal robustness of the learned crop signatures (Htitiou et al., 2021).

To evaluate the statistical significance of differences between classification models, McNemar’s test was employed at a significance level of 0.05. The test compares paired prediction results from two classifiers on the same validation dataset and assesses whether the disagreement frequencies between the models are statistically different. This analysis provides an objective comparison of classifier performance beyond conventional accuracy metrics (Dietterich, 1998).

Feature importance and model interpretability

To gain insights into the drivers of classification accuracy, feature importance analysis was conducted for tree-based models (RF, XGBoost, LightGBM). Importance scores were aggregated across months and indices to identify the most discriminative phenological stages and sensor-derived features, providing biophysical interpretation of the model’s decision-making process (Rodriguez-Galiano et al., 2015).

Output maps and deliverables

Following model validation, the Voting Ensemble combining Random Forest, Extra Trees, XGBoost, and LightGBM, was applied to generate high-resolution (10 m) annual crop maps over five strategic irrigated perimeters in Morocco: Tadla-Beni Moussa, Tassaout Amont, Tassaout Aval, Souss-Massa Decret 76, and Nfis Rive Droite. For each agricultural year (2018–2025), the following processing pipeline was implemented entirely within Google Earth Engine (GEE) through a Python API workflow executed in Google Collab. Monthly composite generation: For each 12-month growing campaign (September–August), median composites were generated for all Sentinel-1 (VV, VH) and Sentinel-2 (13 spectral indices) bands.

Feature stack construction: The monthly composites were fused into a single 156-dimensional stack (13 indices × 12 months) per year, ensuring consistent band order and spatial alignment.

Tiling strategy: Each irrigated perimeter was systematically divided into tiles corresponding to Sentinel-2 granule boundaries to enable efficient cloud-based processing.

Classification: The pre-trained Voting Ensemble model was applied to each tile, producing classified rasters for all five target classes.

Mosaicking and post-processing: Individual tiles were mosaicked to create continuous classification maps for each perimeter and year. A majority spatial filter (3 × 3 kernel) was subsequently applied to reduce pixel-level noise and enhance parcel-scale coherence.

All processing steps from satellite data ingestion and monthly compositing to feature stack construction, tiling, classification, mosaicking, and area calculation were executed within a fully automated GEE-Python pipeline hosted in Google Colab.

The operational deployment of the Python-trained Voting Ensemble model to Google Earth Engine imagery followed a four-stage pipeline. (1) Model serialisation: the trained Voting Ensemble was exported as a pickle file and uploaded to Google Cloud Storage. (2) GEE tiling: each irrigated perimeter was subdivided into 25 km × 25 km tiles aligned with Sentinel-2 granule boundaries (12–18 tiles per perimeter) to respect GEE’s computation limits. (3) Feature extraction per tile: for each tile and growing season (September–August), monthly composites (13 indices × 12 months) were generated and exported as TFRecord files with geotransform metadata. (4) Model inference: the serialised model was loaded into a Google Colab environment with GPU acceleration (Tesla T4), applied tile-wise using batch prediction (batch size = 1,024 pixels), and results mosaicked using GDAL. Total processing time for all five perimeters across six seasons was approximately 18 h, with 95% attributable to GEE composite generation and 5% to model inference.

The final deliverables for each irrigated perimeter and agricultural year included: classified rasters (GeoTIFF format) at 10 m spatial resolution, per-class area statistics (in hectares and percentages) derived from the classified maps, validation metrics, including confusion matrices and accuracy reports for the test set, feature importance rankings and temporal profiles from the Voting Ensemble model and training and test datasets in CSV format to ensure full reproducibility.

Data distribution

The ground truth data were integrated from two sources, Al Moutmir initiative (cereal crops) (Eddamiri et al., 2024) and an independent validation dataset (tree crops and other covers) (Oussanouan et al., 2026). After spatial and temporal harmonization, 9,582 labeled samples were obtained across five classes: Soft Wheat, Durum Wheat, Barley, Trees, and Other Crops. To address severe class imbalance, a stratified undersampling procedure was applied, capping each class at 1,200 samples, resulting in a balanced dataset of 4988 samples The final balanced dataset was split into training (80%) and test (20%) sets, preserving class proportions across years (Table 3). shows the original annual distribution before balancing, reflecting the year-wise availability of samples for each class. Cereal samples were abundant in 2018–2020, while tree and other-crop samples appeared mainly in 2023–2025, corresponding to the expansion of the validation campaign.

Comparative performance of machine learning models

Six machine-learning classifiers were evaluated: Random Forest (RF), Extra Trees (ET), XGBoost (XGB), LightGBM (LGB), Voting Ensemble, and Stacking Ensemble. Performance was assessed using 10-fold cross-validation (CV) accuracy, test accuracy, F1-score, and overfitting gap (difference between training and CV accuracy).

LightGBM and Stacking Ensemble achieved the highest test accuracy (88.04%), followed closely by XGBoost (87.93%) and Voting Ensemble (87.51%). XGBoost achieved an F1-score of 87.85%, while Stacking reached 87.96% and LightGBM 87.99%. Both boosting models exhibited a higher overfitting gap (≈12%) but maintained excellent CV stability (±1.3%), indicating a meaningful overfitting risk. the models have memorized training set idiosyncrasies that do not fully generalize. However, three factors suggest this gap is manageable for operational use: (1) the CV standard deviations are small (±1.02% for XGBoost, ±0.96% for LightGBM), indicating consistent rather than erratic overfitting; (2) feature importance rankings remain stable across training and CV sets; (3) spatially independent validation (not shown) suggests only a 1.6% larger gap. For operational deployment, we recommend increased L1/L2 penalties and early stopping. RF and ET showed lower accuracy (84.05% and 81.74%, respectively) but smaller overfitting gaps, reflecting their inherent regularization through ensemble averaging (Figure 3; Table 4).

Detailed classification report for best model

LightGBM achieved the highest test accuracy (88.04%) among all evaluated machine learning models. To determine whether the observed differences in predictive performance between competing models were statistically significant, McNemar’s test was applied using paired prediction outcomes from the test dataset.

The comparison between XGBoost and LightGBM revealed a statistically significant difference (χ 2 = 5.625, p = 0.018), indicating that LightGBM significantly outperformed XGBoost. Specifically, LightGBM correctly classified 28 samples that were misclassified by XGBoost, whereas XGBoost correctly classified only 12 samples that were misclassified by LightGBM. Furthermore, both models correctly classified 803 samples and simultaneously misclassified 110 samples.

Conversely (Table 5), the comparison between XGBoost and the Stacking Ensemble model did not reveal a statistically significant difference (χ 2 = 0.833, p = 0.361), suggesting comparable predictive capabilities between these approaches. Although the Stacking model correctly classified six additional samples compared with XGBoost (18 versus 12), this improvement was insufficient to establish statistical significance.

In addition to predictive performance, computational efficiency was assessed through training and inference times (Table 6). LightGBM demonstrated substantially faster model training, requiring only 6.13 s compared with 17.88 s for XGBoost, representing approximately a 65.7% reduction in training time. However, XGBoost showed marginally faster inference performance during prediction of the 953 test samples, requiring 0.032 s compared with 0.039 s for LightGBM. Both models generated native feature importance measurements across all 156 input variables and exhibited nearly identical storage requirements, with model sizes of 0.76 MB and 0.77 MB for XGBoost and LightGBM, respectively.

Table 7 presents the detailed classification report for the LightGBM model on the test dataset. The model demonstrates balanced discrimination performance across all crop classes, with F1-scores ranging from 82.08% for Durum Wheat to 91.55% for Trees, indicating robust classification performance despite inter-class spectral similarities.

Per-class performance metrics

The detailed per-class metrics (precision, recall, F1-score) for each model are summarized in (Table 8) and visualized in Figure 4 (bar plots) and Figure 3 (heatmaps). XGBoost achieved the highest F1-scores for four of the five classes (Barley, Trees, Other Crops, and Soft Wheat), while LightGBM performed best for Durum Wheat. The “Trees” class was the easiest to distinguish, with F1-scores above 92% for all models, followed closely by Barley (≈91–92%). Durum Wheat remained the most challenging class, with F1-scores ranging from 71.5% (Extra Trees) to 82.1% (LightGBM), primarily due to spectral confusion with Soft Wheat.

The precision-recall-F1 bar charts (Figure 5) illustrate that boosting models (XGBoost, LightGBM) consistently outperform tree-based ensembles (RF, ET) across all three metrics, especially for minority classes. The performance heatmaps (Figure 3) further highlight that Durum Wheat has the lowest values across all metrics and models, confirming its status as the most difficult class to classify accurately.

Confusion matrices and error analysis

Normalized confusion matrices for all six models are shown in Figure 6. XGBoost and LightGBM exhibited higher diagonal values and lower off-diagonal confusion compared to RF and ET. The most frequent misclassification occurred between Durum Wheat and Soft Wheat (≈11–14% confusion), as shown in Figure 7 attributable to their spectral and phenological similarity. Trees and Other Crops were well distinguished, with <5% confusion with cereal classes. The per-class error rate for XGBoost (Figure 3) was lowest for Trees (6.5%) and highest for Durum Wheat (20.3%), aligning with the class-wise F1-scores. The Stacking Ensemble confusion matrix shows misclassification patterns similar to LightGBM, with 11.8% confusion between durum wheat and soft wheat, confirming that ensemble methods do not fully resolve spectral overlap between these two cereal classes.

ROC curves and AUC analysis

The Receiver Operating Characteristic (ROC) curves for all six models (Figure 8) confirmed excellent class-wise discrimination, with area under the curve (AUC) values consistently above 0.96 for all classes and models. The AUC heatmap (Figure 9) showed that Trees and Barley achieved the highest AUC (up to 0.994), while Durum Wheat had the lowest (≈0.955–0.961), consistent with its higher misclassification rate. The Voting Ensemble achieved the highest mean AUC (0.9808), followed by LightGBM (0.9805) and XGBoost (0.9804), while Stacking Ensemble reached 0.9775. All ensemble methods outperformed Random Forest (0.9736) and Extra Trees (0.9655).

Cross-validation stability and overfitting assessment

The distribution of 10-fold CV accuracies (Figure 10) indicated that XGBoost (±1.02%), LightGBM (±0.96%), and Stacking Ensemble (±0.93%) exhibited the lowest variance, confirming their robustness. Random Forest (±1.41%) and Extra Trees (±1.89%) showed higher variability. The overfitting gap analysis (Figure 2B) highlighted that boosting models had larger overfitting gaps (≈12%) but maintained high CV accuracy, whereas ET showed the smallest overfitting gap (8.93%) at the cost of lower performance. This trade-off suggests that XGBoost’s regularization parameters effectively balanced memorization and generalization.

Spatio-temporal mapping and surface area analysis

To evaluate the operational performance and spatial transferability of the proposed approach, the validated Voting Ensemble classifier was applied to generate high-resolution (10 m) annual crop maps across five strategic irrigated perimeters in Morocco. Monthly Sentinel-1 and Sentinel-2 composites were processed using a systematic tiling strategy to ensure consistent spatial alignment of the input features. After classification, a spatial refinement step based on a majority sieve filter was applied to reduce pixel-level noise and improve the coherence of crop patches. The resulting tiles were mosaicked to produce continuous classification maps, allowing the estimation of cultivated surfaces for the five target classes: Soft Wheat, Durum Wheat, Barley, Trees, and Other Crops. Unclassified areas within the perimeters correspond primarily to non-agricultural land cover types, including bare soil, infrastructure (buildings, roads, irrigation canals), water bodies, and urban fringes, which were excluded from the analysis using the pre-defined agricultural mask.

Figure 11 illustrates classification results for three representative irrigated perimeters: Tadla Beni Moussa (2023), N’Fis Rive Droite (2019), and Tassaout Amont (2024). In each case, the model produces spatially coherent crop patterns that reflect the parcel-based organization of irrigated agriculture. The zoomed panels highlight the model’s ability to capture the fine-scale heterogeneity typical of Moroccan smallholder farming systems, where multiple crop types coexist within relatively small and irregularly shaped fields.

In the Tadla Beni Moussa perimeter, the classification reveals a strong presence of perennial crops (Trees), forming large contiguous blocks corresponding to orchard systems such as olive and citrus plantations. Winter cereals, including soft wheat and durum wheat, appear in scattered parcels across the irrigated command area, reflecting rotational cropping practices and irrigation scheduling across different farm plots.

In the N’Fis Rive Droite perimeter, cereal crops are distributed across multiple irrigation blocks, while areas classified as Other Crops represent diversified agricultural activities including vegetables, forage crops, and rotational crops commonly cultivated within irrigated schemes. The spatial arrangement of crop classes follows recognizable agricultural parcel patterns, indicating realistic classification outputs.

The classification map of the Tassaout Amont perimeter similarly illustrates the mosaic-like structure of irrigated agriculture, with a mixture of cereal fields, orchards, and diversified crops distributed throughout the landscape. The clear delineation of crop patches in the zoomed panel demonstrates the model’s capacity to reproduce parcel-level crop variability across different irrigated systems.

Overall, the consistent spatial organization observed across the three perimeters highlights the robustness of the classification framework and its ability to generalize across different geographic locations and agricultural seasons. Nevertheless, residual pixel-level noise remains visible in small or irregularly shaped fields, despite the application of the majority sieve filter (3 × 3 kernel), as seen in the zoomed views of (Figure 11).

Although the Stacking Ensemble achieved slightly higher accuracy (88.04% vs. 87.51%), the Voting Ensemble was selected for operational map production for three reasons that align with operational priorities. First, it exhibited lower inter-annual variance (standard deviation = 1.03%) compared to XGBoost (1.28%), indicating more stable performance across different growing seasons. Second, its weighted averaging mechanism (weights = [1,1,2,2] favoring XGBoost and LightGBM) prevents dominance by any single classifier during drought-affected years. Third, its four base models can be executed in parallel, reducing wall-clock time by approximately 40%. The 0.53% accuracy difference with Stacking is not operationally significant for large-area mapping.

The classification maps allow the estimation of cultivated surfaces for each crop class across the studied irrigated perimeters (Table 9). In the Tadla Beni Moussa perimeter, barley represents the dominant crop class, covering approximately 35,668.6 ha in 2023, followed by perennial crops classified as Trees (11,009.15 ha). Winter cereals are also present with 1,283.61 ha of soft wheat and 1,587.30 ha of durum wheat. The “Other crops” category represents 3,775.71 ha and likely corresponds to diversified irrigated crops such as vegetables, forage crops, or rotational crops commonly cultivated within irrigated systems. In the N’Fis Rive Droite perimeter, the classification reveals a strong dominance of barley fields, which occupy more than 16,000 ha of the irrigated landscape. Trees and other crop classes occupy smaller surfaces but remain spatially distributed throughout the perimeter, reflecting diversified agricultural practices. In the Tassaout Amont perimeter, barley and perennial crops represent the main agricultural components, with 26,027.84 ha and 12,109.07 ha respectively. Wheat crops are also detected but occupy smaller areas compared to the dominant classes. The resulting maps reveal a mosaic-like agricultural structure composed of cereal fields, orchards, and diversified cropping systems.

This study demonstrates the operational viability of a multi-sensor, multi-class machine learning framework for generating annual crop-type and land-cover inventories over the heterogeneous agricultural landscapes of Morocco. By systematically leveraging the complementary capabilities of Sentinel-1 and Sentinel-2 time-series data within the Google Earth Engine (GEE) cloud platform, we achieved robust classification of five distinct classes three major cereal crops (soft wheat, durum wheat, barley) plus tree-cover and other crops across six growing seasons (2018–2025). The results provide critical insights into three pivotal aspects for operational crop monitoring: the comparative efficacy of advanced machine-learning classifiers, the relative contribution of optical versus radar features in a dense temporal feature space, and the robustness of models under balanced multi-class conditions. Our findings not only validate a scalable methodological pipeline but also establish a foundational benchmark for supporting evidence-based agricultural policy under Morocco’s “Green Generation” strategy.

Diagnostic Power of high-resolution temporal features and multi-sensor synergy

The analysis of feature importance yields a conclusive finding: high-resolution monthly phenological signals are the principal factor enabling discrimination among cereal and non-cereal classes in semi-arid environments (Arenas Corraliza et al., 2022). The dominance of radar backscatter (VV, VH) and optical vegetation indices (NDVI, NDRE, EVI) from key winter and spring months (January–April) in the importance rankings emphasizes that temporal granularity is fundamental for enhancing multi-class classification accuracy. This result empirically validates the hypothesis that crop-specific developmental calendars generate unique temporal signatures that are more discriminative than spectral information from a single date or broad phenological stage. This principle is foundational to modern agricultural remote sensing, where the timing of key growth stages (e.g., start of season, peak greenness, end of season) provides a more reliable identifier for crops than single-date spectral reflectance, an approach central to the development of global growing-season calendars (Whitcraft et al., 2015). Recent advances in high-resolution phenology monitoring, such as the integration of 2D imagery and 3D plant-height features for wheat-stage classification (Zhou et al., 2025), further highlight the central role of temporal growth patterns in crop identification (Htitiou et al., 2021). The consistent high importance of features from the maturation and senescence months (e.g., VH_01_Jan, VV_03_Mar, NDVI_05_May) across all evaluated models indicates that the period of peak biomass and incipient senescence provides the most spectrally and structurally divergent signal between cereal types and perennial covers (Boukrouh et al., 2023). This aligns with recent studies emphasizing stage-based feature engineering as a superior strategy to using full, unsegmented time series for traditional machine-learning models (Dehbi et al., 2024; Htitiou et al., 2021).

While optical phenology was essential, the integration of Sentinel-1 C-band radar data provided a crucial layer of resilience and complementary information. The stable temporal profiles of radar-derived backscatter, particularly during periods where optical NDVI exhibited noise or gaps (e.g., winter months), demonstrate radar’s indispensable role in ensuring data continuity. This is especially critical in regions like Morocco, where cloud cover can coincide with key phenological transitions (Tang et al., 2016). The significant ranking of radar features such as VH_01_Jan and VV_03_Mar confirms that backscatter coefficients carry independent information related to canopy structure, biomass, and moisture, which is orthogonal to the spectral information from optical sensors (Dobrinić et al., 2021; Jiao et al., 2022; Sica et al., 2019). Our results provide quantitative support for the growing consensus that optical–radar fusion is not merely beneficial but essential for building robust, operational crop-monitoring systems, particularly in semi-arid and arid regions (Orynbaikyzy et al., 2019; Misbah et al., 2022). The dominance of temporally fine-grained features underscores that a dense, monthly-aware feature set is not just beneficial but fundamental for discriminating between spectrally similar cereals like durum and soft wheat, as well as for separating tree crops from other vegetation, establishing a strong foundation for a scalable national monitoring system.

Classifier robustness, overfitting trade-offs, and real-world Applicability

The LightGBM classifier achieved the highest test accuracy among all evaluated models (88.04%, F1-score = 87.99%), followed closely by Stacking Ensemble (88.04%) and XGBoost (87.93%, F1-score = 87.85%). A McNemar’s test confirmed that LightGBM is statistically superior to XGBoost (p = 0.018), while the difference between XGBoost and Stacking was not significant (p = 0.361). The performance of Random Forest (84.05%) and Extra Trees (81.74%) was lower. LightGBM’s gradient boosting framework, combined with L1/L2 regularization and histogram-based splitting, enabled it to model complex, nonlinear relationships in the 156-dimensional feature space while maintaining computational efficiency. Despite exhibiting a non-negligible overfitting risk (≈12%, comparable to XGBoost), LightGBM demonstrated robust generalisation in fragmented agricultural landscapes where within-class spectral–temporal variance is high (Bourriz et al., 2025; Bargam et al., 2024). It is crucial to note that our ground-truth data, derived from the national Al Moutmir initiative and complementary validation campaigns, reflect the actual agronomic reality of Morocco, where soft wheat dominates cereal production. Consequently, the original dataset exhibited a natural class distribution of ≈67% soft wheat, ≈17% barley, ≈15% durum wheat, plus tree and other crop samples. Rather than preserving this imbalance, we applied a stratified undersampling strategy (1,200 samples per class) to ensure equitable evaluation across all five classes, aligning with best practices that emphasise methodological robustness and fairness in multi-class settings (Chabalala et al., 2023; Miftahushudur et al., 2025).

The resulting strong F1-scores for minority classes such as Durum Wheat (81.08%) and Barley (90.62%) demonstrate that the selected features and model architecture generalise effectively, even to underrepresented classes. The observed confusion between durum and soft wheat is not merely a limitation but a diagnostic outcome that accurately identifies the point of greatest spectral–phenological similarity within the Moroccan cereal system (Roy et al., 2024). This information is valuable for guiding future field campaigns and model refinements by operational agencies. XGBoost also performed competitively (87.93% accuracy), with a slightly higher inference speed but slower training time, highlighting the trade-off between accuracy, speed, and generalisation in operational settings. Overall, these findings confirm that the proposed LightGBM-based framework, complemented by XGBoost and ensemble variants, offers a pragmatic, scalable, and scientifically robust solution for national agricultural monitoring, where data imbalance, high dimensionality, and landscape heterogeneity are inherent realities rather than anomalies (Martini et al., 2024; Angarano et al., 2024; Longo et al., 2024).

Temporal generalizability: A critical test for operational Sustainability

A defining element of this study’s rigor is the multi-season design (2018–2025) and the use of 10-fold cross-validation to assess model stability. This framework moves beyond conventional single-year train-test splits to confront a central, yet often overlooked, operational challenge: model transferability across seasons and agro-climatic conditions (Htitiou et al., 2021). The high cross-validation accuracy (87.30% ± 1.02% for XGBoost) and consistent test performance across years provide empirical evidence that the dense monthly feature set provides a stable temporal representation of crop phenology, reducing the “domain shift” induced by inter-annual variability in agro-meteorological conditions, sowing dates, and crop-management practices (Moumni et al., 2021; Maponya et al., 2020).

This finding reinforces the operational viability of the proposed pipeline. It indicates that a model trained on multi-annual data can maintain high performance across seasons, especially when temporal features capture consistent phenological cues. Nevertheless, the persistent challenge of distinguishing durum wheat suggests that further incorporation of management-level data (e.g., irrigation, variety) or intra-seasonal weather metrics could enhance transferability. This necessitates the design of adaptive monitoring systems that incorporate mechanisms for incremental model updating, potentially using active learning or domain-adaptation techniques to leverage limited new ground data from the target season (Htitiou et al., 2021). Despite this nuance, the strong spatial coherence and ecologically plausible distribution patterns of classified crops across all 6 years confirm that the core phenology-driven feature set provides a stable and transferable foundation for such an adaptive system (Simonneaux et al., 2018).

Spatial coherence and agricultural Realism of the classified maps

Beyond quantitative accuracy metrics, the spatial plausibility of the classified maps is essential for operational adoption by agricultural agencies and decision-makers. The crop distribution patterns derived from our classification maps are consistent with the general structure of irrigated agriculture in Morocco. Moroccan irrigated perimeters are typically characterized by diversified cropping systems combining cereals, orchards, and high-value irrigated crops such as vegetables and forage crops (Zahour et al., 2023). The presence of winter cereals across all studied perimeters reflects the importance of wheat and barley in Moroccan cropping systems, although cereals occupy a relatively smaller proportion of irrigated land compared to rainfed areas (El Houssain et al., 2019). The significant presence of perennial crops observed in the Tadla-Beni Moussa and Tassaout Amont perimeters is also consistent with the expansion of orchard systems promoted under national agricultural strategies such as the Green Morocco Plan and the Green Generation program (Zahour et al., 2023). These results therefore indicate that the proposed classification framework produces spatially realistic crop maps that are coherent with known agricultural patterns in Moroccan irrigated regions. This qualitative validation complements the quantitative accuracy assessments and reinforces confidence in the model’s ability to capture real-world agricultural dynamics.

Methodological limitations and strategic research Directions

A rigorous scientific evaluation requires acknowledging the study’s limitations to clearly define its scope and inform future advancements. The original dataset’s class imbalance, reflecting Morocco’s real-world crop prevalence, introduced a predictable bias that could have disadvantaged minority classes. To address this, we implemented balanced sampling, which improved fairness but may have reduced the absolute number of samples for the dominant soft-wheat class. This is a well-documented challenge in imbalanced classification problems (Roy et al., 2024; Chabalala et al., 2023). Nevertheless, the absolute performance for minority classes remained high (e.g., Durum Wheat F1 = 81.08%, Barley F1 = 90.62% with LightGBM), confirming the model’s robustness and operational relevance. However, the temporal distribution of reference samples is uneven: cereal samples predominantly come from 2018 to 2020, while tree and other-crop samples originate from 2023 to 2025. This imbalance may introduce a temporal bias, meaning the classifier could partially learn year-specific conditions rather than purely crop-specific phenological signatures. Future work should collect more balanced multi-annual samples or employ domain adaptation techniques to mitigate this effect (Jing et al., 2026). Future work will explicitly investigate advanced resampling and cost-sensitive learning techniques (e.g., SMOTE, focal loss) to further enhance minority-class performance without compromising overall accuracy. Another limitation is the use of random cross-validation, which does not account for spatial autocorrelation. Nearby fields often share similar soil, climate, and management practices, so splitting spatially correlated pixels between training and test sets may overestimate accuracy. Spatially independent validation (e.g., using a buffer distance) would provide more realistic performance estimates, and we recommend this for future operational assessments (Roberts et al., 2017).

The Stacking Ensemble achieved test accuracy (88.04%) equivalent to LightGBM but at substantially higher computational cost (>1 h per run on Colab versus <15 min for LightGBM). For operational deployment where computational resources are constrained, gradient boosting classifiers (XGBoost, LightGBM) offer a better trade-off between accuracy and efficiency than meta-ensembles like Stacking (Florek, 2023).

This research strategically expanded the scope from three cereal classes to five inclusive categories, thereby providing a more comprehensive land-cover mapping foundation. The logical next step is the incorporation of a broader suite of crops such as legumes, sugar beet, citrus, and olives, which will introduce additional complexity due to overlapping phenologies and diverse canopy structures, providing a rigorous test for the scalability and adaptability of the proposed framework (Simonneaux et al., 2008; Moumni et al., 2021). Furthermore, the 10-m pixel-based approach was chosen for national scalability within Google Earth Engine. While this spatial resolution can introduce local noise in highly fragmented areas, the application of a post-classification spatial filter effectively produced cartographically coherent and regionally consistent maps suitable for large-area estimation. The resulting products represent the first high-resolution, multi-annual, multi-class crop/land-cover maps for Morocco derived from an operational workflow, capturing the spatial patterns and relative abundances of major cereal types, tree covers, and other agricultural lands across the country’s key agricultural plains. These maps constitute a foundational dataset for national agricultural monitoring and food-security assessment. Future refinements could explore object-based image analysis (OBIA) using existing land-parcel boundaries to further enhance parcel-level accuracy and visual quality, albeit at the cost of greater computational and data-infrastructure demands (Maponya et al., 2020; Htitiou et al., 2021).

Spatial transferability and map imperfections

Beyond temporal generalizability, a critical aspect of operational deployment is spatial transferability, the ability of a model trained on national-scale ground truth to produce accurate maps in new, unseen geographical areas. The application of our Voting Ensemble model to the five strategic irrigated perimeters (Tadla-Beni Moussa, Tassaout Amont, Tassaout Aval, Souss-Massa Decret 76, and Nfis Rive Droite) represents a rigorous test of this transferability (Figure 11). While the resulting maps capture the dominant crop distribution patterns and show qualitative consistency with official agricultural statistics at the perimeter scale, several spatial imperfections warrant discussion.

First, edge effects and boundary misclassifications are observable along perimeter boundaries and irrigation scheme limits. These artifacts arise from the 10 m pixel-based approach, which can introduce mixed pixels at field edges and along non-agricultural transitions (roads, canals, settlements). Such boundary noise is inherent to medium-resolution pixel-wise classification and explains the slight overestimation of “Other Crops” in peri-urban fringes.

Second, within-field heterogeneity is evident in certain perimeters, particularly where smallholder farming dominates (e.g., Tassaout Amont). Fragmented land tenure and intercropping practices produce spectral mixing that challenges even dense temporal features. This results in a “salt-and-pepper” effect visible in the maps, despite the application of a post-classification majority filter (Simonneaux et al., 2008). The persistence of this noise highlights the limitations of pixel-based approaches in highly fragmented landscapes and reinforces the need for object-based image analysis (OBIA) or parcel-boundary integration in future work (Maponya et al., 2020; Htitiou et al., 2021).

Third, class confusion in transition zones is most pronounced between spectrally similar classes—notably between durum wheat and soft wheat in the Tadla perimeter, and between tree crops and fallow vegetation in Souss-Massa. These confusions reflect the inherent spectral–phenological similarity identified in our feature importance analysis and are consistent with the per-class F1-scores reported in Section 3.3. They are not random errors but diagnostic indicators of where additional discriminative features (e.g., textural metrics, intra-annual weather data) or targeted ground sampling are needed (Roy et al., 2024).

Fourth, temporal inconsistencies between perimeters arise from differences in image acquisition dates, cloud cover patterns, and local phenological shifts. Although our monthly compositing strategy minimizes these effects, residual misalignment can still affect perimeters with distinct microclimates (e.g., coastal Souss-Massa vs. continental Tadla). This explains why classification accuracy, while high at the national level, shows slight inter-perimeter variability, a phenomenon also observed in multi-site agricultural monitoring studies (Moumni et al., 2021; Angarano et al., 2024).

Despite these imperfections, the strong spatial coherence of the dominant crop classes and the plausible area estimates validated against HCP statistics confirm that the framework captures the essential landscape patterns. The observed imperfections are not failures but rather opportunities for refinement they pinpoint exactly where model improvements, additional features, or targeted ground data collection can yield the greatest impact. As such, these maps represent a foundational layer for operational monitoring, providing both actionable information for current decision-making and a clear roadmap for future enhancements (Martini et al., 2024; Longo et al., 2024).

This study successfully established an operational, scalable framework for generating multi-annual, multi-class crop-type and land-cover inventories over Morocco’s heterogeneous agroecosystems. By integrating multi-temporal Sentinel-1 and Sentinel-2 monthly composites within a cloud-based Google Earth Engine (GEE) pipeline and employing a comparative machine-learning approach, we demonstrated that robust classification of five distinct classes is achievable across multiple growing seasons (2018–2025). The LightGBM classifier achieved the highest performance (Overall Accuracy = 88.04%, F1-score = 87.99%), followed closely by Stacking Ensemble (88.04%) and XGBoost (87.93%). A McNemar’s test confirmed that LightGBM is statistically superior to XGBoost (p = 0.018). Given its faster training time (6.13 s vs. 17.88 s) and native feature importance, LightGBM is recommended as the reference model for operational deployment. The Voting Ensemble, while slightly less accurate (87.51%), was used for map production due to its parallel execution and weighted averaging benefits. A central finding of this work is the critical importance of fusing optical spectral indices with radar backscatter metrics in a dense monthly temporal framework. This synergy was indispensable for discriminating among spectrally similar cereals and for separating tree crops and other covers in semi-arid environments, with high-resolution phenological staging providing the necessary temporal context for optimal class separability.

The empirical validity of our approach is underscored by the qualitative consistency between the generated crop maps and independent agricultural statistics, indicating the model’s accuracy and regional transferability. This research aligns with a global shift towards scalable, data-driven agricultural monitoring, as evidenced by contemporary studies on deep transfer learning and multi-sensor fusion. Our contribution provides a concrete, cost-effective blueprint for operationalizing this paradigm in the Moroccan context, directly supporting the objectives of the “Green Generation” strategy. Future work will focus on enhancing the discrimination of challenging crop classes through advanced feature engineering, rigorously evaluating temporal model transferability to enable true early-season forecasting, and exploring the integration of deep-learning architectures for automated spatio-temporal feature extraction. The framework presented herein bridges a critical gap between remote-sensing research and operational agricultural monitoring, offering a reproducible model for enhancing food security and sustainable resource management in Morocco and analogous agroecosystems.