## 1 Introduction  
Although shrublands, savannas and grasslands account for 37% of the world's terrestrial area, many studies have been conducted to estimate how forests contribute to the global carbon cycle [1, 2, 3, 4, 5, 6, 7, 8], but not many monitor the implications of grasslands/shrublands on the carbon cycle at a large scale [9, 10]. Monitoring these non-forested ecosystems in dry areas is important since they are very sensitive to changes in water availability [11] and because they have been suggested as potential candidate areas for major carbon storage efforts [12]. Some of these areas were recorded as carbon sinks, while during years of drought the same areas were carbon sources [13], so long-term measurements at the regional scale are essential for examining the seasonal and inter-annual variability of C fluxes [14]. The net ecosystem exchange (NEE) is the CO2 exchange between the atmosphere and the whole terrestrial ecosystem [15] and it can be decomposed into two terms: gross primary production (GPP) and total ecosystem respiration (Re) [16]. While NEE is directly measurable using Eddy covariance (EC) flux towers [17], it is not possible to obtain direct, integrated observations of either GPP or Re, because these processes represent a multitude of responses by a combination of autotrophic and heterotrophic organisms [15]. However, the GPP estimates from EC towers (GPPec) only represent fluxes at the tower footprint scale (between one hundred metres and several kilometres) [18], and therefore, the estimation of GPP at regional scale can only be made by using ecosystem models (e.g., [19, 20] and/or remotely sensed data (e.g., [2, 21, 22, 23]). Although in the past satellite sensors such as the Advanced Very High Resolution Radiometer (AVHRR) aboard the NOAA satellite and the VEGETATION (VGT) sensor aboard the SPOT satellite have been used to estimate carbon fluxes [24, 25, 26], the most frequently used remote sensing data source to monitor the C cycle at regional or global scale is the Moderate Resolution Imaging Spectroradiometer (MODIS) (e.g., [18, 21, 22, 27]). MODIS estimates of GPP (GPPm) are derived from the MOD17 algorithm, which uses input data from three sources: the eight class 1 km MODIS land cover classification [27], the global scale meteorology obtained from NASA's Data Assimilation Office GEOS-4 global climate model [28], and the fAPAR (fraction of Absorbed Photosynthetically Active Radiation). The utility and applicability of GPPm is linked to its validation, which is challenging because of scaling issues and also the logistical constraints of measuring NPP in the field [16]. Site specific studies such as [2, 29, 30, 31] as well as the two most comprehensive validation projects (the Bigfoot [6] and the AmeriFlux [8, 13]) focused on the inter-annual variation of GPPm in forested biomes in non-water deficient areas [13, 16, 18]. The lack of research focused on GPPm in grasslands was partially filled by [32], which analysed the agreement with GPPec on an eight-day interval on five grasslands in three dry and two continental areas. They found that GPPm agreed well with tower data for three sites (r \(>\) 0.80), but did not capture the seasonal dynamics for the two other sites, and they recommended to extend this work to other grassland ecoregions and ecosystems. Therefore, there is a deficit of studies where GPPm is validated at an annual and at an intra-annual scale on non-forested biomes in water limited areas, since the previous similar works [32] only considered the performance of GPP in one biome (grasslands) and at one temporal scale (weekly). There is also a need for operational methods that can increase the accuracy of the GPPm estimations in these biomes, such as identifying the variables which can reduce the uncertainty in the prediction models. The aim of this work is to cover that gap.