1

Density-dependent immigration promotes population stability in a long-

2

distance migratory bird

3

S. Wilson1, A.E. McKellar2, M.W. Reudink3, P.P. Marra4, and L.M. Ratcliffe5

4
5
6

1

7

Research Centre, 1125 Colonel by Drive, Ottawa, ON, K1A 0H3, Canada

Wildlife Research Division, Environment and Climate Change Canada, National Wildlife

8
9

2

10

Wildlife Research Centre, 115 Perimeter Road, Saskatoon, SK, S7N 0X4, Canada,

Canadian Wildlife Service, Environment and Climate Change Canada, Prairie and Northern

11
12

3

13

Canada

Department of Biological Sciences, Thompson Rivers University, Kamloops, BC V2C 0C8,

14
15

4

16

7012, USA

Migratory Bird Center, Smithsonian Conservation Biology Institute, Washington, DC 20013-

17
18

5

Department of Biology, Queen’s University, Kingston, ON K7L 3N6, Canada

19
20
21
22

Corresponding Author: Scott Wilson, email: scott.wilson@canada.ca, phone:613-990-9740

23

Abstract

24

The spatial structure of populations determines the relative importance of reproduction, survival

25

and movement on temporal population dynamics. However, the mechanisms by which local

26

individuals and immigrants interact and the subsequent effects of immigrants on population

27

productivity are poorly known, especially in large, open populations. We developed an

28

integrated population model (IPM) to study the extent and consequences of immigration on the

29

dynamics of a Neotropical migrant (American redstart, Setophaga ruticilla) over an 11-year

30

period in Ontario, Canada. New immigrants represented the majority of the study population

31

each year with higher rates of immigration for males than females and for first-year breeders

32

than breeders in their second year or older. Immigration was negatively density dependent, with

33

immigrants replacing previously established breeders in a compensatory manner following their

34

death or emigration. Because of the tradeoff between immigration and apparent survival, neither

35

rate had a strong influence on population growth and reproductive output was most strongly

36

correlated with a change in abundance between years. However, if immigration ceased, the study

37

population would become locally extinct within 5 years and thus immigrants were essential for

38

local population persistence. Moreover, we found no evidence for reduced breeding success

39

when immigrants represented a higher proportion of the study population. Our research

40

highlights the importance of movement in the stability of open populations and the strong

41

correlation between the fates of local breeders and the number of immigrants entering the

42

population. We recommend the use of IPMs to address the spatial scale over which immigration

43

occurs and how different scales influence its contribution to temporal population dynamics.

44

Keywords: emigration, integrated population model, movement, population structure,

45

reproductive success, survival

[Type text]

46

Introduction

47

The extent to which population dynamics are determined by within-site processes (reproduction,

48

survival) versus among-site processes (immigration, emigration) should differ depending on

49

where a population lies along a continuum of spatial structure (Hanski and Gilpin 1997; Thomas

50

and Kunin 1999; Bowne and Bowers 2004). This continuum ranges from isolated populations

51

with low immigration through a series of patchy populations with varying connectivity to a

52

highly interconnected network with frequent movement among patches within a larger

53

population (Harrison and Taylor 1997; Thomas and Kunin 1999). The ability of a species to

54

disperse relative to the heterogeneity of the habitat in the ecosystem will largely influence where

55

a population falls along this continuum (Addicott et al.1987; Bowne and Bowers 2004). For a

56

given species and ecosystem, the relative influence of movement on population dynamics will

57

also depend on the spatial scale at which a study is conducted (Camus and Lima 2002; Bowler

58

and Benton 2005). In particular, as the spatial scale increases, immigrants will typically comprise

59

a lower proportion of the individuals in a patch and thus within-patch processes become

60

relatively more influential.

61

Early studies on the importance of movement for populations focused mainly on the

62

influence of immigration in systems with a clear distinction between habitat patches and an

63

inhospitable intervening matrix (Brown and Kodric-Brown 1977, Fahrig and Merriam 1985,

64

Hanski 1998). Numerous animal species inhabit environments with such a spatial structure either

65

naturally (e.g., freshwater ponds) or due to anthropogenic causes (e.g., forest patches in

66

agricultural landscapes). In these environments, immigration can provide an essential rescue

67

effect allowing for population persistence even when immigration rates are relatively low

68

(Stacey and Taper 1992, Vilá et al. 2003, Sanz-Aguilar et al. 2016). Less well known is how

[Type text]

69

local demographic processes (i.e., reproduction and survival) and movement interact to drive

70

patch-level population dynamics when there is a high degree of connectivity among patches of a

71

larger population. Many species occupy landscapes where it can be difficult to distinguish

72

between habitat and non-habitat, and individuals may be able to easily move among preferred

73

patches (Thomas and Kunin 1999; Franklin et al. 2002). Empirical studies based in patches of

74

large, open populations support the high extent and influence of movement on patch population

75

dynamics (Schaub et al. 2012, 2013) but in general few studies have examined the role of within

76

vs. among patch processes in these situations.

77

There are several potential mechanisms by which local individuals and immigrants might

78

interact to influence local abundance within a patch. One possibility is that immigration is

79

negatively density-dependent (Massot et al. 1992, Sæther et al. 1999) with higher immigration

80

following a greater loss of local residents due to mortality or emigration. Negative density-

81

dependent immigration has been observed in a range of terrestrial and aquatic taxa (Massot et al.

82

1992; Gundersson et al. 2002; Greathouse et al. 2005, Wilson and Arcese 2008; Schaub et al.

83

2010; Turgeon and Kramer 2012). Because immigrants are at least partially compensating for the

84

loss of local individuals, this process limits declines in abundance within the patch (Brown and

85

Kodric-Brown 1977; Stacey and Taper 1992; Püttker et al. 2011) and can help to stabilize the

86

broader population system (Sæther et al. 1999; Bowler and Benton 2005). Immigration may also

87

be positively density dependent and increase following years of higher local productivity or

88

survival. Positive density-dependent immigration may occur if individuals use conspecific

89

presence as an indicator of patch quality (Stamps and Krishnan 2005; Betts et al. 2008) but can

90

be destabilizing if immigration declines with decreasing densities leading to local patch

91

extinction (Sæther et al. 1999). Immigration may also be unrelated to local demographic

[Type text]

92

processes, for instance with populations that occupy ephemeral habitats where there is a high

93

degree of inter-annual movement in the system and little correlation between within and among

94

patch processes.

95

Our limited knowledge on how movement affects the dynamics of open populations with

96

high patch connectivity is due in part to the difficulty in estimating immigration rate, defined as

97

the proportion of the recipient population represented by immigrants. A recently developed

98

technique, integrated population modeling (IPM, Brooks et al. 2004; Schaub et al. 2007, Abadi et

99

al. 2010), allows for the estimation of immigration rate by combining data on abundance and

100

demographic rates under a single modeling framework. The estimate of apparent survival from

101

live encounter mark-recapture techniques incorporates mortality and permanent emigration

102

(Lebreton et al. 1992). With additional data on breeding productivity, the only unknown

103

demographic parameter that remains is immigration, but this can be estimated because the

104

abundance data in an IPM represents the sum of all demographic processes. This technique not

105

only allows for an unbiased estimate of immigration rate but also the potential to test hypotheses

106

on how immigration influences population dynamics (Abadi et al. 2010; Schaub et al. 2012).

107

We developed an IPM to examine the extent of immigration and how it influenced the

108

population dynamics and population-level productivity of a Neotropical migratory bird, the

109

American redstart (Setophaga ruticilla). This species has served as an important model for the

110

study of population limitation and regulation in migratory animals, with previous research

111

focusing primarily on factors that influence reproduction and survival (McKellar et al. 2014;

112

2015; Marra et al. 2015; Sherry et al. 2015). American redstarts are distributed across boreal and

113

eastern hardwood ecosystems of North America and select early to mid-successional deciduous

114

and mixed forests (Sherry and Holmes 1997). Like many other forest songbirds in these

[Type text]

115

ecosystems, American redstarts are not localized in discrete habitat patches but are more broadly

116

distributed across forest ecosystems containing different successional stages. Therefore, they

117

provide an excellent model species to understand the extent and influence of movement on the

118

dynamics of large, open populations.

119

We first used the IPM to quantify annual variability in the rates of immigration by age and

120

sex, and examined the contributions of immigration, reproduction and apparent survival to

121

annual population growth. We then tested whether immigration patterns were negatively or

122

positively density dependent as described by the hypotheses above. Finally, we examined

123

whether average breeding success declined when immigrants represented a higher proportion of

124

the study population. While the addition of immigrants can provide numerical stability to a

125

population, it is possible that average breeding success declines when immigrants represent a

126

larger fraction of the population. Several factors might contribute to lower productivity in this

127

situation, including increased agonistic interactions among immigrants and residents as new

128

territories are established (Smith and Ivins 1983; Hayes et al. 2004), lower familiarity of

129

immigrants to the local environment (Stamps 1987), and the possibility that immigrants are

130

lower quality individuals that were only able to secure territories after vacancies arose (McKellar

131

et al. 2013).

132

Methods

133

Study area and sampling methods

134

We conducted field work for this project at the Queen’s University Biological Station (QUBS),

135

near Chaffey’s Lock, ON, Canada (44°34’ N, 76°19’ W) from May–July 2001–2011. American

136

redstarts were monitored in a mixed-deciduous forest consisting primarily of sugar maple (Acer

137

accharum) and Eastern hop hornbeam (Ostrya virginiana). The ~100 ha study area consists of a

[Type text]

138

25 ha campground and a largely undisturbed 75 ha forest separated by a two lane country road

139

(for a detailed description see McKellar et al. 2015). The study site is nested within a broader

140

area of heterogeneous habitat consisting of forest, wetlands, lakes, agriculture and urban areas.

141

Hereafter we refer to the individuals in the 100 ha study area as the ‘study population’. American

142

redstarts exhibit protandry, with males arriving on the breeding grounds several days before

143

females. Upon arrival, males begin singing and establish territory boundaries. From May 1-31,

144

we surveyed the study population from 0600 to 1200, mapping territory boundaries of males,

145

recording whether females were present on territories and identifying any previously banded

146

individuals that returned. The abundance of males and females provided our annual count

147

estimates that were used in the IPM. We also recorded pairing dates and the location of the nest,

148

which we monitored to record laying date, number of eggs laid, hatching success, and fledging

149

success (see McKellar et al. 2014).

150

Birds were captured as soon as possible upon arrival using mist nets and simulated

151

territorial intrusions consisting of song playback and a decoy or fledgling distress calls. All

152

captured individuals were banded with a single Canadian Wildlife Service-issued aluminum leg

153

band and a unique combination of 2–3 colour bands. Nestlings were banded with a single

154

aluminum band, generally between days 5 and 8 after hatching. Age classes for males were

155

determined based on plumage coloration. Males exhibit delayed plumage maturation, with

156

second-year (SY) males exhibiting female-like grey and yellow plumage and after-second-year

157

(ASY) males exhibiting black plumage with orange patches on the wings, tail, and flanks, and a

158

white underside. Female age can only be reliably determined in the hand and was based on retrix

159

coloration and wear, and molt limits when individuals were captured (Pyle 1997).

160

Integrated population model

[Type text]

161

The IPM combined three data sets of demographic information: 1) annual breeding counts of

162

females and males in the study population, 2) annual number of fledged young produced by a

163

sample of breeding-age females and 3) an encounter history of colour-banded individuals to

164

estimate apparent survival (see below). The core of the IPM is a state-space model with a state

165

process that projected the true but unknown development of the population over time as a

166

function of the demographic rates and an observation process that linked the observed counts to

167

the true size assuming observation error.

168

We used a pre-breeding count with separate projections of population size for females

169

and males. We assumed that females bred for the first time at one year of age (i.e., as an SY)

170

(Sherry and Holmes 1997). As noted above, females are not easily distinguished by age, and we

171

previously found similar estimates of apparent survival among SY and ASY females at this site

172

(McKellar et al. 2015). Therefore, our annual population counts included breeding-age females

173

(Fcount) and the demographic information allowed us to estimate the annual number of female

174

fledglings (Ffl), local SY female recruits (Frec), local ASY returning females (Fret) and female

175

immigrants in at least their second year (Fimm). Because male American redstarts exhibit delayed

176

plumage maturation, we were able to distinguish annual male counts between second-year

177

(SYMcount) and after-second-year (ASYMcount) individuals (Sherry and Holmes 1997). This

178

distinction allowed us to estimate additional male age classes: annual number of male fledglings

179

(Mfl), local SY male recruits (SYMrec), SY male immigrants (SYMimm), local third year returning

180

males (TYMret), local after third year (ATY) returning males (ATYMret) and ASY male

181

immigrants in at least their second year (ASYimm).

182

Fledgling production was estimated as Nfl,t ~ Poisson (fratet · Ftot, t) where Nfl,t is the total

183

number of fledglings in year t, fratet is the number of fledglings per breeding female in year t and

[Type text]

184

Ftot, t is the number of breeding females in year t. The number of female fledglings (Ffl,t) in year t

185

was then estimated as a binomial process assuming a 50% probability of being female. Male

186

fledgling number (Mfl,t) was estimated as the difference between total and female fledglings.

187

For each sex, the abundance of local SY recruits was estimated as a binomial process

188

with the number of recruits in year t+1 dependent on the number of fledglings produced in year t,

189

a first year apparent survival probability of juveniles that was assumed to be the same for both

190

sexes due to low return rates (ϕjuv, t) and the probability that a returning juvenile was a male

191

(pMale). The encounter history for juveniles was based only on marked nestlings that were

192

known to have fledged. The number of returning individuals of each sex in year t+1 was also

193

estimated as a binomial process based on the number of SY males, ASY males and females in

194

year t and the apparent survival probabilities for each age-sex class from year t to t+1 (ϕm.sy,t,

195

ϕm.asy,t, and ϕf,t respectively). The total number of individuals in year t also included immigrants

196

that were assumed to enter the population just prior to breeding in year t and contributed to

197

reproduction in that year. Immigration was assumed to follow a Poisson process with Ia,s,t ~

198

Poisson (ωa,s,t) with ωa,s,t equal to the expected number of immigrants entering in year t for age

199

class a and sex s. We chose to estimate the number of immigrants rather than an immigration rate

200

because precision tends to be higher for the former (Schaub and Fletcher 2015). Immigration rate

201

was subsequently measured as a derived parameter (see below). The total number of individuals

202

for each sex at the start of each year t was then the sum of local recruits, returning local breeders

203

and immigrants, which in the case of males could be further distinguished between immigrants

204

born in the prior year (i.e., SY) and those in at least their second breeding year (i.e., ASY). The

205

observation model describes the relationship between the count of females, SY males, and ASY

206

males, and the true population size. We assumed a Poisson distribution Ca,s,t ~ Pois (Tota,s,t)

[Type text]

207

where Tot refers to the estimated annual population size of each age and sex class from the state

208

process model.

209

We also used the IPM to measure two derived parameters 1) population growth rate

210

based on breeding females, Fλ,t = Ftot, t+1 / Ftot, t and breeding males Mλ,t = Mtot, t+1 / Mtot,t and 2) the

211

rates of female, SY male, and ASY male immigration, where the rate is estimated as the number

212

of immigrants for the respective age-sex class divided by the total number of individuals in that

213

class at time t. This definition is based on our expectation that immigrants will enter the study

214

population just prior to breeding at the start of time interval t and is different from similar studies

215

that estimated immigration rate based on the number of individuals entering the population in

216

year t+1 relative to the number that were present in year t (e.g. Schaub et al. 2013, Duarte et al.

217

2016). Note that λt refers to population change between years t and t+1. We also estimated

218

female λ without immigration using the following estimation as a derived parameter within the

219

IPM:
Fapp.λ,t = fratet * 0.5 * ϕjuv,t + ϕf,t

220
221

Estimation of reproductive output and apparent survival

222

The time-varying demographic parameters were specified with a hierarchical structure where

223

annual estimates of each parameter were assumed to originate from a random process with a

224

common mean and temporal variance. For these estimates we used logit link functions for

225

apparent survival and a log link for fledge rate. Each year we had complete reproductive fates for

226

a sample of all territories in the study area. The nesting attempts from this sample were used to

227

determine the annual fledging rate (fledget) ~ Poisson (fratet · Frep,t) where the number of fledged

228

young in the sample of territories (fledget) was the product of the fledge rate (fratet) and the

229

number of females monitored for reproduction (Frep,t).

[Type text]

230

Age and sex-specific survival was estimated with a state-space likelihood and a multi-

231

state representation (Lebreton et al. 1999, 2009) for juveniles (ϕjuv), adult females (ϕf), SY males

232

(ϕm.sy) and ASY males (ϕm.asy). With the multi-state representation for apparent survival, the state

233

equation represents the state (i.e. age-sex class) of the individual:

234

Zi,t|Zi,t-1 ~ categorical (ΩZi,t-1, i, t, 1…S)

235

where i is an individual, t is a time period and (ΩZi,t-1, i, t, 1…S) is a matrix of survival and

236

transition probabilities from each state at time t-1 to each state at time t with S representing the

237

number of true states. Initially, an individual can move from any state to any other state but

238

constraints are applied to the matrix probabilities such that biologically unreasonable transitions

239

are not allowed (e.g. probability SY male to female = 0). The observation equation represents the

240

observation of individual i at time t given its state at time t:

241

yi,t|Zi,t ~ categorical (ΘZi,t, i, t,1…O)

242

with O equal to the number of observed states. If an individual survives and returns to the study

243

area it may be detected with re-sighting probability p. We used the earlier findings from

244

McKellar et al. (2015) to inform variation in re-sighting probability and apparent survival by age

245

and sex in this analysis. Re-sight probability does not differ between age classes but does differ

246

between the sexes as males are more easily detected. Therefore, annual estimates of re-sight

247

probability were estimated with a mean (μ) and temporal residual (ε) for each sex:
logit(ps,t) = μp,s + εp,s,t

248
249

We previously found that apparent survival differs between the sexes and between SY

250

males and ASY males but not between the two female age classes (McKellar et al. 2015).

251

Therefore, apparent survival for the IPM was estimated as:
logit(ϕa,s,t) = μϕ,a,s + ε ϕ,a,s,t

252

[Type text]

253

where annual survival for each age class a (males only) and sex s is estimated with a mean and

254

temporal residual.

255

Model implementation

256

We fit the IPM using a Bayesian analysis with Markov Chain Monte Carlo Sampling in

257

OpenBUGS version 3.2.3 (Lunn et al. 2000) implemented through R (version 3.2.3) using the

258

package R2OpenBUGS (R Development Core Team 2004, Sturtz et al. 2005). With this

259

approach, inference is based on sampling from the posterior distribution, which is proportional to

260

the prior distribution and the likelihood from the data. We used non-informative priors for all

261

model parameters except starting population size where the priors were weakly informative. We

262

ran two chains for 25,000 iterations, discarded the first 15,000 as a burn-in and used a thinning

263

rate of 2 to give 10,000 samples from the posterior distribution for inference. We assessed model

264

convergence through the parameter history plots and R-hat convergence diagnostics (Sturtz et al.

265

2005). There is currently no goodness-of-fit tests for integrated population models (Schaub and

266

Abadi 2011). To evaluate the fit of the mark-recapture model we used a bootstrap goodness of fit

267

test in program MARK on a model with apparent survival varying by age, sex and time, and re-

268

sighting probability varying by sex and time. The estimated c-hat from this model was 1.08

269

indicating no evidence for a lack of fit (Lebreton et al. 1992).

270

Assessing relationships between abundance, demographic rates and population growth

271

We used an approach similar to that employed by Schaub et al. (2012, 2013) to assess

272

relationships among demographic rates and between demographic rates and population growth.

273

For all samples from the posterior distributions we estimate the correlation between the two

274

variables of interest thus incorporating the sampling uncertainty in the estimation of the

275

variables. This approach was used for two specific sets of tests. First, we examined the

[Type text]

276

relationship between each rate at time t and population growth from time t to t+1. Second, we

277

tested for the influence of density dependence on immigration with a correlation between

278

immigration rate in year t+ 1 and the apparent survival probability of local males and females

279

between years t and t+1. Because a change in apparent survival would also influence the number

280

of previous breeders returning to the study population in the following year, we also tested the

281

relationship between the immigration rate in year t+1 and the number of local individuals

282

returning in year t+1. Posterior distributions for the correlation coefficients are typically skewed

283

and we present the posterior median of the correlation coefficient and the percent of the

284

coefficient posterior mass above 0. For correlations of demographic rates and population growth,

285

values near 1 reflect a high probability that a positive change in the demographic rate is

286

correlated with positive population growth. For density dependent relationships, values near 1

287

and 0 indicate a high probability of positively density-dependent and negatively density-

288

dependent relationships respectively.

289

Finally, we tested whether the proportion of immigrant males and females in year t had a

290

negative influence on population level productivity (mean fledge rates) in year t. Because

291

reproductive output at the population level is also influenced by density (McKellar et al. 2014),

292

we included male abundance from the IPM as a covariate in this analysis. We used a similar

293

approach as described above and used the simulations from the posterior distribution in a linear

294

model with fledging rate as the response variable and, immigration rate and male abundance as

295

predictors. This effect was tested for both female and male immigration rate and we report the

296

95% CI for the beta coefficient in each case.

297

Results

298

Population growth and demographic rates

[Type text]

299

Mean population growth (λ) in the QUBS study population was 1.04 for both females (95% CI:

300

1.00, 1.08) and males (95% CI: 0.99, 1.10) with much of the increase at the start of the study

301

period (Figure 1). The estimated total population size varied from a low of 65 individuals to a

302

high of 94 although for 9 of 11 years it fluctuated between 73 and 84 individuals. The fledging

303

rate averaged 1.69 (95% CI: 1.23, 2.35) fledglings per female and varied from a high of 2.94 in

304

2001 to a low of 0.68 in 2011. The apparent survival probability of juveniles was low, averaging

305

only 0.06 (95% CI: 0.03, 0.10), while the probability that a returning juvenile that was detected

306

in the study area was a male was 0.68 (95% CI: 0.42, 0.90). Of 464 banded nestlings that

307

survived to fledge, only 16 returned to the study population and were detected in a subsequent

308

year with an average of 4.2 local female recruits (range: 2 - 7) and 2.3 local male recruits (range:

309

1 - 4) per year. Average apparent annual survival for SY and ASY females combined was 0.38

310

(95% CI: 0.25, 0.53) and averaged higher than probabilities for males where the mean was 0.34

311

(95% CI: 0.21, 0.48) for ASY males and only 0.20 for SY males (95% CI: 0.11, 0.37). Re-

312

sighting probabilities were higher for males at 0.74 (95% CI: 0.47, 0.92) than for females at 0.34

313

(95% CI: 0.17, 0.56). Estimates of temporal process variance for all demographic parameters are

314

provided in Appendix A.

315

For females, the estimated average immigration rate (immigrants/total females) was 0.52

316

(95% CI: 0.37, 0.65, Figure 1). Across years, the female immigration rate varied from a low of

317

0.46 to a high of 0.59 with rates >0.5 in all but two years indicating that in most years the

318

majority of breeding-age females came from outside the study population (Figure 1, Appendix

319

A). Average immigration rates were higher for males at 0.78 (95% CI: 0.59, 0.91) for SYs and

320

0.60 (95% CI: 0.49, 0.69) for ASYs. The number of immigrant males always exceeded the

321

number of local recruits and returning breeders (Figure 1). The estimated mean lambda for

[Type text]

322

females based only on within population rates of reproduction and survival was 0.44 (95% CI:

323

0.33, 0.56) and with this estimate, an initial population size of 40 females would drop to less than

324

one in only seven years if immigration ceased. Moreover, the highest estimate of lambda without

325

immigration for any year was 0.50 (95% CI: 0.32, 0.74) in 2005 and thus the study population

326

required immigrants to prevent a decline (i.e. to maintain λ > 1) in all years.

327

Contributions of demographic rates to population growth

328

Annual fledge rate in year t had the strongest contribution to population growth of females

329

(r=0.43) and males (r=0.44) with probabilities of a positive correlation for 94 and 95 percent of

330

simulations (Table 1). A weaker correlation was observed between population growth and adult

331

apparent survival of females (r=0.34) and ASY males (r=0.35) with probabilities of a positive

332

correlation for 81 and 86 percent of simulations respectively. There was no evidence for

333

correlations with population growth for apparent survival of juveniles, SY males or immigration

334

rate for any age and sex class (all r < 0.2 and probability of a positive correlation < 0.7). Thus,

335

despite being essential for maintaining the size of the QUBS study population, immigration did

336

not consistently contribute to higher population growth for either sex (Table 1).

337

Density-dependent effects on immigration

338

To test our hypotheses for density-dependent immigration, we examined whether the

339

immigration rate for each age and sex class at the start of year t+1 was correlated with the

340

survival of local breeders from year t to t+1 and whether immigration rate was correlated with

341

the number of local breeders at the start of year t+1. We found evidence for negative density-

342

dependent effects in both cases. The median correlation between apparent survival from year t to

343

t+1 and immigration rate in year t+1 for ASY males was r = -0.59 (Pr > 0 = 0.06, Figure 2) and

344

the correlation between ASY male immigration rate and the number of local males returning was

[Type text]

345

r = -0.84 (Pr > 0 = 0.0007). For females the correlation between apparent survival and

346

immigration rate was -0.49 (Pr > 0 = 0.14) and between the immigration rate and the number of

347

local females was -0.86 (Pr > 0 = 0.0004). Because immigration rate is defined as a proportion of

348

total individuals, the strong correlation between immigration rate and returning local individuals

349

is to be expected if the number of immigrants remains constant or increases as the number of

350

returning local individuals declines. In contrast to these patterns, the immigration rate of SY

351

males was not correlated with apparent survival of SY (r = -0.02, Pr > 0 = 0.47) or ASY males (r

352

= -0.04, Pr > 0 = 0.45) and only weakly correlated with the number of returning local males (r = -

353

0.34, Pr > 0 = 0.16).

354

Effects of immigration rate on population level productivity

355

We hypothesized that productivity would be lower in years when a higher proportion of the

356

population was comprised of immigrants but found no support for this hypothesis. The beta

357

coefficient for the effect of female immigration rate on fledging rate was 𝛽" = 0.38 (95% CI: -

358

3.43, 4.44) and for the effect of male immigration rate was 𝛽" = 0.20 (95% CI: -3.06, 3.49).

359

Discussion

360

Immigration has been one of the most difficult demographic rates to study, particularly in open

361

populations where there is a high degree of connectivity among patches of suitable habitat. Using

362

an integrated population model, we estimated annual rates of immigration and tested hypotheses

363

on how immigration and within-patch processes (reproduction, apparent survival) interact to

364

influence patch dynamics for a Neotropical migratory bird breeding in eastern North America.

365

Immigration rates were high with immigrants representing on average 52% of females, 60% of

366

after-second-year males and 78% of second-year males annually. We found strong support for

367

the hypothesis that immigration was negatively density dependent with immigrants replacing

[Type text]

368

previously established breeders in a compensatory manner following their death or emigration.

369

Because of this relationship, immigration rate was not correlated with growth of the study

370

population. However, if immigration were to cease, the study population would become locally

371

extinct in seven years and thus immigration provided an essential rescue effect (Brown and

372

Kodric-Brown 1977) allowing for stability around a carrying capacity despite low juvenile

373

apparent survival and strong annual variability in adult apparent survival. American redstarts

374

exhibit conspecific attraction (Hahn and Silverman 2006) and we have shown previously that

375

small scale movements within the study area are common with individuals more likely to move

376

towards areas with higher density (McKellar et al. 2015). However, at the entire patch scale there

377

was no support for positive density-dependent immigration. While individuals may select

378

breeding sites in closer proximity to one another, our results seem to indicate that there is still a

379

broader patch scale carrying capacity that regulates immigration rate even if breeders are not

380

uniformly distributed throughout the patch.

381

Our findings are similar to those from studies on other territorial species showing

382

negative density-dependent immigration into a population. Schaub et al. (2013) found high rates

383

of immigration for females (0.56 + 0.02) and males (0.43 + 0.02) in an open population of red-

384

backed shrikes (Lanius collurio) with immigration rates negatively correlated with population

385

density. Wilson and Arcese (2008) showed that male immigration rate in a song sparrow

386

(Melospiza melodia) island metapopulation declined with increasing population size on the

387

recipient islands. Experimental approaches have also been used to test the influence of density on

388

immigration rate for non-avian taxa. Turgeon et al. (2012) found strong evidence for

389

compensatory immigration in response to localized mortality in two damselfish species

390

(Stegastes spp.) in Barbados. They also found that for one species, compensatory immigration

[Type text]

391

was related to a perceived increase in habitat quality for immigrants and that it had an effect on

392

the size of immigrant source populations. Gundersson et al. (2002) manipulated patch densities

393

of root voles (Microtus oeconomus) and showed that immigration was negatively density-

394

dependent and that individuals most likely to survive after immigration were those that had

395

previously dispersed and successfully established in other patches. Other experimental examples

396

of an increase in immigration following declines in local density include common lizards

397

(Lacerta vivipara) in France (Massot et al. 1992) and freshwater fish in streams in Virginia, USA

398

(Albanese et al. 2009).

399

Our study population is a strong sink that is dependent on immigration for persistence

400

(Pulliam 1988) but how we interpret the reliance on immigration depends on the extent to which

401

low apparent survival was driven by mortality or permanent emigration of local breeders (Runge

402

et al. 2006). Our estimates of annual apparent survival are low relative to previous estimates for

403

this species from other populations. Male apparent survival in our study population averaged

404

only 0.20 and 0.34 for SY and ASY males respectively with substantial annual variation from

405

0.20 to 0.55 for the latter. Annual apparent survival for male American redstarts is more

406

frequently in the 0.50 to 0.60 range (Johnson et al. 2006, Marra et al. 2015). Our estimate of

407

female apparent survival was 0.38, which is more similar to, yet still lower than, estimates from

408

these other studies at 0.40 to 0.50. Although poorly studied, average apparent juvenile survival is

409

often in the 15-35% range for migratory passerines (e.g., Gardali et al. 2003; Duarte et al. 2014),

410

but averaged only 6% in our study population. The low estimates of apparent survival for all age

411

and sex classes suggest that emigration was likely high and especially so for males and juveniles.

412

The causes and extent of inter-annual variability in emigration are poorly known but other

413

studies have shown that emigration rates vary with prior reproductive success (Cline et al. 2013),

[Type text]

414

environmental conditions on the breeding grounds (Rushing et al. 2015, Wilson et al. 2016) and

415

population density (Kim et al. 2009).

416

We do not know the source of immigrants that filled vacancies in the habitat following

417

mortality or permanent emigration of local breeders and recruits. Within the study population,

418

the majority of inter-annual movements for established breeders that switched territories was less

419

than 100 m but for movements greater than 100 m there was more variability in distance moved

420

out to the limit of monitoring at 1.4 km (McKellar et al. 2015). This finding might suggest that

421

most immigrants are from sites near the study area. Experimental studies of this species on the

422

breeding and wintering grounds have shown that individuals will upgrade to higher quality

423

habitat following the disappearance of territorial birds (Studds and Marra 2005; McKellar et al.

424

2013). After controlling for density-dependent effects on reproductive output we also found no

425

evidence that years with a higher proportion of these immigrants led to a reduction in

426

productivity at the population level. This appears to indicate that immigrants are able to quickly

427

adapt to the new habitat and/or that they are not necessarily lower quality individuals than the

428

local individuals they replaced. Additional study on reproductive success of immigrants and

429

local birds at the individual level would be a useful complement to our study here at the

430

population level.

431

Conclusions

432

Understanding the extent and scale of immigration and emigration has important implications for

433

conservation because it affects the spatial scale at which monitoring and habitat protection needs

434

to be considered for a population (Schaub et al. 2012, Duarte et al. 2016). If regular natal and

435

breeding movement occurs among patches then we must not only consider management needs

436

for a local patch where individuals breed but also for the broader landscape over which dispersal

[Type text]

437

of juveniles and adults regularly occurs. The extent to which populations are synchronized over

438

broader scales is related to the immigration rate per generation. In general, populations within

439

patches are not independent if, on average, individuals occupy more than one habitat patch over

440

their lifetime (Harrison 1991, Harrison and Taylor 1997). The movement rates in our study

441

exceeded this level and thus not only would immigration have a strong influence on our study

442

population but would also be expected to synchronize the dynamics of the broader population

443

system (Hastings 1993).

444

Because our study took place in a single patch covering ~ 1km2 we were not able to

445

estimate how the strength and influence of within vs among-patch processes changes with spatial

446

scale. This would be an interesting and important area for future study and integrated population

447

models provide a methodology for this approach without having to track the movements of

448

individuals across a large area. Collecting demographic data over large scales is difficult but it

449

would be practical to use methods that account for detection probability to estimate abundance

450

(Buckland et al. 2001, Royle 2004) along with a sample of data on reproduction and apparent

451

survival over several years and with increasing spatial scale increments (e.g. 0.5 km2, 1.5 km2,

452

2.5 km2). This approach would only require a study area covering the largest scale but with

453

sufficient sample sizes to estimate rates separately at the smaller scales. An IPM could then be

454

used to estimate immigration rate and the contributions of immigration, reproduction and

455

apparent survival to population growth at each scale, and thus allow us to estimate how

456

movement and the relative importance of within vs among-patch processes changes as the spatial

457

scale increases.

458

Acknowledgements

[Type text]

459

We thank the many field assistants, technicians, and students who worked on this project over

460

the years and an anonymous reviewer for providing valuable comments on an earlier version of

461

this manuscript. Funding was provided by the National Science Foundation (P.P.M) and Natural

462

Sciences and Engineering Research Council of Canada (L.M.R. and students).

463

Literature

464

Abadi F, Gimenez O, Ullrich B, Arlettaz R, Schaub M (2010). Estimation of immigration rate

465
466
467
468
469

using integrated üopulation models. J Appl Ecol 47:393–400
Addicott JF, Aho JM, Antolin MF, Padilla DK, Richardson JS, Soluk DA (1987) Ecological
neighborhoods: scaling environmental patterns. Oikos 49:340–346
Albanese B, Angermeier PL, Peterson JT (2009) Does mobility explain variation in colonisation
and population recovery among stream fishes? – Freshwater Biol 54:1444–1460

470

Betts MG, Hadley AS, Rodenhouse N, Nocera JJ (2008) Social information trumps vegetation

471

structure in breeding-site selection by a migrant songbird. Proc R Soc B 275:2257–2263

472

Bowler DE, Benton TG (2005) Causes and consequences of animal dispersal strategies: relating

473
474
475
476
477
478
479
480
481

individual behaviour to spatial dynamics. Biol Rev 80:205–225
Bowne DR, Bowers MA (2004) Interpatch movements in spatially structured populations: a
literature review. Landscape Ecol 19:1–20
Brooks SP, King R, Morgan BJT (2004) A Bayesian approach to combining animal abundance
and demographic data. Anim Biodivers Conserv 27:515–529
Brown JH, Kodric-Brown A (1977) Turnover rates in insular biogeography: effect of
immigration on extinction. Ecology 58:445–449
Buckland ST, Anderson DR, Burnham KP, Laake JL, Borchers DL, Thomas L (2001)
Introduction to Distance Sampling. Oxford University Press, Oxford.

[Type text]

482
483

Camus PA, Lima M (2002) Populations, metapopulations and the open-closed dilemma: the
conflict between operational and natural concepts. Oikos 97:433–438

484

Cline MH, Strong AM, Sillett TS, Rodenhouse NL, Holmes RT (2013) Correlates and

485

consequences of breeding dispersal in a migratory songbird. Auk 130:742–752

486

Duarte A, Hines JE, Nichols JD, Hatfield JS, Weckerly FW (2014) Age-specific survival of male

487

Golden-cheeked Warblers on the Fort Hood Military Reservation, Texas. Avian Cons Ecol

488

9:4.

489

Duarte A, Weckerly FW, Schaub M, Hatfield JS (2016) Estimating golden-cheeked warbler

490

immigration: implications for the spatial scale of conservation. Anim Cons 19:65-74

491

Fahrig L, Merriam G (1985) Habitat patch connectivity and population survival. Ecology

492
493
494

66:1762–1768
Franklin AB, Noon BR, George TL (2002) What is habitat fragmentation? Stud Avian Biol
25:20–29

495

Gardali T, Barton DC, White JD, Geupel GR (2003) Juvenile and adult survival of Swainson’s

496

Thrush (Catharus ustulatus) in coastal California: annual estimates using capture-recapture

497

analyses. Auk 120:1188–1194.

498
499

Greathouse EA, March JG, Pringle CM (2005) Recovery of a tropical stream after a harvestrelated chlorine poisoning event. Freshwater Biol 50:603–615

500

Gunderssen G, Andreassen HP, Ims RA (2002) Individual and population level determinants of

501

immigration success on local habitat patches: an experimental approach. Ecol Lett 5:294–

502

301

503
504

Hahn BA, Silverman ED (2006) Social cues facilitate habitat selection: American redstarts
establish breeding territories in response to song. Biol Lett 11:337–340

[Type text]

505
506

Hanski I, Gilpin ME (1997) Metapopulation Biology: Ecology, Genetics and Evolution.
Academic Press, San Diego.

507

Hanski I (1998) Metapopulation dynamics. Nature 396:41–49

508

Harrison S (1991) Local extinction in a metapopulation context: an empirical evaluation. Biol J

509

Linnean Soc 42:73–88

510

Harrison S, Taylor AD (1997) Empirical evidence for metapopulation dynamics. In: Hanski I,

511

Gilpin ME (eds) Metapopulation Biology: Ecology, Genetics and Evolution. Academic

512

Press, San Diego, pp 27–42

513
514
515

Hastings A (1993) Coupled interactions between dispersal and dynamics: lessons from coupled
logistic equations. Ecology 72:896–903
Hayes LD, Lin YK, Solomon NG (2004) The effect of female prairie vole (Microtus

516

ochrogaster) immigrants on space use of conspecific female residents. Am Midl Nat 151:88–

517

92

518
519
520
521

Johnson MD, Sherry TW, Holmes RT, Marra PP (2006) Assessing habitat quality for a wintering
songbird in natural and agricultural areas. Conserv Biol 20:1433–1444
Kim SY, Torres R, Drummond H (2009) Simultaneous positive and negative density-dependent
dispersal in a colonial bird species. Ecology 90:230–239

522

Lebreton JD, Burnham KP, Clobert J, Anderson DR (1992) Modelling survival and testing

523

biological hypotheses using marked animals: a unified approach with case studies. Ecol

524

Monogr 62:67–118

525
526

Lebreton JD, Almeras T, Pradel R (1999) Competing events, mixtures of information and
multistratum recapture models. Bird Stud 46:S39–S46

[Type text]

527

Lebreton JD, Nichols JD, Barker RJ, Pradel R, Spendelow JA (2009) Modeling individual

528

animal histories with multistate capture-recapture models. In Caswell, H (ed): Adv. Ecol.

529

Res.Vol. 41, Burlington: Academic Press, 2009, pp. 87–173.

530
531

Lunn DJ, Thomas A, Best N, Spiegelhalter D (2000) WinBUGS – a Bayesian modeling
framework: concepts, structure, and extensibility. Stat Comput 10:25–37

532

Marra PP, Studds CE, Wilson S, Sillett TS, Sherry TW, Holmes RT (2015) Non-breeding season

533

habitat quality mediates the strength of density dependence for a migratory bird. – Proc R

534

Soc B 282 DOI: 10.1098/rspb.2015.0624

535

Massot M, Clobert J, Pilorge T, Lecomte J, Barbault R (1992) Density dependence in the

536

common lizard: demographic consequences of a density manipulation. Ecology 73:1742–

537

1756

538

McKellar AE, Marra PP, Ratcliffe LM (2013) Starting over: experimental effects of breeding

539

delay on reproductive success in early-arriving male American redstarts. J Avian Biol

540

44:495–503

541
542
543
544

McKellar AE, Marra PP, Boag PT, Ratcliffe LM (2014) Form, function and consequences of
density dependence in a long-distance migratory bird. Oikos 123:356–364.
McKellar AE, Reudink MW, Marra PP, Ratcliffe LM, Wilson S (2015) Climate and density
influence annual survival and movement in a migratory songbird. Ecol Evol 5:5892–5904

545

Pulliam HR (1988) Sources, sinks and population regulation. Am Nat 132:652–661

546

Püttker T, Bueno AA, Dos Santos de Barros C, Sommer S, Pardini R (2011) Immigration rates in

547

fragmented landscapes – empirical evidence for the importance of habitat amount for species

548

persistence. PLOS One 6(11): e27963

549

Pyle P (1997) Identification guide to North American birds. Part 1. Slate Creek Press, Bolinas.

[Type text]

550
551
552
553
554
555
556

R Development Core Team (2004) R: a language and environment for statistical computing. R
Foundation for Statistical Computing, Vienna, Austria. www.r-project.org.
Royle JA (2004) N-mixture models for estimating population size from spatially replicated
counts. Biometrics 60:108–115
Runge JP, Runge MC, Nichols JD (2006) The role of local populations within a landscape
context: defining and classifying sources and sinks. Am Nat 167:925–938
Rushing CS, Dudash MR, Studds CE, Marra PP (2015) Annual variation in long-distance

557

dispersal driven by breeding and non-breeding season climatic conditions in a migratory

558

bird. Ecography 38:1006–1014

559
560
561

Sæther B-E, Engen S, Lande R (1999) Finite metapopulation models with density-dependent
migration and stochastic local dynamics. Proc R Soc B Lond 266:113–118
Sanz-Aguilar A, Igual J, Tavecchia G, Genovart M, Oro D (2016) When immigration masks

562

threats: The rescue effect of a Scopoli’s shearwater colony in the Western Mediterranean as

563

a case study. Biol Conserv 198:33–36

564

Schaub M, Gimenez O, Sierro A, Arlettaz R (2007) Use of integrated modeling to enhance

565

estimates of population dynamics obtained from limited data. Conserv Biol 21:945–955

566

Schaub M, Aebischer A, Gimenez O, Berger S, Arlettaz R (2010). Massive immigration balances

567

high human induced mortality in a stable eagle owl population. Biol Conserv 143:1911–

568

1918

569
570
571
572

Schaub M, Reichlin TS, Abadi F, Kéry M, Jenni L, Arlettaz R (2012) The demographic drivers
of local population dynamics in two rare migratory birds. Oecologia 168:97–108
Schaub M, Jakober H, Stauber W (2013) Strong contribution of immigration to local population
regulation: evidence from a migratory passerine. Ecology 94:1828–1838

[Type text]

573
574

Schaub M, Fletcher D (2015) Estimating immigration using a Bayesian integrated population
model: choice of parameterization and priors. Environ Ecol Stat 22:535–549

575

Sherry TW, Holmes RT (1997) American Redstart (Setophaga ruticilla), The Birds of North

576

America Online (A. Poole, Ed.). Ithaca: Cornell Lab of Ornithology; Retrieved from the

577

Birds of North America Online: http://bna.birds.cornell.edu/bna/species/277.

578

Sherry TW, Wilson S, Hunter S, Holmes RT (2015) Impacts of nest predators and weather on

579

reproductive success and population limitation in a long-distance migratory songbird. J

580

Avian Biol 46:559–569

581
582
583
584
585
586

Smith AT, Ivins BL (1983) Colonization in a pika population: Dispersal vs. philopatry. Behav
Ecol Sociobiol 13:37–47
Stacey PB, Taper M (1992) Environmental variation and the persistence of small populations.
Ecol Appl 2:18–29
Stamps JA (1987) The effect of familiarity with a neighborhood on territory acquisition. Behav
Ecol Sociobiol 21:273–277

587

Stamps J, Krishnan VV (2005) Nonintuitive cue use in habitat selection. Ecology 86:2860–2867

588

Studds CE, Marra PP (2005) Nonbreeding habitat occupancy and population processes: an

589
590
591

upgrade experiment with a migratory bird. Ecology 86:2380–2385
Sturtz S, Ligges U, Gelman A (2005) R2WinBUGS: A Package for Running WinBUGS from R.
J Stat Softw 12:1–16.

592

Thomas CD, Kunin WE (1999) The spatial structure of populations. J Anim Ecol 68:647–657

593

Turgeon K, Kramer DL (2012) Compensatory immigration depends on adjacent population size

594

and habitat quality but not on landscape connectivity. J Anim Ecol 81:1161–1170

[Type text]

595

Vilá C, Sundqvist AK, Flagstad Ø, Seddon J, Björnerfeldt S, Kojola I, Casuli A, Sand H,

596

Wabakken P, Ellegren H (2003) Rescue of a severely bottlenecked wolf (Canis lupus)

597

population by a single immigrant. Proc R Soc B 270:91–97

598
599
600

Wilson AG, Arcese P (2008) Influential factors for natal dispersal in an avian island
metapopulation. – J Avian Biol 39:341–347
Wilson S, Alisauskas RT, Kellett DK (2016). Environmental and population effects on

601

permanent emigration of Ross’s and Snow Geese from an Arctic breeding area. J Wildlife

602

Manage 80:117–126

[Type text]

1

Table 1. Median correlation and probability of a positive correlation (Pr > 0) between

2

demographic rates and population growth of American redstarts based on 10,000 replicate

3

simulations from the posterior distribution of each rate. Correlations for adult survival and

4

immigration are only included for growth rate of the respective sex.
Rate

Correlation

Pr > 0

Correlations with female population growth (Fλ)
Fledge rate

0.43

0.94

Juvenile survival

0.08

0.59

Female survival

0.34

0.81

Female immigration rate

0.17

0.65

Correlations with male population growth (Mλ)

5
6
7
8
9
10

Fledge rate

0.44

0.95

Juvenile survival

0.03

0.54

SY male survival

-0.07

0.41

ASY male survival

0.35

0.86

SY male immigration rate

-0.10

0.39

ASY male immigration rate

0.17

0.67

11

Figure Legends

12
13

Figure 1. Estimated mean population size and 95% credible intervals for female, after-second-

14

year male and second-year male American redstarts at Queen’s University Biological

15

Station, 2002-2011. Local recruits and returning breeders are shown with solid circles,

16

immigrants with open squares.

17
18

Figure 2. Posterior distribution histogram plots for the correlation of immigration rate with

19

annual apparent survival and local abundance for ASY male and female American redstarts

20

at Queen’s University Biological Station, 2001-2011. Estimates of the correlation are based

21

on 10,000 replicates from the posterior distributions. The dashed line in each panel

22

represents the median correlation coefficient.

23
24
25

[Type text]

Figure 1.

[Type text]

Figure 2.

[Type text]