Implementing spatially explicit wind-driven seed and pollen dispersal in the individual-based larch simulation model: LAVESI-WIND 1.0
It is of major interest to estimate the feedback of arctic ecosystems to the global warming we expect in upcoming decades. The speed of this response is driven by the potential of species to migrate, tracking their climate optimum. For this, sessile plants have to produce and disperse seeds to newly available habitats, and pollination of ovules is needed for the seeds to be viable. These two processes are also the vectors that pass genetic information through a population. A restricted exchange among subpopulations might lead to a maladapted population due to diversity losses. Hence, a realistic implementation of these dispersal processes into a simulation model would allow an assessment of the importance of diversity for the migration of plant species in various environments worldwide. To date, dynamic global vegetation models have been optimized for a global application and overestimate the migration of biome shifts in currently warming temperatures. We hypothesize that this is caused by neglecting important fine-scale processes, which are necessary to estimate realistic vegetation trajectories. Recently, we built and parameterized a simulation model LAVESI for larches that dominate the latitudinal treelines in the northernmost areas of Siberia. In this study, we updated the vegetation model by including seed and pollen dispersal driven by wind speed and direction. The seed dispersal is modelled as a ballistic flight, and for the pollination of ovules of seeds produced, we implemented a wind-determined and distance-dependent probability distribution function using a von Mises distribution to select the pollen donor. A local sensitivity analysis of both processes supported the robustness of the model's results to the parameterization, although it highlighted the importance of recruitment and seed dispersal traits for migration rates. This individual-based and spatially explicit implementation of both dispersal processes makes it easily feasible to inherit plant traits and genetic information to assess the impact of migration processes on the genetics. Finally, we suggest how the final model can be applied to substantially help in unveiling the important drivers of migration dynamics and, with this, guide the improvement of recent global vegetation models.
How fast vegetation communities can follow their shifting climate envelope in a changing environment is determined by their ability to migrate. This is exceptionally challenging under current global change and plants might strongly lag behind their moving climate envelope (Harsch et al., 2009; Loarie et al., 2009; Moran and Clark, 2012). Temperatures are increasing most strongly in the Arctic. Accordingly, forests in the tundra–taiga transition zone are expected to respond by migration into the tundra (Bader, 2014; Holtmeier and Broll, 2005; MacDonald et al., 2008). However, empirical studies show diverse responses to the warming, including treelines being stable, advancing or even retreating (Harsch et al., 2009). A taiga range expansion though, might positively feedback to a global temperature increase due to albedo reduction (Bonan, 2008; Piao et al., 2007; Shuman et al., 2011).
To predict forest responses to climate, computer models were designed with different scopes of complexity, between highly general to very specific (Grimm and Railsback, 2005; Thuiller et al., 2008). Among these, simulation studies with dynamic global vegetation models (DGVMs) tend to overestimate the turnover of treeless tundra into forests (Brazhnik and Shugart, 2015, 2016; Frost and Epstein, 2014; Kaplan and New, 2006; Roberts and Hamann, 2016; Sitch et al., 2008; Snell, 2014; Yu et al., 2009; Zhang et al., 2013). On the other hand, forest landscape models (e.g. Snell et al., 2014; Shifley et al., 2017; Epstein et al., 2007) and small-scale models (forest-gap or individual-based) provide sufficient detail to realistically represent the responses at a stand level, but need a lot of effort for parameterization, have higher computational expenses, and are therefore typically not applied over large areas (Martínez et al., 2011; Pacala et al., 1996; Pacala and Deutschman, 1995; Zhang et al., 2011) or lack the implementation of wind-driven seed and pollen dispersal (e.g. Epstein et al., 2007). Further problems of DGVMs arise from the use of plant functional types as they consist of species with a wide variety of traits (e.g. Lee, 2011; Snell et al., 2014; Svenning et al., 2014). Nonetheless, the ability to form a closed canopy forest depends mainly on species traits acting at a fine-scale level such as (1) time needed to mature (life cycle, high generation time) and produce viable seeds, (2) dispersal distance and the chance for long-distance seed dispersal, and (3) germination and establishment of new individuals (Svenning et al., 2014). One source of the overestimation of migration rates of DGVMs is the unconstrained seed availability when climate variables allow a vegetation type to establish, which was recently pointed out by using a dispersal function between the grid points in simulations with a DGVM (Snell, 2014; Snell and Cowling, 2015). However, connecting grid cells to allow dispersal among them increases the computational complexity of such models (e.g. Nabel, 2015), but would be necessary to simulate realistic large-scale vegetation responses. In addition, the structure of a tree stand, and its response to changes in external forcing, is determined by further local processes, such as spatially explicit competition among individuals of all ages and their interactions. Of special interest is the adaptation of the traits of individuals of local populations, which are influenced by gene flow through seed or pollen distributions across populations. High exchange can lead to outcrossing that hinders local adaptation, but also prevents negative consequences from diversity losses caused by inbreeding within isolated populations due to founder effects in the process of colonization over large distances (Austerlitz et al., 1997; Burczyk et al., 2004; Fayard et al., 2009; Nishimura and Setoguchi, 2011; Ray and Excoffier, 2010). These processes have, so far, not been implemented continuously over a large scale in simulation models.
During the past decades treeline stands in the Siberian Arctic were densifying, but only rather slowly colonizing the tundra (Frost et al., 2014; Kharuk et al., 2006; Montesano et al., 2016), which could be attributed to seed limitation (Wieczorek et al., 2017). We developed the Larix vegetation simulator, LAVESI, to simulate tree stand dynamics at the Siberian treeline on the southern Taymyr Peninsula and use it as a framework to explore impacts of climate change on larch forests (Kruse et al., 2016). In the first version, the dispersal function randomly dispersed seeds by a probability density function describing a Gaussian term with a fat-tail. This could be parameterized to fit observed stand patterns. The model simulates tree stands on plots, representing a homogeneous forest, which can easily be enlarged to simulate wider areas. However, for simulations on larger transects passing from forests to treeless areas, wind direction, and strength become more important for seed dispersal and needed to be included in the model. Seed dispersal processes are well studied (Nathan et al., 2011a; Nathan and Muller-Landau, 2000) and are sometimes implemented in vegetation models but rarely coupled with wind speed and direction (e.g. Lee, 2011; Levin et al., 2003; Snell, 2014). Also wind patterns might change over time, as the pressure levels vary in a changing climate (Trenberth, 1990), or are directed (Lisitzin, 2012) so that an implementation of wind-dependent dispersal would enable a more realistic simulation of migration (cf. Nathan et al., 2011b).
The new spatially explicit pollination function tracks the full genealogy of a simulated tree stand and furthermore allows the inheritance of individually varying traits of each tree, rather than randomly drawing the actual trait value from the pool of available traits (cf. Scheiter et al., 2013). Additionally, the implementation of spatially explicit seed dispersal and pollination would enable us to align the model to detailed biogeographical knowledge gained from molecular methods (e.g. Navascués et al., 2010; Polezhaeva et al., 2010; Semerikov et al., 2007, 2013; Sjögren et al., 2017). We started with a very detailed small-scale model that can later be used to inform large-scale models especially about plot connectivity through seed dispersal and pollination and subsequent gene flow in landscapes.
We aim with this study to enable the simulation of spatially explicit and wind-dependent seed dispersal and pollination in the individual-based model LAVESI. After the coupling and verification of the seed dispersal kernel to prevailing winds and the incorporation of the pollination we test the model's sensitivity to its parameterization in local sensitivity analyses and the influence on stand development, migration rates, and pollination distances.
2.1 General model description of the Larix vegetation simulator LAVESI
LAVESI is an individual-based spatially explicit model that currently simulates the life cycle of larch species as completely as possible from seeds to mature trees (Kruse et al., 2016). It was set up to improve our understanding of past and future treeline displacements under changing climates, focusing on the open larch forest ecosystem in northern Siberia, which is underlain by permafrost. The relevant processes (growth, seed production and dispersal, establishment and mortality) are incorporated as submodules, which were parameterized on the basis of field evidence and complemented with data from literature. Simulation runs proceed in yearly time steps and are forced by monthly temperature and precipitation time series. The area simulated represents spatially homogeneous forest plots of variable size with the use of an environment grid (e.g. competition) with 20 cm tiles and where the handling of seeds dispersed beyond the plot borders can be set to deletion or reintroduction from the other side to simulate a forest patch. The model is programmed in C using standard template libraries. This and its modular structure allow a straightforward implementation of further extensions.
The model was successfully applied to conduct temperature-forcing experiments, where simulations revealed that the responses of the larch tree stands in Siberia – densification and northwards migration – could lag the applied hypothetical warming by several decades, until the end of the 21st century (Kruse et al., 2016; Wieczorek et al., 2017).
Here we present the implementation of wind-dependent seed dispersal as well as the newly introduced pollination. The absorbing boundary condition had to be revised to allow the simulation of larger areas. Hence, we introduce a new mode of periodic boundary conditions that allows seeds leaving the simulated area (100×100 m) to reenter on the opposite side, so that the borders of a simulation plot are connected along all borders. This mimics a tree stand within a homogeneous forest, similar to forest gap models (e.g. Brazhnik and Shugart, 2016; Pacala et al., 1996; Pacala and Deutschman, 1995; Zhang et al., 2011) and we used it in the simulations used for verification and parameterization for this manuscript. A second mode was implemented for simulations of hypothetical north–south transects (100×1000 m), which were used in the sensitivity analyses, allowing seed dispersal only on the meridional borders but not the latitudinal limits.
2.2 Implementing dispersal processes coupled to wind speed and direction
2.2.1 Pollination probability
Pollen was not represented in the former LAVESI version, but is needed to independently track gene flow by seeds and pollen through time. Accordingly, Figure 1 illustrates how we implemented an individual based pollination for each seed's ovule using a wind-determined and distance-dependent probability distribution function for pollen dispersal (similar to Gregory, 1961). It makes use of the von Mises distribution, which is an angular equivalent to the Gaussian normal distribution, for the two-dimensional representation (Abramowitz and Stegun, 2012).
A pollen dispersal function was newly implemented as a distance-dependent probability function for pollination of each individual seed's ovule, rather than simulating the large amount of pollen released by each tree (Gregory, 1961; Kuparinen et al., 2007). For each seed-bearing tree, the probability of pollen donating trees is calculated and out of the list of potential fathers for each seed one tree is randomly determined according to this probability. The pollination probability of each seed's ovule on a tree is proportional to the amount of pollen in the air column around it, which is, for simplification in the current implementation, not additionally dependent on the performance of the tree so that every tree that bears cones is taken into account. This aspect might be included in future versions. The following function is used here as the distance-dependent probability distribution of arriving pollen:
where r is the distance in m, pe is the ratio of pollen descending velocity Vd, pollen estimated for Larix gmelinii (Eisenhut, 1961) and wind speed Vw and Gregory's parameters C and m are set to C=0.6 cm and m=1.25 (p. 167 in Gregory, 1961).
The probability distribution pr described in Eq. (1) is multiplied by the von Mises distribution (Eq. 2), a continuous probability distribution on the circle, to include pollen distribution over a certain area and couple the process to the wind direction (illustration in Fig. 1; Abramowitz and Stegun, 2012).
where κ is the inverse of the von Mises distribution's variance, and I0(κ) is the modified Bessel function of order 0 as a function of κ, θ is the angle between trees and the actual wind direction. The modified Bessel function in the von Mises distribution is programmed in its integral representation using the Simpson integration scheme (Abramowitz and Stegun, 2012).
Consequently, following Gregory (1961) the pollination probability of a seed's ovule is as follows:
2.2.2 Seed dispersal
In the initial version of LAVESI, seeds are dispersed in random directions and at a distance r in m, estimated by a Gaussian and negative exponential (fat-tailed) dispersal function (Eq. 5, Kruse et al., 2016):
where E0, originally named “width”, is the Gaussian distribution's standard deviation in m, “rand” stands for a random number and “distance ratio” is a weighing factor for the fat tail in m2. Parameter estimates were based on a sensitivity analysis in Kruse et al. (2016) and numerical experiments.
The wind-dependent distance estimation was implemented as a ballistic flight following the assumptions of Matlack (1987). Accordingly, seed dispersal distances depend on the height of the releasing tree top Ht in m, currently estimated as 75 % of Ht (factor ), and are modified by wind speed VW in m s−1 and a species-specific fall speed of propagules (seed plus wing) Vd=0.86 m s −1 for L. gmelinii:
Finally, the direction for the seed dispersal is determined by wind direction, which was randomly selected from a set of observations (see Sect. 2.2.5 for details).
2.2.3 Parameterization to fit field data
The model's parameters had to be revised after implementing the model extensions to achieve simulated tree densities comparable to field data. Forest inventory data were recorded for each larch individual with explicit positions on plots of a minimum area of 20×20 m for several locations along a density gradient from single-tree stands in the north to dense forest tundra stands in the south visited on summer expeditions in the years 2011 and 2013 in northern central Siberia, Russia (Wieczorek et al., 2017). We conducted simulations on 100×100 m areas with closed boundaries initialized by introducing 1000 seeds in the first 100 years of a stabilization period of 1000 years, with forcing climate data randomly sampled from the available data. For the final 80 years of each simulation we used the climate series from the corresponding field site (TY04, see Sect. 2.2.4 for details). We visually compared the number of trees at year 2011 from the central 20×20 m area to the field survey data, which was the first year of fieldwork. The parameters were manually tuned and we iteratively performed simulation runs to improve the simulation results until finally achieving similar stand densities (numbers of trees) as observed (data not shown; parameter values in Table 1).
2.2.4 Temperature and precipitation
Simulations are forced with monthly mean temperature and precipitation sum series from the CRU TS 3.22 database (Harris et al., 2014). These are used to estimate long-term responses and derive the auxiliary climate variables active air temperature (sum of temperatures above 10 ∘C, AAT10) and vegetation length (number of days exceeding the freezing point, net degree days, NDD0) to calculate tree growth, estimate individual tree mortality and establishment from seeds (details in Kruse et al., 2016). We selected a grid box intersecting a location with a known northern taiga tree stand (CH06 at 70.66∘ N; 97.71∘ E, site CF in Wieczorek et al., 2017) and a northern forest tundra stand (TY04, 72.41∘ N; 105.45∘ E, site FTe in Wieczorek et al., 2017). From the available data we excluded years before 1934, because of missing climate station data and hence unreliable extrapolations in the data set (Mitchell et al., 2004). Furthermore, the final year was set to 2013, which is the latest year of fieldwork. The climate at these sites either allows strong tree growth with mean July temperatures of 13.50 ∘C, coldest temperatures during January of −33.24 ∘C and a precipitation sum of ∼328 mm year−1 or only sparse stands to emerge with temperatures of 13.11 and −36.07 ∘C in July and January, respectively, and ∼247 mm annual precipitation (cf. Kruse et al., 2016).
2.2.5 Wind speed and direction
The model is driven with pairs of wind speed in m s−1 and wind direction in degrees [∘]. The winds at 10 m above the surface for the years 1979–2012 at 6-hourly resolution were extracted from the ERA-Interim reanalysis data set (Fig. 3; Balsamo et al., 2015). Because of the coarse spatial resolution (80×80 km), we considered only the grid box over the climate station Khatanga, which is situated roughly in the centre of the treeline ecotone on the southern Taymyr Peninsula (71.9∘ N; 102.5∘ E; Wieczorek et al., 2017). During simulation runs, values are randomly drawn from the year's vegetation period (May–August; Abaimov, 2010) for each seed dispersal event and for the determination of pollination. For simulated years in which climate data are available but no corresponding wind data, a year is randomly selected.
2.3 Sensitivity analyses for dispersal processes
To test the influence of the parameterization of the variables from the newly introduced functions on the model's results, we ran local sensitivity analyses (Grimm and Railsback, 2005; Cariboni et al., 2007). For each simulation repeat, the input parameters (Table 2) were changed by 5 % and 50 % and a sensitivity value was calculated by comparing the results with the reference run:
where V is the variable of interest derived from each simulation run and P is the parameter of interest, both plus (+) and minus (−) 5 % of the estimated parameter, or with the reference value (Kruse et al., 2016).
The simulations were carried out on hypothetical north–south transects with a width of 100 m and length of 1000 m using the new model version and allowing seeds to be dispersed along the meridional borders. Populations were initiated on empty areas only in the lowermost 100 m wide and 100 m long area by randomly distributing 1000 seeds during the first 10 years of a 1000 year long stabilization period. During this phase, seeds exceeding the lowermost 100×100 m area were removed from the simulation. In the following simulation period seeds could enter the area above 100 m and colonize this empty area. The simulation model randomly drew weather conditions for each year from the complete available period 1934–2013 during the stabilization and simulation period. These simulations were repeated 30 times and the positions of each individual tree were recorded at the end of the simulation (500 years). To directly compare results from simulations with changed parameters to reference runs the simulation period was repeated for each parameter variation starting with an identical state of the simulation at the end of the stabilization period and using the same climate series.
For the evaluation of migration rates we selected three target output variables for the area ahead of the 100 m initialization area: (1) stemcount is the total number of stems (trees with a height above 130 cm), (2) forested area is the area covered with >100 stems ha−1, and (3) peak recruit position is the position of the maximum number of stems on the basis of a running mean with a 50 m window. Additionally, the variable stand density, which is the number of stems in the 20×20 m plot in the centre of the lowermost area, was selected to assess impacts on plot level. Furthermore, the pollination distance expressed as the mean distance between the pollen-donating and seed-producing trees was calculated for the evaluation of the pollination function. The resulting sensitivity values were tested for significant changes from the reference results (mean of 0) with a t test with a confidence level of 95 %.
2.4 Model-performance experiments
The memory load was estimated by adding up the size of all data types within each handled structure simulating a plot of one hectare (Table S1). These were multiplied by the actual number of elements in each of the structures. We calculated mean values of the number of handled items of the final 80 years of the simulations for the evaluation of dispersal processes to estimate the total memory needed for the arrays of trees and seeds and the grid representing the environment (Kruse et al., 2016).
To reduce the computation time, we parallelized the code for estimating pollination probabilities, seed dispersal, and tree density computation of the model using the OMP-library and conducted simulations using 1, 4, 8, and 16 CPUs. The performance of the model was evaluated by recording the computation time of each single simulation year for complete simulation runs (1080 years). We conducted four different runs, one with only wind dispersal of seeds (SEED), one with seed and pollination (+ POLL), and two different parallelized pollination computations. First, we tried to simply compute equally sized parts of the complete list of tree individuals including trees that have not produced seeds on the selected number of CPUs (+POLL_PAR-A). In a second variant (+POLL_PAR-B), we attempted to decrease the potential computational overheads of idle CPUs that had finished their job faster because of fewer individuals that needed to estimate pollination for produced seed's ovules, by cutting the list to only trees that produce seeds. The computation time increases with the actual number of trees and seeds present in simulations. In consequence, we analysed the dependency between the time needed for each simulated year and the number of trees and, additionally, the number of produced seeds by generalized nonparametric regression (using the “gam”-function in R-package “gam”; Hastie, 2017). The dependent variable time t was log-transformed prior to analysis. The explanatory variables – number of trees Nt and seeds Ns – were non-parametrically fitted and tested for non-linearity by comparing the deviance of a model that fits the terms linearly with a chi-squared test. In the initial model formula, we also included the interaction between the explanatory variables and excluded non-significant terms from the linear model (p>0.05) until yielding the final best model.
3.1 Verification of wind-dependency
The simulated seeds were solely dispersed in a north or south direction in coherence to the forcing winds (Fig. 2, Tables S2 and S3). The median seed dispersal distances were ∼12.2 m with a north wind and ∼12.0 m with a south wind with a majority of 95 % falling within ∼43 m of the seed tree, but with rare (∼0.1 %) dispersal events >1000 m (Fig. 2). The distance is equally highly correlated with the release height for both wind directions (ρ= 0.63, p<0.0001; Fig. 2).
The pollination events were mainly coming from the direction of the forcing winds: however, ∼18 % deviated from the forcing wind direction (Table S4). This variance is introduced by the formulae used for calculating the pollination probability for each seed's ovules on a tree and is further increased by the random selection of a father out of a subset of all possible mature trees based on the probability density function. The median distance along forcing winds of ∼38 m is, in general, shorter by ∼3–5 m than in other directions (Table S4).
* The integer variable maximum age of
seeds was excluded from the 5 % change sensitivity analysis as only
50 % changes had valid values. Highly significant with p<0.01;
Highly significant with p<0.05; remaining values are
non-significantly different from the reference run.
Values are the mean over 30 simulations.
In northern central Siberia, the main wind directions observed during the vegetation period are a combination of both west and east (Fig. 3a). In some years, one of these directions predominates, and is also characterized by stronger wind events. Accordingly, simulated seeds are dispersed into the general direction of the forcing wind data (Fig. 3b). Dispersal distances can reach up to a maximum of several thousand kilometres, yet the majority of seeds fall within a few hundreds of metres, and these are dispersed over distances depicting the wind speeds as well.
The median pollen flight distances are generally larger than the seed's, with a technically fixed maximum of about the distance from the central plot to the borders (Fig. 3c). Similar to seed dispersal, pollination follows the wind directions and fathers are positioned in the upwind direction of the main occurring winds.
3.2 Sensitivity analyses for implemented dispersal processes
The sensitivity analyses for the implemented seed dispersal function was extended for further model parameters that have an influence on the migration rate. In general, the four target variables have the same response direction towards changes in the parameters (Table 3). The stronger the changes, the more apparent becomes the change in the result so that the significance increases strongly from only 25 % to 79 %. The sensitivity values were of the same order of magnitude with the extreme values of −1.89 and 3.26 for each percent change in the input parameter. Most sensitive is the position of the peak recruitment for the observed migration rate (mean absolute sensitivity of 1.09 and 0.92 for 5 % and 50 %), whereas the impact on the stand level is of minor importance (with sensitivities of only 0.28 and 0.19). The factor of seed productivity fS and the influence factor of weather on germination rate fWeather Germination led to strongest advances of the peak recruit position if increased, which the seed mortality rate on trees PSeed Mortality, Cones caused when lowered.
The sensitivity values for resulting pollination distances for varied parameters were an absolute mean of change of 0.11 for 5 % and 0.02 for 50 % with extremes of −0.08 and 0.30 (Table 4). The stronger the change, the more apparent is a change in the results (40 % to 70 % significant values), although the direction of the changes was similar. However, the change is a magnitude smaller when changed by 50 % but the directions were consistent with those expected, and increasing Gregory's m led to farther pollination distances and vice versa for pollen descending velocity Vd, Pollen. The maximum sensitivities increased from −0.07 to 0.09 on the southernmost plot to about −0.11 to 0.14 for the northernmost plot where a higher proportion of significant values could also be observed.
3.3 Model performance
3.3.1 Memory consumption
The dynamic arrays need 120 bytes for each tree and 98 bytes for each seed. A further 54 bytes are needed for each of the environmental map tiles and another 117 bytes for the storage of output variables for each simulated year (Table S1). The constant containers use 390 bytes for the weather list and the parameter structures contain 642 bytes. On the basis of a simulated typical dense forest with ∼92 000 seeds and ∼25 000 tree individuals stored in the structures for each hectare, a simulation will need roughly ∼15 MB of RAM in a setup of a 1000 year initializing phase and a subsequent 80 year simulation phase.
3.3.2 Computation time
The simulation time increased with the number of trees in the simulation and for the contrasting simulation setups – either only wind-dependent seed dispersal SEED, or also with the calculation of pollination +POLL (Fig. 4). The generalized additive model, including the number of seeds, best explained the increase in computation time and had the lowest AIC value among all simulation types (Table S6). All incorporated variables, namely number of trees and number of seeds, significantly explained the computation time. The number of trees is the most important explanatory variable at ∼79.0 %, followed by an interaction term of the number of trees and seeds at ∼14.6 %, and number of seeds at ∼4.4 % and a residual of ∼2.0 % unexplained variation.
Without including the pollination events, the computation takes ∼0.6 s to calculate a year of a simulated plot on which 30 000 to 40 000 tree individuals are present (Fig. 4). In contrast, this increases to ∼120 s year−1 for a similar stand when calculating the pollen donor for each produced seed (+POLL). The first implemented parallelization of the pollination process (+POLL_PAR A) shortened the computation time by roughly half to ∼65 s year−1 when using eight cores. The second variant (+POLL_PAR B) outruns the first when using one to four cores by a factor of ∼4, but did not decrease the computation time significantly when using more cores than four. The increase to 16 CPUs led to a further decrease of computation time only for the first variant.
The assumption of unlimited seedbeds – allowing species in models to grow as soon as climate space allows them – causes high uncertainty in future predictions with dynamic global vegetation models (e.g. Midgley et al., 2007; Neilson et al., 2005; Sato and Ise, 2012). Implementing time-lagged responses in such models highlighted the need for a proper understanding and implementation of processes that limit species' migrations (Snell, 2014; Snell and Cowling, 2015). To reveal and understand the underlying processes that cause time lags, we designed the model LAVESI that represents all life-cycle stages of larches in high detail from seeds to mature trees, producing seeds themselves, which are then distributed in the environment (Kruse et al., 2016). We built this model to simulate responses of the Siberian treeline ecotone, which is solely covered over vast areas by a single tree species of the genus Larix. Here we describe the model enhancements to achieve, for the first time, a coupled implementation of wind-driven seed dispersal and pollination in the larch forest simulator LAVESI.
4.1 Wind-dependent seed dispersal
The simulated seed dispersal strictly followed the wind forcing and seeds settled in a downwind direction as expected, and not, as in the original model, in a purely ballistic manner (Kruse et al., 2016). We tested, in a local sensitivity analysis, the influence of different parameters on the stand level and the migration process. Sensitivity values were generally low with mean values between ∼0.2 and 1.1, respectively, for the stand level and for the migration rate. They are smaller compared to other parameters found in the sensitivity analysis of the first version of LAVESI by Kruse et al. (2016). In accordance with these findings, the new model is more sensitive to changes in parameters at transient stages such that higher values are found for the peak recruit position. Furthermore, only strong changes by at least 50 %, led in many cases to significant changes in the results, strengthening the robustness of the model to the parameterization. Those parameters leading to more available seeds and higher proportions of recruits (seed production rate, germination) had the highest sensitivity values and if increased they led to a faster migration and stand infilling. As expected, the parameters seed release height, wind speed, distance ratio, and fall speed of propagules became significant for the migration rate but not for the local stand development. Seed dispersal is dependent on the release height of the seed (Matlack, 1987), which is low in the focus region (Wieczorek et al., 2017) and thus leads to low dispersal distances, compared to other taxa (Pinus, González-Martínez et al., 2002, 2006; Picea, Piotti et al., 2009). When winds are constantly blowing at 2.78 m s−1 (10 km h−1), simulated distances seldom reach more than ∼43 m and only on very rare occasions are they observed with distances exceeding 1000 m. The dispersal kernel can thus be described as a combination of a Gaussian distribution, with its maximum fraction reaching ∼12 m from the releasing tree, and a long tail, best described by an exponential function. This aligns well with the implemented function in the model LAVESI (see details in Kruse et al., 2016). The model results are similar when driven with quasi-real wind data from the reanalysis data set ERA-Interim. The short seed dispersal distances depict well the generally observed values of other larch species. For example, Duncan (1954) found for Larix laricina in the northern USA that 94 % of seeds fell within 18 m of the releasing trees. Furthermore, Pluess (2011) found in dense forests of Larix decidua in the Swiss Alps an effective seed dispersal distance of 2–48 m. Moreover, the directions are now more realistically represented and follow the predominant west and east winds as expected (Fig. 3).
The use of winds from only the vegetation period might have introduced a bias, but it is based on the observation that this is the time when seeds are primarily dispersed (Abaimov, 2010). However, secondary dispersal by winds, due to uplift in strong winds, travel in winter on frozen surfaces over long distances (Nathan et al., 2011a; cf. Pluess, 2011), or due to wind-independent animal-mediated zoochory (Evstigneev et al., 2017), is currently not represented but could further facilitate the migration into tundra. When applying this model over historical periods, which are not covered by observations, one must be careful as the wind regimes could have shifted their main wind direction from the past to the current setting and might even change in the future (Lisitzin, 2012; Trenberth, 1990). A change, for example, from north–south to the current east–west wind directions could have limited the recent potential migration rate. This could explain the slow response of the treeline in northern Siberia to global warming, in addition to the long life cycle of larches, as well as prevailing seed limitation in the north (Kruse et al., 2016; Wieczorek et al., 2017).
4.2 Pollination coupled to prevailing wind conditions
Pollen dispersal functions are frequently used to reconstruct vegetation composition from palaeo-archives, for example in the Landscape Reconstruction Algorithm by Sugita et al. (2010), whereas other models have been used to track pollen clouds in tree stands (review in Jackson and Lyford, 1999; Prentice, 1985). Calculating every pollen dispersal event for each tree and seed is computationally challenging, but it can be simplified following the assumptions of Kuparinen et al. (2007). Hence, we implemented a density-dependent probability function and found in the sensitivity analysis that the pollination process was less affected by changing the input parameters than by the seed dispersal process. Values only reached a mean of ∼0.02 when changed by 50 % and increased from south to north, which covers a density gradient. Pollen influx from farther distances is more apparent in the more open stands, which is further supported by the findings in the sensitivity analysis of the original model (Kruse et al., 2016). The pollination distance increases linearly with Gregory's m, which increases the probability for farther standing trees, and decreases as expected for higher pollen descending velocities Vd, Pollen. The density-dependent probability function assigns pollen donors mostly in an upwind direction, but also has a small angular scattering, which was introduced by the use of the von Mises distribution to capture the stochasticity of this process (Gregory, 1961; Kuparinen et al., 2007). This uncertainty could lead to an overestimation of pollination distances, but this seems unlikely because simulated pollen are travelling distances of ∼38 m, which is longer compared to seeds and which is in concordance with observations (e.g. Pluess, 2011). However, the pollen amount and thus the probability of a distant tree to reach a seed-producing tree is dependent on the available resources and could further influence the resulting pollination distances. This relationship was not explicitly included but is partly covered by the use of the tree top height and from this, better performing trees have a higher pollination probability. The pollination distance often reaches the maximal possible distance between two trees in our simulation setup, which was the diagonal of the simulated area. These hypothetical pollen grains could probably reach distances of several hundred metres to kilometres, which would be in the range of the general dispersal distance observed in larches (Dow and Ashley, 1996; Hall, 1986).
4.3 Model performance
The individual-based approach of the model LAVESI-WIND, with the extension of wind-dependent seed dispersal and pollination, bears a high potential of knowledge gain, but this comes with some challenges: (1) repeated calculations for millions of individuals (seeds and trees) are computationally intense (e.g. Snell et al., 2014; Svenning et al., 2014; Nabel, 2015), and (2) they require a certain amount of memory during the simulation runs. Whereas the memory during each simulation run could be minimized to the needs of the simulation setup, the computational power was historically the limiting factor (Grimm and Railsback, 2005). But with the development of recent computer clusters with hundreds of CPUs, it seems very likely that one can overcome this, allowing us to use detailed and spatially explicit models at a regional scale (e.g. Paik et al., 2006; Zhao et al., 2013).
4.3.1 Memory consumption
We estimated the requirements for a hectare of a dense simulated forest as 15 MB of RAM. This means, on typical computer servers, even broad-scale simulation runs are easily feasible for 5000×5000 m, which would need ∼38 GB RAM and take approximately 40 h. The current LAVESI-WIND version was not fully optimized to lower the needs of memory and many variables that might not be needed for a specialized simulation experiment could be excluded. Although, the original simulation program was not intended to be run over continuous square kilometres of forests (Kruse et al., 2016), this is already possible with the current version. The programming language C and the process-based structure of the code support an easy and fast forward development of this model.
4.3.2 Computation time needed for millions of trees, seeds, and pollination
The computational effort of pollination for each seed's ovule increases with the number of mature trees present on a simulated plot. Therefore, to allow simulations to be run on standard computers in manageable time, it was a major goal to minimize the time needed for each simulated year. To meet this requirement, we parallelized parts of the program code that are computationally intensive, namely the processes of pollination and seed dispersal. With our approach, we have been able to decrease the time so far by a factor of 2 when using 8 CPUs, in comparison to using only one. Still, overheads from using a standard template library (STL)-list container lead to a negative exponential progression of the computation time needed per year rather than linear improvements (Fig. 4). Additional gains for other not yet parallelized processes are much smaller than these, but there is further potential to reduce the computation time by using different implementations of the parallelization.
4.4 Potential model applications
The new model version LAVESI-WIND allows for the evaluation of the importance of driving processes, which determine the response speed of tree stands growing at the treeline in Siberia. It can therefore be used for a very detailed evaluation of intra-stand processes determining migration speeds and help to improve abstract dynamic global vegetation models (e.g. Sato et al., 2007; Sitch et al., 2003), forest landscape models (e.g. Seidl et al., 2012), or regional forest gap models (e.g. Brazhnik and Shugart, 2015, 2016). Such a detailed representation of forest stands, as in the model presented here, is unlikely to be able to simulate forest dynamics on a continental to global scale (cf. Neilson et al., 2005). Nonetheless, the model can be used to parameterize dispersal kernels constraining inter-grid cell migration in DGVMs (Snell, 2014; Snell and Cowling, 2015). This could be achieved by comparing the migration rate in a continuous landscape in LAVESI-WIND, which covers grid cells of the DGVM to achieve a better representation of processes constraining or enhancing the spread of a plant species (cf. Lehsten et al., 2018). With this new model version, we can approach novel research questions, such as “do wind regime shifts explain faster or slower migration rates in past climate changes?” Furthermore, one could test how different treeline types determine the migration behaviour in changing environments. These can vary widely, based on the treeline type, being abrupt or with stand densities decreasing with the abiotic gradient and might further be influenced by shrubs that respond faster to current climate warming (e.g. Frost and Epstein, 2014), but which are not represented in the model yet. In addition, this may be influenced by single-tree stands growing ahead of the migration front (Holtmeier and Broll, 2005). Further interesting questions could be addressed, such as the role of refugia during past glacial periods and their influence on present-day tundra colonization by trees (Wagner et al., 2015), with a simplified and thus computationally effective approach. This is a necessary step because the current model version is computationally to demanding to track the full genealogy over simulated areas and time periods. Upscaling approaches could decrease generally the computation time and allow to expand the simulation over larger areas (e.g. Nabel, 2015; Epstein et al., 2007), however, the individual genetic information that passes thorough the landscape would be lost, which might be of interest. By connecting the borders of a simulation plot along the meridional borders we already implemented boundary conditions that allow the simulation of south-to-north transects, which are representative of the treeline area where highest tree densities occur in the south and treeless areas in the north. Thus, with this model, past migration corridors and timings can be revealed by a landscape-scale simulation, potentially answering important questions of the past biogeography of larch species in Siberia.
Before applying this new model version, however, a proper parameterization is necessary. Because pollen productivity and pollination distances as well as seed dispersal distances are not yet available for forests of the northernmost treeline area, the next important step would be to evaluate the modelled seed dispersal and pollination processes with field-based data, and finally, to apply this model to achieve realistic predictions of a future treeline. Molecular methods can help to improve the seed dispersal function, especially microsatellite markers, which can uncover connections among subpopulations and even kinships by parentage analyses at the stand level, which would make the effective seed dispersal distances directly inferable (Ashley, 2010; Dow and Ashley, 1996; Piotti et al., 2009; Pluess, 2011). Additionally, these methods can be used to estimate the fat tail of the dispersal function indirectly (Piotti et al., 2009).
Another interesting application would be to use this model to estimate the pollen influx in lakes (cf. Sugita, 2007). Pollen influx rates are widely used for vegetation reconstructions at the tundra–taiga transition zone (e.g. Klemm et al., 2016) and could now be used either to tune the dispersal parameters for a more precise population dynamics prediction, or inversely, to reconstruct ancient tree stands by simulations. Before the genetics or the influx rates are included in the model, however, a revision of boundary conditions for pollen in the model is necessary. This must include a relevant source area for the pollen (cf. Sugita, 2007) to determine to what extent genetic traits are delivered by pollen from beyond the borders of the simulated area. If this can be efficiently parameterized, the model could further be used to track genetic lineages in time.
We conclude that it is feasible to implement wind-driven seed dispersal and pollination in an individual-based model, which is then able to run across broader areas. However, the simulated area and duration of the simulation are constrained by available computer power and memory, and thus further effort is needed to minimize the computational load of this model in order to allow landscape-scale simulations on a standard computer. With the new model setup, further applications in combination with the genetics of the represented species are now feasible and can bring us detailed knowledge about the behaviour of the treeline and the biogeography of larch species through time.
The source code of the host model is available from GitHub, https://github.com/StefanKruse/LAVESI/releases/tag/v1.01 (last access: 26 October 2018), and stored in the zenodo database: https://doi.org/10.5281/zenodo.1155486 (Kruse and Wieczorek, 2016). The updated version presented here is named LAVESI-WIND and the first version 1.0 is accessible at GitHub at https://github.com/StefanKruse/LAVESI/tree/v1.0 (last access: 26 October 2018) and stored at https://doi.org/10.5281/zenodo.1165383 (Kruse et al., 2018).
The supplement related to this article is available online at: https://doi.org/10.5194/gmd-11-4451-2018-supplement.
SK and AG planned the study, NK, AG, and SK updated the model and implemented new functions. AG and SK performed the simulations and the statistical analysis. SK and AG wrote the manuscript. UH provided substantial advice in the process of data analysis and paper writing.
The authors declare that they have no conflict of interest.
We acknowledge Sven Willner for valuable advice in the process of
parallelizing the program code and Cathy Jenks for proofreading and improving
the manuscript. Furthermore, we particularly thank the handling topic editor
Hisashi Sato as well as Julia Nabel and one anonymous reviewer for their
valuable comments on the previous version of the manuscript. The position of
Stefan Kruse is funded by the Helmholtz Initiative Fund.
The article processing charges for this open-access
publication were covered by a Research
Centre of the Helmholtz Association.
Edited by: Hisashi Sato
Reviewed by: Julia Nabel and one anonymous referee
Abaimov, A. P.: Geographical distribution and genetics of Siberian larch species, in: Permafrost Ecosystems – Siberian Larch Forests, vol. 209, edited by: Osawa, A., Zyryanova, O. A., Matsuura, Y., Kajimoto, T., and Wein, R. W., Springer, Netherlands, Dordrecht, 41–58, 2010.
Abramowitz, M. and Stegun, I. A.: Handbook of mathematical functions: with formulas, graphs, and mathematical tables, Dover Books on Mathematics, Dover Publications, 2012.
Ackerly, D. D.: Community assembly, niche conservatism, and adaptive evolution in changing environments, Int. J. Plant Sci., 164, S165–S184, https://doi.org/10.1086/368401, 2003.
Ashley, M. V.: Plant parentage, pollination, and dispersal: How DNA microsatellites have altered the landscape, CRC. Crit. Rev. Plant Sci., 29, 148–161, https://doi.org/10.1080/07352689.2010.481167, 2010.
Austerlitz, F., Jung-Muller, B., Godelle, B., and Gouyon, P.-H.: Evolution of coalescence times, genetic diversity and structure during colonization, Theor. Popul. Biol., 51, 148–164, https://doi.org/10.1006/tpbi.1997.1302, 1997.
Bader, J.: The origin of regional Arctic warming, Nature, 509, 167–168, https://doi.org/10.1038/509167a, 2014.
Balsamo, G., Albergel, C., Beljaars, A., Boussetta, S., Brun, E., Cloke, H., Dee, D., Dutra, E., Muñoz-Sabater, J., Pappenberger, F., de Rosnay, P., Stockdale, T., and Vitart, F.: ERA-Interim/Land: a global land surface reanalysis data set, Hydrol. Earth Syst. Sci., 19, 389–407, https://doi.org/10.5194/hess-19-389-2015, 2015.
Ban, Y., Xu, H., Bergeron, Y., and Kneeshaw, D. D.: Gap regeneration of shade-intolerant Larix gmelini in old-growth boreal forests of northeastern China, J. Veg. Sci., 9, 529–536, https://doi.org/10.2307/3237268, 1998.
Bonan, G. B.: Forests and climate change: forcings, feedbacks, and the climate benefits of forests, Science, 320, 1444–1449, https://doi.org/10.1126/science.1155121, 2008.
Brazhnik, K. and Shugart, H. H.: 3D simulation of boreal forests: structure and dynamics in complex terrain and in a changing climate, Environ. Res. Lett., 10, 105006, https://doi.org/10.1088/1748-9326/10/10/105006, 2015.
Brazhnik, K. and Shugart, H. H.: SIBBORK: A new spatially-explicit gap model for boreal forest, Ecol. Model., 320, 182–196, https://doi.org/10.1016/j.ecolmodel.2015.09.016, 2016.
Burczyk, J., DiFazio, S. P., and Adams, W. T.: Gene flow in forest trees: How far do genes really travel?, For. Genet., 11, 179–192, https://doi.org/10.1016/j.foreco.2004.05.003, 2004.
Cariboni, J., Gatelli, D., Liska, R., and Saltelli, A.: The role of sensitivity analysis in ecological modelling, Ecol. Model., 203, 167–182, https://doi.org/10.1016/j.ecolmodel.2005.10.045, 2007.
Dow, B. D. and Ashley, M. V.: Microsatellite analysis of seed dispersal and parentage of saplings in bur oak, Quercus macrocarpa, Mol. Ecol., 5, 615–627, https://doi.org/10.1111/j.1365-294X.1996.tb00357.x, 1996.
Duncan, D. P.: A study of some of the factors affecting the natural regeneration of tamarack (Larix laricina) in Minnesota, Ecology, 35, 498–521, https://doi.org/10.2307/1931040, 1954.
Eisenhut, G.: Untersuchungen über die Morphologie und Ökologie der Pollenkörner heimischer und fremdländischer Waldbäume, P. Parey, Berlin, 1961.
Epstein, H. E., Yu, Q., Kaplan, J. O., and Lischke, H.: Simulating future changes in arctic and subarctic vegetation, Comput. Sci. Eng., 9, 12–23, https://doi.org/10.1109/MCSE.2007.84, 2007.
Evstigneev, O. I., Korotkov, V. N., Murashev, I. A., and Voevodin, P. V.: Zoochory and peculiarities of forest community formation: A review, Russ. J. Ecosyst. Ecol., 2, 1–16, https://doi.org/10.21685/2500-0578-2017-1-2, 2017.
Fayard, J., Klein, E. K., and Lefèvre, F.: Long distance dispersal and the fate of a gene from the colonization front, J. Evol. Biol., 22, 2171–2182, https://doi.org/10.1111/j.1420-9101.2009.01832.x, 2009.
Frost, G. V. and Epstein, H. E.: Tall shrub and tree expansion in Siberian tundra ecotones since the 1960s, Glob. Change Biol., 20, 1264–1277, https://doi.org/10.1111/gcb.12406, 2014.
Frost, G. V, Epstein, H. E., and Walker, D. A.: Regional and landscape-scale variability of Landsat-observed vegetation dynamics in northwest Siberian tundra, Environ. Res. Lett., 9, 1–11, https://doi.org/10.1088/1748-9326/9/2/025004, 2014.
Fyllas, N. M., Politi, P. I., Galanidis, A., Dimitrakopoulos, P. G., and Arianoutsou, M.: Simulating regeneration and vegetation dynamics in Mediterranean coniferous forests, Ecol. Model., 221, 1494–1504, https://doi.org/10.1016/j.ecolmodel.2010.03.003, 2010.
González-Martínez, S. C., Gerber, S., Cervera, M. T., Martínez-Zapater, J. M., Gil, L., and Alía, R.: Seed gene flow and fine-scale structure in a Mediterranean pine (Pinus pinaster Ait.) using nuclear microsatellite markers, Theor. Appl. Genet., 104, 1290–1297, https://doi.org/10.1007/s00122-002-0894-4, 2002.
González-Martínez, S. C., Burczyk, J., Nathan, R., Nanos, N., Gil, L., and Alía, R.: Effective gene dispersal and female reproductive success in Mediterranean maritime pine (Pinus pinaster Aiton), Mol. Ecol., 15, 4577–4588, https://doi.org/10.1111/j.1365-294X.2006.03118.x, 2006.
Gregory, P. H.: The microbiology of the atmosphere, Leonard Hill, London & Interscience Publishers, New York, 1961.
Grimm, V. and Railsback, S. F.: Individual-based Modeling and Ecology, Princeton University Press, Princeton, NJ, USA, 2005.
Hall, J. P.: Growth and development of larch in Newfoundland, in Proceedings of the Sixth International Workshop on Forest Regeneration, U.S. Department of Agriculture, Forest Service, General Technical Reports, 21–26, 1986.
Harsch, M. A., Hulme, P. E., McGlone, M. S., and Duncan, R. P.: Are treelines advancing? A global meta-analysis of treeline response to climate warming, Ecol. Lett., 12, 1040–1049, https://doi.org/10.1111/j.1461-0248.2009.01355.x, 2009.
Hastie, T.: gam: Generalized Additive Models, R Package, version 1.14-3, available at: http://cran.r-project.org/package=gam (last access: 18 January 2018), 2017.
Holtmeier, F.-K. and Broll, G.: Sensitivity and response of Northern Hemisphere altitudinal and polar treelines to environmental change at landscape and local scales, Glob. Ecol. Biogeogr., 14, 395–410, https://doi.org/10.1111/j.1466-822x.2005.00168.x, 2005.
Jackson, S. T. and Lyford, M. E.: Pollen dispersal models in Quaternary plant ecology: Assumptions, parameters, and prescriptions, Bot. Rev., 65, 39–75, https://doi.org/10.1007/BF02856557, 1999.
Kaplan, J. O. and New, M.: Arctic climate change with a 2 ∘C global warming: Timing, climate patterns and vegetation change, Clim. Change, 79, 213–241, https://doi.org/10.1007/s10584-006-9113-7, 2006.
Kharuk, V. I., Ranson, K. J., Im, S. T., and Naurzbaev, M. M.: Forest-tundra larch forests and climatic trends, Russ. J. Ecol., 37, 291–298, https://doi.org/10.1134/S1067413606050018, 2006.
Klemm, J., Herzschuh, U., and Pestryakova, L. A.: Vegetation, climate and lake changes over the last 7000 years at the boreal treeline in north-central Siberia, Quat. Sci. Rev., 147, 422–434, https://doi.org/10.1016/j.quascirev.2015.08.015, 2016.
Kruse, S. and Wieczorek, M.: LAVESI – Larix Vegetation Simulator (Version 1.01), Zenodo, https://doi.org/10.5281/zenodo.1155486, 2016.
Kruse, S., Wieczorek, M., Jeltsch, F., and Herzschuh, U.: Treeline dynamics in Siberia under changing climates as inferred from an individual-based model for Larix, Ecol. Model., 338, 101–121, https://doi.org/10.1016/j.ecolmodel.2016.08.003, 2016.
Kruse, S., Gerdes, A., and Kath, J. N.: LAVESI-WIND 1.0 (Version v1.0), Zenodo, https://doi.org/10.5281/zenodo.1165383, 2018.
Kuparinen, A., Markkanen, T., Riikonen, H., and Vesala, T.: Modeling air-mediated dispersal of spores, pollen and seeds in forested areas, Ecol. Model., 208, 177–188, https://doi.org/10.1016/j.ecolmodel.2007.05.023, 2007.
Lee, E.: Impacts of Meteorology-Driven Seed Dispersal on Plant Migration?: Implications for Future Vegetation Structure under Changing Climates, Massachusetts Institute of Technology. Doctoral dissertation, available at: https://globalchange.mit.edu/sites/default/files/Lee_PhD_2011.pdf (last access: 27 May 2018), 2011.
Lehsten, V., Mischurow, M., Lindström, E., Lehsten, D., and Lischke, H.: Simulating migration in dynamic vegetation models efficiently with LPJ-GM, Geosci. Model Dev. Discuss., https://doi.org/10.5194/gmd-2018-161, in review, 2018.
Levin, S. A., Muller-Landau, H. C., Nathan, R., and Chave, J.: The ecology and evolution of seed dispersal: A theoretical perspective, Annu. Rev. Ecol. Evol. Syst., 34, 575–604, https://doi.org/10.1146/annurev.ecolsys.34.011802.132428, 2003.
Lisitzin, A. P.: Sea-ice and iceberg sedimentation in the ocean: recent and past, Springer Science & Business Media, 2012.
Loarie, S. R., Duffy, P. B., Hamilton, H., Asner, G. P., Field, C. B., and Ackerly, D. D.: The velocity of climate change, Nature, 462, 1052–1055, https://doi.org/10.1038/nature08649, 2009.
MacDonald, G. M., Kremenetski, K. V., and Beilman, D. W.: Climate change and the northern Russian treeline zone, Philos. Trans. R. Soc. B Biol. Sci., 363, 2283–2299, https://doi.org/10.1098/rstb.2007.2200, 2008.
Martínez, I., Wiegand, T., Camarero, J. J., Batllori, E., and Gutiérrez, E.: Disentangling the formation of contrasting tree-line physiognomies combining model selection and Bayesian parameterization for simulation models, Am. Nat., 177, E136–E152, https://doi.org/10.1086/659623, 2011.
Matlack, G. R.: Diaspore size, shape, and fall behavior in wind-dispersed plant species, Am. J. Bot., 74, 1150–1160, https://doi.org/10.2307/2444151, 1987.
Midgley, G. F., Thuiller, W., and Higgins, S. I.: Plant Species Migration as a Key Uncertainty in Predicting Future Impacts of Climate Change on Ecosystems: Progress and Challenges, in: Terrestrial Ecosystems in a Changing World, edited by: Canadell, J. G., Pataki, D. E., and Pitelka, L. F., Springer, Berlin, Heidelberg, 129–137, 2007.
Mitchell, T. D., Carter, T. R., Jones, P. D., Hulme, M., and New, M.: A comprehensive set of high-resolution grids of monthly climate for Europe and the globe: the observed record (1901–2000) and 16 scenarios (2001–2100), Tyndall Centre for Climate Change Research, Working Paper 55, 2004.
Montesano, P. M., Sun, G., Dubayah, R. O., and Ranson, K. J.: Spaceborne potential for examining taiga–tundra ecotone form and vulnerability, Biogeosciences, 13, 3847–3861, https://doi.org/10.5194/bg-13-3847-2016, 2016.
Moran, E. V. and Clark, J. S.: Between-site differences in the scale of dispersal and gene flow in red oak, PLoS One, 7, e36492, https://doi.org/10.1371/journal.pone.0036492, 2012.
Nabel, J. E. M. S.: Upscaling with the dynamic two-layer classification concept (D2C): TreeMig-2L, an efficient implementation of the forest-landscape model TreeMig, Geosci. Model Dev., 8, 3563–3577, https://doi.org/10.5194/gmd-8-3563-2015, 2015.
Nathan, R. and Muller-Landau, H. C.: Spatial patterns of seed dispersal, their determinants and consequences for recruitment, Trends Ecol. Evol., 15, 278–285, https://doi.org/10.1016/S0169-5347(00)01874-7, 2000.
Nathan, R., Katul, G. G., Bohrer, G., Kuparinen, A., Soons, M. B., Thompson, S. E., Trakhtenbrot, A., and Horn, H. S.: Mechanistic models of seed dispersal by wind, Theor. Ecol., 4, 113–132, https://doi.org/10.1007/s12080-011-0115-3, 2011a.
Nathan, R., Horvitz, N., He, Y., Kuparinen, A., Schurr, F. M., and Katul, G. G.: Spread of North American wind-dispersed trees in future environments, Ecol. Lett., 14, 211–219, https://doi.org/10.1111/j.1461-0248.2010.01573.x, 2011b.
Navascués, M., Depaulis, F., and Emerson, B. C.: Combining contemporary and ancient DNA in population genetic and phylogeographical studies, Mol. Ecol. Resour., 10, 760–772, https://doi.org/10.1111/j.1755-0998.2010.02895.x, 2010.
Neilson, R. P., Pitelka, L. F., Solomon, A. M., Nathan, R., Midgley, G. F., Fragoso, J. M. V., Lischke, H., and Thompson, K.: Forecasting regional to global plant migration in response to climate change, Bioscience, 55, 749–759, https://doi.org/10.1641/0006-3568(2005)055[0749:FRTGPM]2.0.CO;2, 2005.
Nishimura, M. and Setoguchi, H.: Homogeneous genetic structure and variation in tree architecture of Larix kaempferi along altitudinal gradients on Mt. Fuji, J. Plant Res., 124, 253–263, https://doi.org/10.1007/s10265-010-0370-1, 2011.
Pacala, S. W. and Deutschman, D. H.: Details that matter: The spatial distribution of individual trees maintains forest ecosystem function, Oikos, 74, 357–365, https://doi.org/10.2307/3545980, 1995.
Pacala, S. W., Canham, C. D., Saponara, J., John, A. S. J., Kobe, R. K., and Ribbens, E.: Forest models defined by field measurements: estimation, error analysis and dynamics, Ecol. Monogr., 66, 1–43, 1996.
Paik, S. H., Moon, J. J., Kim, S. J., and Lee, M.: Parallel performance of large scale impact simulations on Linux cluster super computer, Comput. Struct., 84, 732–741, https://doi.org/10.1016/j.compstruc.2005.11.013, 2006.
Piao, S., Friedlingstein, P., Ciais, P., Viovy, N., and Demarty, J.: Growing season extension and its impact on terrestrial carbon cycle in the Northern Hemisphere over the past 2 decades, Global Biogeochem. Cy., 21, GB3018, https://doi.org/10.1029/2006GB002888, 2007.
Piotti, A., Leonardi, S., Piovani, P., Scalfi, M., and Menozzi, P.: Spruce colonization at treeline: where do those seeds come from?, Heredity, 103, 136–45, https://doi.org/10.1038/hdy.2009.42, 2009.
Pluess, A. R.: Pursuing glacier retreat: genetic structure of a rapidly expanding Larix decidua population, Mol. Ecol., 20, 473–485, https://doi.org/10.1111/j.1365-294X.2010.04972.x, 2011.
Polezhaeva, M. A., Lascoux, M., and Semerikov, V. L.: Cytoplasmic DNA variation and biogeography of Larix Mill. in northeast Asia, Mol. Ecol., 19, 1239–1252, https://doi.org/10.1111/j.1365-294X.2010.04552.x, 2010.
Prentice, I. C.: Pollen representation, source area, and basin size: Toward a unified theory of pollen analysis, Quat. Res., 23, 76–86, https://doi.org/10.1016/0033-5894(85)90073-0, 1985.
Ray, N. and Excoffier, L.: A first step towards inferring levels of long-distance dispersal during past expansions, Mol. Ecol. Resour., 10, 902–914, https://doi.org/10.1111/j.1755-0998.2010.02881.x, 2010.
Roberts, D. R. and Hamann, A.: Climate refugia and migration requirements in complex landscapes, Ecography, 39, 1238–1246, https://doi.org/10.1111/ecog.01998, 2016.
Ronce, O.: How does it feel to be like a rolling stone? Ten questions about dispersal evolution, Annu. Rev. Ecol. Evol. Syst., 38, 231–253, https://doi.org/10.1146/annurev.ecolsys.38.091206.095611, 2007.
Sato, H. and Ise, T.: Effect of plant dynamic processes on African vegetation responses to climate change: Analysis using the spatially explicit individual-based dynamic global vegetation model (SEIB-DGVM), J. Geophys. Res.-Biogeo., 117, 002056, https://doi.org/10.1029/2012JG002056, 2012.
Sato, H., Itoh, A., and Kohyama, T.: SEIB–DGVM: A new dynamic global vegetation model using a spatially explicit individual-based approach, Ecol. Model., 200, 279–307, https://doi.org/10.1016/j.ecolmodel.2006.09.006, 2007.
Scheiter, S., Langan, L., and Higgins, S. I.: Next-generation dynamic global vegetation models: learning from community ecology, New Phytol., 198, 957–69, https://doi.org/10.1111/nph.12210, 2013.
Seidl, R., Rammer, W., Scheller, R. M., and Spies, T. A.: An individual-based process model to simulate landscape-scale forest ecosystem dynamics, Ecol. Model., 231, 87–100, https://doi.org/10.1016/j.ecolmodel.2012.02.015, 2012.
Semerikov, V. L., Iroshnikov, A. I., and Lascoux, M.: Mitochondrial DNA variation pattern and postglacial history of the Siberian larch (Larix sibirica Ledeb.), Russ. J. Ecol., 38, 147–154, https://doi.org/10.1134/S1067413607030010, 2007.
Semerikov, V. L., Semerikova, S. A., Polezhaeva, M. A., Kosintsev, P. A., and Lascoux, M.: Southern montane populations did not contribute to the recolonization of West Siberian Plain by Siberian larch (Larix sibirica): a range-wide analysis of cytoplasmic markers, Mol. Ecol., 22, 4958–4971, https://doi.org/10.1111/mec.12433, 2013.
Shifley, S. R., He, H. S., Lischke, H., Wang, W. J., Jin, W., Gustafson, E. J., Thompson, J. R., Thompson, F. R., Dijak, W. D., and Yang, J.: The past and future of modeling forest dynamics: from growth and yield curves to forest landscape models, Landsc. Ecol., 32, 1307–1325, https://doi.org/10.1007/s10980-017-0540-9, 2017.
Shugart, H. H., Leemans, R., and Bonan, G. B. (Eds.): A Systems Analysis of the Global Boreal Forest, Cambridge University Press, Cambridge, UK, 1992.
Shuman, J. K., Shugart, H. H., and O'Halloran, T. L.: Sensitivity of Siberian larch forests to climate change, Glob. Change Biol., 17, 2370–2384, https://doi.org/10.1111/j.1365-2486.2011.02417.x, 2011.
Sitch, S., Smith, B., Prentice, I. C., Arneth, A., Bondeau, A., Cramer, W., Kaplan, J. O., Levis, S., Lucht, W., Sykes, M. T., Thonicke, K., and Venevsky, S.: Evaluation of ecosystem dynamics, plant geography and terrestrial carbon cycling in the LPJ dynamic global vegetation model, Glob. Chang. Biol., 9, 161–185, https://doi.org/10.1046/j.1365-2486.2003.00569.x, 2003.
Sitch, S., Huntingford, C., Gedney, N., Levy, P. E., Lomas, M., Piao, S. L., Betts, R., Ciais, P., Cox, P., Friedlingstein, P., Jones, C. D., Prentice, I. C., and Woodward, F. I.: Evaluation of the terrestrial carbon cycle, future plant geography and climate-carbon cycle feedbacks using five dynamic global vegetation models (DGVMs), Glob. Change Biol., 14, 2015–2039, https://doi.org/10.1111/j.1365-2486.2008.01626.x, 2008.
Sjögren, P., Edwards, M. E., Gielly, L., Langdon, C. T., Croudace, I. W., Merkel, M. K. F., Fonville, T., and Alsos, I. G.: Lake sedimentary DNA accurately records 20th Century introductions of exotic conifers in Scotland, New Phytol., 213, 929–941, https://doi.org/10.1111/nph.14199, 2017.
Snell, R. S.: Simulating long-distance seed dispersal in a dynamic vegetation model, Glob. Ecol. Biogeogr., 23, 89–98, https://doi.org/10.1111/geb.12106, 2014.
Snell, R. S. and Cowling, S. A.: Consideration of dispersal processes and northern refugia can improve our understanding of past plant migration rates in North America, J. Biogeogr., 42, 1677–1688, https://doi.org/10.1111/jbi.12544, 2015.
Snell, R. S., Huth, A., Nabel, J. E. M. S., Bocedi, G., Travis, J. M. J., Gravel, D., Bugmann, H., Gutiérrez, A. G., Hickler, T., Higgins, S. I., Reineking, B., Scherstjanoi, M., Zurbriggen, N., and Lischke, H.: Using dynamic vegetation models to simulate plant range shifts, Ecography, 37, 1184–1197, https://doi.org/10.1111/ecog.00580, 2014.
Sugita, S.: Theory of quantitative reconstruction of vegetation I: Pollen from large sites REVEALS regional vegetation composition, Holocene, 17, 229–241, https://doi.org/10.1177/0959683607075837, 2007.
Sugita, S., Parshall, T., Calcote, R., and Walker, K.: Testing the Landscape Reconstruction Algorithm for spatially explicit reconstruction of vegetation in northern Michigan and Wisconsin, Quat. Res., 74, 289–300, https://doi.org/10.1016/j.yqres.2010.07.008, 2010.
Svenning, J.-C., Gravel, D., Holt, R. D., Schurr, F. M., Thuiller, W., Münkemüller, T., Schiffers, K. H., Dullinger, S., Edwards, T. C., Hickler, T., Higgins, S. I., Nabel, J. E. M. S., Pagel, J., and Normand, S.: The influence of interspecific interactions on species range expansion rates, Ecography, 37, 1198–1209, https://doi.org/10.1111/j.1600-0587.2013.00574.x, 2014.
Thuiller, W., Albert, C., Araújo, M. B., Berry, P. M., Cabeza, M., Guisan, A., Hickler, T., Midgley, G. F., Paterson, J., Schurr, F. M., Sykes, M. T., and Zimmermann, N. E.: Predicting global change impacts on plant species' distributions: Future challenges, Perspect. Plant Ecol. Evol. Syst., 9, 137–152, https://doi.org/10.1016/j.ppees.2007.09.004, 2008.
Trenberth, K. E.: Recent observed interdecadal climate changes in the Northern Hemisphere, B. Am. Meteorol. Soc., 71, 988–993, https://doi.org/10.1175/1520-0477(1990)071<0988:ROICCI>2.0.CO;2, 1990.
Wagner, S., Litt, T., Sánchez-Goñi, M. F., and Petit, R. J. R. J.: History of Larix decidua Mill. (European larch) since 130 ka, Quat. Sci. Rev., 124, 224–247, https://doi.org/10.1016/j.quascirev.2015.07.002, 2015.
Wieczorek, M., Kruse, S., Epp, L. S., Kolmogorov, A., Nikolaev, A. N., Heinrich, I., Jeltsch, F., Pestryakova, L. A., Zibulski, R., and Herzschuh, U.: Dissimilar responses of larch stands in northern Siberia to increasing temperatures-a field and simulation based study, Ecology, 98, 2343–2355, https://doi.org/10.1002/ecy.1887, 2017.
Yu, Q., Epstein, H., and Walker, D.: Simulating the effects of soil organic nitrogen and grazing on arctic tundra vegetation dynamics on the Yamal Peninsula, Russia, Environ. Res. Lett., 4, 45027, https://doi.org/10.1088/1748-9326/4/4/045027, 2009.
Zhang, J., Zhou, Y., Zhou, G., and Xiao, C.: Structure and composition of natural Gmelin larch (Larix gmelinii var. gmelinii) forests in response to spatial climatic changes, PLoS One, 8, e66668, https://doi.org/10.1371/journal.pone.0066668, 2013.
Zhang, N., Yasunari, T., and Ohta, T.: Dynamics of the larch taiga–permafrost coupled system in Siberia under climate change, Environ. Res. Lett., 6, 24003, https://doi.org/10.1088/1748-9326/6/2/024003, 2011.
Zhao, G., Bryan, B. A., King, D., Luo, Z., Wang, E., Bende-Michl, U., Song, X., and Yu, Q.: Large-scale, high-resolution agricultural systems modeling using a hybrid approach combining grid computing and parallel processing, Environ. Model. Softw., 41, 231–238, https://doi.org/10.1016/j.envsoft.2012.08.007, 2013.