Young Birds Switch but Old Birds Lead: How Barnacle Geese Adjust Migratory Habits to Environmental Change

Long-distance migratory animals must contend with global climate change, but they differ greatly in whether and how they adjust. Species that socially learn their migration routes may have an advantage in this process compared to other species, as learned changes that are passed on to the next generation can speed up adjustment. However, evidence from the wild that social learning helps migrants adjust to environmental change is absent. Here, we study the behavioral processes by which barnacle geese (Branta leucopsis) adjust spring-staging site choice along the Norwegian coast, which appears to be a response to climate change and population growth. We compared individual-based models to an empirical description of geese colonizing a new staging site in the 1990s. The data included 43 years of estimated annual food conditions and goose numbers at both staging sites (1975–2017), as well as annual age-dependent switching events between the two staging sites from one year to the next (2000–2017). Using Approximate Bayesian Computation, we assessed the relative likelihood of models with different “decision rules”, which define how individuals choose a staging site. In the best performing model, individuals traveled in groups and staging site choice was made by the oldest group member. Groups normally returned to the same staging site each year, but exhibited a higher probability of switching staging site in years with larger numbers of geese at the staging site. The decision did not depend on food availability in the current year. Switching rates between staging sites decreased with age, which was best explained by a higher probability of switching between groups by younger geese, and not by young geese being more responsive to current conditions. We found no evidence that the experienced foraging conditions in previous years affected staging site choice. Our findings demonstrate that copying behavior and density-dependent group decisions explain how geese adjust their migratory habits rapidly in response to changes in food availability and competition. We conclude that considering social processes can be essential to understand how migratory animals respond to changing environments.

Long-distance migratory animals must contend with global climate change, but they differ greatly in whether and how they adjust. Species that socially learn their migration routes may have an advantage in this process compared to other species, as learned changes that are passed on to the next generation can speed up adjustment. However, evidence from the wild that social learning helps migrants adjust to environmental change is absent. Here, we study the behavioral processes by which barnacle geese (Branta leucopsis) adjust spring-staging site choice along the Norwegian coast, which appears to be a response to climate change and population growth. We compared individual-based models to an empirical description of geese colonizing a new staging site in the 1990s. The data included 43 years of estimated annual food conditions and goose numbers at both staging sites , as well as annual age-dependent switching events between the two staging sites from one year to the next (2000)(2001)(2002)(2003)(2004)(2005)(2006)(2007)(2008)(2009)(2010)(2011)(2012)(2013)(2014)(2015)(2016)(2017). Using Approximate Bayesian Computation, we assessed the relative likelihood of models with different "decision rules", which define how individuals choose a staging site. In the best performing model, individuals traveled in groups and staging site choice was made by the oldest group member. Groups normally returned to the same staging site each year, but exhibited a higher probability of switching staging site in years with larger numbers of geese at the staging site. The decision did not depend on food availability in the current year. Switching rates between staging sites decreased with age, which was best explained by a higher probability of switching between groups by younger geese, and not by young geese being more responsive to current conditions. We found no evidence that the experienced foraging conditions in previous years affected staging site choice. Our findings demonstrate that copying behavior and density-dependent group decisions explain how geese adjust their migratory habits rapidly in response to changes in food availability and competition. We conclude that considering social processes can be essential to understand how migratory animals respond to changing environments.

INTRODUCTION
The choices that animals make in response to their environment have typically been shaped by evolution, and are therefore expected to maximize the animal's survival and reproduction. However, environments can change in ways that are hard to predict (Dall et al., 2005). In those cases, animals must deal with uncertainty in the consequences of their decisions. To understand those decisions, it is necessary to know which environmental factors individuals use to inform their decision, and how they integrate those factors to make the decision (i.e., their "decision rules"; Bauer et al., 2011;Budaev et al., 2019). This is particularly true for long-distance migrants, which must make decisions in anticipation of future and distant conditions (Kölzsch et al., 2015).
Animals use current environmental conditions on which to base their decisions, but also previous experiences may affect decisions (Berbert and Fagan, 2012). Memories allow animals to predict habitat quality by deducing temporal trends in stochastic seasonal environments (Abrahms et al., 2019). Furthermore, exploration of the environment can extend such experiences and thereby contribute to making better decisions in the future (Mettke-Hofmann et al., 2002;Tebbich et al., 2009), for instance, by informing the animal about the spatial distribution of resources. Another mechanism that can help the animal to make better decisions is social learning, which allows animals to exploit the experiences of others (Danchin et al., 2004;Couzin et al., 2005;Guttal and Couzin, 2010). Social learning can be an effective means to solve complex problems (Hoppitt and Laland, 2013), especially when combined with learning from previous individual experiences (Rendell et al., 2010). Recent semi-natural experiments suggest that animal populations can indeed accumulate improvements of migratory routes over several generations by combining individual learning with social learning (Sasaki and Biro, 2017;Jesmer et al., 2018), but evidence from natural populations is lacking. It remains largely unknown how migratory animals combine current and previous individual experiences with social learning to make decisions, and whether this combination helps them to adjust their migrations to environmental change.
A good candidate for further investigation is the barnacle goose, which is a social migratory species that has shown striking changes in migratory behavior in response to population growth and climate change (Eichhorn et al., 2009;Jonker et al., 2013). Barnacle geese follow the green wave of grass growth in spring (van der Graaf et al., 2006), but the sites where they stop along the way to accumulate crucial fat reserves for breeding (Drent et al., 2007) seem to be largely determined by tradition. For example, the barnacle goose population that migrates north along the Norwegian coast to breed on the Svalbard archipelago traditionally staged exclusively in Helgeland ( Figure 1A; Black et al., 2014). Recently, a striking change has occurred in this tradition (Tombre et al., 2019). After a small group of birds in the 1990s colonized a new staging site 250 km further north, Vesterålen, the majority of the population has moved to the new site within a few generations ( Figure 1C). The increasing number of birds in Vesterålen coincided with a strong increase in population size, which increased competition for resources at the traditional staging site. The shift in distribution also fits with an increase over the years of suitable habitat in Vesterålen due to climate change. Spring has advanced at both staging sites by 3 weeks since 1975. Grass growth simulations indicated that this advance has led to a higher grass production during the staging period at both sites, and simultaneously to a strong decrease in grass digestibility in Helgeland, but not in Vesterålen where spring starts ∼4 weeks later. As a result, the total production of digestible biomass per square meter of grass during the staging period has more than doubled in Vesterålen, but remained constant in Helgeland ( Figure 1B).
Tombre et al. estimated from ring resightings of individually marked birds that ∼62% of the increasing use of Vesterålen can be attributed to birds that switched from the traditional to the new staging site in subsequent years, suggesting that the choice of staging site might be partly determined by geese responding to the changes in resource availability. However, in a year-toyear comparison switching rates did not correlate with foraging conditions, neither in the current nor in the previous year. This suggests a lack of direct response to changes in food availability, and implies that optimal foraging models (e.g., Bauer et al., 2006;Klaassen et al., 2006) are unlikely to explain the observed dynamics in staging site choice. Furthermore, young birds exhibited higher switching rates than older birds ( Figure 1D). This implies that age-dependent changes in decision-making, which may (partly) result from social processes, affected the observed changes in migratory behavior.
We reason that the observed dynamics in staging site choice may be better understood when explicitly taking into account the ecological and social information that is available to individual animals, and the "decision rules" by which they integrate this information. To this end, we designed a set of simulation models, in which we implemented different potential sets of decision rules by which each individual in a simulated population of barnacle geese decides whether it stages in Helgeland or in Vesterålen. We used individual-based models, which are particularly suitable when the decisions by individuals and interactions between individuals are expected to affect the dynamics of the population (Bauer and Klaassen, 2013). Specifically, we analyzed which set of decision rules best explains the observed changes in stagingsite use, by comparing the performance of different models. Using Approximate Bayesian Computation (Beaumont, 2010), we simultaneously test which model is the most plausible given the empirical data, and estimate the values of the parameters in the selected model(s). Each model contains a different combination of the following five components: (i) adjusting choice to the expected quality of the current staging site, obtained by memorized individual experiences in previous year(s), (ii) comparing expected quality of the current staging site with expected quality of the alternative staging site, which is obtained through explorative behavior in previous year(s), (iii) leaving the choice to others by traveling in a group, (iv) reconsidering staging site choice at arrival in Helgeland, dependent on the current number of geese and/or grass cover, and (v) impact of age on any of the previous four processes. FIGURE 1 | Barnacle goose spring-staging sites. All panels are reproduced from Tombre et al. (2019). (A) is a map of the migratory route (green arrows), and the two staging sites in red and blue. The geese winter at the Solway Firth, and breed on Svalbard. (B) shows the annual estimated staging site quality at both staging sites, estimated as the sum of the daily digestible biomass growth of grass leaves during the staging period. The lines are linear regressions and the shaded areas delineate the 95% confidence interval of local regressions. Panel (C) shows the number of spring staging geese at the two sites as found by the same study. Lines are the trends estimated by local regression, the colored areas depict confidence intervals. (D) shows the probability of geese of particular ages (y-axis) in each year (x-axis) to switch from staging in Helgeland to staging in Vesterålen in the subsequent year, as obtained from resightings of individually marked geese.

Individual-Based Models
We simulated barnacle goose population dynamics in individualbased population models with discrete time steps of one year (see Figure 2 for a visual description). In each model, the simulation runs started in 1970 with a population of 3,000 individuals with randomly assigned sex (50% chance of either male or female) and age (the initial age distribution was derived from a pilot simulation). Each individual was also assigned an age at which to become available as a partner, determined by drawing randomly from the Poisson distribution + 1 and λ = 1.5. This specific distribution with a mean of 2.5 and a standard deviation of 1.2 matches the empirically observed distribution (mean = 2.5, SD = 1.1; Choudhury et al., 1996). At the start of each time step, partnerships were determined, with pairs randomly assigned between available individuals (i.e., at or above the age of first reproduction and unpaired) of the opposite sex.
Individuals remained with the same partner in subsequent years, only becoming available again as a partner when the partner died (Black et al., 2014; in reality the annual chance of a pair to separate is 2%, which we chose to ignore). All unpaired birds and a randomly assigned bird within each pair then chose a staging site: either Helgeland or Vesterålen. During the first time step, all individuals were set to choose Helgeland. In later time steps, individuals could instead decide to visit Vesterålen (see section Staging Site Decision Rules). Subsequently, each paired female reproduced with probability b s,t , where s is the visited staging site and t is the calendar year.
Previous simulation studies of goose behavior have focused on energetics Klaassen et al., 2006). While explicitly modeling density-dependent energy gain at staging sites and the consequences for reproductive success, they simplified the process of decision-making by assuming optimal behavior. We focused on the decision-making process and instead simplified the energetic part of the model. We  Table 2, differ only in the presence or absence of five components: "memory" (yellow), "reconsider" (purple), "exploration" (green), "groups" (red) and "aging" (blue). Mortality occurs with a probability of 0.11, reproduction occurs with a probability that depends on the conditions of the staging site (grass production as well as the number of geese at the staging site, see Equation 1) that the female has visited, being either Helgeland (H) or Vesterålen (V).
assumed that b s,t depends linearly on the annually estimated grass production at the staging site that she visited, and also decreases linearly with an increasing number of birds at that staging site, depending on the surface of foraging area: where N s,t is the abundance of birds at the visited staging site and K s is the total surface area of suitable foraging habitat at that staging site (m 2 ). The probability of reproduction in absence of competition, r s,t , is a linear function of the digestible grass production per m 2 during the staging period in year t at staging site s, q s,t (measured in g/m 2 , see next section): where a is a conversion factor (m 2 /g). The lower boundary of 0.1 reflects the low probability of reproduction observed for geese with very low body condition before departing Helgeland (Prop et al., 2003). Instead of deriving K H (carrying capacity in Helgeland) and conversion factor a mechanistically, we fitted them by performing model simulations without staging site choice, assuming all individuals to stage in Helgeland. The simulated population sizes were compared to the population count data between 1970 and 1997, when virtually all individuals visited Helgeland (see Figure 1). K H and a were estimated as 44,000 and 0.0082, respectively, by selecting the values that minimized the distance between the simulated population sizes and the empirically derived values (see section Calculating the Distance of Each Simulation to Empirical Data). Based on the ratio of agricultural land in the two areas (summed surface of agricultural land in 2017 was estimated at 27.6 and 88.5 km 2 for the main goose areas in Helgeland and Vesterålen; data downloaded from www.ssb.no), and given that barnacle geese in Helgeland also make use of natural salt marshes and that barnacle geese in Vesterålen face competition for food with pink-footed geese (Tombre et al., 2019), we estimated conservatively that K V was two times K H . A higher value of K V had no strong effect on the model selection results, as the population in Vesterålen remained far below carrying capacity in all simulations (see Appendix I and Table S1).
The number of offspring produced by a reproducing female was drawn from a Poisson distribution + 1 with λ = 1, resulting in a mean of two offspring, which equals the distribution in the number of juveniles associated with successful breeders in the wintering area (Black et al., 2014). At the end of each time step, individuals had a probability of dying, d, estimated at 0.11 (Black et al., 2014). Each simulation consisted of 48 time steps, representing the period from 1970 to 2017.

Grass Production at the Staging Sites
The digestible grass production per m 2 during each spring staging period t at staging site s, q s,t (g/m 2 ), was taken from Tombre et al. (2019). It was estimated as the sum of the daily digestible biomass growth of grass leaves from 30 April to 20 May (Prop and Black, 1998). The daily values were calculated by means of the simulation model CATIMO (Canadian Timothy Model; Bonesmo and Bélanger, 2002a,b). CATIMO simulates the daily growth of cell walls and cell contents in the leaves of timothy, Phleum pratense. Timothy is one of the main agricultural grass species and an important food source for barnacle geese in Norway (Black et al., 1991). Daily grass growth (g/m 2 ) was converted to digestible daily grass growth (g/m 2 ) by taking into account that the digestible proportion for barnacle geese is 0.16 and 0.64 for cell wall and cell content respectively (Prop and Vulink, 1992). The simulations were based on daily local temperature and radiation values. See Tombre et al. (2019) for a full explanation.

Staging Site Decision Rules
We compared 22 models, all with different decision rules determining the choice of staging site. Each set of decision rules is a combination of five components. The first component is memory, which we incorporated as an effect of staging site quality that the focal individual experienced in previous years. The second component is exploration, which we modeled as an effect of staging site quality at the alternative staging site in previous years when the individual was alive. The third component is traveling in groups. This is an effect of the staging site choice of others, in most cases the group leader, and hence a consequence of social learning. The fourth component is to reconsider staging site choice at arrival in Helgeland, with each individual continuing migration to Vesterålen with a probability that depends on the number of geese in Helgeland and/or the grass cover in Helgeland. As the fifth component, we included age-dependent differences between individuals in any of the four previous components (see also Figure 2).
In all models, paired birds stay together and normally return to the staging site of the previous year. In case newly paired birds did not visit the same staging site in the previous year, they make a random choice between both staging sites. Analysis of ring resightings before and after pair formation does not suggest a sex bias (TO and JP, unpublished data). Unpaired birds normally visit the staging site of the previous year. We assumed that each individual has an 18% probability of remaining with its parents during the first spring migration (Black et al., 2014), thereby copying the staging site choice of the parents. The first-year birds that do not stay with their parents follow others, based on one of the following criteria (denoted by parameter c juv , for all parameters see Table 1): (1) follow a random non-first-year bird, (2) follow a parent (i.e., an individual that has produced offspring in the previous year), or (3) follow an individual of at least 10 years old, which is approximately the top 30% of the age-distribution (Black et al., 2014).
On top of this basic scheme, each individual can decide to switch staging site relative to the previous year. In the first model, each individual has a fixed annual probability of switching staging site (parameter n). Subsequent models incorporate different combinations of the five components as described below.

Memory
In each year, the expected probability of reproducing when returning to the current staging site (as opposed to switching to the other staging site), E b c , is given by a weighted average of its past experiences at that site. The weight of each of those experiences is given by the decay function e − y m , where y is the 'age' of the experience (in years) and parameter m determines the rate at which memories fade. We assumed that individuals start switching to the other staging site when E(b c ) falls below a threshold that is given by parameter x a . Below this threshold, the probability of switching increases with decreasing E(b c )with a rate that is determined by parameter x r , where: Exploration Individuals explore the alternative staging site at the end of the staging period with probability (v), enabling them to inform their expectation of the reproduction probability when visiting the alternative staging site, E b a . If the difference between E b a and E b c is larger than x b , then the probability of switching staging site in the next year is given by: where parameter x r determines how fast the probability of switching increases as the difference between E b a and E b c increases. This component only affects the model results when memory is also implemented, with equation 3b replacing 3a.

Groups
Instead of individually deciding where to go, birds may also choose to follow another individual, thereby copying its choice of staging site. We modeled this by assigning each bird to a group, and determining staging site choice per group instead of per individual. In this case, juveniles do not join an individual, but a group. We assumed that 18% of the juveniles joins the group of their parents (Black et al., 2014), and the rest joins a randomly chosen group. Group decisions may be made in different ways, denoted by parameter c group . We assumed that individuals either (1) form groups with a single leader, which may be (1) a random bird, (2) a randomly chosen parent (i.e., an individual that has produced offspring in the previous year), or (3) a randomly chosen bird from among the oldest ones. Alternatively, each group member first makes an individual choice as explained above, after which the group reaches consensus by adopting the "majority vote" (4). Note that simple and plausible behavioral mechanisms allow individuals to follow any of these rules, without having an overview of the process (Couzin et al., 2005). We further assumed that individuals join the same group as in the previous year (but see component v, Aging). Maximum group size is determined by parameter g, with groups splitting into two equally sized groups when larger than g, and merging with a random other group when smaller than 0.25g.

Reconsidering Staging-Site Choice
At arrival in Helgeland, individuals have the possibility to reconsider their choice, and continue to Vesterålen. The probability to continue is either linearly dependent on the number of geese, N H (Reconsider geese ), or on the grass cover at  (Black et al., 2014) Randomly chosen parameters (drawn from transformed uniform distribution) c juv 1, 2 or 3 -Criterion for juveniles choosing whom to follow: 1 = randomly chosen individual, 2 = randomly chosen parent, 3=old individual (>10 years) c group 1, 2, 3 or 4 -Criterion for group decision: 1 = randomly chosen group leader, 2 = randomly chosen parent is group leader, 3 = oldest individual is group leader, 4 = majority vote; each individual has a preference and the majority preference is chosen by all g 20-1,000 Individuals Maximum group size (groups split in two when at size g, and merge with another group at size 0.25g). ). Used in models without 'exploration' Threshold value of the difference between the expected probability of reproducing in the alternative and the current staging site, below which the probability of switching staging site starts to increase linearly (see Equation 3b). The switching probability is set at zero if the individual has no expectation of the alternative staging site. Used in models with 'exploration' x r 1-200 -The slope of increase in probability of switching as the expected probability of reproduction decreases. Aging We explored four different potential effects of age. The first assumed that the influence of previous experiences on the current decision decreases with the age of the individual. We modeled this by multiplying the probability of switching (see Equations 3a and 3b) with an age-factor a that changes with age according to the function a = e 1−age ar , where age is measured in years. Parameter a r (also in years) determines the strength of the age-effect. A second possibility is that the probability of exploring (v) decreases with the individual's age, which is modeled by multiplying v by age-factor a. Thirdly, if the animals make migratory decisions in groups (see component i), then there may be an age-effect in the probability of changing to a randomly chosen new group, w 0 , which is then multiplied by the age-factor a. Fourthly, there could be an ageeffect in the tendency of individual geese to reconsider their staging site choice upon arrival in Helgeland. This is modeled by multiplying the probability to reconsider (Equations 4a and 4b) with the age-factor a.

Empirical Data
To determine which model is most plausible, we compared the simulations to two different sets of empirical data, both published by Tombre et al. (2019). The first set consists of the annual number of spring staging barnacle geese in Helgeland and in Vesterålen. This set contains 86 data points, being the estimated numbers of birds at each site in each year from 1975 to 2017 ( Figure 1C). They were derived from annual counts during the staging period in Helgeland and Vesterålen, and annual counts of the total population size in the wintering area. The second set of data points consists of the probabilities of individual birds switching from staging in Helgeland to staging in Vesterålen in subsequent years (from here on referred to as "switching probabilities"). Each data point is the switching probability for an individual of a given age (age 1 to age 20) in a particular calendar year between 2000 and 2016 ( Figure 1D). These data points were derived from resightings of marked individuals at both staging sites, as well as the wintering and the breeding area. Further details can be found in Tombre et al. (2019). As hardly any geese were observed staging in Vesterålen from 1975 to 1995, we infer that switching probabilities from Helgeland were zero from 1975 to 1995 for all ages. The years 1996 to 1999 were not part of the analysis. This resulted in a total of (21 + 17) × 20 = 760 data points. We did not compare the switching probabilities in the other direction (from Vesterålen to Helgeland), because these could not be estimated in years when the simulated bird numbers were zero in Vesterålen.

Model Selection: Approximate Bayesian Computation
We evaluated the relative strength of the different models by comparing simulations to the empirical data using Approximate Bayesian Computation (ABC; Beaumont, 2010) in R (R Core Team, 2018). This statistical tool has been developed to quantify the fit of different individual-based models to different sets of empirical data simultaneously. The ABC-method allows the fit of different models to be compared (e.g., models with and without memory), as well as comparing the fit of different parameter values within each model (e.g., values of a parameter determining the rate of memory loss). The method is called "Bayesian" because the method updates the degree of belief in each model given the empirical data. It is "Approximate" because it is not an analytical method, which is generally not an option for individual-based models, but instead relies on simulations (van der Vaart et al., 2015). We used rejection-ABC, the simplest and most accessible type of ABC that can be used for ecological models with multiple parameters (van der Vaart et al., 2015(van der Vaart et al., , 2016. Calculations were performed as in the R-package "abc" (Csilléry et al., 2012), except for indicated differences. Below, we explain the method step by step. First, parameter values are defined. Where possible, parameters were estimated from the literature (see Table 1). For the other parameters, distributions were defined such that all possible values are included (see Table 1). These distributions are referred to as "prior distributions." Then, 10,000 simulation runs were performed for each model. For each simulation, the values of all parameters in the model were drawn at random from the prior distributions. After all simulation runs were performed, we calculated the distance between each run and the empirical data (see next section). To give equal weight to both used datasets (bird numbers and switching probabilities), we calculated the distance of each simulation run to the empirical data separately for each dataset, and then took the mean of the two to arrive at a single distance estimate for each simulation run. Finally, the 100 runs with the smallest distance were selected. The evidence for model x relative to model y is expressed by the Bayes factor (B x,y ), which in this context is defined as the ratio of simulations from each model among the selected runs (van der Vaart et al., 2016).
To test whether the result would change with more simulations, we ran a bootstrapping test of the model selection accuracy by repeating the procedure 100 times, each time with a randomly chosen half of all simulation runs. To evaluate the ability of the ABC-method to distinguish between different models, we carried out cross-validation as implemented in the function "cv4postpr" in the "abc" R-package, and described in Csilléry et al. (2012). First, 100 simulation runs are randomly selected from each model. Then, for each of these runs, the complete model selection procedure is repeated after removing this run from the simulation data and replacing the empirical data with this run. The result is a "confusion matrix", where each row represents the number of simulations under each model, and each column represents the number of simulations assigned to that model by the model selection procedure.
The distribution of parameter values among the selected simulation runs ("posterior distributions") can be regarded as a probability distribution for each parameter, and acts as a sensitivity analysis. To test whether the posterior distributions were significantly different from the prior distributions (distribution of parameter values among all runs), we performed a Chi-square test after dividing the data into 10 equally-sized bins with the function "bin" in R-package OneR (von Jouanne-Diedrich, 2017). To correct for multiple testing, we applied a Bonferroni correction to the standard significance level of 0.05.

Calculating the Distance of Each Simulation to Empirical Data
Distance (ρ) is defined as the standardized Euclidian distance between all data points j in simulation i (M i ) and the same data points in the empirically derived data (D): where M i,j is the output of run i for datapoint j, D j is the empirically derived value of data point j, and sd(M j ) is the standard deviation of data point t in all simulation runs. As in van der Vaart et al. (2015), we used standard deviation instead of median absolute deviation (as is done in the "abc" package; Csilléry et al., 2012), because the median was zero for several datapoints and this led to undefined distances. To Selected models are in bold. *Model selection was performed in two steps: first only with models without "Age." The best model and competitive models were tested in a second step, together with a new set of models based on those models that included "Age".
avoid overfitting, we chose to compare the simulations to the statistically estimated trends (Figures 1B,C), rather than to the raw empirical data. We made this decision because an unknown part of the inter-annual variation in the empirical data is caused by non-modeled processes, such as annual conditions in the breeding area and observation errors.

Reducing the Number of Simulations
To reduce the required number of simulation runs, we adopted a two-step model selection procedure. First, we performed a model selection of scenarios without the "age" component (models 1 to 15 in Table 2), and executing 10,000 simulation runs per model. We then composed seven additional models based on the selected models, but also including an age-effect (models 16 to 24 in Table 2), and executed 10,000 simulations per model. We did not consider models with an age-effect on more than one component, to further limit the number of models to be tested. These additional models were tested in a new model selection procedure, also including the selected models from the first model selection. For parameters in the first selection where the posterior distribution was significantly different from the prior distribution ( Figure S1, Table 3), we updated the prior distributions for the simulations in the second model selection procedure (Table 3).

Model Selection
The simulation most similar to the empirical data was produced by model 16, which includes "reconsider geese ", "groups", and an age-effect on "groups." The pattern resulting from this simulation corresponded to the observed annual bird numbers at both staging sites (Figures 3D,E). Moreover, it showed a decrease in switching probability with age ( Figure 3F), which was similar to the pattern in the empirical data ( Figure 1D). This model was also the best represented model among the 100 best simulation runs (60 out of 100 runs, Table 2). The same model but with "memory" (model 21) was represented with 34 runs. With a Bayes factor of B 16,21 = 1.8 there is no evidence that memory does not play a role, but it does not improve the performance of the model in explaining the empirical data. Roughly, a Bayes factor of 3 to 10 is regarded as "substantial evidence" and above 10 as "strong evidence" (Kass and Raftery, 1995;van der Vaart et al., 2016). Apart from models 16 and 21, only model 18 (model 21 but without "reconsider geese ") occurred among the 100 best models, with 6 runs (B 16,18 = 10 and B 21,18 = 5.7), meaning that there is substantial evidence for models 16 and 21 over model 18, and strong evidence over all other models. Hence, the results suggest that staging site choice is made in groups, with a decrease over age in the probability that individuals change groups, and that groups switch to another staging site based on the current number of geese at the staging site. The results are indefinite regarding the role of previous experiences at the alternative staging site.
There is no evidence that exploration of the other staging site in previous years plays a role, nor that there is an effect of current food conditions at the staging site.

Model Validation
Because models 16 and 21 both came out as likely to underlie the empirical data, we focused on these models in the model validation. When repeating the model selection analysis 1,000 times with a randomly chosen half of the data, models 16 and 21 together always made up the majority of the selected runs (range 91-100 out of 100 selected runs, mean 96, Figure 4). Hence, the evidence for models 16 and 21 relative to the others is robust. The only other model that appeared among the selected simulations was model 18 ("memory", "groups" and an age-effect on groups, mean 4, range 0-9). The cross-validation procedure suggested that the model selection performs badly in estimating the underlying model of randomly drawn simulations: of the runs that were produced by model 16 or 21, only 67% were also estimated as such ( Figure S2). This result was to be expected, because simulations were similar between models for a large proportion of parameter combinations. For example, switching did not occur at all in many simulation runs of all models with "groups" (between 6 and 40%). When performing the cross-validation procedure with the 100 best fitting runs of each model instead of randomly drawn runs, then 98.5% of the runs produced by model 16 or 21 were also estimated as such ( Figure S2). Hence, when the data was close to the observed trends, the model selection performed well.

Parameter Estimation
In the 100 selected simulations runs of the first step in the model selection (see Figures 3A-C for simulation results), the distribution of values (posterior distributions) of 10 out of 15 parameters were significantly different from the defined prior distributions, of which five were in models that were represented among the best simulations (Table 3, Figure S1). For those parameters, we defined new prior distributions for use in the simulation runs for the main model selection (Table 3). In the selected simulations after the second step in the model selection, the posterior distributions of five out of ten parameters were significantly different from the defined prior distributions (c group , g, ge 0 , w 0 and a r , Table 3, Figure 5). In all of the selected runs the birds traveled in groups. Smaller groups occurred more often among the selected runs than larger groups (see Figure 5B). In most of the selected runs, the oldest individuals led the group (78 out of 100 runs, Figure 5A). Simulations where group decisions were made by a majority vote always performed badly (see Table 3 and Figure S1; it did not occur in the selected runs in the first step, and was therefore removed from the prior distribution of the main selection procedure). Individuals switched between groups in all selected simulations, with most of the runs having an initial switching probability below 0.4, and a relatively slow decrease with increasing age (Figure 5F). In the selected runs where the probability for a group to switch increased with goose numbers at Helgeland (94 out of 100), birds started to switch when numbers were between 10,000 and 15,000 geese (parameter ge 0 ). The selected runs including "memory" and "reconsider geese " responded less strongly to density (parameter ge r ) than the runs with "memory" but without "reconsider geese " (Figure 5D). There was no pattern in the maximum probability to reconsider staging site (parameter ge m ; Figure 5D). The selected runs with memory (40 out of 100) showed no clear pattern in the rate of memory loss (parameter m; Figure 5F), suggesting that the rate of memory loss is not importantly affecting the dynamics. The same was the case for x r , the rate at which the probability of switching increases when the expected probability of reproducing declines (the slopes in Figure 5C).

DISCUSSION
Simulations resembled the empirical data best when geese were assumed to travel in small groups that are led by the oldest individuals, and when young geese switched more between groups in subsequent years than did older individuals (Table 2, Figure 5). Further, the results suggest that the current food conditions are of minor importance to staging site choice, but that the abundance of geese in Helgeland does increase the probability for groups to reconsider their choice and continue to Vesterålen. The model results are indecisive about whether experiences acquired by the group leaders in previous years, i.e., the "memory" component, influence the decision to switch staging site. We found no evidence that experiences at the alternative staging site in previous years contributes to the decision ( Table 2). Below we discuss the implications of these results in more detail.

Grouping
The well-known fact that geese operate in groups need not inherently imply that each individual's choice of staging site is influenced by other members the group. For example, group-foraging pink-footed geese during spring staging decided individually on their specific daily foraging locations FIGURE 4 | Bootstrapping test of model selection. The ABC-analysis was repeated 100 times, each time using a randomly drawn 50% of all simulation runs. Plotted are the frequency distributions of the representation of each model among the 100 best simulation runs in each repeat, excluding zero. The only models that always occurred among the selected runs are models 18 (green, mean = 4), 21 (blue, mean = 32) and 16 (red, mean = 64). Models 16 and 21 together always represented at least 91 of the 100 selected runs (gray, mean = 96). These are two models that include "reconsider geese ", "grouping", and an age-effect on "grouping." Model 21 additionally includes "memory" (see Table 2 and Figure 2). (Chudzinska et al., 2016). Our results are the first to suggest that group decisions do play a role in the choice of staging site. In all selected simulations (i.e., best fitting with the empirical data), staging site choice was made in groups.
The results further suggest that these decisions are not arrived at by a majority vote. The gradual increase in numbers in Vesterålen in the 1990s is not compatible with this decision rule, which requires a high proportion of all individuals to prefer switching, before the first geese start to switch. This aligns with the idea that strong conformity is generally not a good strategy in changing environments, because innovative behavior is unlikely to spread even when highly adaptive (Eriksson et al., 2007;Kandler and Laland, 2009). The most likely group decision rule was to follow the oldest, and therefore most experienced, bird of the group. This rule performed better than following parents (Chi-squared test, χ 2 1 = 36.6, p < 0.0001), which in turn performed better than following a random individual (Chi-squared test, χ 2 1 = 7.7, p = 0.006). Following experienced birds might be adaptive because the annual food conditions at the staging site vary stochastically (Figure 1B), and longer experience will provide a better prediction of next year's staging site conditions. In contrast, following an individual that produced offspring in the previous year is hardly predictive of the chances to reproduce in this year if annual stochasticity is high (Baldini, 2012). This may explain why the model results indicated that following an individual that raised offspring was less likely than following an experienced leader. That following experienced birds is better than following successful breeders also could explain why in reality most firstyear barnacle geese choose not to follow their parents on their first spring migration; on average, it would provide a higher pay-off to follow old and experienced individuals than to follow the parents. However, inclusive fitness arguments predict that unrelated group members may be more hostile than parents or other related individuals. Indeed, this also holds for barnacle geese (Black et al., 2014). Nonetheless, there are more examples of animals that are more likely to copy old (Amlacher and Dugatkin, 2005) and knowledgeable (Kendal et al., 2015) individuals, and to copy experienced others rather than the parents (Agostini et al., 2017). In bird flocks, leaders have been shown to be the more experienced individuals (Flack et al., 2012;Mueller et al., 2013). Our results imply that following experienced birds is especially advantageous when recent success needs not be a good predictor of subsequent success, but multiple-year averages of success are.

Reconsideration of Staging Site Choice at Arrival in Helgeland
The component "reconsider geese " featured in all selected simulation runs. In models with this component, group leaders are more likely to reconsider their staging site choice after arrival in Helgeland in years when the number of birds in Helgeland is high. Simulations with this density-dependent effect corresponded better to the empirical data, because this effect keeps individuals from switching to Vesterålen before 1990. This also explains why simulation runs with "reconsider grass " do not perform well, not even when combined with "reconsider geese ". In models with "reconsider grass ", the probability of reconsidering staging site choice increases as the grass phenology is more advanced at arrival in Helgeland. In those simulation runs, individuals do often colonize Vesterålen before 1990 because years with an early spring also occurred before 1990 ( Figure 1B). Hence, these results suggest that the choice between Helgeland and Vesterålen is not a direct response to the "green wave" of spring phenology (van der Graaf et al., 2006). Instead, the growing preference for Vesterålen follows from a response to other geese, both positive (grouping) and negative (densitydependent switching).

Memory and Exploration
From an optimal foraging perspective, it is expected that any knowledge about the conditions at the current or alternative staging site should play a strong role in the decision whether or not to return to the current site in the following year Krebs, 1986, Abrahms et al., 2019). This influence was captured in the "memory" and "exploration" components of the model. The "memory" component was part of 40 out of 100 of the selected simulations (models 18 and 21; see Table 2), Although this is not evidence against memory playing a role, we conclude that there is no need to assume that geese memorized foraging conditions at the staging site in the previous year(s). Note that this only concerns memory of foraging conditions. In all models, individuals (or at least group leaders) are assumed to have spatial memory, and remember the migration route and staging site of the previous year (Mettke-Hofmann and Gwinner, 2003).
Adding the "exploration" component also did not improve the fit of simulations to the data, as the best model in the first step of the model selection with exploration (model 12) was less well represented than the same model without exploration (model 11). Hence, the current results are also indecisive with regard to the importance of exploration for decision-making. Geese have only rarely been observed to spend a significant amount of time at both staging sites in one spring, but they  Table 2) are given in dark gray. In the line graphs, each model is shown in a different color. For the bar plots, patterns did not differ between models. (A) gives the following criterion in each simulation of the models with "groups" (parameter c group ). (B) is the frequency distribution of the maximum group size in each simulation run (parameter g). The lines in (C) define how the annual switching probability depends on the individual's expected probability of reproducing at the current staging site, E(b c ). They are determined by parameters x a (threshold value below which the probability becomes non-zero) and x r (the slope of the line below x a ) in the models with memory but without exploration. In (D), the lines are determined by parameters ge 0 , ge r and ge m in the models with "reconsider geese ", which determine how the probability to switch preference after arrival at the staging site, depending on the number of geese there. (E) shows how, resulting from differences in m, the weight of each memory declines over the years. The lines in (F) depend on parameters w 0 and a r , and define how the probability of switching between groups decreases with age in each simulation (only models 16 and 21).
occasionally made a short stop in Helgeland before continuing to Vesterålen (PS, IT and JP, unpublished data). Less frequently, geese staging in Helgeland were also observed in Vesterålen at the end of the staging period, although most geese fly directly north after staging in Helgeland (PS and Larry Griffin, unpublished visual observations of departing geese and satellite tracks). A potential way forward is to add a third set of empirical data to the comparison, for example containing information on individuals that were (or were not) observed at multiple staging sites, in relation to their switching behavior. However, exploring individuals may be easily missed by observers if they land only shortly or not at all, making it hard to determine the rate of occurrence by ring resightings. More information on the rate of exploration and age-dependent changes in exploration could be derived by tracking individuals with gps-tags. Another possibility is to model the effect of exploration in more detail, which might lead to a better fit with the current empirical data. For example, new simulations could allow the probability of exploring Helgeland when staging in Vesterålen to be different from the probability of exploring Vesterålen when staging in Helgeland.

Aging
The finding that migratory decisions are age-dependent confirms a general trend that young birds become more consistent in their migratory decisions as they grow older (Lok et al., 2011;Oppel et al., 2015;Vansteelant et al., 2017). In Eurasian spoonbills (Lok et al., 2011) as well as in pink-footed geese (Clausen et al., 2018), a higher probability for young individuals to switch wintering site between years was attributed to young birds being more explorative. This has also been the main hypothesis to explain the higher probability of staging site switching by barnacle geese (Tombre et al., 2019). However, our results suggest that juveniles do not explore new staging sites deliberately. Instead, they are more likely to travel with different groups in subsequent years, which results in a higher probability of ending up at different staging sites. Also this groupswitching behavior might be understood as being "explorative", but it is social exploration rather than spatial exploration. This is an important distinction because it implies that migratory innovation needs not start with young and naïve individuals, as was suggested before. The modeling exercise indicates that the colonization of Vesterålen is more likely to have been initiated by old and experienced individuals, which were being followed by young animals.

Suggestions for Future Research
By comparing our simulations to the statistical trends in the empirical data (instead of the raw empirical data) uncertainty in the empirical trends were not conveyed to the statistics of our model selection procedure. We think that this method is to be preferred over using the raw data in this case, because the empirical trends were derived from a resighting analysis of individual bird observations. That analysis takes into account the probability of either not observing a bird when it is actually there (resighting probability) or not observing a birds because it is dead (mortality probability; for details see Tombre et al., 2019). To take these considerations into account when comparing individualbased models to the raw data, it would be necessary to also simulate the process of data collection within these models. We propose that incorporating data collection in the simulation exercise could be an interesting venue for future research.
In individual-based models, each decision must be modeled explicitly (Bauer and Klaassen, 2013). The advantage is that all of the underlying assumptions, many of which remain implicit and are often ignored in other types of modeling, become explicit. A disadvantage is that it remains unknown how much the way a process is modeled affects the results. For example, we modeled the process of pair formation and group formation in a basic way, with young individuals choosing a partner or a group at random. There are indications that individuals will be more likely to group with others that they grew up with (Choudhury and Black, 1994;van der Jeugd et al., 2002). Other studies have shown that social structure within groups can have strong effects on group dynamics (e.g., Bateman et al., 2013). Modeling these aspects more precisely could produce further insights into the causes and consequences of group formation by barnacle geese.
We stress that we only investigated the tip of the iceberg when it comes to individual differences. There may well be differences in decision-rules between individuals other than those mediated by age. Research on individual differences (Dingemanse et al., 2010), including those in barnacle geese (Kurvers et al., 2009), has shown that animals within the same population and of the same age can differ greatly in personality characteristics such as dominance, aggression, and exploration. Although beyond the scope of this study, such individual variation could be incorporated as an extension of the current study by assuming that individuals within the same population can act according to different sets of decision rules.

Cultural Evolution of Migratory Behavior
Social learning is an essential part of migratory inheritance and development for many migratory bird species (Sutherland, 1998;Helm et al., 2006;Németh and Moore, 2014), and for barnacle geese in particular (e.g., Eichhorn et al., 2009;Jonker et al., 2013). This study is the first attempt to infer the details of the learning processes in migratory decision-making from empirical data. The results indicate that geese travel in groups led by the oldest individual whose decisions are density-dependent, and the modeling explains how barnacle geese are able to respond so rapidly to long-term trends in competition and climate change at the staging sites. This is in line with the long-held conviction that cultural evolution allows for faster adaptation than genetic evolution (Boyd and Richerson, 1985;Sutherland, 1998). However, copying the behavior of conspecifics can also inhibit behavioral adjustment, and cause sub-optimal traditions to be maintained (Warner, 1988;Day et al., 2001;Németh and Moore, 2014). In order for social learning to lead to rapid adaptation, it typically needs to be combined with some low amount of individual learning, or other processes that introduce variation (Rendell et al., 2010). Intriguingly, the decision process that we identified here as being the most likely for migratory decisions by barnacle geese, does exactly this. Most geese follow others, but some of the experienced geese that lead the groups alter their decisions in response to current conditions.
We have discussed how the adaptive value of the observed decision rules is expected to depend on the amount and nature of the environmental variation. Although we expect that our results will also apply to other decisions, both by barnacle geese and by other social species, we stress that care should be taken when generalizing the results. An interesting venue for future research will be to apply the methods presented here to other published studies of migratory behavior across taxa and across situations. Finding general patterns between decision rules and environmental and social context will help to understand why some populations are more vulnerable to environmental change than others, and allow for better predictions of the ecological consequences of climate change. Currently, most studies of population dynamics do not consider the specific processes by which animals make their decisions. While arguably in some cases this may be a legitimate simplification, in cases like the present one, social and developmental aspects of decisionmaking turn out to be essential for understanding the populationscale response to environmental change.

DATA AVAILABILITY STATEMENT
All datasets generated for this study are included in the article/Supplementary Material.

ETHICS STATEMENT
Ethical review and approval was not required for the animal study because all the data for this study has been published elsewhere, and where necessary were approved by an Animal Ethics Committee.

AUTHOR CONTRIBUTIONS
TO, JP, KL, and GR conceived the study design. JP, IT, and PS provided the data. TO performed the analyses and wrote the manuscript with contributions from all authors.