We integrated diverse datasets into a harmonized framework, disaggregating to crop level (Supplementary Table 1). We applied the Monte Carlo approach to estimate gridded global cropland GHG emissions for 46 crop classes in 2020 following the IPCC guidelines (Extended Data Fig. 1). Emissions were expressed as carbon dioxide equivalents (CO2e) using the 100-year global warming potential values from IPCC AR6. Crop-specific GHG areal and caloric intensities are calculated by normalizing emissions against physical area and calorie production, respectively. For regional analyses, we aggregated results by World Bank regions (Supplementary Fig. 1).
Climate zone, runoff and irrigation
We obtained yearly temperature, number of frost days, precipitation and potential evapotranspiration at 5-arcmin resolution over the period of 2010–2040 from FAO Global Agro-Ecological Zones (GAEZ v.4). The period comprises the historical and predicted climate using representative concentration pathway RCP 4.5. Following the IPCC 2019 refinement guidelines, we generated global gridded climate zones (Supplementary Fig. 2). The NASA Land Data Assimilation System40 estimated daily surface and subsurface runoff at 5-arcmin resolution. We stacked the daily surface and subsurface runoff maps for 2020 to indicate where runoff occurred (Supplementary Fig. 3). By integrating census from national and subnational survey and reports from FAO, World Bank and other international organizations, irrigation statistics including the pixel percentage equipped for irrigation and the percentage actually irrigated were spatialized at 5-arcmin resolution41. We excluded the pixels percentage with <1% as the 99th outlier. We multiplied these two percentage maps and imply the pixels that were irrigated when their product is greater than zero (Supplementary Fig. 4).
Cropland and managed grassland
We obtained and linearly interpolated national crop-specific harvested area during 1984–2021 from FAOSTAT12. For the spatial distribution in 2020, we used the Spatial Production Allocation Model (SPAM 2020 v.2.0)25, which provides gridded crop-specific harvest area, physical area and production at a 5-arcmin resolution. SPAM 2020 ensured consistency with FAOSTAT national totals and integrated satellite-based datasets, crop-specific suitability information and expert judgement with national or subnational statistics to generate global crop distribution maps for 46 crop classes (Supplementary Table 2). We estimated gridded calorie production by multiplying crop production with crop energy conversion factors for each of 46 crop classes (Supplementary Table 3). We aggregated FAOSTAT crop categories to match SPAM crop classes and constructed time series of national crop-specific harvested areas for 1984–2021.
We determined the spatial extent of managed grasslands using the latest History Database of the Global Environment (HYDE v.3.3)26,42. For each country, we averaged managed grasslands for the period 2017–2023, including managed grazing lands for high stocking densities, hay production and silage. Although FAO land use statistics theoretically encompass these categories, the available data are insufficient for consistent global mapping. Furthermore, distinguishing pasture from rangeland remains challenging and uncertain, as this differentiation is not explicitly addressed in FAOSTAT. However, for some countries, FAO reports permanent meadows and pastures as ‘cultivated’ and ‘naturally growing’ categories.
The harvested area for the 46 SPAM 2020 crop classes and managed grassland areas at subnational administrative divisions (state/province) were subsequently extracted from the generated land use maps. These subnational areas were used for developing synthetic fertilizer data.
Synthetic fertilizer
We derived national total nitrogen (N) fertilizer consumption and application rate from 1984 to 2020 from FAOSTAT12. Missing years were gap-filled using a temporal–spatial–economic referenced interpolation approach, accounting for time-series trend, spatial location and economic status. We assumed that countries within the same subcontinent and income level experienced similar interannual fluctuations in N application rates43. For gaps longer than 3 years, missing periods with available data on one (equation (1)) or both sides (equation (2)) were filled using subcontinental and income-level averages; for gaps of 3 years or less, cubic spline interpolation was applied27. Total N consumption for missing years was then computed as the product of interpolated application rate and harvest area.
$${{\rm{Interpolated}}\,{\rm{rate}}}_{i+k}=\frac{{{\rm{Referenced}}\,{\rm{rate}}}_{i+k}}{{{\rm{Referenced}}\,{\rm{rate}}}_{i}}\times {{\rm{Raw}}\,{\rm{rate}}}_{i}$$
(1)
$$\begin{array}{l}{\mathrm{Interpolated}\,\mathrm{rate}}_{i+k}=\displaystyle\frac{{\mathrm{Referenced}\,\mathrm{rate}}_{i+k}\times {\mathrm{Raw}\,\mathrm{rate}}_{i}}{{\mathrm{Referenced}\,\mathrm{rate}}_{i}}\times \displaystyle\frac{k-i}{j-i}\\\qquad\qquad\qquad\qquad\quad+\displaystyle\frac{{\mathrm{Referenced}\,\mathrm{rate}}_{i+k}\times {\mathrm{Raw}\,\mathrm{rate}}_{i}}{{\mathrm{Referenced}\,\mathrm{rate}}_{j}}\times \displaystyle\frac{j-k}{j-i}\end{array}$$
(2)
Where ‘raw rate’ refers to the incomplete data containing missing values, ‘referenced rate’ is the complete data providing interannual variation reference, i and \(j\) denote the beginning and ending years of a gap and i + k represents the kth missing year.
We compiled subnational total N fertilizer consumption for 40 countries from national agencies, fertilizer associations and literature (Supplementary Table 4 and Supplementary Fig. 5). When 2020 subnational consumption was unavailable, we estimated it using the share of each subnational unit relative to the national total from the closest available year. If N-specific consumption was missing, we adopted the share of aggregated N–P–K fertilizer consumption, assuming it reflects the share of N fertilizer.
The share of grassland N fertilizer to national total N consumption was derived from the International Fertilizer Association (IFA)44 and linearly interpolated for missing years during 1984–2020. For countries lacking grassland fertilizer data but with managed grassland, we assigned a minimum value of 0.1% of national N consumption, consistent with IFA lower bounds. We capped application rates at 400 kgN ha−1 yr−1, roughly 50% greater than the IFA maximum rate for grassland and updated the share of grassland N consumption accordingly. The updated grassland N consumption accounts for 5% of global N fertilizer use, comparable to 3% reported by IFA. We assumed negligible fertilizer use in forestry, aquaculture and rangelands45. Country-specific shares were then applied to partition FAOSTAT national total N consumption between cropland and managed grassland for 1984–2020.
IFA periodically reported national crop-specific N, P and K fertilizer application rates from 1984 to 2019. We obtained the compiled fertilizer dataset46, standardized crop names and supplemented missing harvest area from FAOSTAT. For countries reporting only total N–P–K consumption, the total was partitioned into N fertilizer using subcontinental ratios of N to N–P–K consumption. Area-weighted mean rates for integrated crops (for example, other cereals and other roots) were based on SPAM 2020 crop classes.
We then gap-filled national crop-specific N fertilizer application rate for 2020. For countries with at least one available rate, missing periods on one (equation (1)) or both sides (equation (2)) were filled using national cropland average application rate as a reference. For countries with no crop-specific data for the entire period, we substituted the average rate of the same crop within the corresponding subcontinent and income-level group. To ensure consistency with national or subnational total N consumption from FAOSTAT and subnational information, crop-specific rates were adjusted using equation (3).
$${\rm{Adjusted}}\,{\rm{rate}}=\frac{{\rm{Consumption}}}{{\sum }_{i=1}^{46}{{\rm{Rate}}}_{i}\times {{\rm{Area}}}_{i}}\times {\rm{Rate}}$$
(3)
Where ‘consumption’ refers to national or subnational N fertilizer consumption, ‘rate’ is crop-specific N application rate, ‘area’ is the corresponding national or subnational crop-specific harvest area and i denotes 46 crop classes.
The share of individual N fertilizer type relative to national total N consumption was calculated from FAO for the year 2020. For countries not reported by FAO, we used the data reported by IFA instead. These shares were applied to disaggregate crop-specific rate by N fertilizer types for each country (Supplementary Fig. 6).
Animal manure applicationManure production
To estimate the gridded annual nitrogen production rate from manure for 2020, we used the harmonized Gridded Livestock of the World dataset (GLW v.4.0)47, developed by FAO. This dataset provides the spatial distribution of cattle, swine, poultry, sheep, goats and buffalo at a 5-arcmin spatial resolution. GLW database was developed with an asymmetric method that assigns different weights to pixels using environmental and social predictor variables48. Animal census counts were allocated according to these weighted distributions using Random Forest models49,50. A full description can be found in ref. 51.
We characterized animal populations and herd/flock dynamics following the Global Livestock Environmental Assessment Model (GLEAM) equations52, disaggregating animal populations by species, liveweight, production phase and feeding practices at the pixel level. This enabled the estimation of cohort-specific manure production. Regional N excretion rates for ruminants across agro-ecological zones (thermal and moisture regimes) were obtained from ref. 20, while those for poultry and pig were sourced from IPCC 2019 refinement.
Manure application to croplands
Managed manure N was calculated by multiplying N excreted amount with the share treated in manure management system (MMS), excluding the uses as fuel or deposited on pasture, range and paddock. The shares in each MMS by livestock type, commodity and region were derived from IPCC 2019. MMS categories include lagoon, slurry, solid storage, drylot, daily spread, digester, pit (<1 month and >1 month) and other. Volatilization and leaching losses were applied by MMS using IPCC 2019 fractions, while the losses through disposal, waste discharge, collection and transport were not considered due to limited data availability (Supplementary Dataset). The remaining manure N was summed across livestock types within each pixel. Application rates were calculated by dividing the summed manure N by cropland harvested area. Owing to the lack of crop-specific information, we assumed that manure N was evenly distributed across all crops within the pixel20 and combined rate from N fertilizer and manure was less than 700 kg N ha−1 yr−1 following ref. 8.
Crop residuesCrop residue production
Gridded crop-specific yield and harvested area at 5-arcmin resolution were sourced from SPAM 202025. We estimated aboveground crop residue yield by applying crop- and region-specific residue-to-yield ratios, derived from fixed values or empirical yield-dependent equations (Supplementary Tables 5 and 6). Residue yields from all applicable equations were averaged to minimize bias. Aboveground residue mass per pixel was obtained by multiplying harvested area by residue yield and total crop-specific N per pixel was calculated using IPCC 2019 N content values.
On-site burning
We estimated on-site burning of agricultural biomass using monthly estimates from 2019 to 2021 at 0.25-degree resolution from Global Fire Emissions Database (GFED4s)53. These estimates combined satellite observations with biogeochemical modelling. Our gridded residue biomass of all crops was aggregated to match GFED4s resolution to estimate the crop-specific residue burned fraction, which was then downscaled to 5-arcmin pixels. Total burned biomass per pixel was capped at 90% of cereal residue production.
Animal use, feeding and bedding
Cereal residues are an important feed source for ruminants but generally unpalatable to monogastric animals. To quantify residue use for ruminant feed, country-level cattle, sheep and goat meat and milk production (2019–2021, FAO) were spatially allocated to pixels using GLW v.4 (around 2020). In each pixel, meat and milk production were converted into residue demand using feed conversion ratios (FCRs) and crop residue feed fractions (CRFFs). FCRs represent the feed dry matter needed to produce meat or milk protein at herd level, accounting for unproductive animals. CRFFs represent the fraction of an animal’s diet composed of crop residues. Regional averages of FCRs and CRFFs were sourced from ref. 20 and were disaggregated by animal type, production type, livestock production system and climate zone. The dominant livestock production system per pixel was identified using Global Distribution of Ruminant Livestock Production Systems (GRPS v.5). We calculated protein content using standard values: 138 g protein kg−1 for bovine meat, 137 g protein kg−1 for sheep and goat meat and 33 g protein kg−1 for milk.
Data on residue use for bedding are scarce but important in wealthier regions where residues are less critical as feed54. European Union bedding requirements were used as a baseline and applied globally, assuming they reflect high-income country patterns. Daily bedding requirements were 0.375 kg for cattle, 0.1 kg for sheep and goats and 0.0625 kg for pigs, accounting for seasonal and production system differences. Poultry bedding was excluded, as industrial poultry farms typically use sawdust or other materials. Bedding use was considered only for cattle, sheep and goats within mixed, urban or other LPS. These requirements were converted into gridded yearly demands using GLW v.4 animal population data.
Other off-field uses
In low-income countries, crop residues are often burned as domestic fuel, depending on alternative fuel availability and competing animal feed. In China, where monogastric livestock dominate and forest cover is low, about 30% of residues were historically used for fuel in the mid-1990s55. We adopted this 30% as an upper limit for other off-field uses (fuel, construction materials, mushroom cultivation and so on) in low-income countries. Low-income countries are defined as those with gross domestic product (GDP) per capita below US$1,035 (2019–2021 prices, World Bank).
In high-income countries (GDP per capita >US$12,535), industrial applications, biogas, construction and cultivation substrates account for around 10% of crop residues. For middle-income countries, the maximum share was linearly interpolated between 30% (low-income) and 10% (high-income) according to GDP per capita. Over time, changes in GDP adjusted this upper limit dynamically.
Fractional allocation at pixel level is determined by GDP and residue availability. Combined residues removal from burning, animal use and off-field uses were capped at 90% to ensure minimum retention.
Trading of crop residues
In some pixels, estimated residue for animal usage exceeded available residues. Such imbalances are expected since residues are often transported between regions for animal feeding and bedding. To address this, we redistributed residues from surplus to deficit pixels, prioritizing nearby supply within a maximum distance of 80 pixels (approximately 750 km at the equator). A pixel is considered in deficit if combined on-field burning and/or animal usage exceeded 90% of local crop residue production. Remaining deficits after redistribution were assumed to reflect local overestimation of burning or animal demand. If residues remained after burning and animal use, they were allocated to other off-field uses if estimated for that pixel, while maintaining at least 10% residue retention.
Left on field
Crop residues not assigned to burning, feeding, bedding or other off-field uses were assumed to remain on the field. A minimum of 10% residues retention was enforced, reflecting the practical difficulty of fully removal.
Belowground residues
Belowground residues were estimated as a crop-specific ratio of aboveground biomass (including harvested yields and aboveground residues) using the parameters from the IPCC 2019 refinement. Since aboveground residue production varies with crops and regions, belowground biomass estimates also vary spatially across pixels.
Rice cultivationGrowing season
We determined the number and duration of growing seasons using RiceAtlas24, which compiled national and subnational data from 115 countries (58 with subnational details) reviewed by regional experts. We converted these original vector layers to gridded maps at a 5-arcmin resolution.
Water regime
SPAM 2020 distinguishes rice as irrigated or rainfed. In line with IPCC 2019, we further disaggregated rainfed rice into upland, regular rainfed, drought-prone and deepwater categories and irrigated rice into continuously flooded, single-aeration and multiple-aeration systems (Supplementary Table 7). By integrating national statistics, United Nations Framework Convention on Climate Change (UNFCCC) reports, FAO Rice Information and literature sources, we reconstructed national rice water regime data for 85 countries with rainfed rice (Supplementary Table 8) and 35 countries for irrigated rice (Supplementary Table 9). When several years were available, the most recent data were adopted. The compiled rainfed water regime data cover most rice-growing counties (Supplementary Fig. 7). For countries lacking rainfed regime data, we applied subcontinental mean ratios. For irrigated rice, we used the global mean ratio (excluding China, Japan and South Korea) (Supplementary Fig. 8). These ratios were applied to SPAM 2020 rainfed and irrigated areas to generate gridded rice water regime maps. We determined the precultivation water status for different rice water regime systems following ref. 21 (Supplementary Table 10).
Organic amendments
Organic amendments, including straw incorporation, compost, animal manure and green manure, are considered as key factors affecting rice CH4 emission in IPCC 2019. Owing to limited global data and the negligible quantities of compost and green manure applied, only straw incorporation and animal manure were included in this study.
Drained peatland
We delineated global peatland extent by combining four representative datasets. Peat-ML23 served as the base map, as it provides the broad coverage by leveraging 165 predictors (climate, soil, vegetation and terrain) and a machine learning model trained on 14 major peatland and peatland-free ecoregions. To complement Peat-ML, we incorporated three indicative products: (1) Global Forest Watch56, which compiles five expert-rule-based regional peatland layers; (2) Global Peatland Map 2.0 (GPM v.2.0)57 produced by United Nations Environment Programme-led GPA from over 200 international experts; and (3) FAO drained organic soils58, which uses Histosol soil layer as a proxy for peatland distribution.
We merged these products into a hybrid peatland extent map (equation (4) and rescaled pixel-level peatland fractions to match national peatland areas reported by GPA (Supplementary Fig. 9). Drained peatlands under cropland were initially identified by multiplying peatland fractions with SPAM 2020 physical cropland area (46 crops). We then adjusted pixel-level drained-peatland area to align with the corresponding national totals by GPA. After these adjustments, regions with a long history of peatland drainage or recent expansion such as North America, Europe and Southeast Asia show high pixel-level drainage ratio of peatland (>20%), whereas Africa and South America generally exhibit very limited disturbance (<1%) (Supplementary Fig. 10). While this national-level constraint ensures that the magnitude of drained-peatland emissions is realistic at the country scale, substantial uncertainty persists in the fine-scale spatial allocation of drainage within countries. Some regional mismatches remain relative to field-based peatland studies, underscoring the need to incorporate field evidence into future global peatland drainage datasets to better capture regional heterogeneity and refine future assessments (Supplementary Results Section 3).
$${\mathrm{Per}}_{\mathrm{peat}}=\left\{\begin{array}{l}0,\mathrm{At}\,\mathrm{least}\,\mathrm{two}\,{\mathrm{Per}}_{\mathrm{ind}}=0\,\mathrm{and}\,{\mathrm{Per}}_{\mathrm{peatML}} < 10 \% \\ \mathrm{Mean}\,\mathrm{of}\,\mathrm{three}\,{\mathrm{Per}}_{\mathrm{ind}},\,\mathrm{At}\,\mathrm{least}\,\mathrm{two}\,{\mathrm{Per}}_{\mathrm{ind}}\\ > 0\,\mathrm{and}\,{\mathrm{Per}}_{\mathrm{peatML}}=0\\ {\mathrm{Per}}_{\mathrm{peatML}},\mathrm{Otherwise}\end{array}\right.$$
(4)
Where Perpeat, Perind and PerpeatML are pixel-level percentage of peatland, indicative products and Peat-ML, respectively.
Monte Carlo simulations using IPCC guidelines
We followed the latest IPCC guidelines for each emission source (Supplementary Methods Sections 12–15). Specifically, we adopted the IPCC 2019 refinement to estimate N2O emissions from synthetic fertilizer, manure and crop residue incorporation and CH4 emissions from rice paddies. The IPCC 2006 guidelines were used for crop residue burning, and the IPCC 2013 supplement for emissions from cultivated drained peatlands. For rice cultivated on peatlands, only CH4 emissions were estimated following the IPCC 2019 methodology for rice paddies, assuming no on-site peat emissions and including only drainage ditch emissions.
IPCC emission factors were derived from different approaches depending on data availability. For N2O emissions from N inputs and CH4 emissions from rice paddies, where observations are abundant, linear mixed-effect models were used that capture the effects of climate, crop type, nitrogen form, water management and organic amendment. By contrast, for drained peatlands and crop residue burning, where observational data are sparse, mean values stratified by climate zone and crop type were applied.
To represent the local variability and quantify uncertainty, we performed Monte Carlo simulations with 500 iterations. For each emission factor and parameter, we constructed probability density functions using reported means, uncertainty ranges and sample sizes in the IPCC guidelines (Supplementary Table 11). Following the IPCC 2006 guidelines (vol. 1, ch. 3), we used lognormal distributions for positive variables, beta distributions for fractions bounded between zero and one and normal distributions for in situ CH4 emissions from drained peatlands, which can include negative values.