Gridded soil surface nitrogen surplus on grazing and agricultural land: Impact of land use maps

Excess N application on agricultural land greatly impacts the environment in multiple ways, driven by population growth and improving quality of human diets. Therefore, it is essential to quantify the sources of the emissions of N compounds and their determinants (e.g. biological N fixation (BNF), mineral fertilizer, manure N and N deposition) to develop adequate mitigation measures. Here we aim at comprehensively mapping and quantifying N fluxes on agricultural land to analyze these sources on different scales. As underlying grazing land maps used for such calculations are fairly different in terms of methodology and definition and thus spatial extent and pattern, we investigate how this diversity in grazing land maps affects quantification of N indicators. We compared three different global grazing land maps and analyzed the propagation of differences to discrepancies in N indicators calculated from them. We discovered that (i) area differences propagated to high discrepancies in N surplus mostly in Asia, and to a minor extent also in Europe and Northern Africa. (ii) BNF constitutes an important translator for differences on grazing land to N indicators, while also being a source of further uncertainty, which warrants further scrutiny. (iii) A more inclusive definition of grazing land results in overall less N surplus given the larger areas included but allows to provide a more comprehensive estimate of the influence of human activity on the N cycle. This study is the first to provide an in-depth analysis of the effect of grazing land and agricultural land area differences on various N budget terms and N indicator calculation, highlighting opportunities for further research, and the importance of a comprehensive accounting of N surplus when using an inclusive definition of grazing land.


Introduction
Nitrogen (N) compounds as plant nutrients play a key role in food production and the increase in agricultural yields. Many of the current agricultural practices are associated with an excess application of N to agricultural land (Vitousek et al 1997, Galloway et al 2013. In order to increase agricultural output and to replenish any N lost from soils, under industrial settings, mineral N-fertilizers are used abundantly, and livestock manure is released merely as an add-on. As a consequence, human activity results in a perturbation of the global N cycle, with impacts on the local and global environment through eutrophication, acidification and through the formation of N 2 O, a powerful greenhouse gas (Fowler et al 2013, Reis et al 2016. With global N 2 O emissions from agriculture on the rise and the N cycle reported to having overshot a 'safe operating space' within which the functioning and resilience of the Earth system can be sustained (Steffen et al 2015), action to reduce N inputs is pressing (IPCC 2020).
For the development of policies to manage and reduce the globally increased N input it is key to quantify and monitor major N flows (Meisinger et al 2008, Winiwarter andExpert Panel on Nitrogen Budgets 2016). To take into account not only the global but also local effects of excess N application, it is essential to know the spatial distribution of major N flows to and from agricultural land. This can be achieved by calculating N indicators such as N surplus and Nitrogen Use Efficiency (NUE) using global land-use maps. Such maps have been developed for cropland (Liu et al 2010, West et al 2014, Bouwman et al 2017 being complimented with maps provided by several process-based models that also include calculations on pastures and rangeland (Tian et al 2018). Covering also grazing land in the N 2 O calculation is essential as this is the origin of around half of the total agricultural N 2 O emissions (Dangal et al 2019, IPCC 2020.
Previous studies demonstrate that the underlying land-use maps for such calculations show distinct variations. For cropland, Kaltenegger and Winiwarter (2020) show that the variations in crop maps also influence the N indicators calculated from them. Including grazing land in spatially explicit N surplus calculations (from which N 2 O emissions can be estimated), is particularly challenging due to the lack of robust data and the high variability in existing maps (Erb et al 2007, 2016, Fetzel et al 2017, Phelps andKaplan 2017. Maps on the extent of global grazing lands suffer from differences in classification, methodologies, errors in georeferencing or the use of different satellite sensors and the particular difficulties to separate grazing land from other land uses; note that grazing can occur in a wide range of ecosystems and is not bound to the occurrence of grassland ecosystems such as steppes or artificial grasslands (Milchunas andLauenroth 1993, Robinson et al 2011). Thus, mapping of grazing land by means of remote sensing is challenging due to the weak link between the land use (grazing) and land cover (artificial grasslands, natural grasslands with and without trees, semideserts, shrublands, etc) (Kuemmerle et al 2013, Erb et al 2016. In consequence, key differences in grazing land maps result from the use of different grazing area definitions. The most commonly used maps are based on or are adjusted to fit the freely available, global statistics by the Food and Agricultural Organization (FAO), which only include information on permanent meadows and pastures. These grazing lands are defined as land used permanently, i.e. five years or more, for herbaceous forage crops, either cultivated or growing wild. This definition excludes natural grazing lands that are used for non-permanent land uses such as occasional or sporadic grazing, land uses typical for e.g. transhumance. The omission of such areas can result in an area difference as large as 10% of terrestrial ice-free surface (Erb et al 2016, Arneth et al 2019. Furthermore, besides this omission, the quality of the national FAO data is highly variable between countries and inconsistencies become apparent when combined with other land-use information , Stadler et al 2018. Such large differences in land use maps have mostly been neglected in the calculation of global spatially explicit N surplus or N 2 O emissions, with studies assessing N 2 O emissions from grazing land often using the same land-use map (Dangal et al 2019, Tian et al 2018). Even the global N 2 O Model Intercomparison Project (NMIP: Tian et al 2018) that compared ten process-based models and analyzed differences to establish a global N 2 O budget consistently used cropland and grassland extent from just one source, the Historical Database of the global Environment -HYDE (Klein-Goldewijk et al 2017). Understanding the uncertainty, that origins not only from differences in data, but also from differences in definitions, however, is key for deriving robust interpretations, especially for developing or assessing N management policies (Oenema et al 2003).
Considering the methodological and semantic discrepancies in the derivation of grazing land area datasets as described above, combined with the effect crop map discrepancies have on N surplus calculation, we here aim to explore (i) the influence of methodological discrepancies in grazing land maps on N surplus calculations and (ii) the effect of a more inclusive grazing land definition.

Methodology
In order to assess how discrepancies in grazing land maps translate to differences in patterns and magnitude of various N indicators, we calculated N surplus, a measure for the amount of N remaining on the soil surface, and Nitrogen Use Efficiency (NUE), a measure for the efficiency of N uptake by the crop, using three different landuse maps. Following the summary table of book-keeping for soil surface N budgets as described by Oenema et al (2003), we include mineral fertilizer, manure N, biological N fixation and N deposition as N inputs and N in crop harvest as well as N grazed by livestock as N outputs for our N indicator calculations. We exclude mineralization, the decomposition to inorganic material or crop residues. These stocks are assumed to stay within the system. N surplus was calculated by subtracting N outputs from N inputs and dividing it by the respective area, and NUE was calculated by dividing N output by N input.
N budget calculations on cropland were added to our grazing land calculations using the M3 crop map  so that a possible effect of discrepancies on grazing land to discrepancies on agricultural land generally (cropland and grazing land) could be detected. M3 was chosen because it was identified as being the most suitable product to expand cropland calculations to grazing land (Kaltenegger and Winiwarter 2020).
The base year for all calculations is 2010.

Grazing land
Three different grazing land maps were chosen to make differences and their implication for the final calculations visible: (1) Ramankutty et al (2008), (2) HYDE (Klein-Goldewijk et al 2017) and (3) an updated map based on the approach outlined in Erb et al (2007). Ramankutty et al (2008) and HYDE use the FAO definition for pastureland ('Permanent Pastures') but use different approaches to downscale national information from the FAO's Statistical Database (FAOSTAT) to the spatially explicit grid. HYDE data for this category is split into rangeland (less accessible land) and pastureland, whereas Ramankutty et al (2008) do not differentiate grazing land categories, but correct national pasture areas on basis of plausibility checks. The more comprehensive approach by Erb et al (2007) aims at including permanent and non-permanent grazing land, integrated in a 'closed-budget' land-use dataset (comprising cropland, infrastructure area, forests, wilderness and non-productive area and grazing land and covering 100% of the area in each grid cell). A similar approach was followed in the construction of the Global Agroecological Zone (GAEZ) mapping by IIASA and FAO (Fischer et al 2008). More information on the differences between the maps and an analysis of their discrepancies can be found in the publications by Erb et al (2016) and Fetzel et al (2017). All maps are available at the spatial resolution of 5 arc minutes.
HYDE grassland (pastureland & rangeland) is calculated using FAOSTAT data combined with a spatially explicit map of present land cover by the European Space Agency (ESA) based on Medium Resolution Imaging Spectrometer (MERIS) satellite data for the year 2010 (Klein-Goldewijk et al 2017). Each land cover type from the satellite dataset is assigned a probability of cropland and pastureland to occur in a grid cell. Through four allocation steps, FAOSTAT cropland and pastureland data is distributed according to the map's probability of a certain land use type occurring in a grid cell, thereby producing a map with spatially distributed land use data which is representative for the year 2010. Protected areas (from the UN Environment Programme World Conservation Monitoring Center -UNEP/WCMC) and Australian 'non-used areas' (from the National Land and Water Resources Audit -NLWRA) were excluded as well as inaccessible tundra areas in Northern America and Russia. Maximum 90% of total area of a grid cell could be assigned to agricultural land due to the assumption that at least 10% of the area (in a 5' grid cell) would be used for small infrastructure or would be unsuitable for cultivation. Ramankutty et al (2008) combined two satellite datasets to create the spatial distribution pattern for cropland and pastureland in the year 2000. One dataset was from Boston University (BU-MODIS) where 17 different land cover classes were identified from October 2000 to October 2001. This dataset was developed by combining SPOT (satellite type) vegetation data, which is monitoring global vegetation daily or in a ten-day period, and expert assessments using a flexible classification scheme. The other dataset was the Global Land Cover 2000 (GLC2000), where 22 different land cover classes were identified with the help of regional experts and regional classifications between November 1999 and December 2000. The combination of these two maps rendered the best results because the maps counteracted each other's weaknesses. However, to include all possible land use categories (new ones emerged due to the overlap of two different classifications of BU-MODIS and GLC2000) the final map contained 391 land cover types. Multiple linear regression was used to calculate the parameters for the cropland and pastureland distribution. Bootstrap technique was used to estimate the uncertainty for the distribution parameters using inventory data. In a last step, the derived land use data was aggregated to compare it to inventory (FAOSTAT plus census data on subnational levels) cropland and pastureland data for approximately 15,990 different administrative units (an agricultural census usually only takes place every 5-10 years, in which time also administrative units can change) for the year 2000. From this comparison, a correction factor was calculated for the satellite data, which was restricted to a predefined range, trusting satellite data more than inventory data. As the data provided by Ramankutty et al (2008) was only available for the year 2000, for the purpose of this paper we adjusted it to fit FAOSTAT country data for the year 2010 (see supplementary information for more details (available online at stacks.iop.org/ERC/3/055003/ mmedia)). The grassland map provided by Ramankutty et al (2008), will hereafter be referred to as 'Ramankutty'.
In contrast to the two approaches described above, where the extent of permanent meadows and pastures was derived from satellite data and adjusted to FAOSTAT statistics, the approach by Erb et al (2007) maps grazing land which is defined as 'potential grazing area'. Here, this basic approach was followed but updated information and newly available maps for cropland, forests, wilderness and non-productive areas were used. Owing to the lack of data for non-permanent grazing land, the basic approach of mapping the entirety of grazing land is based on subtracting built-up area, cropland area, unproductive areas, wild areas and forest area from total land area for each grid cell. In the updated version used here, non-productive land, built-up land, cropland area and permanent pastures and rangelands are taken from HYDE (version 3.2.1) (Klein-Goldewijk et al 2017). The wilderness mask is generated by combining human footprint data for 1993 and 2009 (Venter et al 2016a(Venter et al , 2016b with data on intact forest landscapes, available for 2000 and 2013 (Potapov et al 2017). Core wilderness areas were defined as areas with no human artefacts and, within the forest zone as defined by Potapov, located in an intact forest landscape. Peripheral wilderness within the forest zone only satisfies one of those criteria. Forest area is derived from the CCI land cover dataset of the European Space Agency (ESA) and distinguishes between closed forests, with tree cover greater than 40%, and open forests, with tree cover between 15 and 40% (ESA 2017).
Total grazing land thus subsumes permanent pastures and rangelands from HYDE 3.2.1, a share of 25% of peripheral wilderness, 50% of open forests, and, following the 'closed-budget' approach, the remainder of land area per grid cell that has not been assigned to any of the other land uses above (other land, may be grazed) (Erb et al 2007). This definition is in line with the definition of the grassland category by the IPCC (Bickel and Köhl 2006), with the exception that it does not include unused (wild) grasslands.
This approach leads to an inclusive grazing land definition, encompassing, for example, browsed land such as shrubs, trees and succulents consumed by livestock (Phelps and Kaplan 2017). It is important to note that all uncertainties in underlying maps accumulate in this map. However, comparisons with remote sensing data and census statistics reveal a good fit (Erb et al 2007(Erb et al , 2016. This grazing land map, based on the approach suggested by Erb et al (2007), with extensions as described, will hereafter be referred to as 'Erb'.
Given the differences in definition between these land-use maps, the term 'grazing land' as it is used throughout this paper includes both, permanent pastures as defined by FAO, as well as areas used for occasional grazing like browsing or transhumance, while the term 'Grassland' only refers to 'permanent pastures'.

Fertilizer
Synthetic fertilizer application on cropland per crop type and country category was taken from Heffer (2013), complimented with additional information on grass crop fertilization from Heffer et al (2017), and distributed on M3 harvested areas per crop type and country (see supplementary information S3 for details). Because synthetic fertilizer application within an IFA country category such as 'ROW -Rest Of World' which includes Latin American countries as well as African countries is only given as total application per crop type, no further differentiation between application rates within such a country category was possible. As Latin American and Asian countries can have an up to ten-fold higher application rate compared to African countries, data was later adjusted to fit FAOSTAT annual country totals (FAOSTAT 2019, The World Bank Group 2020).
Factors provided by Lassaletta et al (2014) were used to calculate synthetic fertilizer application to grazing land. These factors provide the fraction of total synthetic fertilizer applied per country that is diverted to such. This information was only needed for calculations referring specifically to grazing land. Calculations performed for total agricultural land do not require this separation as we assume a different fertilizer distribution affects only the respective grid cell, and cropland and grazing land are added together in such a grid cell for the calculation of N surplus and N use efficiency (NUE) on agricultural land.

Manure N
Spatially explicit data on livestock numbers per livestock system were taken from the Gridded Livestock of the World (GLW) (Gilbert et al 2018) and combined with N excretion rates per animal type and country from the Greenhouse Gas -Air Pollution Interactions and Synergies (GAINS) model (Amann et al 2011). Fractions of manure left on grazing land as well as managed and lost during storage per livestock system were taken from Herrero et al (2013b). Country-specific application rates to cropland were taken from Liu et al (2010) to calculate the fraction of manure N applied to cropland as well as the fraction of manure N recycled to grazing land. All calculations were made on the grid cell level.
As manure N spatial distribution was determined by GLW and did not necessarily match the spatial distribution of the cropland and grazing land combinations used here, all data was filtered to only include cells where manure N and cropland or grazing land area were collocated (see supplementary information for further information).

Biological N fixation
BNF on cropland was calculated using crop specific factors such as harvest index, amount of N content from BNF (Ndfa) and below ground N provided by Herridge et al (2008) for crops and BNF rates for leguminous and non-leguminous grass crops by Smil (1999). Grass crops in this context refer to grass crops found in the M3 crop maps that we interpreted as temporary pastures included in cropland (Kaltenegger and Winiwarter 2020) as defined by Monfreda et al (2008).
BNF on grazing land was calculated using two different methods. The first method made use of spatially explicit data on global net primary production (NPP)the energy provided by a plant minus the energy for its own metabolism ( The basis for the saturation function between NPP and BNF was taken from Cleveland et al (1999) who hypothesized that NPP may be a proxy for carbon potentially available to N fixers.
The second method to calculate BNF was with the use of data on evapotranspiration (ET)the sum of evaporation from soil and plant and transpiration through plant canopy (Jørgensen and Fath 2008). This relationship also follows an estimate given by Cleveland et al (1999), who hypothesized a correlation between ET and BNF. The corresponding equation was taken from Meyerholt et al (2016). Data on ET was taken from two sources, the Inter-Sectoral Impact Model Intercomparison Project (ISIMIP) (https://esg.pik-potsdam.de/ search/isimip/) and the GLEAM model (https://www.gleam.eu/). .31) (IPCC 2019). We rearranged the equation used to calculate N excretion to obtain the amount of N intake per animal: ) K annual N excretion rates of animal of species/category T, kg N animal −1 yr −1 N intake(T) K annual N intake per head of animal of species/category T, kg N animal −1 yr −1 N retention(T) K fraction of daily N intake that is retained by animal of species/category T N retention(T) was taken from IPCC (2019 ) Table 10.20, 10A.1-10A.4 and combined with shares of livestock in a specific growth stage of animal type in herd from the same tables.

N indicator calculation
The N inputs included for each calculation depended on whether grazing land or agricultural area was used as a basis. For N surplus and NUE calculations on grazing land, synthetic fertilizer spread on grazing land only, manure left and spread on grazing land and BNF and N deposition on grazing land were considered as N inputs while N grazed was considered as being the only N output. When N indicators were calculated for total agricultural land, total synthetic fertilizer application, manure N application on cropland and grazing land (with volatilization losses and losses due to storage being subtracted) as well as manure N left on grazing land and BNF and N deposition on agricultural land were included. N contained in crop harvest as well as N grazed by livestock were N outputs.

Results
In order to trace back N indicator differences to grazing land area discrepancies we first take a closer look at allocation and area variations between the three different grazing land sources before arriving at a detailed analysis of the discovered grazing land discrepancies on the N indicator calculation.

Grazing land discrepancies
Further to summing up grazing land areas per region for all three sources (figure 1), we calculated its relative standard deviation in each grid cell to show the agreement or disagreement between the land use maps in more detail (figure 2). While differences between HYDE and Ramankutty only exceed 10% in Western Asia and North Africa, area size and allocation of grazing land in parts differ quite strongly between all three sources in most regions. Discrepancies in global grazing land area between Erb and each of the other sources amount to 1,061 km 2 , an area slightly bigger than the territory of Canada. For instance, Erb assign more than tenfold more grazing land area to India compared to Ramankutty and HYDE, highlighting the potential scale of the conceptual differences between permanent and other forms of grazing. However, there are regions with little differences between all three sources such as the Caribbean, Australia & Oceania, Central America and Eastern Asia. In general, given their more inclusive definition of grazing land, Erb allocate larger areas to most regions. Exceptions are Western Asia and Northern Africa. As seen in figure 2, allocation differences between Erb and the other two sources become visible for Northern areas and Southern and Southeastern Asia which is also reflected in area differences. However, allocation and area discrepancies can be detached. In Australia, differences in allocation become visible, however these do not affect the regional results for total grazing land area to a greater extent, as the areas assigned to Australia only differ in allocation and not size.

N indicator discrepancies on grazing land
Differences in N surplus, calculated from each of the grazing land maps, are shown as as total N surplus per region in figure 3 and relative standard deviation in figure 4. Grid cells containing less than 5% grazing land area were excluded from being displayed in the map to avoid extreme values for N input produced by a difference in allocation of grazing land and N inputs. All N indicators compared in this section use NPP as a basis for BNF calculations as discrepancies between the NPP and ET approach remain small. This does not to a greater extent affect the results, nor the conclusion of this study as is explained in more detail in the supplementary information.
Most differences between N indicators calculated from the different grazing land maps are smaller than the area differences between the sources. While grazing land area for almost all regions is largest in Erb, N surplus does not necessarily reflect this pattern. The biggest discrepancy when calculating N surplus based on the three grazing land maps, remain in Central and Eastern Europe and China (figures 3 and 4). Regardless, regions where grazing land area discrepancies are high, like Southern and Southeastern Asia, also show higher discrepancies in the N indicators calculated from them as can be seen in figure 3. These discrepancies are driven by allocation differences between the three maps as can be observed in figure 4. N surplus difference is high for Western Industrial Europe with N surplus based on Erb being over 20 kg ha −1 yr −1 lower than N surplus based on Ramankutty or HYDE which is rather similar. The same is true for Eastern and South Eastern Europe as well as Southeastern Asia.  For a more detailed analysis per region, we took a closer look at how differences in area and allocation between the three datasets on grazing land area affect differences in N input and N output. Additionally, we investigated the drivers for the discrepancies encountered by inspecting the countries dominating the regional results and the budget terms constituting the biggest share of N inputs. As can be seen in table 1, differences between the budget terms often exceed differences observed for the N indicators ( figure 3). This can be explained by N inputs and N outputs showing discrepancies of similar direction and size which then balance each other as they are put in relation for N indicator calculation.
Looking at Southern Asia, showing the second highest relative standard deviation, India and Pakistan are mainly responsible for the differences in N surplus. These are regions where the Erb approach results in much larger grazing land areas than the approaches by Ramankutty and HYDE. In these regions, the more encompassing grazing land map results in more N input, especially BNF and N deposition, in the N surplus and NUE calculations. Although also higher manure N and synthetic fertilizer input caused by allocation discrepancies can be found in the manure N results based on Erb, NUE is not to a greater extent affected by this, as input discrepancies are balanced by grazing discrepancies due to the large share of manure N in total N inputs (more than 70%).
A similar phenomenon can be found in Western and Eastern Asia as well as North Africa, where low NUE discrepancy and high N surplus discrepancy can be found. In both cases highest N inputs of the respective source are not matched with largest areas thereby heightening the differences in N surplus. In Western as well as Eastern Asia, BNF and N deposition contribute about half to total N input to an equal share. Although Erb allocates less area than both other sources to Western Asia, they allocate more than twice as much area than Ramankutty and HYDE to Turkey, where highest BNF rates in Western Asia can be found. These area discrepancies lead to about three times higher BNF and N deposition values when calculated from the grazing land map from Erb compared to the other sources. In North Africa, where manure N contributes 59%-64% of N inputs, higher values of this budget term are reached when calculated based on Erb leading to total N input based on Erb being second highest while grazing land area taken from this source is smallest.
In Eastern Asia, China is the biggest influence on the regional results, being responsible for 87%-90% of total N inputs to this region. However, it is the area difference in Mongolia that leads to HYDE reporting the largest grazing land area for Eastern Asia, while derived N inputs are only the second highest (highest N inputs have been assigned to the grazing land map derived from Erb). N inputs derived from each grazing land map remain similar because BNF rates are generally low in this region and hence area differences do not result in bigger BNF differences. Table 1. Relative Standard Deviation between area, N inputs (manure N (N man ), synthetic fertilizer N (N syn ), BNF based on NPP (N BNF_NPP,, )BNF based on ET from ISIMIP (N BNF_ISIMIP ), BNF based on ET from GLEAM (N BNF_GLEAM ) and N deposition (N dep )), N output (N grazed ) and N indicators (NUE and N surplus) calculated based on Ramankutty, HYDE and Erb. Values larger than 10% are printed in bold.

Relative Standard Deviation [%]
The opposite is true for Central Africa, South America, Eastern Africa, North America, Southern Africa and Central Asia and Russian Federation, where allocation differences lead to NUE discrepancies. However, discrepancies are not visible when looking at N surplus as area discrepancies lower the differences between the N budget terms due to assigning relatively larger areas to large N inputs and outputs. For Central Africa, BNF is an important input constituting between 53% and 57% to total N input. As Erb allocates around five times as much area to the Democratic Republic of Congo, this leads to five times more BNF in this country, which is the largest contributor in this region. BNF also plays a central role in South America, where larger grazing land areas from Erb in Brazil lead to almost twice as much BNF input in this country as they coincide with high BNF rates due to high NPP rates. All other N inputs are rather similar and as BNF constitutes more than a third of total N inputs, it is a driver for discrepancies. In Eastern Africa, largest area discrepancies can be found in Ethiopia, South Sudan and Kenya whereas in North America large area discrepancies can be found in the US, which contributes over 85% of N inputs in that region. These discrepancies lead to large discrepancies in BNF as well, which contributes more than one third to total N inputs in these regions and shows the highest differences between the different sources. The same is true for Southern Africa and Central Asia and Russian Federation, where BNF is the largest N input. Influential countries in these regions are South Africa and Botswana and Russia respectively.
In Western Industrial and Eastern and South Eastern Europe, allocation differences are not externalized in manure N and synthetic fertilizer differences, as hot spots for these N inputs seem to be covered by all land use maps. However, area and allocation discrepancies causing differences in BNF and N deposition lead to NUE and N surplus discrepancies. Large differences in BNF in Western Europe can be found in Norway. Using the updated Erb dataset results in allocating more than 100 times as much BNF to this country group than the assessment based on Ramankutty or HYDE. This leads to BNF contributing over 50% to the total N inputs when calculated based on Erb as opposed to as little as 10% when calculated using HYDE data. However, Norway is a country where BNF shows high discrepancies between the different calculation modes. When NPP is used for the calculation, BNF is up to twice as high as when calculated using ET due to NPP being relatively high in the South of Norway while ET is relatively low for Norway as a whole. BNF differences in Eastern and South Eastern Europe are driven by Ukraine and Poland. Differences between different BNF calculation modes are rather low.
Due to similar methodologies, N indicators calculated based on Ramankutty and HYDE are more similar than when comparing them to N indicators based on Erb (table 2). Biggest differences between Ramankutty and HYDE based results can be found in all regions of Asia (except Central Asia). On the global level, discrepancies between N surplus results based on Erb and based on Ramankutty or HYDE are up to seven times larger than discrepancies between Ramankutty and HYDE alone, while for NUE results, this number increases to over 28.

N indicator discrepancies on agricultural land
Some of the discovered differences become smaller when expanding to agricultural land by including cropland (table 3). This can be observed when comparing results for Western and Southeastern Asia. For both regions, synthetic fertilizer application can explain this effect. Due to only 18% (HYDE) to 36% (Ramankutty) of agricultural area being cropland area in Western Asia, compared to 89% (HYDE) to 90% (Ramankutty) in Southeastern Asia, the difference between grazing land and agricultural land NUEs in Western Asia is smaller than in Southeastern Asia. In Southeastern Asia, the difference of NUE on grazing land as described before almost disappears due to synthetic fertilizer becoming the primary N input for HYDE and Ramankutty based calculations. However, in Southern Asia NUE discrepancies between HYDE and Ramankutty based calculations arise when including cropland area. This can be explained by the additional manure applied from livestock kept indoors balancing out the beforehand observed discrepancies in manure N. However, as this is only true for the combination of indoor and outdoor excreted manure, discrepancies remain for N grazed leading to differences in NUE.
N surplus and NUE discrepancies between calculations based on Erb versus Ramankutty also become smaller when cropland area is added. Discrepancies in synthetic fertilizer disappear completely and manure N discrepancy only remains in Central Africa, where Ramankutty allocates neither enough cropland nor grazing land area to cover all cells where manure N is located. This can on the one hand be explained by more area being covered by agricultural area, smoothing out allocation discrepancies and by additional manure N and most synthetic fertilizer N being located on cropland. Because cropland N is equal for Erb and Ramankutty based calculations it acts as a moderator to N surplus and NUE discrepancies. Higher discrepancies remain in Northern America, where N input differences decreased more than area differences, leading to a lower NUE based on Erb grazing land. In Southeastern Asia, where manure N and synthetic fertilizer account for more than 70% of N input the decrease of discrepancies between them while area differences remain leads to similar NUE but differing N surplus.

Discussion
Calculating spatially explicit N indicators from land-use maps is common practice, however little attention has been given to the differences in these land-use maps and their effect on the resulting NUE or N surplus. Our results suggest that conceptual differences originating from different grazing land definitions are highly influential on the N indicator calculation, especially on grazing land. Comparing the global N surplus ratio between Erb and each of the other two sources (1.12 & 1.14) to the ratio between HYDE and Ramankutty (0.98), hints towards the impact of conceptual grazing land differences being about six times stronger than impacts of other uncertainties.
We also identified BNF as playing a crucial role in the translation of the discrepancies on grazing land, with uncertainties related to the calculation of this N input being high in some regions.

Differences in grazing and agricultural land
Discrepancies rooted in the methodology used to develop the respective land-use map (observable in the comparison Ramankutty versus HYDE) offer a way to identify areas where data uncertainties prevail. In contrast, discrepancies rooted in the definition of grazing land (Ramankutty & HYDE versus Erb) can be interpreted as an opportunity to evaluate human perturbation of the N cycle in a more encompassing way. Areas, that were not considered as being affected by human activity before, are enclosed in this approach. This could be of importance when looking at countries, where livestock may be kept as an asset or a draught animal and might not graze on permanent pastures but rather on smaller areas with variable land use (Teufel et al 2010, Herrero et al 2013a. Such activity is excluded from the approaches used by HYDE and Ramankutty who refer only to a subset of grazing activities. While, ensuing, the approaches by HYDE and Ramankutty are biased to the lower end, the approach by Erb et al (2007) is probably biased to the higher end. It can be assumed that area derived by the subtractive approach (see method section) results in the inclusion of other land uses that prevail on the areas identified as grazing land by the Erb approach. This indicates the importance and need for further research on grazing land in order to enhance the mapping and assessment capacities aimed at evaluating the effect of human activity in the agricultural sector on the N cycle. While discrepancies in N indicators on grazing land are quite striking, differences on agricultural land are not as prominent. Nevertheless, using a more inclusive definition of grazing land (Erb) leads to a larger account of N inputs to agricultural lands by approximately 29 Tg yr −1 as compared to results based on the maps following the FAO 'permanent grass land' definition. Combined with N outputs, which show a difference of only 11 Tg yr −1 between the exclusive (Ramankutty, HYDE) and inclusive (Erb) grazing land definition, this leads to an increased N balance on agricultural lands of approximately 17 Tg yr −1 . This excess N found in the calculations based on Erb equals about 6% of the total global N input and over 12% of the global N balance. Still, global N surplus based on Erb is lowest for agricultural land and rather low for grazing land due to larger areas being considered in the calculation. Regional and sub-regional differences remain where in some cases N surplus based on Erb is higher and in other cases lower compared to the other two sources.
However, not only the effect of differing area definitions on N surplus is lowered when including cropland but also the overall N surplus calculation. This effect is visualized in figure 5 where we subtracted a 'best estimate' for cropland only (Kaltenegger and Winiwarter 2020), which is based on very similar approaches, from our calculations for agricultural land based on Erb (figure 6). We note that for large parts in the world, only a part of N surplus in cropland area also remains when looking at total agriculture. This indicates the transfer of N material via livestock (extracted from grazing land) to cropland, where it causes excess surplus. Although still being rather high, biggest reduction in N surplus when considering its distribution over agricultural area rather than cropland only, can be observed in China.
For some N poor areas, a similar effect is visible the other way around: areas of significant N deficiency that appear in cropland only, e.g. in Western Africa, become invisible when total agricultural land is considered. This may point to land use rotation (not visible from the land use maps) causing a more or less balanced situation, at least for N, on the longer timescale and the grid size used (not excluding sub-grid and shorter-term effects).
Nevertheless, our results still show high N surplus in China, India, Central Europe and Northeastern North America (figure 6). Still, even in those regions, distribution of N towards grazing land creates a much smaller surplus than otherwise expected (see Kaltenegger and Winiwarter 2020). Moreover, while N surplus is very low in many developing countries, a consistent long-term depletion of nitrogen is not visiblewhich also would be difficult to maintain over extended periods.

Biological N fixation (BNF)
Additional to uncertainties related to grazing land differences, uncertainties related to BNF calculations were identified to influence N indicator calculations.
Area differences in grazing land affect BNF most strongly, leading to discrepancies in N input as BNF contributes more than 50% to total N inputs in some. Its influence on N indicator calculations on grazing land highlights the importance to investigate BNF and uncertainties related to it more closely.
We calculated global BNF on grazing land based on NPP and ET to range from 31 Tg yr −1 using ISIMIP ET data on the Ramankutty grazing land map to 53 Tg yr −1 using GLEAM ET data on grazing land taken from Erb. Due to the correlation of NPP and ET, both being linked to water availability and solar radiation, differences between NUE calculated using BNF based on ET or NPP remain small. Highest discrepancies between BNF calculation based on NPP and ET can be found in Central Asia and Russian Federation (9% difference between NPP and ISIMIP ET). In this region more NPP can be found than ET leading to higher N input and consequently a lower NUE. This is also true for the North Africa region.
While only small discrepancies between the different modes for calculation in this paper were discovered, differences to other sources are slightly bigger. Herridge et al (2008) used FAOSTAT data on yield and area combined with values on BNF taken from literature to calculate global BNF rates per agricultural system. They estimated BNF from pasture and fodder legumes to lie between 12 Tg yr −1 and 25 Tg yr −1 and BNF from savannas to amount to less than 14 Tg yr −1 . Adding these two estimates amounts to a range of 26-39 Tg N yr −1 being fixed biologically, being rather at the low end of our global BNF results on grazing land. However, this can most likely be explained by Herridge et al (2008) only including legumes and no other grazing areas. Meyerholt et al (2016) used the land surface model 'O-CN' to calculate BNF spatially explicit and per plant functional type (PFT). They estimated about 40 Tg N yr −1 globally being fixed by C3 and C4 grasses, which is well within our calculated range. The reason for the estimate by Meyerholt et al (2016) being lower than the estimate based on the inclusive grazing land definition, following Erb, could be that vegetation types such as Mediterranean and arid shrubland' as well as wet savannas and moist tundras are not included in the C3 and C4 PFTs but could, at least partially, be considered as grazing land by Erb. Cleveland et al (1999) derived BNF rates per ecosystem from a literature review and arrived at a global BNF estimate between 100 Tg yr −1 and 290 Tg yr −1 . Selecting all ecosystems where grazing land could be found, such as savannas, tundra, shrubland etc, amounts to an estimate between 52 Tg yr −1 and 193 Tg yr −1 which fits the results we provide here based on the grazing land definition used by Erb. Further research to produce more robust results will be needed to undermine the hypothesis that BNF plays a crucial role in the translation of grazing land discrepancies to differences in N indicators. This is especially true for the calculation method used. Although a link between NPP and BNF as described by Cleveland et al (1999) was also found by Ashworth et al (2018), recent work suggests limitations in the use of NPP and ET as proxies for BNF calculation (Davies-Barnard and Friedlingstein 2020). This is also true when aiming to include effects on BNF from grazing land management such as legume intercropping (Rumpel et al 2015), which comes with the challenge of species distinction on grazing land where no global dataset is available. Additionally, Fetzel et al (2017) find NPP to be responsible for >30% of the variance found in mapping grazing intensity at the global level (ranging from about 10% up to 50% at the world-regional level). This variance is caused by the huge range of available NPP estimates due to different modelling approaches. From this one could assume a strong impact from NPP uncertainty on BNF and subsequently N indicator calculation.
Further research aimed at gathering more information on global species distinction on grazing land and narrowing the uncertainty related to NPP is certainly warranted. However, as this affects all grazing land maps to the same extent, this remains outside the scope of our work where we aim to scrutinize the discrepancies resulting from the use of different maps.

Comparison of N indicators to other studies
While the largest differences in N indicators and their compounds result from differences in grazing land definitions, a comparison between the results for permanent pastures from this study (HYDE & Ramankutty) with the same reference presented by Dangal, reveals that a discrepancy of similar magnitude relates to the role of manure N (49.02 Tg yr −1 Ramankutty, 52.46 Tg yr −1 HYDE and 74.20 Tg yr −1 Dangal). This difference can most likely be explained by differing animal numbers and excretion rates. While Dangal et al (2019) use manure N data from Zhang et al (2017), which combines a previous version of GLW with IPCC 2006 guidelines default excretion values, we use regional and animal specific excretion factors as provided by the GAINS model.
To check for further discrepancies, our global and regional results on NUE and N surplus from agricultural land were compared to Bouwman et al (2009), while results for Europe were additionally compared to data presented by Leip et al (2011). Bouwman et al (2009) calculated global N budgets on agricultural land for the years 1970, 2000, 2030 and 2050. Their budgets for the year 2000 are similar to the ones calculated by us with some discrepancies in the details (see table 4). BNF is calculated from yield and N contents for leguminous crops while an invariant BNF rate is applied to all pastures. This leads to rather low BNF fixation (30 Tg yr −1 ) compared to our results for agricultural land (75-78 Tg yr −1 ). Manure N, on the contrary, is very high in Bouwman et al (2009) with 101 Tg yr −1 compared to our results of 76 Tg yr −1 . This could be explained by our low manure N management and application rates taken from Herrero et al (2013b). Looking at table 4, one can see that the regional N surplus as calculated in this study fits well with the results presented by Bouwman et al (2009), keeping in mind that the base year of our calculations is 2010 as oppose to 2000 as used by Bouwman et al (2009).
The NUE calculated for Europe for 2001-2003 by Leip et al (2011) is 54% which fits well to our calculation of 56%-57%. However, looking at each country in more detail, greater discrepancies become visible, especially for Estonia, Finland, Italy, Ireland and Sweden (figure 7). In Italy, synthetic fertilizer input influences the results the greatest with Leip et al (2011) ascribing nearly 70% higher synthetic fertilizer input. This discrepancy can be explained by the difference in base years between Leip et al (2011Leip et al ( ) (2001Leip et al ( -2003 and our NUE calculations (2010). Mineral fertilizer application in Italy was reduced by more than 60% from 2001 to 2010 (UNFCCC 2020). In Sweden, Leip et al (2011) calculates a higher N withdrawal leading to a higher NUE than our results. This could be explained by a reduction in crop production from 2001 to 2010 (FAOSTAT 2020). The same is true for crop production in Ireland, which decreased by about 25% from 2002 to 2010 (FAOSTAT 2020). In Estonia, BNF highly influences the discrepancies as the approach used by Leip et al (2011) (using a fixed share of N needed for uptake being covered by BNF for all grazing land) leads to a 30-fold lower result.

Conclusion
As was shown throughout this paper, N surplus calculations vary greatly between regions, with discrepancies in grazing land definition being more influential on the results than methodological differences. Highest discrepancies in both cases were found in Asia. We found that BNF constitutes an important translator for these differences, while also being a source of further uncertainty, highlighting the opportunity for further research.
The right choice of a land-use map is fully determined by its objective as each map has its own advantages and disadvantages. While HYDE offers the opportunity to develop a time series for the past 12.000 years, Ramankutty facilitates the combination of a cropland map with a crop map (M3). Both maps also enable an easier comparison with other studies as they follow the most common definition of grazing landpermanent pastures.
When aiming at evaluating the impact of human activities on the global N cycle, we argue that the Erb grazing land map is a better choice. Its more encompassing grazing land definition leads to the inclusion of N surplus in areas that would otherwise have not been considered. While these areas may be low in productivity, their contribution to overall N supply remains relevant hinting towards a more profound anthropogenic influence on natural systems.
Given this extension, that would also redistribute some additional N over a considerably larger area, the overall pattern provides new insights into the anthropogenic alterations on N and specifically on NUE and N surplus.
These insights offer a way to expand the coverage of previous work (e.g. Bouwman et al 2013, Tian et al 2018, Dangal et al 2019 to not only include permanent pastures but also browsing and transhumance land. This expansion is essential for improving the understanding of the role of human activity in terms of planetary boundaries. Furthermore, it paves the way towards a more complete accounting of potential N release from agricultural lands, e.g. when attempting to assess non-linear effects of N availability on N 2 O emissions (see e.g. Shcherbak et al 2014) as well as for forging effective policy strategies aimed at a more sustainable development of the food and land systems.