Desert dunes transformed by end-of-century changes in wind climate

Sand dunes in arid regions are significant mobile landforms which require adaptation and mitigation strategies for protecting human infrastructure and economic assets from encroachment, and which play a significant role in desertification and atmospheric dust emissions. Here we show how the shape, migration speed and direction of active desert dunes around the world are projected to change by the end of this century, in response to shifts in sand-moving wind regime associated with climate change under the SSP5-8.5 scenario. We find significant transformations in dune dynamics for many sand seas and dune fields across the Sahara, The Horn of Africa, the Southern Arabian Peninsula, South-Asia, China, and Australia – as well as increased potential of sand sea expansion and re-activation of dormant dune fields – which we can link to alterations in the Hadley circulation, extra-tropical cyclone activity, and monsoon systems resulting from global warming.

Desert dunes and sand seas cover approximately 20% of the world's arid zones 1 and their morphology and patterning are an important diagnostic of environmental surface conditions not only on Earth but also on other planetary bodies 2 . Crescentic shaped barchan dunes, for example, are indicative of a unidirectional sand-moving wind regime in low sediment supply conditions 3 , while the appearance of linear dunes on the surface of Titan has been interpreted as evidence of a bidirectional atmospheric flow regime shaping granular surface material in its equatorial region 4,5 .
The encroachment of moving desert dunes can pose significant threats to transportation infrastructure, agriculture, industry, and settlements [6][7][8][9] , either by direct burial under migrating dune bodies or by increased incoming wind-blown sand flux off the dunes, hence any change in the mobility and direction of migrating dunes can greatly alter the local risk assessment 10,11 . Migrating sand dunes can be agents of desertification, as their passage leaves long-term detrimental impacts on soil productivity and ecological richness 12,13 . They play an important role in dust emissions into the atmosphere at globally significant dust sources, where the blasting impact of sand flux accompanying mobile dunes releases silt and clay particles from dried-up lake beds into the air 14 , an important climate-change feedback process 15 . Understanding potential future changes in desert dune morphology, mobility, and migration direction due to changes in wind climate therefore has a great socio-economic relevance.
Changing wind climate also plays a key role in the potential expansion of dune fields and sand seas, as well as reactivation of currently dormant fields due to climate change 16 . An increase in wind power driving sand transport may trigger the transformation of a vegetation-controlled sand surface into mobile bare-sand dunes 17 , particularly if accompanied by overgrazing or increasing drought conditions 18 . The resulting increase in surface albedo from dark vegetation cover to bright bare sand can affect the regional radiation balance in a further climate-change feedback mechanism 19 .
Research on future changes in wind climate has predominantly been pursued in the context of the impact on power generation from current and planned wind turbine installations 20 . This work is, however, focussed on regions that are most relevant to the wind power industry, i.e. areas of existing and prospective turbine arrays, mainly North America, Europe, and East Asia, and with an emphasis on oceanic wind conditions for off-shore installations. These regional projections therefore have limited applicability to arid zone desert dunes and sand seas. Wind climate projections for the power industry are furthermore concerned only with scalar magnitudes and do not consider directional variability and change of prevailing wind vectors. Wind directional regime is, however, crucial to understanding and predicting the shape and migration of sand dunes.
A sand-moving wind regime can be quantified by the average scalar magnitude of Drift Potential (DP) vectors, proportional to the wind speed cubed above an initiation threshold, and the vector-resultant of the drift potentials, in magnitude (RDP) and direction (RDDir) 21 . DP reflects the overall capacity of the wind to shift sand, while RDP and RDDir quantify the net resultant displacement vector of sand. The ratio RDP/DP quantifies the directional variability of the sand drift regime, with values near unity indicating unidirectional regimes and small values indicating multiple wind vectors which together yield very little net movement of sand. DP, RDP, and RDDir are typically derived from meteorological records of wind speed and direction at a 10-minute temporal resolution spanning periods of 10 years or more to represent the wind climate 22 . Past comparisons between measured DP and observed sand sea activity 21,23 has suggested that low sand-drift environments typically show DP<27 m 3 s -3 , whereas high sand-drift environments are associated with DP>54 m 3 s -3 (units will be omitted henceforward, per convention).
In this study we analyse wind data from Coupled Model Intercomparison Project Phase 6 (CMIP6) climate simulations to determine projected changes, by the end of this century, in sand-moving wind regime parameters in the world's arid zones under the SSP5-8.5 Shared Socioeconomic Pathway emission scenario. We interpret the projected changes in different desert regions around the globe to infer potential increases as well as decreases in dune field activity, shifts in migration direction of mobile sand dunes, changes in dune shapes and pattern, and impacts on currently dormant dune fields.

Dune fields under current wind regimes
In order to qualify the associations between sand-moving wind regime and sand dunes, we examined contemporary satellite imagery of the Earth's surface in each cell of the CMIP6 climate simulation grid located in the arid zone (N=9498), to record the presence and areal extent of desert dunes and sand seas, their morphological types, and the associated inferred net sand drift direction. We find that active desert dunes are present across 23% of the arid zone land surface, with areal extent ranging from local pockets (<10% of a grid cell) to full coverage (>80%). Of the four major sand dune types, transverse-crescentic (T) dunes are found in 15%, linear-seif (L) dunes in 11%, barchan (B) dunes in 9%, and star (S) dunes in 4% of the arid zone. The bulk of active desert dunes are found in Africa (61%) and Asia (35%), with far fewer in the Americas (4%) and Australia (<1%). Just under half (46%) of desert dune areas are occupied by a single type of dune, while the other 54% sees co-occurrence of different types, most commonly transverse-crescentic dunes with barchans, linear, or star dunes.
Comparison of the global dune field inventory with sand-moving wind regime parameters extracted from the CMIP6 hind-cast simulations for the start of this century (2000)(2001)(2002)(2003)(2004)(2005)(2006)(2007)(2008)(2009)(2010)(2011)(2012)(2013)(2014)(2015) reveals that areas with active sand dunes experienced a median DP of 28, with extreme values of more than 90 at the 95 th percentile of the frequency distribution. While the traditional differentiation between low, moderate, and high sand-drift environments has previously been based on a subjective classification of 13 desert regions 21 , we propose boundaries based on both the DP frequency distribution as well as the physics of sand dune movement. We class low sand-drift environments where DP < 10, a threshold value that equates to the annual turn-over of the minimum size 24 (10 m long) dune as well as the 25 th percentile of the frequency distribution, whereas high sand-drift environments are classed as DP > 50, equating to the 75 th percentile of the distribution and associated with the yearly turnover of a 25 m long dune. The latter threshold is similar to the traditional boundary (DP > 54) and 25 m is close to the mode (most common) length of terrestrial barchan dunes 25 . For the first 15 years of this century areal mean Drift Potential in the low sand drift environment was 4, in the moderate drift environment 28, and 80 in the high drift environment.
The formation of the four major dune types (B, T, L, S) depends on two variables: the wind directional variability (quantified by RDP/DP), and the sediment availability (quantified as an Equivalent Sediment Thickness, EST) 3 . Our arid zone inventory of active desert dunes, mapped against the sand-drift regime parameters for 2000-2015, yields frequency distributions of RDP/DP for each type that suggest a separation between barchans and linear dunes -both known to form at low EST -and a separation between transverse-crescentic and star dunes -both known to form at high EST (Figure 1). Intersections of the percentile distributions indicate that at low EST the boundary between barchan dunes and linear dunes may be set at RDP/DP=0.725, with barchans more common above this value (more uni-directional sand drift) and linear dunes below this value. At high EST the cross-over between transverse-crescentic and star dunes is found at RDP/DP=0.435, with star dunes more common below this value (more multi-directional sand drift) and transverse-crescentic dunes at higher values. Of the arid zone active desert dune area, 29% experienced a multi-directional sand-drift regime (RDP/DP<0.435), 31% a bi-directional (0.435<RDP/DP<0.725), and 39% experienced a uni-directional regime (RDP/DP>0.725).
There is a good agreement between the Resultant Drift Direction (RDDir) obtained from the hind-cast simulations of the 2000-2015 period and the direction of dune migration that was inferred from the satellite imagery. For barchans (B) the median angular error between the two is 17 degrees (N=838 observations), for linear dunes (L) 19 degrees (N=1007), and for transverse dunes (T) 21 degrees (N=1241).

Future changes in sand-movement regime
A comparison of sand-moving wind regime parameters extracted from the 2000-2015 hind-cast simulations and from SSP5-8.5 emission scenario simulations for the end of this century (2085-2100) indicate that 73% of currently active desert dune areas are projected to experience a significantly different Drift Potential (p<0.05) under this scenario. Around 1/3 rd of desert dune area sees an increase in DP, the other 2/3 rd a decrease. Increase and decrease in DP beyond +10 and -10, respectively, are roughly balanced in areal extent (~7% each). Changes in Resultant Drift Potential show a less skewed ratio of 43% areal extent of increase and 57% areal extent of decrease, while the extent experiencing a rise in RDP greater than +10 is markedly higher (at 10%) than the extent experiencing a decline below -10 (2%).
In the aggregate, median DP across these areas declines by 12% (from 27.5 to 24.2), while RDP declines by only 3% (from 11.7 to 11.3). The spatial variability of these sand drift parameters greatly increases, however, showing an 11% (for DP) and 13% (for RDP) rise in standard deviation. The more extreme values of sand drift are evident in higher 95 th -percentiles, increasing from 88 to 95 for DP, and from 72 to 77 for RDP. Transitions between sand-drift environment classes are small, with more than 90% of desert dune areas remaining in their class.
Directional variability in sand drift sees more significant changes, with areas under multi-directional as well as uni-directional RDP/DP ratios shrinking, by 10% and 7% respectively, but areas with bi-directional regimes expanding by 19%. In three-quarters of active desert dune areas in moderate and high sand drift environments (i.e., DP>10) the resultant drift direction, RDDir, changes by less than 10 degrees, while 14% of these areas see a directional change of more than 15 degrees, with 8% seeing a backing RDDir (anticlockwise) and 6% a veering RDDir (clockwise).
The global distribution across the arid zone of changes in sand drift regime parameters by the end of this century (Figure 2) reveals a number of notable regions of interest. The most prominent feature is the increasing DP in a zonal band across the Central Sahara, the Horn of Africa, the Southern Arabian Peninsula, and South-east Pakistan & Rajasthan. Conspicuous decreases in Drift Potential are visible in Patagonia, the Mediterranean Maghreb, Iraq, around the boundary between Iran and Afghanistan, the Tibetan Plateau, Eastern Taklamakan and the Gobi Desert. The distribution of RDP changes in addition to those of DP shows significant increase in central Algeria & Libya, as well as moderate increases in central Western Australia. Large regions that experience little change or slight reductions in DP and RDP are North America, Southern Africa, and the Iran-Turkmenistan-Uzbekistan-Kazakhstan region. The consequences of combined changes in DP and RDP affecting the sand-drift directional variability (RDP/DP) are most notable in the Northern half of the Sahara turning more uni-directional, the Sahel turning more multi-directional, a strong contrast between Yemen & Oman turning more uni-directional while the Rub' al Khali turns more multi-directional, and large parts of the Australian interior experiencing a more uni-directional regime. Changes in Resultant Drift Direction for areas with DP>10 are more dispersed and localised, except for large parts of southern Australia where the direction is backing significantly (turning anti-clockwise). Table 1 lists currently active sand seas and desert dune fields around the globe that are particularly affected by changes in sand drift regime parameters.

Impacts on desert dunes
Changes in sand drift regime can have four main impacts on desert dunes. First, a change in RDP affects their migration rate, particularly for barchan and transverse dunes, being proportional to resultant drift potential and inversely proportional to dune height. Figure 2 shows the Bodélé Depression experiencing the globally largest increase in RDP by the end of this century. The region currently already sees extreme RDP on the order of 300, reflected by the world's fastest barchan dunes 26 with heights of 20-40 m migrating at ~50 m yr -1 , and the projected increase of 15-20% will see a commensurate speed-up (Figure 3a). As the ballistic impact of wind-blown sand is the key mechanism for dust emissions from the dried-up lakebed into the atmosphere and the Bodélé Depression is already the greatest single source of global atmospheric dust 27 , this contribution will likely increase significantly.
Second, a change in directional variability (RDP/DP) leads to dunes transforming in shape. Given their relatively small and self-contained size (hence short turn-over time), the most likely shape change is for barchans transforming into linear dunes under a wind regime with more directional variability (Figure 3b). Roughly 18% of currently uni-directional sand dune areas see a transition to a bi-directional regime, including 39 grid cells where barchans are currently present (5% of the areal extent of barchan dune fields). Nearly 30% of dune areas currently under a multi-directional regime are projected to transition to a less variable (bi-directional) regime. This includes 91 grid cells containing star dunes (23% of this dune type's areal extent), which therefore are potentially subject to transforming into transverse dunes, although the very large size of star dunes (generally several 100s of metres in diameter) means that changes in shape will be very slow.
Third, a change in Resultant Drift Direction (RDDir) alters the direction of migration or extension of dunes, particularly for barchans and linear dunes. Significant changes in drift direction of more than 15 degrees are projected for 119 grid cells (10%) experiencing moderate or high drift environments and currently containing barchans or linear dunes. Barchans can change direction as individual bedforms, while linear dunes may break up into rows of separate dunes if the drift direction shifts more perpendicular to the established crest line. Migrating dunes changing course can significantly impact human infrastructure and settlement developed around historic dune patterns, as illustrated in Figure 3c. Changes in RDDir can also undermine the effectiveness of existing sand drift mitigating measures, such as sand fences, which were optimised based on past sand drift regime parameters.
Fourth, significant future increases in DP or RDP may initiate re-activation of dormant dune fields which are currently held down by vegetation, especially in regions where climatic change leads to an increased soil moisture deficit. Of particular concern in this respect are western Australia and the Thar Desert. In the former, large areas of the Great Victoria and the Gibson Desert are projected to experience RDP increase on the order of 50%, while directional variability is set to reduce from RDP/DP values of around 0.35 (multidirectional) to around 0.60 (bi-directional). Combined with an expected reduction in soil moisture by the end of this century 28 , and the likely associated decrease in vegetation cover, the more focused Resultant Drift Potential may trigger a reactivation of currently dormant linear dunes. The Thar Desert, meanwhile, is expected to see an increase in DP & RDP on the order of 25% under a highly uni-directional wind regime (Figure 3d). While the area is projected to receive increased monsoon rainfall 29 , the relatively high population density and subsistence agriculture may expose currently vegetated long-walled parabolic dunes to reactivation as well as expansion of the existing transverse-crescentic dune fields under the stronger drift potential.

Discussion
The predicted changes in sand drift regime in the various regions around the globe are likely rooted in shifts in the global atmospheric Hadley Cell circulation and related monsoonal systems, altering the seasonality of sand moving winds from contrasting directions. The origins of these changes are linked to the projected future weakening and poleward expansion of the Hadley circulation due to global warming 30 , particularly in the Northern Hemisphere, together with changes in extra-tropical cyclone activity 31 and monsoon systems 32 . In Figure 4 we explore these regionally distinct impacts by comparing the projected changes in the sand rose for each typical month of the year for several select locations. Changes in the West African Monsoon are apparent in Southern Mauretania (Figure 4a), which sees significantly increasing sand drift under stronger monsoonal westerly winds in Jul-Aug-Sep, while the North-easterly Harmattan winds during the dry winter months remain largely unchanged, resulting in a significant shift in RDDir potentially threatening settlements in interdune areas. In Northern Algeria (Figure 4b), by contrast, the predicted decline in overall DP is due to significant decreases in the winter and spring North-westerlies (Nov-Apr), likely due to the projected weakening in extra-tropical cyclone activity over the Mediterranean 33 . The Northward expansion of the South-West Monsoon over the southern Arabian Peninsula is conspicuous in Oman (Figure 4c), showing mild decrease in drift potential under the North-westerly winds during the winter months (Dec-Jan-Feb), but with major change in sand drift regime during the middle of the year (Apr-Sept) from the current weak westerlies to strong South to South-westerly monsoon winds generating significantly higher drift potential toward the North-east, particularly in Jul-Aug. The resulting change in sand drift regime may alter the shape of local barchan dunes into linear dunes migrating in a North-easterly direction (Figure 3b). The poleward movement of the Southern Hemisphere storm track is apparent in the sand drift regime changes for the Great Victoria Desert in Western Australia (Figure 4d), which show a significant weakening of sand drift under westerlies during austral winter (May-Sep), conform with previous predictions 34 , while the austral summer period (Oct-Mar) sees a mild increase in sand drift under the South-easterly trade winds. This yields a higher RDP in a more uni-directional sand drift regime.
Our finding of a global 12% decline in median DP over active desert dune areas is commensurate with a recent prediction of average reduction in saturated sand flux by the end of this century for 45 sand seas 35 , while our projections of relatively small decreases in DP over most of central and western North America align well with regional models indicating a 5% decline in mean annual wind energy density for wind power generation 36 . The patterns of future changes in DP and RDP (Figure 2) furthermore match closely those predicted from a previous generation of climate model for northern Africa and Australia 17 . Our results however reveal the great spatial variability in predicted changes in sand drift regime parameters -both increases and decreases in activity -and the potential impacts on desert dune morphology and migration in arid zones around the world. Sand drift regime is one of the two key controls on dune dynamics, the other being sediment supply and in particular the role of vegetation restricting the movement of sand by wind. Climate change projections for most semi-arid zones generally predict increases in persistent drought conditions 37 , likely leading to a reduction in vegetation cover and increasing desert dune activity in tandem with changes in sand drift regime.

Calculation of sand drift regime parameters
For every day in each dataset a daily average Drift Potential vector was derived for each N216 grid cell from the four wind variables. The scalar magnitude of average Drift Potential is defined as: where U is surface wind speed at 10 m above the ground, and U t is its threshold for initiation of sand transport. We used a threshold of U t = 6.2 m s -1 , which is equivalent to the impact threshold shear velocity for D 50 = 0.3 mm median grain size quartz sand extrapolated via a logarithmic velocity profile from a Nikuradse roughness length z 0 = 5 D 50 /30 to the 10 m height wind speed. This same threshold is typically used in drift potential studies 21-23 . Because sand transport is a cubic function of wind speed involving a threshold, the temporal resolution of wind data has a critical impact on DP estimates, progressively under-estimating DP with longer wind speed averaging times due to the filtering out of high wind speed events that contribute disproportionally to sand drift. The standard time-scale for DP calculations is the 10-minute temporal resolution of WMO wind speed recordings, but CMIP6 climate models only provide 3-hourly wind speed data at best, leading to an underestimation of DP by 10-35% 22 . Instead, we calculated daily average DP from the reported daily mean wind speed and maximum wind speed by deriving a two-parameter Weibull PDF with a shape parameter, k, and a scale parameter, λ, where we equate the mean of the PDF to the daily mean wind speed, and we equate the (1-(1/96)) percentile of the corresponding CDF to the maximum wind speed, which is associated with an internal modelling time-span of 15 minutes (96 such time-spans in a day). The Weibull PDF is the statistical distribution that best represents meteorological wind speed data 41 . We take the thus derived PDF to draw 144 random samples of 10-minute wind speeds, averaged over 25 realisations, to calculate DP magnitude using equation 1 (the difference between Weibull PDFs at 15-minute time scale vs 10-minute time-scale is negligible).
The vector direction of the daily average Drift Potential was determined from the zonal and meridional daily averaged wind velocity reported for the grid cell, linearly interpolated from the Because the Drift Potential equation is the functional core of a traditional sand transport equation 42 , DP units can be converted to a sand flux using a suitable coefficient. For the standard 0.3 mm median grain size quartz sand: DP = 10 is equivalent to a saturated sand flux of 6.23 m 3 m -1 yr -1 (assuming a sand bulk density of 1600 kg m -3 ).

Comparison with ERA5 data
We compared DP values derived from the 2000-2015 hindcast simulations of the HadGEM3 model to DP values calculated from the ECMWF Re-analysis V5 (ERA5) surface (10 m height) wind data for the same period, the variables '10 metre U wind component' (zonal) and '10 metre V wind component' (meridional). ERA5 is a globally complete and consistent gridded reanalysis dataset with a high spatial resolution of ~30 km and temporal resolution of 1 hour 43 and may be considered the best global dataset of past weather conditions. HadGEM3 has previously been found to have small anomalies with ERA5's surface wind data compared to other CMIP6 models 44 . We used an approach for deriving daily average DP similar to our method for HadGEM3 data, except that we fitted a Weibull PDF through the 24 hourly scalar wind speed values from ERA5 and then re-sampled 144 wind speeds from that PDF to generate a 10-minute scale daily average DP. Note that wind speed data in ERA5 represent 'instantaneous' winds, not hourly averages, and so the PDF derived from these data does not filter out high wind speed events. The direction of the daily average DP vector was determined from the daily averaged U and V components.
Because the ERA5 grid is at a much higher spatial resolution than the N216 grid used in HadGEM3, the DP values at ERA5 grid-nodes were linearly interpolated and up-scaled to the areal-weighted average inside each N216 grid cell for comparison with the HadGEM3 results. The regression on DP values from ERA5 vs those from HadGEM3 for all N216 land surface grid cells returned a goodnessof-fit R 2 = 0.81 and a proportionality coefficient of 1.07, which we deemed sufficiently close to unity to avoid applying a scaling factor to the HadGEM3 results.

Dune field inventory
We used current satellite imagery in Google Earth Pro to visually assess the presence and estimate the areal extent of active desert dunes of four generalised types in each grid cell of the HadGEM3 model domain located over land surface in the arid zone between ±60° latitude (N=9498). The arid zone was defined as the B-zones in the Köppen-Geiger climate classification system, taken from Beck et al. 45 , and land surface grid cells were defined as containing 50% or more land fraction (as opposed to water). Active desert dunes were identified as having fully bare sand surfaces, visible crest lines and slip-faces. We used a first-order classification to distinguish between four main desert dune types that are typically recognized 3,46 : Barchan (B), Transverse-crescentic (T), Linear-seif (L), and Star (S) dunes, and we only considered the spatial scale of individual dune bedforms that can respond to the decadal sand drift regime, rather than large scale sand sea structures such as mega-dunes, draas, and transgressive ridges that evolve over centuries. For each dune type present in a grid cell we estimated its areal extent in four classes: localized pockets (<10%), sparse (~30%), partial (~50%), and full (>80%) coverage. Where possible, we inferred overall migration direction of the dunes, from slip-face evidence and spatial structure of dune fields, according to eight cardinal-ordinal compass directions (N, NE, E, …). We find active desert dunes in 2191 arid-zone grid cells, with Barchans present in 843, Transverse-crescentic in 1441, Linear-seif in 1008, and Star dunes in 401. Further statistics are reported in main text.

Sand roses
Sand roses shown in the figures are generated following the conventions of Fryberger & Dean 21 , based on the daily Drift Potential vectors (N=21600 for each grid cell). Each of 32 azimuth sectors (11.25° wide) shows the average daily drift potential from that direction, weighted by its timefraction. The sum of the weighted average daily drift potential over all sectors is equal to the overall DP for the whole 15-year period. The RDP vector is indicated by the red line. Note that drift potential per sector indicates the direction the sand drift is coming from (as per a conventional wind rose), while the RDP vector shows the resultant direction that sand drift is moving toward.

Supplementary Information
The supplementary information to this paper includes Google Earth layers of all drift parameter maps shown in Figure 2, as well as numerical labelling of these parameters and their changes in each HadGEM3 grid cell.

Data availability
All the data used in this study are publicly available at the following sites: HadGEM3-GC31-MM in the CMIP6 repository from the World Climate Research Programme: https://esgf-index1.ceda.ac.uk/search/cmip6-ceda/ ERA5 hourly data on single levels from 1979 to present from the Climate Data Store -Copernicus Climate Change Service: https://cds.climate.copernicus.eu/cdsapp#!/dataset/reanalysis-era5-single-levels?tab=overview