Assessing the Impact of Spatial Resolution and Aggregation Method on Sentinel-2 NDVI Time Series in Grasslands of Mainland Spain

High-resolution satellite imagery has substantially improved the monitoring of vegetation dynamics; however, the influence of spatial resolution and pixel aggregation on NDVI time series consistency remains insufficiently quantified, particularly across multiple native resolutions within a single sensor platform. This study evaluates how Sentinel-2 spatial resolutions (10 m, 20 m, and 60 m) and two pixel aggregation methods (pure-pixel and centroid) affect NDVI time series in 14,031 grassland plots across mainland Spain over the period 2018–2023. High-quality NDVI time series were selected using the Interpolation Efficiency Indicator (IEI), and discrepancies relative to a 10 m pure-pixel baseline were quantified through the Time Series Angle Distance (TSAD) and Root Mean Square Error (RMSE). A sensitivity check confirmed that the radiometric differences between Band 8 (10 m) and Band 8A (20/60 m) introduce negligible bias compared with genuine spatial-resolution effects. Formal non-parametric statistical testing—omnibus Kruskal–Wallis with epsilon-squared (ε2) effect sizes and pairwise Cliff’s Delta comparisons—was applied to assess the magnitude and practical significance of the observed differences across plot area categories and Köppen climate groups (B, Cs, Cf). Results show that coarser resolutions (60 m) substantially reduce NDVI reliability, excluding more than half of the plots under the pure-pixel criterion and smoothing temporal variability, whereas 10 m and 20 m resolutions preserve most spectral and temporal information. The 20 m resolution introduces moderate but non-severe phenological distortion (median TSAD ≈ 0.05 rad, RMSE ≈ 0.026) with a 78% reduction in data volume and 72% reduction in processing time. The choice between pure-pixel and centroid sampling has negligible impact at 10–20 m but becomes relevant at 60 m, where pure-pixel selection reduces errors from spectral mixing at the cost of severe sample attrition. Parcel area strongly conditions the error metrics, with large effect sizes (ε2=0.273) in the smallest plots, while Köppen climate classification decisively shapes TSAD (up to ε2=0.447), indicating that spatial degradation distorts phenological patterns differently across climate classes. These findings support a multi-scale monitoring strategy: 10 m for fragmented, heterogeneous grasslands (<3 ha), 20 m as a computationally efficient alternative for homogeneous areas (>10 ha), and outline potential implications for policy frameworks such as the Common Agricultural Policy (CAP).