Roosting ecology and the evolution of bat landing maneuvers

Biomechanics is poised at the intersection of organismal form, function, and ecology, and forms a practical lens through which to investigate evolutionary linkages among these factors. We conducted the first evolutionary analysis of bat flight dynamics by examining the phylogenetic patterning of landing mechanics. We discovered that bats perform stereotyped maneuvers that are correlated with landing performance quantified as impact force, and that these are linked with roosting ecology, a critical aspect of bat biology. Our findings suggest that bat ancestors performed simple, four-limbed landings, similar to those performed by gliding mammals, and that more complex landings evolved in association with novel roost types. This explicit connection between ecology and biomechanics presents the opportunity to identify traits that are associated with a locomotor behavior of known ecological relevance, thus laying the foundation for a broader understanding of the evolution of flight and wing architecture in this extraordinarily successful mammalian lineage.


Introduction 24
Morphologists and biomechanicians often study organismal evolution as a function of 25 three interrelated factors: structure (morphology), function (mechanics or behavior), and context 26 (ecology). Detecting linkages among traits from these categories and discerning where trait shifts 27 correspond with patterns of diversification not only provides evidence of selection but can also 28 point to specific drivers of adaptive radiations, which are a central phenomenon in evolution 29 whereas low-impact landings (two-point) are associated with roosting on stiff horizontal surfaces, 86 such as cave ceilings or tree hollows. 87 The broad biological importance of roosting ecology and interspecific variation in landing 88 mechanics offers an opportunity to discover how the biomechanical basis of landing performance 89 may underlie how bats take refuge and disperse within their environment. In the present study, we 90 Table 1: Study taxa, landing style observed, peak landing impact forces, and roosting ecology category. Roost categories are: cavity in standing tree (CST), exposed   figure 2. See Figure 1-Source Data 1 for raw impact force measurements for each landing and Figure  1-Source Data 2 for each individual's mean peak impact force, the latter of which was used to generate these figures. See Source Code File 1 for R code to reproduce these plots.

Landing impact force increases with points of contact across landing styles 132
Two-point landings from 32 individuals of 15 species, from 4 families, uniformly resulted 133 in low peak impact forces, with a mean of 0.95 ± 0.54 BW (mean ± s.d.,). Four-point landings 134 Among the taxa we sampled, we detected three independent shifts to from four-to two- We investigated the relationship between roosting ecology and landing style by using phylogenetic 178 logistic regression to compare landing style with roosting ecology using alternative roost 179 classification schemes, which aggregated roost categories with similar physical characteristics (table  180 2). Compared to the null model, our aggregated model had greater explanatory power for 181 predicting landing style from roosting ecology, as indicated by AIC score. Our null model, which 182 tested for association between roosting habits and landing style using the most common roost type 183 for each species, revealed a significant positive association only between cavity-roosting and two-184 point landings (ß NullCST-2pt = 3.862; p Null,CST-2pt = 0.01701). Aggregating roosting categories according 185 to physical properties, such as compliance, orientation, and spatial constraint, allowed us to test 186 the hypothesis that these physical properties are significantly associated with the mechanics of the 187 three known landing styles. We found that two-point landings were positively associated with stiff,  Table 2: Correlations between landing style and roosting ecology from phylogenetic logistic regressions. We provide Firth-corrected coefficient estimates (ß) with 194 bootstrapped 95% confidence intervals (in brackets) and Wald p-values (in italics) to denote significant associations between roost type. Significant p-values are 195 bolded and set within shaded cells. AIC scores provide comparison between our Null and Aggregated model (models with smaller AIC are preferred; differences are 196 meaningful when ≥ 2). P-values are conditional upon phylogenetic signal, a, where values near 0 denote strong phylogenetic signal and values approaching 1 197 indicate weak phylogenetic signal. Roosting ecology categories correspond with those listed in Table 1: cavity in standing tree (CST), exposed on standing tree (EST), 198 rocks and/or caves (R/C), termite or ant nests (TAN), foliage-leaf tent (FOL-LT), unmodified foliage (FOL-UF), and rock crevices (CREV). See Table 2

Discussion 205
Using a combination of field and lab-based measurements, we investigated functional links 206 between landing mechanics and roosting ecology, which is a critical biological factor for bats. Our

Other factors that may influence landing style 303
We focused on associations between landing mechanics and roosting ecology in the 304 present study, but other traits could also influence diversity of landing maneuvers among bats. 305 Here we highlight a couple, including sensory ecology and wing morphology. 306 A bat's ability to sense the location, geometry, and surface characteristics of a potential 307 landing site contributes to its capacity to execute accurate, precise landings. Therefore, variation in 308 sensory ecology, specifically echolocation capacity and call structure, could also influence the 309 landing maneuvers of bats. Most bats navigate their environments and detect prey using laryngeal 310 approaching stationary targets, such as roosts (but see (Tian and Schnitzler, 1997)) Interspecific 318 variation in echolocation behavior during landing could reveal patterns that coincide with 319 differences in impact forces and body rotations as bats call to sense the roost during approach. 320 Additionally, in the Pteropodidae, which do not possess laryngeal echolocation, landing behavior 321 could be constrained due to sensory limitations in their capacity to resolve details of potential 322 roosts with high temporal resolution during an approach flight. 323 Interspecific differences in wing morphology may also relate to variation in landing 324 mechanics. Because aerodynamic forces are highly dependent upon the velocity of airflow over the 325 wings, and landing occurs at low speeds, bats accomplish landing maneuvers using inertial forces 326 almost exclusively (Bergou et al., 2015). The wing's capacity to effect body rotation via inertial 327 torques is therefore related to their mass moment of inertia, which in turn is determined by the 328 distribution of mass within the wings. Studies that characterize interspecific differences in wing 329 mass distribution may therefore reveal a relationship between the wing's body mass-normalized 330 mass moment of inertia and the rotational complexity of landing maneuvers. For example, 331 variation in wing inertia could arise from interspecific differences in wing length or relative mass 332 of the bones, muscle and skin that comprise the wing, particularly in the distal regions. 333

334
Estimated ancestral landing mechanics provide support for a gliding bat ancestor 335 Approaches that eliminate the measurement of impact forces will allow for broader sampling with 394 videography because some species are difficult or impossible to train to land on a small force 395 platform. Furthermore, recording landing videos at known roost locations rather than with 396 captured individuals in a field-based flight arena might also permit broader sample while 397 simultaneously documenting variations in landing style on natural roosts. 398 In addition to increased sampling, future work would also benefit from efforts to measure 399 and reconstruct the evolutionary history of morphological traits related to landing maneuvers, 400 such as wing mass distribution, which determines the inertial torques bats use to execute landings 401 (Bergou et al., 2015), and other skeletal features relating to limb stresses and landing impact forces. 402 Taken together, these efforts would determine whether there are clade-based links among roosting 403 habits, landing style, wing morphology, and diversification rates. If shifts in roosting ecology were 404 associated with speciation in certain lineages (e.g., stenodermatines), and if roosting ecology, 405 landing style, and wing morphology were linked, then one should detect significant shifts in 406 diversification rates for clades that arise following coordinated shifts in roosting habits and landing 407 style.

Materials and Methods 437
Focal taxa, field sites, and animal capture 438 We recorded 665 landings from 96 bats, representing 35 species, and 9 families (table1). 439 We collected all measurements from wild-caught bats except for Rousettus aegyptiacus and taxa from 440 We used a published time-calibrated molecular phylogeny (Shi and Rabosky, 2015), pruned 499 to our focal taxa, for all phylogenetic analyses (excluding A. watsoni, which was not included in the 500 Shi & Rabosky tree), using the Phytools R-package (Revell, 2018(Revell, , 2011. We then assigned one 501 landing style as a discrete character to each taxon according to its most-often observed landing style 502 (table 1). 503 We conducted an ancestral state reconstruction using stochastic character mapping 504 (Huelsenbeck et al., 2003), as implemented in the make.simmap function of the R package 505 phytools, to reconstruct the evolutionary history of landing styles among sampled taxa. We used 506 the fitDiscrete function in the R package Geiger (Harmon et al., 2007) to compare the fit of four 507 different models for the transition matrix of the stochastic character mapping procedure: equal 508 rates, symmetric, all rates different, and meristic. The equal rates model yielded the lowest AICc 509 score, thus we selected this model, which gave all state changes equal probability, and computed 510 the posterior probability for each landing style at internal nodes from 1000 simulated stochastic 511

maps. 512
Next, we used phylogenetic generalized least squares regression (PGLS), implemented in 513 the R function pgls from the package Caper (Orme, 2018) to explore the extent to which landing 514 impact force is predicted by points of contact. Here, we estimated phylogenetic signal using the 515 maximum likelihood value of Pagel's lambda and treated points of contact and peak 3D landing 516 impact force as continuous variables. We log-transformed impact forces to ensure normality. We 517 then computed a phylogenetic ANOVA (10000 iterations) with post-hoc tests using the 518 phylANOVA function in the R package phytools to test for pairwise differences in log-peak 519 landing impact forces among landing styles. Peak impact force was the response variable and 520 landing style was the factor. We omitted two species from these analyses due to an inability to 521 unambiguously designate them as a two-, three-, or four-point landing: M. schreibersii due to its high 522 degree of behavioral variability and T. tricolor because it performs a specialized landing maneuver to 523 alight on a vertical substrate (Boerma et al., 2019), rather than beneath a horizonal roost as in the 524 landing experiments for our other sampled taxa. 525 We used phylogenetic logistic regression with Firth's correction (Ives and Garland, 2009), 526 as implemented in the R package, phylolm (Ho and Ané, 2014), to test the hypothesis that landing 527 styles are associated with the physical properties of roosts. We applied 2000 bootstrap replicates to 528 generate confidence intervals for and test the significance of the model coefficients, ß, which relate 529 to the probability of observing a particular landing style (categorical response variable) given a 530 particular roosting ecology (categorical predictor variables). Positive coefficients indicate a positive 531 association between predictor and response variables, whereas negative coefficients denote a 532 negative relationship. We excluded T. tricolor from this analysis because it is the only sampled 533 taxon to perform its landing maneuver and to roost in tubular furled leaves. We compared two 534 models of roosting habits, the latter of which aggregated multiple roost types according to their 535 physical properties, thereby testing our hypothesis that diverse roost types that share physical 536 properties are correlate with landing style. The models were as follows: 1) a null model in which we 537 assigned each taxon's roosting ecology according to its most-commonly cited roost type (table 1), 538 and 2) an alternate model in which we aggregated roosting ecologies that include stiff, primarily 539 horizontal surfaces (CST, EST, and R/C) into a single category, spatially-constrained, horizontal bodyweight, for each recorded landing. Impact forces were recorded at 1000 Hz, 567 normalized to the individual's body mass, and smoothed using a zero-phase 2 nd order low-568 pass Butterworth filter with a cutoff frequency of 100 Hz, parameters which are identical to 569 those of previous bat landing studies. The total (resultant) force into the ceiling was 570 calculated and the peak extracted for each landing. 571 • Figure  for all phylogenetic analyses in this study. The tree was trimmed to include only our focal 576 species prior to any analyses. 577 • Figure 2-Source Data 2: This .csv file contains the source data for the stochastic character 578 mapping and ancestra state reconstruction for landing style. 579 • Table 2-Source Data 1: This .csv file contains the source data for the phylogenetic logistic 580 regression results summarized in Table 2. 581 • PGLS phylANOVA-Source Data 1: This .csv file contains the source data (Taxon, Landing 582 Style, and mean peak impact force) for the phylogenetic generalized least squares regression 583 and phylogenetic ANOVA with pairwise comparisons. This data file omits Thyroptera 584 tricolor from both analyses because it is the only species to perform its specialized four-point 585 landing and to roost in furled leaf tubes. It also omits Miniopterus schreibersii because its 586 landing style was equivocal. See Supplemental File 3 for a version that includes Thyroptera 587 tricolor. 588 • Supplemental File 1: Posterior probabilities for landing style at each node in the phylogeny 589 shown in Figure 2. Posterior probabilities were estimated using an equal rates model for 590 1000 simulated stochastic character maps. 591 • Supplemental File 2: Node legend for posterior probabilities in the phylogeny shown in 592 Figure 2. 593 • Supplemental File 3: This .csv file summarized the landing style and mean peak impact 594 force for each species for which impacted forces were measured, except Miniopterus 595 schreibersii, for which landing style was equivocal. See Table 1, Figure 1, and Source Data 596 files for Figure 1 for impact forces and landing style in M. schreibersii. 597 • Source Code File 1: This .R file contains the code for plotting impact forces shown in 598 Figure 1. 599 • Source Code File 2: This .R file contains the code for computing the ancestral state 600 reconstruction for landing style using stochastic character mapping. 601 • Source Code File 3: This .R file contains the code for computing the phylogenetic 602 generalized least squares and ANOVA that test for associations between landing style and 603 landing impact force, and the phylogenetic logistic regressions that test for associations 604 between landing style and roosting ecology.