The Imprint of Vegetation on Desert Dune Dynamics

Here we demonstrate the qualitative and quantitative influence that vegetation has on stabilizing desert dunes. We use topographic data to isolate translation and deformation of dune patterns, upwind of and across a sharp gradient of vegetation, at White Sands dune field, New Mexico. Barchanoid dunes are unstable due to an aerodynamic surface wave instability. The dynamics of vegetated parabolic dunes are different; deformation becomes localized, and random, once plant density reaches a critical value associated with the barchanoid‐parabolic transition. Plants stabilize dunes not only by slowing them down but also by shutting off the fundamental mechanism that generates new sand waves and destabilizes dunes. Increasing plant density downwind increases vegetation‐induced form drag and results in decreasing dune migration rate. We suggest that similar biological modulation of pattern‐forming instabilities may also occur in other landscapes.


Introduction
Dunes in primarily unidirectional wind regimes act as downstream propagating waves of sand, but they are fundamentally unstable. Elevation perturbations along the crests of transverse dunes cause them to become sinuous, and eventually break up into barchanoids in sand-starved environments (Parteli et al., 2011). Aerodynamic perturbations generate new surface waves on the backs of dunes at a wavelength of order tens of meters (Elbelrhiti et al., 2005;Ping et al., 2014), causing barchanoids to emit smaller dunes from their horns. This counteracts collision-driven dune growth and leads to barchanoid fields that are stable in aggregate, even though individual dunes are unstable (Elbelrhiti et al., 2005(Elbelrhiti et al., , 2008Hersen & Douady, 2005;Hersen et al., 2004). Vegetation acts to slow dune migration, through three mechanisms: (i) reducing boundary shear stress b due to form drag (Durán & Herrmann, 2006b;Lancaster & Baas, 1998); (ii) binding and consolidation of sand by roots (Waldron, 1977) and (iii) facilitating the formation of sand-stabilizing soil crusts . In order for plants to gain a foothold, however, their growth rate must outpace the rate of burial by sand on the lee faces of dunes (Barchyn & Hugenholtz, 2015;Durán & Herrmann, 2006b;Jerolmack et al., 2012;Reitz et al., 2010) (Figures 1e-1g). This occurs at a threshold dune migration rate, below which the stabilizing effects of plants result in a positive feedback that further slows dune migration (Baas & Nield, 2007;Durán & Herrmann, 2006b;Pelletier et al., 2009). Stabilization typically occurs first at the low-elevation horns because they have smaller deposition rates than the center. This initial pinning of the dune horns, while the rest of the dune continues to migrate, results in barchanoids inverting their shape to a parabolic form (Durán & Herrmann, 2006b;Reitz et al., 2010) (Figures 1a and 1b-1d). Parabolic dunes may continue to evolve until vegetation densities reach a point that the dune is completely stabilized by plants.
The current understanding of how vegetation stabilizes dunes is essentially a kinematic picture; plants decrease erosion rate and ultimately pin the substrate in place . Here we posit that plants disrupt the fundamental dune-forming instability by perturbing the sand flux (q s ) field-effectively shutting down the formation of surface waves, and the emission of smaller dunes from the horns of barchanoids. To test this idea, we examined dune dynamics at the White Sands National Monument (New Mexico, USA), a well-studied unidirectional gypsum dune field (Jerolmack et al., 2012;Kocurek & Ewing, 2005;Kocurek et al., 2007;Langford, 2003;Reitz et al., 2010) that arises abruptly from a line source of sediment (downwind distance, x = 0). While not exact, a unidirectional approximation for both the dune migration direction and wind direction is reasonable. Wind speeds above threshold are strongly unimodal (from the Southwest) and the majority of dune migration occurs along these SW winds (Jerolmack et al., Figure 1. Landscape patterns at White Sands National Monument. (a) Aerial photograph of the dune field; black rectangle marks the region of this study. The three zones of distinct dune dynamics and morphology described in the text are indicated as Regions I-III in the figure. (b-d) Representative transverse, barchanoid, and parabolic dunes, respectively. Varied effects of vegetation across scales include (e) stabilization of soil through roots, (f) sediment deposition due to wake effects, and (g) growth of surface crusts.

Figure 2.
Definition sketch of deformation variables. Left "unshifted" shows a dune planform outline sampled at two different years (Year 1 is blue, and Year 2 is red for the entire figure), and below it is the elevation profile of the dune sampled along the green dashed transect. Right "shifted" shows the outlines and profiles shifted between Years 1 and 2 by the profile travel distance obtained during the calculation of the dune migration rate, V c . Vertical deformation is then measured profile by profile along the dune using ∏ . Planform deformation is computed using the metric D aff and D 2 min . For details see main text and the supporting information.
2012; Pedersen et al., 2015). Previous research has identified three regions associated with differing sediment transport and dune dynamics ( Figure 2). Region I lies between the sediment source and x ≈ 2 km, where initially transverse dunes break up into barchanoids and migration rate drops rapidly due to a topographically induced downwind decline in b (Jerolmack et al., 2012;Pelletier, 2015). Downwind in Region II, b approaches a stable value and there is a ∼5-km-long section where barchanoid dunes increase their spacing, but otherwise maintain a consistent size and migration rate (Baitis et al., 2014;Ewing & Kocurek, 2010). Beginning around the location x ≈ 7 km, the barchanoids rapidly become colonized by vegetation and invert their shape-over a scale of 1 km-to a parabolic morphology (Barchyn & Hugenholtz, 2015;Durán & Herrmann, 2006b;Jerolmack et al., 2012;Reitz et al., 2010) (Figures 1b-1d). Region III is associated with parabolic dunes that become increasingly vegetated, elongated, and slower moving as they progress downwind (Jerolmack et al., 2012;Pelletier, 2015).
The barchanoid-parabolic shape inversion at White Sands has been examined in detail (Reitz et al., 2010). Briefly, plants begin to colonize barchanoid dunes when the rates of vertical surface change (erosion/deposition) cross a threshold value, ∼0.2 m/yr, which corresponds to a dune migration rate V c ∼ 2-3 m/yr that occurs at x ≈ 7 km. This threshold crossing progresses from the dune horns to the center as the dune migrates, resulting in a planform pattern inversion because the horns' migration rate slows first. Reitz et al. (2010) used geometry and mass conservation to infer that there should be a linear decline in the dune migration rate with distance across the barchanoid-parabolic transition; they had insufficient data, however, to test this idea. These findings on the barchanoid-parabolic transition at White Sands provide the context for our new results on dune deformation presented here. Connecting these vegetation and dune-pattern transitions (Reitz et al., 2010), to dune deformation associated with the surface wave instability (Elbelrhiti et al., 2005;Ping et al., 2014), is the primary novel contribution of this paper. We also present refined kinematic and vegetation density data that improve upon previous observations.

Results and Discussion
We make use of repeat aerial LIDAR topographic data collected in September 2009 and June 2010, to construct digital elevation models (DEMs) with horizontal and vertical resolutions of 1 and 0.1 m, respectively (see section 4). The dunes at White Sands have characteristic lengths of ∼100 m, heights of several meters, and migration rates of several meters per year, so structure and dynamics are well resolved with this data set (Barchyn & Hugenholtz, 2015;Pelletier, 2015;Xia & Dong, 2016). We assume that the rates of dune migration for the 2009-2010 windy season are typical of the migration that would occur in any other given year. In addition, we use high-resolution (∼0.3 m per pixel) aerial photography from 2004, which is capable of resolving individual plants, to quantify the spatial distribution of vegetation across the study area ( Figure 1). An implicit assumption is that the vegetation remained in a quasi steady state over the course of 6 years between images and topography data. Since we interpret vegetation density averaged over the width of the dune field in this study-and not at the scale of individual plants or dunes-we consider this a reasonable assumption; if the vegetated boundary migrated at the rate of individual dunes, this would result in a maximum offset of 20 m relative to the 2004 location. Anecdotally, vegetation communities have not changed notably over our 10 years of field work in the dune field. The novelty of our approach is that we examine dune deformation-residual changes in both the profile and planform shape after removing net translation-to probe the dynamics of dune-pattern transitions and the influence of vegetation ( Figure 2). Profile deformation is quantified in the downwind (x) direction using a technique first developed for subaqueous bed forms (McElroy & Mohrig, 2009). Dunes are separated from the substrate, and the phase lag associated with maximizing correlation is used to determine a dune's migration rate from the two DEMs (see section 4). The residual elevation differences, ∏ (m/yr), represent deformation at each point on the dune: where is elevation and t indicates profiles separated by time t. We use the magnitude of ∏ throughout the paper instead of the rate itself. A different technique, though similar in spirit, is used to quantify deformation in the planform shape of dunes. A cloud of points that represent the footprint of a dune are tracked between the two surveys. The residual differences in the positions of these points, after subtracting the net translation of the dune (obtained from V c ), are used to compute a quantity called D 2 min (units of m 2 ) and a measure of homogeneous strain we call D aff (units of m 2 ) (Figure 2; see supporting Information). These metrics are commonly employed to study the deformation of granular materials (Falk & Langer, 1998;Griffa et al., 2011;Utter & Behringer, 2008). D 2 min quantifies the amount of non-affine deformation that a body under shear has undergone and is defined as follows: where N(r c ) is the number of dune sample points belonging to a cluster centered in r c , r ′ i is the current point being sampled, and r ′A i is the location that r ′ i would be in if all deformation occurred along the affine deformation matrix A (see supporting information for further explanation). D aff is computed from A in the course of finding D 2 min . In performing this analysis, we are essentially treating the dune field as a thin layer of amorphous (soft and deformable) materials on the field scale, and investigating its planar deformation behavior in response to wind shear and a changing environmental condition (vegetation). We first examine qualitative patterns in the maps of profile-deformation ( | | ∏ | | ) and planform deformation (D 2 min ), in the context of the dune-pattern transitions (transverse to barchanoid, and barchanoid to parabolic) at White Sands. The most striking result is the presence of laterally coherent stripes of large deformation rate in dune profiles, approximately 15-20 m in wavelength, which are pervasive on transverse and barchanoid dunes but nearly absent on parabolic dunes (Figures 3 and S11). Spectral analysis of the ( | | ∏ | | ) map confirms the presence of a characteristic wavelength ( 15.1 m) for deformation on unvegetated dunes (see Figure  S10) that is consistent with the surface wave instability proposed by Elbelrhiti (Elbelrhiti et al., 2005). While field evidence for this instability has been documented in the structure of dune profiles (Elbelrhiti et al., 2005;Ping et al., 2014), it has not to our knowledge been observed from dune dynamics. In the parabolic region these stripes largely disappear; not only is deformation rate lower, but it is also visibly less coherent (Figure 3). These findings provide a first confirmation that plants suppress the dune-forming instability. Turning to planform deformation, the most significant pattern is a zone of large D 2 min and D aff in Region I compared to the rest of the dune field. Large values of planform deformation are associated with the initial upwind margin of the dune field, where protodunes are combining to form transverse dunes. This initial peak is followed by a consistent decline in D 2 min that follows the breakup of transverse dunes into barchanoids (see Figure 4b). This is consistent with the planform instability proposed by Parteli (Parteli et al., 2011). There is also a second, more subtle change in D 2 min and D aff associated with changes in planform pattern across the barchanoid-parabolic transition. We suggest that the dynamics of the barchanoid-parabolic transition can be described as an abrupt shutdown of the surface wave instability and a gradual planform inversion, which are manifest as notable declines in profile deformation and modest increases in planform deformation.
To examine quantitative downwind trends in dune dynamics and vegetation, we compute width-averaged quantities of dune migration, affine deformation D aff , and vegetation density veg , which represents the area fraction of dunes that are covered in plants (Figure 4). Here width averaging refers to averaging in the cross-wind direction; we performed error analysis of all quantities involved to confirm that associated error is second order or smaller (see supporting information for more details). Dune migration data confirm previous findings that V c declines rapidly in the transverse Region I, is roughly constant in the barchanoid Region II, and then gradually decreases toward 0 in the parabolic Region III. The latter spatial decrease in dune migration rate is roughly linear and is mirrored by an approximately linear increase in vegetation density. Such a linear decrease in V c was proposed by Reitz et al. (2010) as a necessary condition for inversion of barchanoid dunes to a parabolic form. The inverse correlation between V c and veg also holds at the individual dune scale ( Figure S11). Because V c ∝ q s by mass conservation, and q s ∝ ( b − c ) where c is the threshold entrainment stress (Durán & Herrmann, 2006a;Martin & Kok, 2017), the decline in dune migration rate is due to either a decrease in the boundary stress or an increase in the threshold entrainment stress (or both). In the unvegetated dunes of White Sands, it has been shown that the dominant winds produce a boundary stress approximately 1.6 times the threshold stress, that is, b ∕ c ≈ 1.6 (Jerolmack et al., 2011). Previous work has shown an inverse relation between vegetation density and boundary shear stress due to the form drag effect (Durán & Herrmann, 2006b) where 0 is the boundary stress in the absence of plants, m = 0.16 is a model parameter, is the ratio of plant to surface drag coefficients, and is the ratio of basal to frontal area of a plant. Considering that plant density increases from veg ≈ 0 to veg ≈ 0.02 over the parabolic region, and utilizing parameter values found by Wyatt to be applicable to creosote desert shrubs ( = 200; = 1.45) (Wyatt & Nickiing, 1997), we would expect the excess boundary shear stress ( b − c ) to decrease by roughly a factor of 5. This is consistent with the magnitude and the trend of decreasing V c , suggesting that vegetation-induced form drag may exert the greatest influence on declining migration rates of parabolic dunes. The observation that the aerodynamic instability shuts down at low vegetation coverage ( veg ≈ 0.01) also supports the notion that form drag, rather than sand stabilization by roots, is the dominant effect of plants.
In order to quantify coherent deformation patterns ( Figure 3) and relate them to plant density, we compute a "coherent deformation density" cd . We perform blob detection on the | | ∏ | | map of the entire study area, compute the ratio of blob to dune area for each dune in the data set, and then perform width averaging to obtain a downwind profile of cd (see section 4). This quantity detects spatial coherence of dune profile deformation, even when magnitudes are low. There are two major features in the profile of cd : First, coherent deformation is highest at the upwind margin and gradually decreases downwind to the end of the barchanoid Region II; and second, there is a marked drop in coherent deformation at the barchanoid-parabolic transition, which corresponds to the location where vegetation density begins to increase. This provides quantitative support for the qualitatively distinct behaviors observed upwind and downwind of the barchanoid-parabolic transition ( Figure 3). The onset of vegetation in the parabolic region leads to a rapid loss of coherent deformation.

Conclusion
The current work has used the lens of deformation to examine two dune-pattern transitions at White Sands that are common in desert environments. This new approach demonstrates that a proposed planform instability (Parteli et al., 2011) is responsible for the breakup of transverse dunes into barchanoids. After this breakup, migration and deformation data indicate that the barchanoids evolve to a state of relative pattern stability that is maintained by the aerodynamic surface wave instability and dune collisions, as proposed previously (Elbelrhiti et al., 2005;Hersen & Douady, 2005;Hersen et al., 2004). The surface wave instability is detected as coherent stripes of profile deformation, which are completely disrupted with the onset of vegetation above a remarkably low areal density, approximately 1%. The disruption of coherent deformation is a local process that occurs in tandem with the previously well-documented morphological transitions associated with vegetation (see Figure 3). While plants may help bring about the barchanoid to parabolic transition by stabilizing dunes with their roots, our analysis at White Sands indicates that the dominant effect of plants is to disrupt and slow down the aerodynamics associated with a barchanoid morphology. In the parabolic region, deformation becomes smaller in magnitude and more localized, while migration rate gradually slows due to the increase in plant density and its associated form drag. While it is well known that increasing plant density results in decreasing erosion (Pelletier et al., 2009), the vegetation gradient at White Sands provides an unusually clear and quantitative demonstration of this effect. More fundamentally, this work shows how plants can drive a qualitative shift in the form and dynamics of a landscape pattern, by modulating an abiotic pattern-forming instability. A related example is the vegetation-induced transition from braided to single-thread river channels; by slowing the rate of bank erosion, plants prevent the formation of

Methods
DEMs of the study area outlined in Figure 1 were collected by the National Center for Airborne Laser Mapping in September 2009 and June 2010. For each DEM, 2-D elevation profiles were sampled along the width of the study area along the time-averaged wind direction. For each elevation profile, individual dune profiles were identified for both DEMs. A unique V c for each dune profile was then obtained by finding the spatial lag that maximized correlation between the sample years. The V c for every dune profile in the study area was then used to subtract the distance the profile had translated from 2009 to 2010 in the June 2010 DEM. This variably shifted DEM was then compared to the 2009 data to obtain the values used to quantify both the vertical and planform deformation, ∏ and D aff , respectively, and the coherent deformation density cd . The vegetation data used to calculate veg were produced from circa 2004 high-resolution aerial photography that covers the same area as the DEMs. For more information on methodology, including the error associated with the various quantities calculated in the paper, see supporting information.

Code Availability
All codes used to analyze the data are available online (https://github.com/dylanlee/WSScripts).