Recent global land cover dynamics and implications for soil erosion and carbon losses from deforestation

(1)

Recent global land cover dynamics and implications for soil erosion and carbon losses from deforestation

Xiangping Hu

^a,

*, Jan Sandstad Næss

^a

, Cristina Maria Iordan

^a

, Bo Huang

^a

, Wenwu Zhao

^b

, Francesco Cherubini

^a

aIndustrialEcologyProgram,DepartmentofEnergyandProcessEngineering,NorwegianUniversityofScienceandTechnology,N-7491,Trondheim,Norway

bStateKeyLaboratoryofEarthSurfaceProcessesandResourceEcology,FacultyofGeographicalScience,BeijingNormalUniversity,Beijing,100875,China

ARTICLE INFO

Articlehistory:

Received3November2020

Receivedinrevisedform12March2021 Accepted24March2021

Availableonline26March2021

Keywords:

Landcoverchange Climatechange Deforestation Landuseemissions Soilerosion

ABSTRACT

Changesinlandcoverareincreasinglyaffectinglandsurfacepropertiesandprovisionofecosystem services. Understanding recent historical land cover changes and their interlinkages with key environmentalprocessesisinstrumentaltobettersupportstrategiesforland-usemanagement.The recentlyreleasedproductsfromtheEuropeanSpaceAgencyandCopernicusClimateChangeService containhigh-resolution(300m)timeseriesofgloballandcovermapsfrom1992to2018.Thisstudy investigatesthelandtransitionsintheseproductsandexplorestheeffectsontwokeyenvironmental aspects,namely,carbonlossesfromdeforestationandsoilerosionrates.Weusedapowerfulserverfor bigdataanalysistoretaintheoriginalspatialresolutionofthedatasets.Wefound722Mha(5.5%ofthe totalice-freelandsurface)ofgrosslandcoverchanges,whichmainlyinvolvedtransitionstoandfrom forest/agriculture. Cropland gains are 205 Mha and losses 126 Mha (net expansion of 79 Mha).

Deforestationoccurringin242Mhawasmainlycausedbyagriculturalexpansion,whereas196Mhawere afforested.Settlementsshowthelargestrelativeexpansion(44Mha,+210%),ofwhich67%(29Mha) occurredatexpensesofagriculturalland.Deforestationcaused12.3(7.6/14.2)GtCarbonlossesfrom below-andabovegroundbiomassfrom2010to2018,correspondingto1.5(1.0/1.8)GtCarbonperyear.

Globalagricultureactivitieshaveincreasedtotalsoilerosionof3.2Gtandsoilerosionratesof0.22Mg ha¹yr¹intheperiod2001–2012,especiallyintropicalregions.Theidentiﬁedlandtransitionsand changesinkeyenvironmentalprocessesreﬂectahuman-dominatedEarthsystemandtheindirecteffects ofclimatechangeonlandcover,especiallyinborealecosystems.

1.Introduction

Landisthebasisforhumanlivelihoodsandwell-being.Itisthe mainsupplieroffood,freshwaterandmultipleotherecosystem services, including biodiversity. Human activities are already affectingmorethan70%(likely69–76%)oftheglobal,ice-free land surface (IPCC, 2019).Anthropogenicactivitieshave caused extensivelandcoverchanges(KleinGoldewijketal.,2017),which alsoaffectedprotectedareaswithinbiodiversityhotspots(Bailey etal.,2016;Huetal.,2020).Changesinlandcoverarebothadriver anda consequenceofglobalenvironmental change(Foleyetal., 2005; Bonan,2008;Alkama and Cescatti,2016), andtheyhave

greatimpactonoursocietyandecosystemsatregionalandglobal scales(Verburgetal.,2011;Turneretal.,2007).

Assessingchangesinlandcoverandtheirenvironmentaleffects isanon-goingandfastdevelopingresearcharea(Nowosadetal., 2019;Huaetal.,2018;MousivandandArsanjani,2019;Liuetal., 2018b;IPCC,2019).Majorapplicationsaimtounderstandspatial and temporal patterns of land cover dynamics (IPCC, 2019;

Mousivand and Arsanjani, 2019; Liu et al., 2018b), explain the underlyingmechanismscausingthechanges(Lepersetal.,2005;

Ceccherinietal.,2020;Jaimesetal.,2010),constructmodelsfor assessing challenges to society and environment (Hurtt et al., 2020;Kuemmerleetal.,2016;Chenetal.,2020),quantifyimpacts toclimatechangeataglobal(AlkamaandCescatti,2016;Bonan, 2008;Presteleetal.,2017)orregionalscale(Huangetal.,2020a;

Huetal.,2019;Lejeuneetal.,2018;Cherubinietal.,2018a),assess effectsonecosystemservices(Tolessaetal.,2017;Chenetal.,2019;

Bayeretal.,2020;Venäläinenetal.,2020),andhelpstakeholders todesignmoresustainablelandusepolices(EllisandRamankutty, 2008;Duveilleretal.,2020;Verburgetal.,2011;Seneviratneetal.,

*Correspondingauthorat:IndustrialEcologyProgram,DepartmentofEnergy andProcessEngineering,NorwegianUniversityofScienceandTechnology,N-7491, Trondheim,Norway.

E-mailaddress:[email protected](X.Hu).

http://dx.doi.org/10.1016/j.ancene.2021.100291

ContentslistsavailableatScienceDirect

Anthropocene

j o u r n al h o m e p a g e : w w w . el s e v i e r . c o m / l o c at e / a n c e n e

(2)

2018;Englundetal.,2020).Differenttypesoflandcoverorland usechangeshavebeenstudiedintheliterature,suchasagriculture expansionandcontraction(Grauetal.,2005;Spawnetal.,2019;

Krause et al., 2009; Næss et al., 2021; Leirpoll et al., 2021), deforestation and afforestation(Setoet al., 2012; Jaimeset al., 2010; Li etal., 2016; Hansenet al.,2013), forestdynamics and management(Ceccherinietal.,2020;Panetal.,2011;Hansenetal., 2013;Cherubinietal.,2018b),urbanization(Arsanjanietal.,2013;

Cakiretal.,2008;Arsanjanietal.,2018),wetlandshrinkage(Song etal.,2012;Ghoshetal.,2018;DebanshiandPal,2020;Xuetal., 2019),anddesertification(Veronetal.,2006;Bestelmeyeretal., 2015; Lamchin et al., 2016). A common approach used is the temporal analysis of variability per pixel to detect land cover transitions.Resultsarefrequentlyshownasaggregatedaveragesat relativelylowresolutionorfocusonspecificregionsandspatial scales. Thisis due totwo main issues that analysts face when assessing changes in land cover at high spatial resolution and globalscale.Thefirstisrelatedtothechallengesinprocessingbig datasets, because timeseries of high-resolution global datasets withdetailedlandcoverclassificationareusuallyhighlycompu- tationallydemandingandtimeconsuming.Thesecondisdueto thelackofboth spatiallyandtemporarilyconsistentlandcover datasetsthatcovertheglobeforalongtimeperiodofobservations (Nowosad et al., 2019). Recently, theEuropean Space Agency’s ClimateChangeInitiative(ESACCI)providedannualgloballand covermapsfrom1992to2015,andtheCopernicusClimateChange Serviceclimatedatastore(C3S-CDS)releasedthegloballandcover annualmapsfrom2016to2018(C3S,2019).Thesetwodatasetsare highlyconsistentandhavehighspatialresolution(300matthe equator),andtheyofferthepossibilitytoanalyzegloballandcover transitionsindetail.

SomestudieshaveinvestigatedtheESACCIlandcoverproduct fromdifferentperspectives.Lietal.(2016)usedanearlyversionof theESACCIlandcoverproducttoanalyzelandtransitionsatthree speciﬁc points in time, 2000, 2005 and 2010.However,due to computational challenges,they disregardedsome transitionsto save resources and time. Liu et al. (2018b) identiﬁed land transitionsandhotspotsoflandusechangesusingthetrajectory analysis method. However, again due to limited computation capacity, they had to aggregate the datasets into 0.5 0.5 resolution and then compute the area proportions in the aggregated datasets.Liet al.(2018)analyzed thegrossandnet

changesinlandcoversforsomeplantfunctionaltypescommonly used in land surface models. They investigated the spatial distributionofchangesinlandcoverforcropland,grasslandand forestsbetween1992and2015.TheyalsocomparedtheESACCI land cover product to the LUH dataset (Hurtt et al., 2011), concludingthattheESACCIdatacanimprovetherepresentation of land cover dynamics in land surface models and the characterizationoflocalorglobalcarboncycledynamics.Nowosad etal.(2019)assessed andvisualizedmainlandcovergainsand lossesbetween1992and2015afterpre-processingthelandcover datatoresizethematalandscapelevel.TheESACCIlandcover product was also used by Mousivand and Arsanjani (2019) to quantifyhistoricalchangesinlandcoversandthenpredictfuture land transitions using Markov chain. Some of these studies apparently did not consider the curvature of the Earth while translatingthegriddedlandcoverdatasetsintoarealextensions, andtherebyintroducedabias(overestimation)oftheareasforthe grid cells at high latitudes (Liu et al., 2018a; Mousivand and Arsanjani,2019;Liuetal.,2018b).Overall,thesestudiesdidnot includetherecentC3S-CDSproduct(whichalsoincludesannual landcovermapsfrom2016to2018),sinceitwasnotavailableat that time, and they did not quantify the speciﬁc land cover transitionsor theirspatialdistributions. Preliminarystudiesare alsointegratingtheESACCIlandcover productswithadjacent ﬁelds, such as remotely sensed data or other ground-level observations tostudy the effects of changes in land covers on regionalclimate(Huangetal.,2020a;Duveilleretal.,2018),and, althoughrelativelylittleexploredsofar,theyofferpossibilitiesfor integrations with multiple datasets to discern the effects of changesinlandcoversonenvironmentalareasofconcerns,likefor examplesoilerosion, vegetation carbonstorage,and ecosystem services.

In this study,we quantitativelyanalyzeand visualizerecent historicallandcoverdynamicsatagloballevelfrom1992to2018 asrepresentedbytheintegrationoftheESACCIandC3S-CDSland coverproducts.Theanalysisretainstheoriginalspatialresolution (300m)ofthemapsthankstotheuseofapowerfulserverforbig dataanalysis,anditshowsgainsandlossesinlandcover,aswellas trendsandspatialpatternsoflandcovertransitions(i.e.,fromone landcoverclasstoanother).The37originallandcoverclassesof thedatasetsaretranslatedintothelandcoverclassesspeciﬁedby theIntergovernmentalPanelonClimateChange(IPCC),whichare

Fig.1.Asimpliﬁedﬂowchartoftheapproachandmainmethodologicalstepsconsideredinthisstudy,fromtheoriginaldatasetstotheestimatesofcarbonemissionsfrom vegetationlossesandsoilerosion.

(3)

morereliablefortheidentiﬁcationoflandcovertransitionsthan theoriginalclasses.

Tothebestofourknowledge,thisistheﬁrstattempttomerge andanalyzethetwolandcoverproducts,toshowglobalmapsand timeseriesofbothchangesinlandcovers(lossesandgains)and landtransitions betweeneach pairoflandcoverclasses,andto identifyareasofpotentialhotspotsofchangesinlandcoversata globallevel.Inaddition,weapplythisintegrateddatasettoassess theeffectsofhistoricalrecentchangesinlandcoversontwomajor global environmental issues: total carbon emissions and their spatial patterns from losses in aboveground and belowground biomasscarbonduetodeforestationtoagricultureandsettlement;

contributiontoglobalsoilerosionfromlanduse(i.e.,agriculture) andchangesinlandcovers(e.g.,fromforesttoagricultureandvice versa).Fig.1 showsa simpliﬁedﬂowchartof theapproach and mainmethodologicalstepsconsideredinthestudy.

Insummary,theaimofourworkistoaddressthefollowing researchquestions:

Whatarethemaingloballandtransitions from1992to2018 accordingtotheESACCIandC3SCDSproducts?

What are the spatial distributions of the main land cover transitions?

Howdoesgloballandcoverandlandusechangesinﬂuencekey environmentalaspects,suchascarbonlossesfromdeforestation andsoilerosionrates?

2.Methodology 2.1.Landcoverdatasets

ThelandcoverproductsfromESACCIandC3S-CDSprovide global annualmapsfrom1992 to2015and from2016to2018, respectively(ESA,2017;C3S,2019),thatarehighlyconsistentover time and space. Theywere speciﬁcally developed to provide a consistent evaluation of the historical temporal evolution of changesinlandcoversatagloballevel.Toensurethecontinuity, thesetwodatasetshavethesamespatialresolution(300matthe equator,or10arc-sec),andeachmapisintheformof129600 64800gridsLon/Latraster.BothproductsusetheWorldGeodetic System84(WGS84)referenceellipsoidasthecoordinatereference system,andtheyprovideglobalmapsdescribingthelandsurface usingthesamecategoriesandcodingscheme.Thesedatasetshave beendevelopedtoincreasetherobustnessoflandcoverdatasets forcarbonandclimatechangemodels,landsurfacestudies,and landscape dynamic assessment (ESA,2017; C3S, 2019; Liet al., 2018,2016).

Thedatasetsarebasedon37landcoverclassesaccordingtothe UnitedNationsLandCoverClassiﬁcationSystemtodescribethe Earth’s terrestrial surface (Di Gregorio, 2016). These maps are obtainedbycombiningseveralearthobservationproductsandby usingtheGlobCoverunsupervisedclassiﬁcationchain(ESA,2017;

Defournyetal.,2009;Poulteretal.,2015).Unliketheothersingle- yearproducts,thelandcoverproductsfromESACCIandC3S-CDS maintaingoodtimeconsistency,andtheiroverallglobalaccuracy is about 71 % (ESA, 2017; C3S, 2019). Some classes, such as cropland,forests,urbanandbareareashavehigheraccuracy,and others, such as mosaic classes, have lower accuracy. The high certainty of cropland classes makes the land cover products especially useful for cropland monitoring. Global user and produceraccuraciesofcropland classesrangebetween85–94% and76–92%acrossbothdatasets,withmediansof89%and82% (ESA,2017;C3S,2019).Thesehighuserandproduceraccuracies indicateahighspatialcroplandmatchandhighprecisionintotal croplandareaextent.Additionally,broadleavedevergreenforests

havehighuserandproduceraccuracies(86%and96%forESACCI, and89%and86%forC3SC-CDS,respectively),therebymakingthe products suitable for monitoring tropical deforestation and afforestation.Forexample,ESACCIhasaregionaloverallaccuracy of85%whenappliedasaforestproductintheBrazilianAmazon (Qinetal.,2019).

The original 37 land cover classes in these two land cover productshavebeentranslatedintothegenericIPCClandcover classeswithacross-walkingtable(Supplementary TableS1) to obtain more reliableresults of land cover transition dynamics (ESA, 2017). This conversion avoids false change detection betweenlandcoverclasseswhicharesemanticallyclosein the original37landcoverclasses(ESA,2017;Liuetal.,2018b).The ESACCIlandcoverproducthasalreadybeenusedtoinvestigate the pattern of changes in land covers (Nowosad et al., 2019;

Tschora and Cherubini, 2020; Liu et al., 2018b), to assess consistency with other products (Hua et al., 2018; Liu et al., 2018a) and to investigate climate-land interactions (Duveiller etal.,2018;Huangetal.,2020a;Lietal.,2018).However,tothe bestofourknowledge,thesetwolandcover productshavenot beenintegratedintoasingledatasetandthespatialpatternsata globallevelhavenotbeendepicted.

2.2.Datasetsofglobalbiomasscarbonandsoilerosion

A harmonized global dataset of biomass carbon density, aboveground(ABC) and belowground(BBC),for 2010at 300m resolutionhasbeenrecentlymadeavailable (Spawnand Gibbs, 2020;Spawn etal.,2020), andit ishereusedtoquantifygross vegetationcarbonlossesduetothemajorlandusetransitionsthat happenedbetween2010and2018.Thiscarbondensitymapwas createdwithanovelmethodthatintegratesremotelysensedmaps ofspeciﬁcvegetationcharacteristicswithancillarymapsoftree cover, land covers and rule-based decision tree (Spawn et al., 2020).Thedatasetcombinespublishedestimatesforvegetation speciﬁcdensitiesandtheirresultshavebeenrigorouslyvalidated.

ThedatasetwasrecentlybuiltfromtheESACCIproduct,andit hasthesameresolutionofthetwolandcoverdatasetsusedinour study.

Landcovertransitionsandanthropogenicactivitiesareamong the main causes of soil erosion (SE), which lead to serious consequencestohumansociety,suchasriskstofoodsecurityand depletionofecosystemservices(Borrellietal.,2017;Tarolliand Straffelini,2020;Huangetal.,2020b).SEisusuallydeﬁnedasthe massofnetsoillossperunitareaandtime(Nearingetal.,2017b).

Borrellietal. (2017)providedaresampledglobaldataset ofSE rates in2001 and2012, witha spatial resolutionof 25km, to study temporal changesin soil erosion rates over time. Their estimatesofSEratesarebasedonthereviseduniversalsoilloss equation(RUSLE) modellingapproach,a well-known empirical method for predicting SE (Risse et al., 1993; Stolpe, 2005;

Bagarelloetal.,2017).IntheRUSLEapproach,SEisexpressedas themass of soilloss perunit areaandtimeand computed by combining the contribution factors, such as land cover and management, rainfall runoff erosivity, soil erodibility, slope steepnessandlength,andsoilconservationpractice.Thisdataset was produced to investigate the nexus between land cover transitionsandSE,especiallytomonitortherisksofincreasingSE ratesduetocroplandexpansion.Forexample,Borrellietal.(2017) showedaglobalincreaseinsoilerosiondrivenbylandusechange of2.5%between2001and2012.

Allthesedatasetsareconvertedtothesamespatialresolution ofthelandcoverdata,i.e.,300m.Iftheoriginalmapsdonotcover thewholeglobe,theyareextendedtocover-90N-90Nand-180

E-180E,andvaluesofpixelsintheextendedpartsaremarkedas missingvalues.

(4)

2.3.Analysisoflandcoverdynamics

Toanswertheﬁrstandsecondresearchquestionsofthisstudy, wecombinedtheESACCIandtheC3S-CDSlandcoverproductsto investigatethehistoricallandcoverdynamicsfrom1992to2018.

Our mainfocusisontotalchangeandtrendof eachlandcover class,landtransitionsbetweendifferentlandcovers,andspatial distributionsoftheobservedlanddynamics.Weﬁrstlycomputed thetotalareaofeachgridforeachspeciﬁcIPCClandcoverclass, consideringthecurvatureoftheearth.Thetotalareaforeachland coverclassm(denotedasAream)iscomputedasfollows, Aream¼X^N

j¼1

LCj¼m

Aj ð1Þ

foreveryyearfrom1992to2018.jistheindexofthegridwitha spatialresolutionof300mattheequator,andLCjisthelandcover classingridj.Ndenotesthetotalnumberofgridsonthewhole Earth.Ajistheareaofgridj,whichvariesbylatitude:0.09km²at theequatorandgraduallysmallerwhileprogressingtowardsthe polebecauseofthecurvatureoftheEarth.Thelatteristakeninto accountasfollows,

Aj¼

p

180R²jsinðlat2Þsinðlat1Þj jlon2lon1j ð2Þ where,RistheEarthradius,lat1,lat2,lon1andlon2,arethevaluesof latitudeandlongitudeofgridAj,and||standsfortheabsolute value.Tobetterillustratetheglobaltrendsofchangesineachland coverclass,resultsfromEq.(1)arenormalizedtothetotalareain 1992ofeachlandcovertype.

We then investigate the pairwise land cover transitions (denoted as LCm,Yearp→n,Yearq)betweentwoyears toseeifa grid changedthelandcoverclass(e.g.,frommton),

LCm;Yearp!n;Yearq¼LCj¼m

YearpLCj¼n

Yearq ð3Þ

where, LCm,Yearp→n,Yearqis a matrix with0and 1 toindicatethe transitionoflandcoverclassminyearp(denotedasYearp)toland coverclassninyearq(denotedasYearq)forallpairsofmandn,m6¼ n.Thetotalareaofthetransitionbetweenthetwoyears(denoted asAream,Yearp→n,Yearq)iscomputedbycombiningEqs.(2)and(3)as follows

Aream;Yearp!n;Yearq¼X^N

j¼1

LCm;Yearp!n;YearqAj ð4Þ Eq.(4)canbeusedtocomputethetotalareaofpairwiseland covertransitionsforanygiventwoyearsbetween1992and2018.

The global land cover transitions are identiﬁed and then aggregated to9km spatial resolutionat theequatorbya local windowof30-by-30gridsformapvisualizationpurposes.Thearea proportion(from0to1)ofeachlandcovertransition(denotedas AreaPm,Yearp→n,Year) betweenany chosen two years in each local windowiscomputedas

AreaPm;Yearp!n;Yearq¼ P^P

j¼1

LCj¼m

YearpLCj¼n

YearqAj

!

P^P

j¼1

Aj

ð5Þ BothsumsinnumeratoranddenominatorinEq.(5)arewithin thelocalwindow,whichcontains900grids,i.e.,P=900.Different low spatial resolutions can be obtained through Eq. (5) if necessary.

Thechangesinlandcoverscanbeobtainedsimilarly.First,we identiﬁedthelandcoverincreaseordecrease(denotedasLC_j,m",

YearpandLC_j,m#,Yearq,respectively)forgridjandeachlandcoverclass mbetweenanygiventwoyearsbetween1992and2018,

LCj;m";Yearp!Yearq¼LCj6¼m

YearpLCj¼m

Yearq

LCj;m#;Yearp!Yearq¼LCj¼m

YearpLCj6¼m

Yearq

ð6Þ UsingEq.(6),theincreaseanddecreaseofeachlandcoverare codedas1and-1.Ifthereisnochangebetweenthetwoyears,we setthe valueof the grid to0. Using this approach, the spatial distributionpatternsofthechangesinlandcoverscanbeobtained.

Wethenaggregatedtheresultto9kmspatialresolutionusinga localwindowof30-by-30gridsforvisualizationpurposes,

AreaPm¼ P^P

j¼1LCj;m";Yearp!YearqþLCj;m#;Yearp!Yearq Aj

!

P^P

j¼1

A_j

ð7Þ

SimilartoEq.(5),bothsumsinnumeratoranddenominatorin Eq.(7)arewithinthelocalwindow.

2.4.Estimatesofcarbonlossesandsoilerosion

To answer the third research question of this study, we estimatedthe biomass carbonlosses and SEdue to landcover transitions.Thelossofbiomasscarbonduetolossinforestcover for the transition from forest to agriculture or settlement is estimatedfromtheglobalABCandBBCdensitydatasetsrecently madeavailable(dataarerepresentativeof2010only)(Spawnetal., 2020).The ABCloss due todeforestation (foresttoagriculture, denotedas

D

^ABCFOR→AGR)iscomputedbyidentifyingthegridsof forestareas in2010that havebeen convertedtoagriculturein 2018,assumingthedifferenceinabovegroundbiomassbetween the two land covers as an instantaneous emission to the atmosphere.Forthesegrids, carbonemissionsareestimatedby thedifferencebetweentheassociated ABCcontentofforestsin eachgridtransitionedtoagricultureandtheglobalaveragecarbon contentofagriculturalland,asfollows

D

^ABC^FOR!AGR¼X

g

ABCg;FOR;2010ABCAGR

Ag ð8Þ

where,ABCAGRistheestimatedaverageABCdensitypergridcellg ofallagriculturalgrids.Thisisdonetoconsideraveragecarbon densityforagriculturalland,becausetheABCdensityfor2010is not available (as the grid was classiﬁed as forest). ABCAGR is assumedtobeequaltotheglobalmean,the10thpercentileorthe 90thpercentileofallABCofagriculturalgridsin2010.Thesumin Eq.(8)computesthetotalcarbonemissionsfromallthegridsthat wereclassiﬁedasforestin2010andagriculturein2018.Similarly, theestimatedABClossforforesttosettlement(

D

^ABCFOR→SET)is computedas

D

^ABC^FOR!SET ¼X

h

ABC_h;FOR;2010A_h ð9Þ

undertheassumptionthatdeforestationtosettlementleadsto a total loss of ABC. The sum in Eq. (9) estimates total carbon emissionsfromallthegridshthatwereclassiﬁedasforestin2010 but settlement in 2018. Eq. (3) can be used to identify the transitionsinEqs.(8)and(9).Itwasthenpossibletoproducemaps of the loss of ABC and show the main emission regions, after aggregationtolocalwindowsof30-by-30gridsthatincludethe sumofthethereinlossofABCtoimprovevisualization.Thesame procedure is usedtoestimate carbon lossesfrom belowground biomasscarbon,butinthiscasetheBBCdensitymapsareused.

Similarly, we study SEassociated with themain land cover transitions by integrating the SE dataset and the land cover

(5)

products.TheSEdatasetprovidedbyBorrellietal.(2017)contains theSEratesintwoyears,2001and2012(denotedasSER2001and SER2012).UndertheassumptionthattheSEratechangeslinearly duringthesetwoyears,therateofchangeoftheSEratesineach gridkcanbecomputedas

bk;m!n¼SERk;m;2012SERk;n;2001

D

^t ^ð¹⁰^Þ

with

D

^t⁼^11.^With²⁰⁰¹^as^the^reference^year^(t⁼^0),^the^total

increasedSEfrom2001to2012(denotedas

D

^SEm,2001→n,2012)can becomputedas

D

^SE^m;2001!n;2012¼P

k

Rt¼11

t¼0 bk;m!ntdtAk

¼11 2

X

k

A_kSER_k;m;2012SER_k;n;2001 ð11Þ

ThesuminEq.(11)isovereachgridkwhichisclassiﬁedasland coverclassmin2001andclassnin2012,andm=nisallowed.Such asettingenablesustocomputetheSEnotonlyduetolandcoveror landusechange,butalsoforthesamelandcoverovertime.Inthis study we focus on SE due to land cover transitions between agriculture and forest, agriculture remaining agriculture, and forestremainingforest.

3.Results

3.1.Globaltrendsoflandcoverchanges

Thetrendsinrecenthistorical changesforeachtype ofland coverclassareshowninFig.2,theunderlyinglosses,andgainsin SupplementaryFigureS1,andthemaintransitionsbetweenpairs oflandcoverclassesinFig.3.Resultsarenormalizedtotherelative areaofeachlandcoverclassin1992.

Inrelativeterms,settlementisbyfarthelandcoverclassthat increasedthemostduringthestudyperiod,asithasexpandedby morethan200%(Fig.2b).Thesecondlargestincrease(3%)isfor

agriculture,whosenetareagainspeakedaround2006andthen remainedrelativelyconstant(Fig.2a).However,thisstabletrend does not imply that there were no transitions to and from agricultureandotherlandcoverclasses(mainlyforest)after2006, becauseitistheresultofabalancebetweenagriculturegainsand losses(Fig.S1).Until2006,conversionofforeststoagriculturewas largerthantheoppositetransition(fromagriculturetoforests), resultinginnetareagainsforagriculture,butafter2006thetwo transitions are approximately of thesame magnitude (Fig. 3a).

Grasslandareasexpandofabout1%(Fig.2a),andthisislargelyat expenses of forest areas (Fig. 3b). All other land covers have decreasedwithvaryingtrends.Wetlandhadthemostsigniﬁcant decreasingtrend(-6%),mainlyasaresultofdryingduetoclimate change,unsustainablewaterwithdrawals,andexcessivenutrient loadsthatfavorplantgrowth(Werneretal.,2013;Greenet al., 2017),whichmakevegetationtograduallytakeoverwetlandareas (Fig.3c).Shrublandandsparsevegetationdecreasedsteadilyfrom 1992 to 2008 (-2%), and thereafter the net trends remained relativelyconstant.Changesinareaextensionsofshrublandand sparsevegetationaredifﬁculttointerpret.Shrublandintropical areascanbetypicallyassociatedwithsavannahorCerrado,andthe shrubland-to-foresttransitioncanindicateadevelopmentofearly- stagetreesin1992thatareprogressivelybecominglargerandwith closercanopy.Theoppositetransition(forest-to-shrubland)canbe indicative offorest degradation(Fig.3d).Net changes of forest areasshowadeclineof1%,whichmainlyoccurredbefore1999,but declinesinforestareasarestillon-going(althoughcompensated by forest expansion). Bareareas are also declining in favor of grassland,a transitionthatcanbeexplainedbytheprogressive observedgreeningoftheEarthasafeedbacktorisingatmospheric carbonconcentrationsandtemperatures(Zhuetal.,2016).

3.2.Totallandcoverchangesandspatialdistributions

Thetotalchangesinlandcoversfromthebeginning(1992)to the end (2018) of ourstudy period are quantiﬁed in terms of millionhectares(Mha)andshowninTable1(numbersinthetable refertothetransitionsfromrowtocolumn).Therearetransitions amongalllandclasses,exceptfromsettlementtootherclasses.The global spatialdistribution of the patternsof relative gains and lossesforeachlandcoverandthemainlandcovertransitionsare showninFigs.4and5,respectively.

Atotalof722Mhaofgrosschangesinlandcoversoccurredin ourplanetbetween1992and2018accordingtotheESACCIand C3S-CDS datasets. Agriculture is the land cover class with the largestareagains(205Mha)andforesttheonewiththelargest decrease(241Mha).Atthesametime,agriculturalhaslost126 Mhaandforesthasgained196Mha,resultinginnetareachanges of +79 Mha and –45 Mha, respectively. Agriculture mainly translated to forests (79 Mha), a transition that is typically associated with either cropland abandonment due to socio- economicreasons,suchasinEasternEuropeafterthefallofthe Soviet Union (Lesivet al.,2018), or afforestation for ecological restoration, suchas in coastal Brazil (Rezende et al., 2018), or measuretocontrastlanddegradation,suchasinChina(Pengetal., 2014).Theseregionsoftheworldareamongthosethatshowthe maintrendsinagriculturaldeclinesandincreasesinforestareas (Fig.4a,band Fig.5a).Settlementisthesecondmaindriverof agriculturearealoss(29Mha,about23%ofthetotalagriculture loss),andmostlyoccurredinEasternChina,India,US,andWestern Europe (Fig. 5b). Agriculture expansion mainly happened at expensesofforest(123Mha),shrubland(30Mha)andgrassland (25.7Mha),andessentiallytookplaceinmanyworldregions,such asCentralandEasternAsia,US,SouthAmerica,andCentralAfrica (Fig.4a).Thelargestpresenceofforest-to-agriculturetransitions occurredatthebordersoftheAmazonbasinandinSoutheastAsia Fig.2. Normalizedtrendsforeachlandcoverfrom1992to2018.Differentcolors

indicatedifferentclasses.Trendsaresmoothedusingaﬁve-yearmovingaverage.(a) Normalizedtrendsforagriculture,forest,grassland,wetland,shrubland,sparse vegetation,andbareareas.(b)Normalizedtrendsforsettlement.AGR:Agriculture, FOR:Forest,GRA:Grassland,WET:Wetland,SET:Settlement,SHR:Shrubland,SPA:

Sparsevegetation,BAR:Barearea.

(6)

(Fig.5c).Globaldeforestationresultedinalossofforestareaequal to241Mha.Afteragriculture,themainforestlossesarecausedby shrubland (52Mha) andgrassland (27Mha),detected in South Americaandtheborealclimate(Fig.5d).Agriculture,shrubland, and grassland together are responsible for 84 % of forest loss.

Expansionofforestareasmainlycomefromagriculture(79Mha), shrubland(54Mha)andwetland(25Mha).Thesetransitionsare mostly due to the effects discussed above, such as cropland abandonmentandafforestationprogramsinthecaseofagricul- ture-to-forest (Fig. 5a), forest degradation or potential incon- sistenciesintheclassiﬁcationsystemfortheshrubland-to-forest, andglobalwarmingand/oranthropogenicactivitiesforwetland- to-forest(Fig.5e).Therearemultipleotherpossibleconsiderations forthechangesinlandcoversandtransitionsamongthedifferent classes. Grassland major gains are fromsparse vegetation (34 Mha),forest(27Mha)andbareareas(20Mha).Thelattermostly occurredinaridandsemi-aridareas(Fig.5g).Shrublandlargely increasedatexpensesofforests(52Mha),andsettlementmostly tookoveragriculturalland(29Mha,67%ofthetotalexpansion) (Fig.5b).Further,settlementappearstobeanon-reversibleland coverclassinthedatasetsconsideredinouranalysis,becauseno

conversionofsettlementtootherlandcoverclassesisregistered throughout the whole study period. Sparse vegetation mainly originated from bare areas (30 Mha) and grassland (20 Mha), whichcanbecorrelatedtoafeedbackofwarmerconditionsand higherCO2atmosphericconcentrationsthatextendthegrowing season and stimulate vegetation activities and tree growth, especiallyathigherlatitudes(Fig.5h).Similarly,sparsevegetation inborealclimatesarealsotransitioningtoforests,asaresultof generally improved conditions for tree growth and woody encroachment at forest-tundra ecotones (Fig. 5f). As expected, wetland losses (29 Mha) mainly occurred at high latitudes, especiallyin theNordicregionandRussia (Fig.5d),and mainly resultedinexpansionofforests.Overall,bareareadecreasesof57 Mha,withanetlossof29Mha.Barearealoses30Mhaand20Mha tosparsevegetationandgrassland,respectively(Figs.5gandh).

3.3.Globalvegetationcarbonlossesandsoilerosion

TheglobalbiomassdensitymapsfromSpawnetal.(2020)and SpawnandGibbs(2020)havethesamespatialresolutionofthe landcoverproductsfromESACCIandC3S-CDS,andtheycanbe Fig.3.Landcovertransitions,inmillionhectares(Mha).Trendsaresmoothedusingaﬁve-yearmovingaverage.(a)Agriculture/Foresttransitions.Thebluelineshowsthe transitionfromagriculturetoforest,andtheredlineshowsthetransitionfromforesttoagriculture.Theyellowlineisthenetgainorlossofforestinforest/agriculture transitions.(b)Forest/Grasslandtransitions.Thebluelineshowsthetransitionfromforesttograssland,andtheredlineshowsthetransitionfromgrasslandtoforest.The yellowlineisthenetgainorlossofgrasslandinforest/grasslandtransitions.(c)Forest/Wetlandtransitions.Thebluelineshowsthetransitionfromforesttowetland,andthe redlineshowsthetransitionfromwetlandtoforest.Theyellowlineisthenetgainorlossofwetlandinforest/wetlandtransitions.(d)Forest/Shrublandtransitions.Theblue lineshowsthetransitionfromforesttoshrubland,andtheredlineshowsthetransitionfromshrublandtoforest.Theyellowlineisthenetgainorlossofforestinforest/

shrublandtransitions.AGR:Agriculture,FOR:Forest,GRA:Grassland,WET:Wetland,SET:Settlement,SHR:Shrubland,SPA:Sparsevegetation,BAR:Barearea.

Table1

Globallandcoverchangesduringtheperiod1992-2018.Thenumbersinthetablerefertothetransitionsfromrowtocolumn.AGR:Agriculture,FOR:Forest,GRA:Grassland, WET:Wetland,SET:Settlement,SHR:Shrubland,SPA:Sparsevegetation,BAR:Barearea,WAT:water(Unit:Mha).

AGR FOR GRA WET SET SHR SPA BAR WAT Sum(loss)

AGR 0 79.0 9.2 0.2 29.1 4.0 1.7 1.2 1.8 126.2

FOR 123.7 0 26.6 14.4 4.3 52.1 9.0 2.6 8.9 241.5

GRA 25.7 19.4 0 0.3 4.7 1.5 20.1 6.9 1.0 79.5

WET 1.1 24.7 0.6 0 0.3 0.4 0.3 0.0 1.4 29.0

SET 0 0 0 0 0 0 0 0 0 0

SHR 30.1 54.4 3.8 0.2 2.2 0 1.2 0.5 0.5 93.0

SPA 18.7 13.6 33.6 0.3 0.7 0.7 0 12.4 0.6 80.7

BAR 4.7 0.4 19.6 0.0 2.0 0.1 29.7 0 0.9 57.4

WAT 1.4 4.4 0.8 2.3 0.4 0.5 0.7 4.6 0 15.0

Sum(gain) 205.2 196.0 94.3 17.7 43.7 59.3 62.7 28.2 15.1 722.2

Net(gain-loss) 79.1 45.5 14.8 11.3 43.7 33.7 18.0 29.2 0.2

(7)

integratedtoinvestigatethebiomasscarbonlossesaboveground (ABC) and belowground(BBC) due tothe recenthistorical land cover transitions. Inthis paper, wefocused onthegrossglobal losses due to the transitions from forest to agriculture and settlement.Between2010and2018,weestimatedthat9GtC(5.7/

10.7GtC,whenthe90th/10thpercentileisusedinEq.(8))arelost from aboveground vegetation clearance driven by agriculture expansion, and 0.23 GtC are lost from conversion of forest to settlements. The spatial distribution of the ABC loss due to agricultureexpansioninforestareasisshowninFig.6.Themajor lossesareinthetropicalband,mainlySouthAmerica,Africa,and Indonesia,andareestimatedat6.5 (4.0/8.1)GtC.Inthetropics, deforestation rates are high and trees typically store a larger amountofcarbonthaninotherbiomes.In termsofBBC,losses followthesamespatialpatternofABC(seeSupplementaryFig.S2), although quantities are smaller. A total of 2.7 (1.8/3.3) GtC is estimated ata global level,and 1.7 (1.1/2.1) GtC in thetropics.

SupplementaryTableS2showsspeciﬁcvalues.

AverageSEratesincreasedfrom5.0Mgha¹yr¹in2001to5.2 Mg ha¹ yr¹ in 2012 in the grid cells that wereclassified as agricultureinbothyears(AGR-AGR).TheSEratethusincreased 0.23 Mg ha¹ yr¹(Table2), avaluethatis morethan 5times highertherateoftheforestareasthatremainedforested(0.041Mg ha¹yr¹),i.e.,gridcellsthatwereclassifiedasforestbothin2001 and2012(FOR-FOR).Thetransitionfromforesttoagriculture(FOR- AGR,gridcellsthatwereclassifiedasforestin2001butagriculture in2012)acceleratedSEwithanaverageincreaseof0.71Mgha¹ yr¹,whichwas17timeshigherthantherateofgridcellswhere forestsremained forests.These changes in ratesof soil erosion resultedinatotalincreaseofSEof3.2Gtand1.1GtintheAGR-AGR andFOR-FORcase,respectively.Transitionfromforesttoagricul- turehascaused a total soilerosionlossof 0.14Gt. Despitethe smallerrates,thetotalSEislargerinFOR-FORthaninFOR-AGR becausetheformer occursoveramuch largerareaoftheglobe thanthelatter.TheaverageSErateafterafforestationofagriculture land(AGR-FOR,i.e.,gridcellsthatwereclassifiedasagriculturein 2001butforestin2012)is0.095Mgha¹yr¹,whichismuchlower thanthecasewhereagricultureremainedagriculture(AGR-AGR).

Theseresultsconﬁrmtheeffectivenessofforestestablishmentasa measure to contrast soil erosion. Large-scale afforestation

programshavebeenimplementedinthelast20yearstoprevent (oreven reverse) land degradation,such as those in the Loess PlateauinChina(Fengetal.,2020,2016).

GlobalSEcaused byagriculturalactivities(AGR-AGR)mainly occurredinthetropicalband,suchasBrazil,Sub-SaharanAfrica andSoutheastAsia(Fig.7).DecreasedratesinSEwerefoundin areascorrespondingtocountrieswithtransitioningoradvanced economies,suchasNorthAmerica,EuropeandEastChina(Fig.7), where agriculture conservation practices are more common (Borrelli et al., 2017). While declines of SEs in forests can be attributed to climatic variables (and to some extent to forest management),agricultural landareas are moreexposed tosoil erosionfromanthropogenicfactors(e.g.,tillage).Anexpansionof practices of conservationagriculture in developing countriesis seenasanoptiontomitigatemostofthemajornegativetrendsof SE,andultimatelypreventfurtherdeforestationfromexpansionof agriculturalareas(Borrellietal., 2020;IPCC,2019;Smithet al., 2020).

4.Discussion

This paper performed an extensive quantitative analysis of globallandcoverdynamicsusingthe27years’timeserieshigh resolution mapsfromtheESACCIand C3S-CDS products. This overviewdepictsthemajortrendsinchangesinlandcoversinour planet.Thespatialdistributionpatternsoflandcover dynamics shown in this studyare broadlyin linewith those previously reported. For example, Li et al. (2018) found net forest losses between1992and2015ofabout60Mha(45.5Mhainourstudy), mainlyoccurringinSouthAmerica,CentralAmericaandIndonesia.

A study that used satellite data to investigate global land transitionsbetween1982–2016forhighlyaggregatedlandcover classes,suchastallvegetation(i.e.,trees5minheight),short vegetationandbareground,founddeclinesinbaregroundof116 Mha(57Mhainourstudy,butforashortertimeperiod),increases in treecoverand reductionsin short vegetation(estimatesnot directly comparable with our study due to different land classiﬁcation systems) (Song et al., 2018). They also found clear regional patterns in land use changes, such as tropical deforestation,agriculturalexpansion, temperatereforestationor Fig.4. Globallandcoverchanges(increaseordecrease)forselectedlandcoverclassesasafractionofagridcell.Allthesubﬁguresareaggregatedto9kmspatialresolutionfor visualizationpurposesonlyandthenumbersinthecolorbarsindicatethefractionofthegridcellaffectedbythechange(positivevaluesindicateexpansion,negativevalues contraction).Landcoverchangeof(a)agriculture,(b)forest,(c)shrubland,(d)wetland.Notedifferentscalesoncolorbars.

(8)

afforestation,urbanization,treecoverexpansioninnorthernand montanesystems,andvegetationlossesinmanyaridandsemi- aridecosystems.

Previousstudieshavevalidatedandassessedtheaccuracyofthe landcoverdatausedinouranalysis(Lietal.,2018;Huaetal.,2018;

Liuetal.,2018a;Defournyetal.,2009).Therearelimitationsinthe ESACCIandC3SCCIlandcoverproducts,mainlyduetopotential misclassiﬁcations.Globaloverallaccuraciesofthetwolandcover productsare71%andvarywithspatiallocationandbetweenclasses.

Regionaloverallaccuraciesarefoundtobe70%inSouthAmerica (Pérez-Hoyosetal.,2017),62%inAfrica(Pérez-Hoyosetal.,2017),72

%inChina(Yangetal.,2017),84%incoastalEurasia(HouandHou, 2019),and64%intheArctic(Liangetal.,2019).Whiletheproducts haverelativelyhighglobalaccuraciesforcropland,forests,andurban classes, there are relatively larger uncertainties related to the mappingofmosaicnaturalvegetation,mossesandlichenclasses.

Thismeansthatouranalysisoffersmorereliabledatafortheland cover transitions involving forests,agriculture, and settlements,

whileithaslargeruncertaintiesforthelandclassesshrublandand sparsevegetation.Additionally,accuraciesofindividualclassesvary slightlyacrossthetwolandcoverproducts.Forexample, inthe ESACCIproduct,theurbanclassgloballyhasahigherthanaverage useraccuracy(88%)andalowerthanaverageproduceraccuracy(51

%),leadingtoaglobalunderestimationofurbanareas.However,this hasbeenimprovedinC3S-CDSwherebothurbanuserandproducer accuraciesareaboveaverage(75%),providingimprovedprecisionin totalareaextent.

Wefoundthatmosttransitionsoccurredbetweenforestand agriculture,andagricultureextensionwasthemaindrivingforce ofdeforestation.Agriculturehadthelargestnetgainsduringthe studyperiod(Table1),andthisismainlybecauseoftheexpansion ofagricultureactivitiesinBrazil,Africa,CentralAsia,EasternChina andSoutheastAsia(Fig.5).Thiscanindicateanincreasingdemand forfoodproductsfromdevelopingcountries,whichisassociated withhigherrisksofconversiontocroplandorpastureofforestor shrublandareas.Onthecontrary,agriculturelandisdecliningin Fig.5.Globalspatialdistributionofmainlandtransitions.Allthesubﬁguresareaggregatedto9kmspatialresolutionforvisualizationpurposesonlyandthenumbersinthe colorbarsindicatethefractionofthegridcellaffectedbythetransition.Transitionfromagriculturetoforest(a),fromagriculturetosettlement(b),fromforesttoagriculture (c),fromforesttoshrubland(d),fromwetlandtoforest(e),fromsparsevegetationtoforest(f),frombareareatograssland(g),frombarearetosparsevegetation(h).Note differentscalesoncolorbars.

(9)

otherlocations,suchasIndia,Europe,EastChinaand Southeast Brazil,where naturalrevegetation islikely ongoing(Mousivand andArsanjani,2019).Settlementexpandedalmostlinearlyduring the entirestudyperiod, mainlyatthecost ofagriculturalland.

Further,thecontinuingtrendsinlossesofforestsandshrublandfor agriculturalexpansionhashighrisksforirreversiblespeciesloss, especially inareaswithhighlyfragmentedhabitats,suchasthe biodiversityhotspots(Bettsetal.,2017;Barlowetal.,2016).Dueto climatewarming,wetlandathighlatitudeisdecreasing,anditis mainlyreplacedbyforests(Zhuetal.,2016).Wefoundgreening trendsforbareareas,withincreasingpresenceofgrassandsparse vegetationbecauseofCO2fertilizationandotherclimatefeedbacks (Zhuetal.,2016).Ontheotherhand,humaneffectsofconservation activitiesandforestrycanleadtolargescaleexpansionofforests, suchasinSouthernChinaandcoastalBrazil(Brandtetal.,2018).

Our analysis also showed an example of theimportance of changes inland coversfor keyglobal environmentalprocesses, suchassoilerosionandcarbonemissions,thankstotheintegration of the land cover products with other datasets of terrestrial ecosystemproperties.WeestimatedglobalannualABClossdueto deforestationfromexpansionofagricultureandsettlementof1.2 (0.7/1.4) GtC/yr, which is of 1.5 (1.0/1.8) GtC/yr when BBC is included.Inthetropics,valuesareof1.0(0.6/1.3)GtC/yrwhenboth ABCandBBCareconsidered.Thesearegrossestimatesofcarbon losses fromvegetation due tomajor drivers (i.e., expansion of agriculture and urban areas), which do not rely on speciﬁc assumptionsaboutthefateofconvertedforestlandorharvested woodproducts.Acomparisonwithexistingstudiesischallenging, especiallyatgloballevels,owingtodifferentapproaches,spatial scales,temporalperiods,andtypesofprocessesconsidered(Zarin, 2012;Ramankuttyetal.,2007).Harrisetal.(2012)estimatedan averagegrosscarbonlossof0.81GtC/yrusingsatelliteobserva- tionsofdeforestationinthetropicsfrom2000to2005,including bothabove-andbelow-groundbiomasscarbon.Theiroutcomes are25–50%ofrecentlypublishedestimates.Baccinietal.(2017) estimatedchangesofcarboninwoodyvegetationfromdeforesta- tionandforestdegradationinthetropicalregionandreportedan averagelossof0.86GtC/yrfrom2003to2014usingtheMODIS Fig.6.SpatialdistributionoftheABCloss(unit:tC)from2010to2018duetodeforestationtoagriculture.Themapisaggregatedto9kmspatialresolutionforvisualization purposes.

Table2

GlobalSEchangebetween2001and2012.AGR:Agriculture,FOR:Forest.

SE AGR-AGR FOR-FOR FOR-AGR AGR-FOR

Rate(Mgha¹yr¹) 0.23 0.041 0.71 0.095

Total(Unit:Gt) 3.2 1.1 0.14 0.018

Fig.7.Spatialpatternofdifferencesinglobalsoilerosionratesbetween2001and2012inagriculturalland(AGR-AGR).NegativevaluesindicateSEalleviation(i.e.,areduction inSErate),andpositivevaluesSEaggravation.AGR-AGRreferstothegridcellsclassiﬁedasagriculturebothin2001and2012(unit:t).Thismapisaggregatedto75kmspatial resolutionforvisualizationpurposesonly.

(10)

pantropical satellite data, and argued that emissions from degradationareprobablylargerthanthosefromdeforestation.

Ouranalysisonlyaccountedforcarbonlossesfromvegetation carbonresultingfromtransitionsfromforestcovertoagriculture and settlements. Most of the existing global studies produced estimatesthatincludeavarietyofotherprocesses,suchaschanges insoilorganiccarbon,forestgrowthordegradation,andexpansion ofnewforests,soreportingnetemissionsfromchangesinlanduse.

Forexample,HoughtonandNassikas(2017)estimatedaglobalnet carbonemissionfromlandusechangesduring2006–2015of1.1 GtC/yrusingthebookkeepingmodeltoprescribevegetationand soilcarbondensitychanges.LeQuéréetal.(2016)reportedaglobal carbonemissionof1.3GtC/yrduring2006–2015fromanensemble oflandsurfacemodels,andFriedlingsteinetal.(2019)provideda recentestimateof1.5GtC/yrfortheperiod20092018.Theyboth computed carbon ﬂuxes from land use changes and included estimates for bothvegetation (deforestationand forestgrowth/

expansion)andsoilcarbon.Alltheseanalysesusedthebookkeep- ingmodeltoprescribevegetationandsoilcarbondensitychanges andlandcoverdatafromtheFoodandAgricultureOrganizationof the United Nations (Keenan et al., 2015). In the bookkeeping method, there are many broad assumptions about the fate of clearedlandsandtheirrespectivecarbonstockstoestimatethe associated net carbon ﬂuxes, which include CO2 ﬂuxes from deforestation, afforestation, logging (forest degradation and harvest activity), shifting cultivation, and regrowth of forests followingwoodharvestorabandonmentofagriculture.

Ouranalysisdidnotattempttoestimatesoilemissions,butthe soilcarbonpoolismuchlargerthanthecarboninlivingbiomass.

However,theeffectsoflandusechangesonsoilcarbonemissions are difficult to quantify and still poorly understood (Guo and Gifford,2002;Donetal.,2011).Existingestimatestypicallyrelyon default factors that approximatesoil carbon fromaboveground biomassdensities(Panetal.,2011),butthesefactorshighlyvary acrossecosystemsandclimates(Piaoetal.,2009;Duarte-Guardia etal.,2019).Forexample,globalandregionalmeta-analysisoffield measurementsfindthatsoilorganiccarboncaneitherincrease, decrease,orremainconstantafterafforestation,becausethesoil carbonresponseissensitivetotheinteractionsofmultiplelocal factors,suchastreespecies,previouslanduse,initialsoilcarbon content,standage,climateandsoiltype(Hongetal.,2020).Our analysis alsotreated allcarbon lossesfrom forest clearance as instantaneousemissions,withoutconsideringpotentialstorageof carboninharvestedwoodproductsorforestremainingforestafter thedisturbance(i.e.,forestmanagement,whichiscommoninthe Northernhemisphere)(Iordanetal.,2018;Ceccherinietal.,2020;

Cherubini et al., 2016). Alternative emission estimates require assumptions about the fate of converted lands and biomass harvested,whichcurrentlycannotbeascertainedwithstatistical conﬁdenceacrosstheglobe.Ifarobustmethodwouldallowforthe considerationofthesedata,theresultwouldbeareductioninour estimatesofcarbonlosses.

Changes in land cover is the one of the main driversof SE (Borrelli etal.,2017),and a betterunderstandingof trendsand spatial patterns of SEfor different land uses can supportlocal policies aimingatmaintaining soilqualityand preventing land degradation. WefoundhighSErates onagricultureland inthe tropics, where currentagriculturalpractices areacceleratingSE due tounsustainablemanagementpractices,suchastillageand over-grazing (Borrelli et al., 2017; Montgomery, 2007). More conservationpracticesareneededinthesecountriestopreserve soillayers.Ourﬁndingsof5.2Mgha¹yr¹in2012ofSErateson agriculturalland are in linewithNearing et al. (2017a), which reported a value of global average SE rate on cropland with traditionaltillageof5.7Mgha¹yr¹.However,wefoundsmaller rates of SE than the ﬁndings of Borrelli et al. (2017), which

estimatedSEratesoncroplandtobemuchhigherthantheglobal averagesoilerosionrate.Thismightbeduetodifferentresolution ofthedataused(25kmvs.250m),andtopossibleinconsistencies inclassiﬁcationoflandusetypesacrossthelanduseandtheSE datasets.

Landcoverchangeis anindicatorforassessingthethreatto biodiversity,eventhoughtheimplicationsarehighlydependenton thelocalcontext(Uchidaand Ushimaru,2014;Strassburget al., 2020). Agriculture expansion on forest and other land covers challenges biodiversity conservation, especially in biodiversity hotspotsandintactlandscape(Bettsetal.,2017).Ourresultsreveal thatsigniﬁcantchangesoccurredatthebordersofpristineareas, suchastheAmazonrainforest,thetropicalforestsofCentral-east Africa,andtheBrazilianCerrado(shrubland),withhighpotential risksforspeciesloss(Huetal.,2020).Thistrendcallsforglobal attentiontomakethetoday’sagri-foodsectormoreefﬁcientandto adopt more sustainable diets, so to alleviate pressure from deforestation and prevent further habitat conversion (Folberth etal.,2020;McElweeetal.,2020).Otherwise,theexpectedrising demandsforfoodandfeedforanalwaysincreasingpopulationwill continue to drive expansion of agriculture land at expenses of naturalareas(IPCC,2019).

5.Conclusions

WecombinedtheESACCIlandcoverproductwiththenewly releasedC3S-CDSlandcoverproducttostudytheirrepresentation ofrecenthistoricallandcoverdynamicsbyretainingtheoriginal high spatial resolution (300 m) of the datasets. We further integratedtheseproductswithglobaldatasetsofbiomasscarbon densityandSEratetoanswertheresearchquestionsraisedinthe introductionsection.

Wefoundthat5.5%ofthetotalice-freelandsurfacetransitioned toadifferentlandcoverfrom1992to2018.Themajorexpansions are registered for settlements (+44 Mha, the largest relative increase)andagriculture(+79Mha),whiledeclinesareobserved forwetland(11Mha,thelargestrelativedecrease),forests(-45 Mha),shrubland(33Mha)andbareland(-29Mha).Declinesin forest areas mostly occurred before 1999, but they are still ongoing today (although net changes are small because of contrasting trends of forest expansion). The main transitions involved areexpansion of agriculturalland intoforestlandor shrubland, afforestation of agricultural land or shrubland, greeningof barelandandwetlanddryingtoforestsorsparse vegetation. Most of these transitions are the results of anthropogenicactivities,eitherdirectlythroughlanduse(i.e., deforestationforagriculturalexpansion)orindirectlyviawell- knownclimatechangefeedbacks(i.e.,vegetationenhancement orwetlanddryinginborealecosystems).

Clear spatial patterns emerged for the major historical land transitions. Forexample,agricultureexpansion atexpensesof forests or shrubland mainly happened at the borders of the Amazonbasin,inSoutheastAsia,EasternAsia,US,andCentral Africa.Afforestationismainlyfoundinassociationwitheither croplandabandonmentduetosocio-economicreasons,suchas inEasternEurope,orforecologicalrestorationtocontrastland degradation, such as in coastal Brazil and China. Settlement expansion mostly occurred in Eastern China, India, US, and Western Europe. Innorthern ormountainous ecosystems,we foundthatclimatefeedbacksarefavoringthetransitionofsparse vegetationtoforests,greeningofbareland,andwetlanddrying.

Wefoundimportantcontributionsfromrecentlanduseandland coverchanges tokeyglobalenvironmental processes,suchas grosscarbonlossesfromvegetationclearanceandeffectsonsoil erosion rates.Between 2010and 2018, we estimated that an