# Is Beak Morphology in Darwin's Finches Tuned to Loading Demands?

Joris Soons, Annelies Genbrugge, Jeffrey Podos, Dominique Adriaens, Peter Aerts, Joris Dirckx, Anthony Herrel

► To cite this version:

Joris Soons, Annelies Genbrugge, Jeffrey Podos, Dominique Adriaens, Peter Aerts, et al., Is Beak Morphology in Darwin's Finches Tuned to Loading Demands?. PLoS ONE, 2015, 10 (6), pp.e0129479. 10.1371/journal.pone.0129479. hal-02145307

HAL Id: hal-02145307

https://hal.science/hal-02145307v1

Submitted on 6 Jan 2025

HAL is a multi-disciplinary open access
archive for the deposit and dissemination of sci-
entific research documents, whether they are pub-
lished or not. The documents may come from
teaching and research institutions in France or
abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est
destinée au dépôt et à la diffusion de documents
scientifiques de niveau recherche, publiés ou non,
émanant des établissements d'enseignement et de
recherche français ou étrangers, des laboratoires
publics ou privés.

<figure><img src="image_1.png" /></figure>
<figure><img src="image_2.png" /></figure>

RESEARCH ARTICLE

# Is Beak Morphology in Darwin's Finches Tuned to Loading Demands?

Joris Soons¹*, Annelies Genbrugge¹,², Jeffrey Podos³, Dominique Adriaens², Peter Aerts⁴,⁵, Joris Dirckx¹, Anthony Herrel⁴,⁶

¹ Laboratory of Biomedical Physics, University of Antwerp, Groenenborgerlaan 171, 2020, Antwerpen, Belgium, ² Evolutionary Morphology of Vertebrates, Ghent University, K.L. Ledeganckstraat 35, 9000, Gent, Belgium, ³ Department of Biology, 221 Morrill Science Center, University of Massachusetts, Amherst, Massachusetts, 01003, United States of America, ⁴ Department of Biology, University of Antwerp, Universiteitsplein 1, 2610, Antwerpen, Belgium, ⁵ Department of movement and sports sciences, University of Ghent, Watersportlaan 2, 9000, Gent, Belgium, ⁶ Département d'Ecologie et de Gestion de la Biodiversité, Museum National d'Histoire Naturelle, 57 rue Cuvier, Case postale 55, 75231, Paris, Cedex 5, France

* joris.soons@ua.ac.be

## Abstract

One of nature's premier illustrations of adaptive evolution concerns the tight correspondence in birds between beak morphology and feeding behavior. In seed-crushing birds, beaks have been suggested to evolve at least in part to avoid fracture. Yet, we know little about mechanical relationships between beak shape, stress dissipation, and fracture avoidance. This study tests these relationships for Darwin's finches, a clade of birds renowned for their diversity in beak form and function. We obtained anatomical data from micro-CT scans and dissections, which in turn informed the construction of finite element models of the bony beak and rhamphotheca. Our models offer two new insights. First, engineering safety factors are found to range between 1 and 2.5 under natural loading conditions, with the lowest safety factors being observed in species with the highest bite forces. Second, size-scaled finite element (FE) models reveal a correspondence between inferred beak loading profiles and observed feeding strategies (e.g. edge-crushing versus tip-biting), with safety factors decreasing for base-crushers biting at the beak tip. Additionally, we identify significant correlations between safety factors, keratin thickness at bite locations, and beak aspect ratio (depth versus length). These lines of evidence together suggest that beak shape indeed evolves to resist feeding forces.

## Introduction

The often tight correspondence between bird beaks and plant morphology well-illustrates the power and precision of natural selection [1]. Within bird populations, subtle variations in beak morphology can affect foraging efficiency [2], and in some cases tip the balance between survival and starvation [3]. For seed-crushing birds, however, analyses of feeding capacity require a broader view than just beak morphology *per se*. This is because a bird's ability to crush seeds is determined mainly by bite force capacity, which in turn depends primarily on the orientation

<figure><img src="image_3.png" /></figure>

### OPEN ACCESS

**Citation:** Soons J, Genbrugge A, Podos J, Adriaens D, Aerts P, Dirckx J, et al. (2015) Is Beak Morphology in Darwin's Finches Tuned to Loading Demands? PLoS ONE 10(6): e0129479. doi:10.1371/journal.pone.0129479

**Academic Editor:** Kornelius Kupczik, Max Planck Institute for Evolutionary Anthropology, GERMANY

**Received:** November 29, 2013

**Accepted:** May 10, 2015

**Published:** June 12, 2015

**Copyright:** © 2015 Soons et al. This is an open access article distributed under the terms of the Creative Commons Attribution License, which permits unrestricted use, distribution, and reproduction in any medium, provided the original author and source are credited.

**Funding:** This work was performed in part at the Center for Nanoscale Systems (CNS), a member of the National Nanotechnology Infrastructure Network (NNIN), which is supported by the National Science Foundation under NSF award no. ECS-0335765. CNS is part of the Faculty of Arts and Sciences at Harvard University. Financial support to this project was provided by the Research Foundation - Flanders (FWO) to JS, an interdisciplinary research grant of the special research fund of the University of Antwerp to PA, JD, AH, a PHC Tournesol project to AH, AG and DA, and JP was supported by NSF grant IBN-0347291. The funders had no role in study design,
data collection and analysis, decision to publish, or
preparation of the manuscript.

**Competing Interests:** The authors have declared
that no competing interests exist

and cross sectional area of the jaw closer muscles [4–6]. These muscles, situated at the back of
the head, generate crushing forces that are transferred to food by means of the upper and lower
beak through a complex kinetic beak apparatus [6–11]. Beak morphology, by contrast, likely
evolves to facilitate successful food manipulation and song production. However, beaks also
need to avoid structural failure during loading [4]. The loading regime itself presumably varies
with food manipulation strategies (e.g. tip biting versus base crushing) and may thus drive the
evolution of beak shape [12].

We here test, in Darwin's finches of the Galápagos Islands, the hypothesis that beak morphology provides structural integrity, i.e., that it evolves to resist feeding forces and avoid fracture [13]. Beaks in Darwin's finches are known to vary broadly across and within species, and to evolve via natural selection as a response to variation in food type, food availability, and interspecific competition [3,10,14–18]. Our main approach is to develop and apply finite element (FE) models of the upper beak. FE models enable exploration of the effect of complex shape variation on stress magnitude and distribution [19–21], and allows us to draw inferences about patterns of loading during biting. We include here 13 species and incorporate data into our models for the bones, the jaw closer muscles, and the keratinous rhamphotheca that encapsulates the beak. Including the rhamphotheca in the model allows us to provide a more realistic estimate of stresses incurred during biting [22,23] that we achieved in previous studies [24,25].

The goal of this paper is twofold. First we examine how natural loadings are reflected in
safety factors for a subset of species for which physiological FE models could be created (i.e. for
which data on muscles could be obtained). Second, we created models of all 13 species scaled to
a common surface area to muscle ratio [26]. This allows us to test whether variation in beak
shape matches differences in food handling behavior. We predict that species will have beak
shapes that allow them to minimize peak stresses during those feeding behaviors most com-
monly employed. In particular, we predict that base-crushing species should show the lowest
stresses when biting at the base of the beak, that tip-biting species should show the lowest
stresses during tip-biting, and that probing beaks should show higher stresses during both tip-
and base crushing compared to those species specialized for these behaviors.

**Materials and Methods**

Specimen collection and CT scanning

No animals were killed for the purpose of this study. Road-killed specimens were collected dur-
ing February-March of 2005 and 2006 on the main road connecting the airport on Baltra with
Puerto Ayora, Santa Cruz Island. Specimens were collected under a salvage permit provided by
the Galápagos National Park Service. Intact specimens were collected and preserved in a 10%
aqueous formaldehyde solution for 24 h, rinsed and transferred to a 70% aqueous ethanol solu-
tion. One individual of *Geospizia fortis*, *Geospizia fuliginosa*, *Geospizia scandens*, *Platyspiza
crassirostris*, *Certhidea olivacea*, and *Camarhynchus parvulus* were scanned at the UGCT scan-
ning facility (www.ugct.ugent.be), using a micro-focus directional type X-ray tube, set at a volt-
age of 80 kV and a spot size of 10 mm. Specimens were mounted on a controllable rotating
table (MICOS, UPR160F-AIR). For each specimen, a series of 1000 projections of 940 by 748
pixels were recorded covering 360 degrees.

Museum specimens of *Geospizia magnirostris*, *Geospizia difficilis*, *Pinaroloxias inornata*,
*Cactospiza pallida*, *Camarhynchus psittacula*, *Geospiza conirostris*, and *Camarhynchus pauper*
from the collections of the California Academy of Sciences maintained in a 70% aqueous etha-
nol solution were scanned at the Harvard CNS facility using an X-Tek XRA-002 micro-CT im-
aging system set at 75 kV. Specimens were mounted on a rotating table and a series of 3142
projections of 2000 by 2000 pixels covering 360 degrees was recorded. The voxel size of the
scans was dependent on the specimen size and ranged from 2.96 µm for the smallest species to 47.00 µm for the larger species. As such, the bone shell of the bony core was represented by multiple voxels, and details within the bony beak were visible.

## Muscle data

Muscle data were available for five of the species collected in the field as road-killed specimens: *G. fortis*, *G. fuliginosa*, *G. scandens*, *C. olivacea* and *C. parvulus*; data for three additional species were derived through dissection of museum specimens from the California Academy of Sciences (*G. magnirostris*, *G. conirostris*, and *G. difficilis*). These species were dissected and all muscle bundles of the jaw removed individually [27]. Muscles were blotted dry and weighed on a Mettler microbalance (±0.01mg). Next, muscles were transferred individually to Petri dishes and submerged in a 30% aqueous nitric acid solution for 18 h to dissolve all connective tissue [28]. After removal of nitric acid, muscles were transferred to a 50% aqueous glycerol solution and fibers were teased apart using blunt-tipped glass needles. Thirty fibers were selected from each muscle bundle (Table 1) and drawn using a binocular scope with attached *camera lucida*. A background grid was also drawn in each image to provide an object for scaling. Drawings were scanned and fiber lengths determined using imageJ (freely available at http://rsb.info.nih.gov/ij/).

Based on muscle mass and fiber length, the physiological cross-sectional area of each muscle bundle was estimated assuming a muscle density of 1036kg/m³ [29]. Since pennate muscles were separated into their individual bundles, no additional correction for pennation angle was included. Corresponding force-generation capacities for each muscle were calculated assuming a muscle stress of 30N/cm² [30]. As the external adductor and pseudotemporalis muscles act only indirectly on the upper mandible [6,9,10,31], the component of the muscle force transferred to the upper mandible was calculated taking into account the position of the muscles and their angles relative to the jugal bone (which transfers the forces from the lower jaw and the quadrate to the upper beak by pulling the beak downward when quadrate is rotated backward during jaw closing). The pterygoid muscle bundles act directly on the upper mandible [6,9,10,31], and muscle forces are directly transmitted through the pterygoid/palatine complex. The combined muscle forces applied to the model are shown in Table 2.

## Segmentation and FE modelling

The segmentation of the bony core and the keratin layer (Amira 4.1 64-bit version, TGS systems) was done for all species. The tetrahedral grid construction [32], the assigned material properties (homogeneous and isotropic linear elastic with $E_{keratin} = 1.7$ GPa and $E_{bone} = 7.3$ GPa) [22] and the boundary conditions are similar to those reported previously [22–24, 33]. Convergence testing on the results, namely the von Mises (vM) stresses and the bite forces, was conducted so they change less than 5% when doubling mesh size. The final FE grids of approximately 2–3 million 4-noded linear tetrahedral elements were obtained. Consequently, multiple layers of elements were available to model the thin bony walls.

For every specimen, a FE model for base and tip loading was established. In this study we were only interested in the stresses on the upper beak, and as such elements were fully constrained (in all 6 degrees of freedom) at the back of the beak at the neurocranium. Hinging of the upper beak is possible due to bending of the thin bone (BA in Fig 1). In order to avoid singularities around point constraints, we constrained elements (in x, y and z direction, no rotation) at the ventral surface of the rhamphotheca simulating the bite constraint (T and Bs in Fig 1); the calculated reaction force of these elements is the bite force which is of equal magnitude yet of opposite sign as the seed reaction force. Bite positions were estimated using in vivo
<table><caption>Table 1. Table summarizing muscle mass and fiber length data used to calculate the external forces acting on the upper beak of the different species of Darwin's finch.</caption><thead><tr><th>Species</th><th>MDM</th><th>MAMER</th><th>MAMEV</th><th>MAMEP</th><th>MAMOQ</th><th>MPsTSI</th><th>MPsTSm</th><th>MPsTP</th><th>MPtVI</th><th>MPtVm</th><th>MPtDI</th><th>MPtDm</th><th>MRPal</th><th>MPPIQ</th></tr></thead><tbody><tr><td>G. magnirostris</td><td>mass</td><td>0.0522</td><td>0.1813</td><td>0.0451</td><td>0.0595</td><td>0.0084</td><td>0.0137</td><td>0.0624</td><td>0.0447</td><td>0.0752</td><td>0.0758</td><td>0.0568</td><td>0.022</td><td>0.0104</td><td>0.0085</td></tr><tr><td></td><td>fl</td><td>3.99</td><td>2.31</td><td>2.01</td><td>2.02</td><td>2.58</td><td>1.36</td><td>1.39</td><td>4.18</td><td>2.58</td><td>2.56</td><td>2.56</td><td>2.63</td><td>2.61</td><td>2.65</td></tr><tr><td>G. cornirostris</td><td>mass</td><td>0.0455</td><td>0.0389</td><td>0.0178</td><td>0.0262</td><td>0.0035</td><td>0.0033</td><td>0.0157</td><td>0.0269</td><td>0.0079</td><td>0.0253</td><td>0.0240</td><td>0.0093</td><td>0.0044</td><td>0.0086</td></tr><tr><td></td><td>fl</td><td>4.12</td><td>2.55</td><td>2.93</td><td>2.37</td><td>2.21</td><td>1.81</td><td>2.62</td><td>4.52</td><td>3.10</td><td>3.20</td><td>3.42</td><td>2.65</td><td>2.73</td><td>2.63</td></tr><tr><td>G. difficilis</td><td>mass</td><td>0.005</td><td>0.0102</td><td>0.0025</td><td>0.0044</td><td>0.0002</td><td>0.0006</td><td>0.0012</td><td>0.0064</td><td>0.0020</td><td>0.0055</td><td>0.0071</td><td>0.003</td><td>0.0014</td><td>0.0004</td></tr><tr><td></td><td>fl</td><td>2.91</td><td>1.67</td><td>2.01</td><td>1.72</td><td>1.64</td><td>1.53</td><td>2.15</td><td>2.76</td><td>2.36</td><td>4.32</td><td>2.54</td><td>2.21</td><td>1.98</td><td>1.62</td></tr><tr><td>G. fuliginosa</td><td>mass</td><td>0.0070</td><td>0.0163</td><td>0.0035</td><td>0.0069</td><td>0.0010</td><td>0.0020</td><td></td><td>0.0020</td><td>0.0027</td><td>0.0036</td><td>0.0052</td><td>0.0032</td><td>0.0014</td><td>0.0020</td></tr><tr><td></td><td>fl</td><td>2.85</td><td>2.13</td><td>1.41</td><td>1.46</td><td>2.02</td><td>1.12</td><td></td><td>1.46</td><td>2.09</td><td>1.99</td><td>1.69</td><td>2.10</td><td>2.60</td><td>2.22</td></tr><tr><td>C. parvulus</td><td>mass</td><td>0.0056</td><td>0.0130</td><td>0.0029</td><td>0.0065</td><td>0.0005</td><td>0.0013</td><td></td><td>0.0022</td><td>0.0044</td><td>0.0016</td><td>0.0027</td><td>0.0029</td><td>0.0012</td><td>0.0040</td></tr><tr><td></td><td>fl</td><td>3.45</td><td>3.07</td><td>2.11</td><td>1.92</td><td>1.44</td><td>1.67</td><td></td><td>1.58</td><td>2.81</td><td>2.81</td><td>3.19</td><td>2.75</td><td>3.18</td><td>2.61</td></tr><tr><td>C. olivacea</td><td>mass</td><td>0.0033</td><td>0.0041</td><td>0.0011</td><td>0.0022</td><td></td><td>0.0007</td><td></td><td>0.0005</td><td>0.0021</td><td>0.0014</td><td>0.0008</td><td>0.0016</td><td>0.0004</td><td>0.0019</td></tr><tr><td></td><td>fl</td><td>2.05</td><td>1.74</td><td>1.14</td><td>1.33</td><td></td><td>1.03</td><td></td><td>1.30</td><td>2.66</td><td>2.27</td><td>2.47</td><td>1.86</td><td>1.80</td><td>1.86</td></tr><tr><td>G. scandens</td><td>mass</td><td>0.015</td><td>0.0329</td><td>0.009</td><td>0.0126</td><td>0.0016</td><td>0.0071</td><td>0.0052</td><td>0.0168</td><td>0.0072</td><td>0.0062</td><td>0.0247</td><td>0.0108</td><td>0.0057</td><td>0.0042</td></tr><tr><td></td><td>fl</td><td>3.89</td><td>2.52</td><td>1.89</td><td>1.98</td><td>2.64</td><td>1.22</td><td>1.22</td><td>1.54</td><td>2.99</td><td>2.96</td><td>2.70</td><td>2.75</td><td>2.60</td><td>2.30</td></tr><tr><td>G. fortis</td><td>mass</td><td>0.02</td><td>0.0651</td><td>0.0159</td><td>0.0217</td><td>0.0026</td><td>0.0094</td><td>0.0124</td><td>0.0283</td><td>0.0199</td><td>0.0119</td><td>0.0284</td><td>0.0333</td><td>0.0041</td><td>0.0059</td></tr><tr><td></td><td>fl</td><td>3.91</td><td>1.87</td><td>1.85</td><td>1.70</td><td>2.02</td><td>1.31</td><td>2.06</td><td>3.58</td><td>2.38</td><td>2.52</td><td>2.07</td><td>1.73</td><td>2.88</td><td>3.40</td></tr></tbody></table>

Mass in gram, fiber lengths in mm. N = 1 for all species except G. fortis where N = 4. For the smaller species C. olivacea, C. parvulus and G. fuliginosa the distinction between the medial and lateral superficial m. pseudotemporalis was extremely difficult to make. As such both layers were grouped. Similarly, for the smallest species, C. olivacea, we were unable to extract the m. adductor mandibulae ossi quadrati (MAMOQ), MDM, m. depressor mandibulae; MAMEP, m. adductor mandibulae externis rostralis; MAMEV, m. adductor mandibulae externus ventralis; MAMEP, m. adductor mandibulae externus profundus; MAMOQ, m. adductor mandibulae ossi quadrati; MPSTSI, m. pseudotemporalis superficialis pars lateralis; MPsTSm, m. pseudotemporalis superficialis pars medialis; MPsTP, m. pterygoideus ventralis pars lateralis; MPtVm, m. pterygoideus ventralis pars medialis; MPtDI, m. pterygoideus dorsalis pars lateralis; MPtDm, m. pterygoideus dorsalis pars medialis; MRPal, m. retractor palatine; MPPIQ, m. protractor pterygoidei et quadrati.

doi:10.1371/journal.pone.0129479.t001
<table><caption>Table 2. Von Mises stress for different loading conditions (LC) for 13 Darwin's finches (PB: physiological base biting; PT: physiological tip biting; FB: base biting scaled to *fortis*; FT: tip biting scaled to *fortis*).</caption><thead><tr><th>Species</th><th>LC</th><th>Fj</th><th>Fp</th><th>vM1</th><th>vM2</th><th>vM3</th><th>GV1</th><th>GV2</th><th>GV3</th><th>F</th><th>V<sub>bone</sub></th><th>V<sub>ker</sub></th></tr></thead><tbody><tr><td>G. fortis</td><td>PB</td><td>11.9</td><td>10.9</td><td>21</td><td>19</td><td>29</td><td>0.83</td><td>0.98</td><td>0.81</td><td>30.4</td><td>104</td><td>128</td></tr><tr><td>G. fortis</td><td>PT</td><td>11.9</td><td>10.9</td><td>47</td><td>36</td><td>41</td><td>0.78</td><td>0.95</td><td>0.82</td><td>21.5</td><td>104</td><td>128</td></tr><tr><td>G. magnirostris</td><td>PB</td><td>23.0</td><td>21.9</td><td>35</td><td>24</td><td>40</td><td>0.86</td><td>0.97</td><td>0.77</td><td>57.9</td><td>232</td><td>195</td></tr><tr><td>G. magnirostris</td><td>PT</td><td>23.0</td><td>21.9</td><td>66</td><td>35</td><td>46</td><td>0.62</td><td>0.98</td><td>0.81</td><td>41.6</td><td>232</td><td>195</td></tr><tr><td>G. magnirostris</td><td>FB</td><td>11.9</td><td>10.9</td><td>22</td><td>16</td><td>26</td><td>0.86</td><td>0.97</td><td>0.77</td><td>29.3</td><td>157</td><td>131</td></tr><tr><td>G. magnirostris</td><td>FT</td><td>11.9</td><td>10.9</td><td>43</td><td>23</td><td>30</td><td>0.62</td><td>0.98</td><td>0.81</td><td>20.9</td><td>157</td><td>131</td></tr><tr><td>G. fuliginosa</td><td>PB</td><td>2.2</td><td>3.2</td><td>18</td><td>15</td><td>18</td><td>0.74</td><td>0.93</td><td>0.85</td><td>7.1</td><td>24</td><td>35</td></tr><tr><td>G. fuliginosa</td><td>PT</td><td>2.2</td><td>3.2</td><td>24</td><td>28</td><td>27</td><td>0.60</td><td>0.88</td><td>0.84</td><td>4.6</td><td>24</td><td>35</td></tr><tr><td>G. fuliginosa</td><td>FB</td><td>11.9</td><td>10.9</td><td>29</td><td>23</td><td>28</td><td>0.74</td><td>0.93</td><td>0.85</td><td>29.2</td><td>115</td><td>164</td></tr><tr><td>G. fuliginosa</td><td>FT</td><td>11.9</td><td>10.9</td><td>35</td><td>41</td><td>42</td><td>0.60</td><td>0.88</td><td>0.84</td><td>19.0</td><td>115</td><td>164</td></tr><tr><td>G. conirostris</td><td>PB</td><td>5.4</td><td>5.3</td><td>16</td><td>25</td><td>42</td><td>0.86</td><td>0.85</td><td>0.76</td><td>15.3</td><td>76</td><td>82</td></tr><tr><td>G. conirostris</td><td>PT</td><td>5.4</td><td>5.3</td><td>25</td><td>40</td><td>47</td><td>0.98</td><td>0.83</td><td>0.70</td><td>10.4</td><td>76</td><td>82</td></tr><tr><td>G. conirostris</td><td>FB</td><td>11.9</td><td>10.9</td><td>25</td><td>33</td><td>67</td><td>0.86</td><td>0.85</td><td>0.76</td><td>32.3</td><td>112</td><td>121</td></tr><tr><td>G. conirostris</td><td>FT</td><td>11.9</td><td>10.9</td><td>45</td><td>63</td><td>67</td><td>0.98</td><td>0.83</td><td>0.70</td><td>18.8</td><td>112</td><td>121</td></tr><tr><td>G. scandens</td><td>PB</td><td>3.2</td><td>4.5</td><td>20</td><td>15</td><td>23</td><td>0.67</td><td>0.92</td><td>0.74</td><td>8.9</td><td>59</td><td>56</td></tr><tr><td>G. scandens</td><td>PT</td><td>3.2</td><td>4.5</td><td>17</td><td>33</td><td>35</td><td>0.88</td><td>0.95</td><td>0.82</td><td>6.1</td><td>59</td><td>56</td></tr><tr><td>G. scandens</td><td>FB</td><td>11.9</td><td>10.9</td><td>40</td><td>28</td><td>43</td><td>0.67</td><td>0.92</td><td>0.74</td><td>26.2</td><td>126</td><td>119</td></tr><tr><td>G. scandens</td><td>FT</td><td>11.9</td><td>10.9</td><td>30</td><td>60</td><td>63</td><td>0.88</td><td>0.95</td><td>0.82</td><td>17.9</td><td>126</td><td>119</td></tr><tr><td>G. difficilis</td><td>PB</td><td>1.5</td><td>1.6</td><td>14</td><td>12</td><td>21</td><td>0.87</td><td>1.00</td><td>0.97</td><td>3.6</td><td>22</td><td>30</td></tr><tr><td>G. difficilis</td><td>PT</td><td>1.5</td><td>1.6</td><td>18</td><td>19</td><td>18</td><td>0.73</td><td>0.98</td><td>0.98</td><td>2.4</td><td>22</td><td>30</td></tr><tr><td>G. difficilis</td><td>FB</td><td>11.9</td><td>10.9</td><td>38</td><td>24</td><td>51</td><td>0.87</td><td>1.00</td><td>0.97</td><td>26.1</td><td>86</td><td>116</td></tr><tr><td>G. difficilis</td><td>FT</td><td>11.9</td><td>10.9</td><td>50</td><td>56</td><td>53</td><td>0.73</td><td>0.98</td><td>0.98</td><td>16.7</td><td>86</td><td>116</td></tr><tr><td>P. inornata</td><td>FB</td><td>11.9</td><td>10.9</td><td>55</td><td>32</td><td>59</td><td>0.66</td><td>0.95</td><td>0.85</td><td>26.2</td><td>109</td><td>96</td></tr><tr><td>P. inornata</td><td>FT</td><td>11.9</td><td>10.9</td><td>76</td><td>67</td><td>61</td><td>0.71</td><td>1.00</td><td>0.83</td><td>13.8</td><td>109</td><td>96</td></tr><tr><td>C. olivacea</td><td>PB</td><td>0.9</td><td>1.0</td><td>23</td><td>31</td><td>45</td><td>0.50</td><td>0.84</td><td>0.81</td><td>2.0</td><td>9</td><td>9</td></tr><tr><td>C. olivacea</td><td>PT</td><td>0.9</td><td>1.0</td><td>23</td><td>32</td><td>22</td><td>1.00</td><td>1.00</td><td>0.94</td><td>1.0</td><td>9</td><td>9</td></tr><tr><td>C. olivacea</td><td>FB</td><td>11.9</td><td>10.9</td><td>58</td><td>66</td><td>100</td><td>0.50</td><td>0.84</td><td>0.81</td><td>24.4</td><td>99</td><td>100</td></tr><tr><td>C. olivacea</td><td>FT</td><td>11.9</td><td>10.9</td><td>58</td><td>78</td><td>59</td><td>1.00</td><td>1.00</td><td>0.94</td><td>12.3</td><td>99</td><td>100</td></tr><tr><td>C. pallida</td><td>FB</td><td>11.9</td><td>10.9</td><td>59</td><td>33</td><td>80</td><td>0.77</td><td>0.71</td><td>0.88</td><td>25.4</td><td>68</td><td>101</td></tr><tr><td>C. pallida</td><td>FT</td><td>11.9</td><td>10.9</td><td>38</td><td>46</td><td>64</td><td>0.67</td><td>0.87</td><td>0.89</td><td>15.1</td><td>68</td><td>101</td></tr><tr><td>C. parvulus</td><td>PB</td><td>1.3</td><td>1.3</td><td>22</td><td>12</td><td>20</td><td>0.57</td><td>0.94</td><td>0.61</td><td>3.3</td><td>11</td><td>18</td></tr><tr><td>C. parvulus</td><td>PT</td><td>1.3</td><td>1.3</td><td>13</td><td>19</td><td>25</td><td>0.66</td><td>0.82</td><td>0.58</td><td>2.2</td><td>11</td><td>18</td></tr><tr><td>C. parvulus</td><td>FB</td><td>11.9</td><td>10.9</td><td>47</td><td>26</td><td>43</td><td>0.57</td><td>0.94</td><td>0.61</td><td>29.6</td><td>88</td><td>143</td></tr><tr><td>C. parvulus</td><td>FT</td><td>11.9</td><td>10.9</td><td>29</td><td>41</td><td>55</td><td>0.66</td><td>0.82</td><td>0.58</td><td>19.5</td><td>88</td><td>143</td></tr><tr><td>C. pauper</td><td>FB</td><td>11.9</td><td>10.9</td><td>23</td><td>17</td><td>37</td><td>0.60</td><td>0.94</td><td>0.52</td><td>31.2</td><td>118</td><td>163</td></tr><tr><td>C. pauper</td><td>FT</td><td>11.9</td><td>10.9</td><td>19</td><td>25</td><td>35</td><td>0.82</td><td>0.99</td><td>0.75</td><td>19.4</td><td>118</td><td>163</td></tr><tr><td>C. psittacula</td><td>FB</td><td>11.9</td><td>10.9</td><td>39</td><td>24</td><td>33</td><td>0.46</td><td>0.88</td><td>0.90</td><td>30.5</td><td>93</td><td>127</td></tr><tr><td>C. psittacula</td><td>FT</td><td>11.9</td><td>10.9</td><td>25</td><td>28</td><td>35</td><td>0.71</td><td>0.63</td><td>0.90</td><td>21.9</td><td>93</td><td>127</td></tr><tr><td>P. crassirostris</td><td>FB</td><td>11.9</td><td>10.9</td><td>21</td><td>21</td><td>31</td><td>0.84</td><td>0.93</td><td>0.81</td><td>32.1</td><td>138</td><td>152</td></tr><tr><td>P. crassirostris</td><td>FT</td><td>11.9</td><td>10.9</td><td>25</td><td>32</td><td>43</td><td>0.82</td><td>0.97</td><td>0.90</td><td>21.7</td><td>138</td><td>152</td></tr></tbody></table>

Model input forces (= muscle forces) for jugal (Fj) and palatine (Fp) are given in N. von Mises stresses are given in MPa for three positions indicated in Fig 1 (vM1: on top of bone, near bite position; vM2: on top of the nasal hinge; vM3: nasal bone ipsilateral side for base biting, both sides for tip biting). Gray values in CT-stack for the same positions are given (0 black, 1 white). The resulting (model) biting force (F) is given in N; the volume of keratin and bone ($V_{ker}$ and $V_{bone}$) are given in mm³.

doi:10.1371/journal.pone.0129479.t002
<figure><img src="image_4.png" /><figcaption>**Fig 1.** Schematic representation of our multi-layered (bone: B, keratin: K) finite element modeling approach, for the medium ground finch *Geospiza fortis*. Bending area (Ba) and bite position (base: Bs, or tip: T) were constrained in our models for translation and rotation, and muscle forces were applied in our models via the jugal (J) and palatine (P) jaw bones (black elements are constrained). Locations of vM Stress recordings are indicated with transparant ellipses (1b: on top of bone, near base bite position; 1t: on top of bone near tip bite position; 2: on top of the beak near the nasal hinge; 3: nasal bone).</figcaption></figure>

doi:10.1371/journal.pone.0129479.g001

observations of animals cracking or manipulating seeds in the field [13]. Elements were con-
strained for a central tip bite close to the tip of the bony core, as birds are never seen to bite at
the very tip in the wild (T in Fig 1). For an unilateral base bite we used a point at one fourth of
the total length of the rhamphotheca, close to the edge of the bony core (Bs in Fig 1), hence cre-
ating homologous bite points across species. The available muscle forces (Table 2) estimated
through dissection were applied to rigid elements at the end of the jugal and pterygoid bones in
the direction of these bones (J and P in Fig 1). This resulted in physiological FE models, for
both base and tip loading, for eight species (16 models). In addition we developed scaled FE
models for all 13 species. We therefore measured the beak surface of all beaks and scaled them
in such a way that they have the same surface area as *G. fortis*. Next we applied the muscle
force of *G. fortis* to these size-scaled beaks. *G. fortis* was chosen because we had muscle data
and scans available for multiple individuals of this species [27], but we could have picked one
of the other 7 finches. As such, the effects of size are eliminated and the effect of beak shape
variation on the mechanical behavior of the beak can be isolated for the 13 species of Darwin's
finches included here [26].

Bite force and maximum vM stress were obtained for all 40 FE models (Table 2): 16 physio-
logical models (tip and base loading for 8 beaks with available muscle data) and 24 models (tip
and base loading for all 13 beaks size-scaled to G. fortis, minus G. fortis, see Table 2). Maximum
vM stress was taken according to Saint-Venant's Principle over three different locations dis-
playing the highest stresses on the bony core, away from areas influenced by model constraints
(Figs 1 and 2, Table 2). We used 98th percentile to select maximum vM stresses presented in
Fig 3 and Table 2. In these locations, the gray values of the corresponding CT images were also
recorded. It has been shown [34–36] that the strength of bone increases with density. Scan pa-
rameters were optimized for each CT scan and no calibration was performed, making it impos-
sible to obtain the density directly from these gray values. Consequently, gray values were
transferred to a qualitative linear scale, where 0 represents the least dense bone in the sample
and 1 the densest. vM stresses on the keratin layer were not taken into account and are not dis-
played in Figs 2 and 4. Nevertheless, a correct modeling of the rhamphotheca is important
since it has a significant influence on the stress distribution and magnitude in the bony core
[23].

During base loading vM stresses are not recorded for the contralateral side of the nasal
bone, because a bird could decrease this stress by reducing the bite force at that side. As a con-
sequence, the resulting reaction force on the contralateral side will be lower. It is, however, very
complex to introduce this bite force tuning into our models and as such is ignored here. To
<figure><img src="image_5.png" /><figcaption>**Fig 2.** Physiological Finite Element Model results, side view, for 4 selected Darwin's finch species known to use their beaks in different ways during feeding [10, 49]: *G. fortis* (base crushing beak), *G. scandens* (probing and crushing beak), *C. olivacea* (probing beak), and *C. parvulus* (tip biting beak). Results are shown for both base (1st column) and tip (2nd column) biting simulations. Arrows indicate the location and magnitude of the calculated seed reaction forces. Warmer colors represent higher von Mises stresses. Additional FE models for the eight finches for which muscle data were available are presented in S2 and S3 Figs.</figcaption></figure>

doi:10.1371/journal.pone.0129479.g002

quantify a birds' risk of beak fracture during biting, we calculated safety factors for each species' beak, by dividing bone strength by peak stress. Values for bone strength were obtained using the linear relationship of Fyhrie and Vashishth, bone strength = 0.0061 E<sub>Bone</sub> [37]. In this study we used 45 MPa for a experimentally obtained Young's modulus of 7.3 GPa for finch beak bone [22,23]. Finally keratin thickness on top and bottom of the beak, and beak dimensions

<figure><img src="image_6.png" /><figcaption>**Fig 3.** Peak vM stress for physiological FE models (tip and base loading, eight species with available muscle data) and for scaled FE models (tip and base loading, scaled to same size as *G. fortis* and with the same muscle forces) for thirteen species of Darwin's finch (* = juvenile).</figcaption></figure>

doi:10.1371/journal.pone.0129479.g003
<figure><img src="image_7.png" /><figcaption>**Fig 4.** Top view of scaled FE models of the upper beaks of 13 species of Darwin's finches. All beaks were scaled to same size and muscle force as G. fortis, with stresses calculated for both base (left) and tip (right) biting. Warmer colors represent higher vM stresses. Note how stresses are lower during the behaviors typically employed by each species, with base crushers showing the lowest stress values during base-loading and tip crushers during tip-loading. Species with probing beaks show generally high stresses under both loading conditions. Top and side views of scaled FE models are presented in S4 and S5 Figs.</figcaption></figure>

doi:10.1371/journal.pone.0129479.g004

(depth and width relative to length) were measured and are reported for four functional groups (Table 3), and Pearson's correlations (with Bonferroni correction) of these values compared to the SFs are shown in Table 4.

## Results

The calculated von Mises stresses for the physiological boundary conditions in different Darwin's finches are presented in Figs 2, 3, S2 and S3 (see also Table 2). The force vectors representing bite forces (Figs 2, S2 and S3) correspond well with field-measured *in vivo* bite forces (Table 5) [4, 5, 25]. Note, however, that *in vivo* bite forces reported here differ proportionally from those reported previously [25], as in that study the bite forces were not corrected for the

<table><caption>**Table 3.** Characterization of keratin thickness and beak dimensions (relative to beak length) in the different functional groups.</caption><thead><tr><th></th><th>bottom keratin thickness / beak length</th><th>top keratin thickness / beak length</th><th>beak depth / beak length</th><th>beak width / beak length</th></tr></thead><tbody><tr><td>Crush (N = 3)</td><td>0.095 ± 0.015</td><td>0.023 ± 0.005</td><td>0.22 ± 0.03</td><td>0.406 ± 0.010</td></tr><tr><td>Probe and Base (N = 3)</td><td>0.067 ± 0.007</td><td>0.0205 ± 0.0013</td><td>0.161 ± 0.013</td><td>0.35 ± 0.02</td></tr><tr><td>Probe (N = 2)</td><td>0.029 ± 0.005</td><td>0.018 ± 0.007</td><td>0.120 ± 0.006</td><td>0.304 ± 0.003</td></tr><tr><td>Tip (N = 5)</td><td>0.08 ± 0.03</td><td>0.035 ± 0.011</td><td>0.22 ± 0.05</td><td>0.44 ± 0.07</td></tr></tbody></table>

Table entries are means ± standard deviations.

doi:10.1371/journal.pone.0129479.t003
<table><caption>Table 4. Pearson correlation and *p*-values between beak dimensions corrected for beak length and safety factors extracted from the different models.</caption><thead><tr><th></th><th>bottom keratin thickness / beak length</th><th>top keratin thickness / beak length</th><th>beak depth / beak length</th><th>beak width / beak length</th></tr></thead><tbody><tr><td>safety factor base biting</td><td><strong>0.80</strong> (p = 0.001)</td><td>0.25 (p = 0.40)</td><td><strong>0.71</strong> (p = 0.006)</td><td><strong>0.68</strong> (p = 0.011)</td></tr><tr><td>safety factor tip biting</td><td><strong>0.90</strong> (p < 0.001)</td><td>0.65 (p = 0.015)</td><td><strong>0.88</strong> (p < 0.001)</td><td><strong>0.88</strong> (p < 0.001)</td></tr></tbody></table>

Bold values illustrate significant correlations (p<0.05/4, Bonferroni correction) between safety factors and anatomical features.

doi:10.1371/journal.pone.0129479.t004

lever arms of the bite force set-up. Among finches that crush seeds at the base of their beaks, safety factors varied between 0.7 and 2.5 (Table 5). Moreover, safety factors during natural loading regimes (i.e. base biting) generally decrease as bite force increases. Importantly, safety factors are lower during tip loading behaviors suggesting that beaks in these species are not optimized for tip biting. For other species, safety factors varied between 1 and 2.3 and were generally more similar between the two loading conditions.

Density variation in bone can be expected and might result in a locally higher strength at high stress locations as suggested by our preliminary data. Indeed, the gray values (Table 2) at high vM stress locations were all higher than 0.5 (0 indicates the least dense bone in the sample, 1 indicates the densest bone). For these natural loading conditions, higher relative bone densities are thus observed at locations with higher stresses which can lead to increased safety factors.

The von Mises' (vM) stresses for the 13 size-scaled models are presented in Figs 3, 4, S4 and S5 (also see Table 2). The maximum vM stress in the size-scaled models differs broadly across different species (Table 2). For example, the beaks of *C. olivacea* and *P. inornata* in our scaled models show highest peak vM stresses of up to 100 MPa, indicating that these shapes are not suited to withstand loading regimes similar to those observed for *G. fortis*. The deep and wide seed crushing beaks of *G. magnirostris*, *G. fortis* and *G. fulignosa*, show lower stresses under base loading (below 30 MPa) than under tip loading, consistent with these birds' emphasis on base biting strategies in nature and again in accordance with our predictions. One can also observe moderate peak stresses for the tip biting beaks (*P. crassirostris*, *C. pauper*, *C. psittacula*,

<table><caption>Table 5. *In vivo* measured bite force (6 species, means ± standard deviations, N = number of specimens) compared with the model bite force and the model safety factors (SF).</caption><thead><tr><th></th><th>behavior</th><th>Measured force at base (N)</th><th>Measured force at tip (N)</th><th>Model force at base (N)</th><th>Model force at tip (N)</th><th>Model base biting SF</th><th>Model tip biting SF</th></tr></thead><tbody><tr><td><strong>G. fortis</strong></td><td>Base crush</td><td>23 ± 9 (N = 382)</td><td>19 ± 7 (N = 382)</td><td>30</td><td>22</td><td>1.6</td><td>1.0</td></tr><tr><td><strong>G. magnirostris</strong><br>(N = 29)</td><td>Base crush</td><td>65 ± 17 (N = 29)</td><td>44 ± 10 (N = 29)</td><td>58</td><td>42</td><td>1.1</td><td>0.7</td></tr><tr><td><strong>G. fuliginosa</strong><br>(N = 115)</td><td>Base crush</td><td>5.5 ± 1.9 (N = 115)</td><td>4.6 ± 1.6 (N = 115)</td><td>7.1</td><td>4.6</td><td>2.5</td><td>1.6</td></tr><tr><td><strong>G. conirostris (*)</strong></td><td>Probe &amp; crush</td><td></td><td></td><td>15</td><td>10</td><td>1.1</td><td>1.0</td></tr><tr><td><strong>G. scandens</strong><br>(N = 64)</td><td>Probe &amp; crush</td><td>10 ± 3 (N = 64)</td><td>7 ± 3 (N = 64)</td><td>8.9</td><td>6.1</td><td>2.0</td><td>1.3</td></tr><tr><td><strong>G. difficilis</strong></td><td>Probe &amp; crush</td><td></td><td></td><td>3.6</td><td>2.4</td><td>2.1</td><td>2.3</td></tr><tr><td><strong>C. olivacea</strong><br>(N = 18)</td><td>Probe</td><td>2.0 ± 0.5 (N = 18)</td><td>1.2 ± 0.4 (N = 18)</td><td>2.0</td><td>1.0</td><td>1.0</td><td>1.4</td></tr><tr><td><strong>C. parvulus</strong><br>(N = 29)</td><td>Tip crush</td><td>5.6 ± 1.3 (N = 29)</td><td>4.2 ± 1.2 (N = 29)</td><td>3.4</td><td>2.2</td><td>2.1</td><td>1.8</td></tr></tbody></table>

doi:10.1371/journal.pone.0129479.t005

<footnote>(* Juvenile data)</footnote>
except *C. parvulus*), ranging from 31 MPa to 39 MPa. The probe and crush beaks (*G. scandens*, *G. difficilis*) show higher stresses, however, ranging from 43 to 51 MPa. Finally, comparatively low maximum stresses are found during tip loading for species that use the tips of their beak during foraging and have beaks with high curvature, as predicted (e.g. *Camarhynchus pauper* and *C. psittacula*, 35 MPa). The crushing beaks (*G. fortis*, *G. fuliginosa*, *G. magnirostris*) on the other hand have higher stresses (42 to 47 MPa, SF = ±1) under these loading conditions. The main exceptions to the otherwise rather precise tuning of beak strength and feeding mode concerns relatively high stresses observed for *G. conirostris* (67 MPa), *C. parvulus* (55 MPa), and *Cactospiza pallida* (80 MPa) under different loading regimes.

The rhamphotheca was thinnest in the probers (Table 3). Base crushers have a thicker keratin layer on the bottom of the beak than tip biters, and tip biters have thicker keratin on the top of the beak. The size scaled models show a strong correlation between the beak safety factor and the thickness of the keratin at the base of the beak (Table 4). It is also interesting to note that base and tip biters seem to have similar aspect ratios of the beak (depth and width relative to length), and that these ratios are lower in the probers (Table 3). Moreover, we observe a strong and positive correlation between aspect ratio and beak safety factor (Table 4).

## Discussion

For all species, the estimated maximum vM stresses, particularly near the nasal hinge, were 4–44% lower than maximum vM stresses calculated in prior models that did not take into account the keratinous rhamphotheca [24]. An earlier study on *Padda oryzivora* [23] showed that changing the bone modulus value also had a linear effect, albeit small, on safety factor values (SFs changed from 2.5 to 3.0 if bone modulus changed from 6.7 to 7.9 GPa). Also, variation in the keratin modulus within its measured interval (Eker = 1.3 till 2.1 GPa) had a little effect on SF. However, if we decreased the keratin modulus below this, a profound effect was observed (for Eker of 1.7 GPa the SF = 2.8, while for Eker = 0.5 GPa the SF = 1.8), indicating the importance of incorporating the keratin layer in the models. Models that integrate information about multiple layers ideally consider well-defined material properties for each constituent material. These parameters were obtained on a different species of finch (*Padda oryzivora*) since experiments cannot be conducted with Darwin's finches due to the protected status of these birds. Stress regimes for multi-layered beaks using these experimentally-obtained elastic moduli of keratin and bone had been modeled with success previously for *Padda oryzivora* and then validated using digital speckle pattern interferometry [22,23].

Keratin itself is also a multi-layered structure [38,39] that typically shows an anisotropic mechanical behavior, with preferential directions of failure depending on cell orientation [38]. The mechanical behavior of keratin also depends on its hydration state [22,40]. However, in living animals the hydration state of the keratin is controlled to the underlying living tissues (epidermis and dermis). In birds specifically, the beak consists of four distinct layers: a layer of dead epidermal cells (stratum corneum), a series of living epidermal cells, the dermis, and finally the bony beak [41]. Keratin cells typically orient such that they are aligned with the principal deformation [41]. In contrast to bone, which is a stiff mineralized tissue prone to fracture, keratin often shows local failure which does not penetrate the entire structure due to the organization of the cells in the different layers of the structure. Moreover, keratin is abrasion resistant and worn layers are shed and renewed by the addition of new layers of the stratum corneum through deposition from the living epidermal cells. As such the keratinous rhamphotheca is a continuously growing structure where damaged cells are shed and replaced. Consequently, we here decided to focus on the mechanical behavior of the underlying bony beak while taking into account the mechanical behavior of the keratinous rhamphoteca.
We found that the beaks of different finch species operate within a range of safety factors be-
tween 0.7 and 2.5 (Table 5). From an engineering perspective this range cannot be regarded as
narrow, although from a biological perspective it is rather narrow if one takes into account the
large variation finches express in bite forces (1N to 65N), beak sizes (1cm to 2cm), shapes, and
biological uses. The *in vivo* realistic range of applicable safety factors is likely even narrower
than this, given that the 0.7 value was calculated for a behavior rarely observed in nature: tip-
biting in the largest species, *G. magnirostris*. A safety factor lower than 1 implies that the beak
would break if the bird would use its full muscle strength while biting at the tip. This is, howev-
er, a biting strategy rarely observed in this base-biting bird (note that this species displays a
safety factor of 1.1 for base loading). Moreover, bone remodeling could further increase safety
factors by decreasing vM stresses in highly loaded areas, a statement supported by the higher
grayscale values in these regions. In future modeling, density dependent elastic moduli and
strength could be used to provide further insights in this material tuning. Nevertheless, calcu-
lated safety factors imply that beaks can withstand loading under normal conditions, although
the exceptional jaw muscle hypertrophy of the ground finches [27], particularly *G. magnirostris*
(S1 Fig), introduces risks to the beak's structural integrity during tip biting.

As a further caveat, the magnitude of the safety factors given here should be interpreted as
relative only when comparing different species, since our measures for the Young's modulus
were based on data for a different species (*P. oryzivora*). Moreover, we used the linear approxi-
mation of Fyhrie and Vashishth [37]. Indeed, a wide range of strengths for denser bone,
ranging from 106 to 224 MPa [42–44] is available in literature and Darwin's finches could
potentially have evolved denser bone with a higher strength. Moreover, stress magnitudes are
sensitive to changes in material properties [23,45–48]. Another assumption we make in calcu-
lating the safety factors is that the beaks of Darwin's finches have similar material properties as
those measured for *P. oryzivora*. We are, however, unable to test this directly since the preser-
vation of the tissues can have a major impact on the mechanical properties of bone and keratin.
In addition, freshly killed animals cannot be obtained.

A unique advantage of FE models, applied in our study, is that they can be size-scaled, thus
allowing a size-independent view on the structural merits of different beak shapes. In our set of
FE models, for which we scaled the beaks of all Darwin's finches to the same size and applied
identical muscle forces [26], our primary result is that maximum vM stresses differ broadly
across different species, in ways that align with the different species feeding strategies (Table 2).
For example, the beaks of *C. olivacea* and *P. inornata* show highest peak vM stresses in our
models (up to 100 MPa). This suggests that their beak shape is not suited to crack very hard
seeds. Indeed, these birds use their beaks almost exclusively to capture insects, rather than to
crush or manipulate hard objects [49].

The high stresses for *G. conirostris* observed in our models might be explained by the fact that the geometry included in our study was that from a juvenile, while the applied muscle forces on this juvenile geometry were those from an adult. Likely, ontogenetic changes in shape and ossification take place to optimize the beak for adult loading conditions [50] as previously demonstrated for *G. fortis*. Results for *C. parvulus* and *C. pallida* remain puzzling, but suggest that these species may be functionally constrained in the use of their beak when applying large muscle forces.

We observed a thicker rhamphotheca in base crushing and tip-biting birds compared to
birds that probe and bite, as well as dedicated probers (Table 3). Moreover, the keratin thick-
ness at the bottom of the beak was strongly correlated to the beak safety factor indicating that
thicker keratin helps protect the beak and increase its safety factor (Table 4). Interestingly, the
thicker top keratin for tip biters matches the patterns of stress generated during tip loading,
suggesting that in both ecotypes keratin is an important part of the stress mitigation strategy.
In addition, positive and high correlations between beak aspect ratio and beak safety factor
indicate that birds with relatively deeper and wider beaks have higher safety factors (Table 4) and are thus better equipped to withstand loading due to biting. This provides evidence that beak shape affects the risk of failure.

In summary, our FE models demonstrate that beak shapes are generally well-suited for mitigating risk of fracture in accordance with a species' predominant feeding habitat. Understanding how beaks evolve to reduce risk of fracture may help to explain patterns of selection on beak size and shape in natural populations, and ultimately should be considered as one of the axes of adaptation and specialization in the Darwin's finch radiation.

## Data availability

Data is available from biomesh.org.

## Supporting Information

S1 Fig. Muscle mass and tip bite force compared to the body mass of finches. Note how *G. magnirostris* and *G. fortis* have an exceptionally high muscle mass and bite force for their size.
(DOCX)

S2 Fig. Top and side view for physiological FE models of upper beak during base biting for 8 Darwin finches.
(DOCX)

S3 Fig. Top and side view for physiological FE models of upper beak during tip biting for 8 Darwin finches.
(DOCX)

S4 Fig. Top and side view for scaled FE models of upper beak during base biting for 13 Darwin finches.
(DOCX)

S5 Fig. Top and side view for scaled FE models of upper beak during tip biting for 13 Darwin finches.
(DOCX)

## Acknowledgments

We thank S. Maas and J. Weis (University of Utah) for allowing us to use the FEBio software package. We also thank three anonymous reviewers for constructive comments on an earlier version of the manuscript. Field work was coordinated through the Charles Darwin Research Station and the Galápagos National Park Service (GNPS). We are particularly grateful for the generosity of the GNPS in granting a salvage permit, and to the California of Sciences for the loan of material. The authors thank Luis De Leon, Ana Gabela, Andrew Hendry, Mike Hendry, Eric Hilton, Sarah Huber, Katleen Huyghe, and Bieke Vanhooydonck for their assistance in the field. This work was performed in part at the Center for Nanoscale Systems (CNS), a member of the National Nanotechnology Infrastructure Network (NNIN).

## Author Contributions

Conceived and designed the experiments: JS, PA, JD, AH, DA. Performed the experiments: JS, AG, PA, AH, JP. Analyzed the data: JS, AG, PA, AH. Wrote the paper: JS, AG, JP, PA, JD, AH, DA.
## References

1. Temeles EJ, Kress WJ. Adaptation in a Plant-Hummingbird Association. Science. 2003; 300:630-3. PMID: 12714743

2. Temeles EJ, Koulouris CR, Sander SE, Kress WJ. Effect of flower shape and size on foraging performance and trade-offs in a tropical hummingbird. Ecology. 2009; 90:1147-61. PMID: 19537537

3. Boag PT, Grant PR. Intense Natural Selection in a Population of Darwin's Finches (Geospizinae) in the Galapagos. Science. 1981; 214:82-5. PMID: 17802577

4. Herrel A, Podos J, Huber SK, Hendry AP. Evolution of bite force in Darwin's finches: a key role for head width. J. Evol. Biol. 2005; 18:669-75. PMID: 15842496

5. Herrel A, Podos J, Huber SK, Hendry AP. Bite performance and morphology in a population of Darwin's finches: implications for the evolution of beak shape. Funct Ecol. 2005; 19:43-48.

6. van der Meij M, Bout R. Scaling of jaw muscle size and maximal bite force in finches. J Exp Biol. 2004; 207:2745-53. PMID: 15235003

7. Bock WJ. Kinetics of the avian skull. J Morphol. 1964; 114:1-42.

8. Bout RG. Biomechanics of the avian skull. In: Gasc JP, Casinos A, Bels VL, editors. Vertebrate Biomechanics and Evolution II. Oxford: BIOS Scientific Publishers Ltd; 2002. pp. 49-62.

9. van der Meij MAA, Bout RG. The relationship between shape of the skull and bite force in finches. J Exp Biol. 2008; 211:1668-80. doi: 10.1242/jeb.015289 PMID: 18456894

10. Bowman RI. Morphological Differentiation and Adaptation in the Galapagos Finches. University of California press: Berkeley and Los Angeles; 1961.

11. Nuijens F, Zweers F. Characters discriminating two seed husking mechanisms in finches (Fringillidae: Carduelinae) and Estrildids (Passeridae: Estrildinae). J Morphol. 1997; 232:1-33.

12. Campàs O, Mallarino R, Herrel A, Abzhanov A, Brenner MP. Scaling and shear transformations capture beak shape variation in Darwin's finches. PNAS. 2010; 107(8):3356-60. doi: 10.1073/pnas.0911575107 PMID: 20160106

13. Herrel A, Podos J, Vanhooydonck B, Hendry AP. Force-velocity trade-off in Darwin's finch jaw function: a biomechanical basis for ecological speciation? Funct Ecol. 2009; 23:119-125.

14. Grant PR. Inheritance of size and shape in a population of Darwin's finches. Proc R Soc Lond B. 1983; 212:403-432.

15. Grant PR. Ecology and Evolution of Darwin's Finches. Princeton Univ: Princeton; 1999.

16. Grant PR, Grant RB. Evolution of character displacement in Darwin's finches. Science. 2006; 313:224-6. PMID: 16840700

17. Lack D. Darwin's Finches. Cambridge Univ Press: Cambridge; 1947.

18. Foster D, Podos J, Hendry AP. A geometric morphometric appraisal of beak shape in Darwin's finches. J Evol Biol. 2008; 21:263-75. PMID: 18021202

19. Rayfield EJ. Finite element analysis and understanding the biomechanics and evolution of living and fossil organisms. Annu. Rev. Earth Planet. Sci. 2007; 35: 541-576.

20. Ross CF. Finite element analysis in vertebrate biomechanics. Anat Rec. 2005; 283A:253-8.

21. Richmond BG, Wright BW, Grosse I, Dechow PC, Ross CF, Spencer MA, et al. Finite element analysis in functional morphology. Anat Rec Part A. 2005; The Anatomical Record Part A: 283(2), 259-74. PMID: 15747355

22. Soons J, Herrel A, Aerts P, Dirckx JJJ. Determination and validation of the elastic moduli of small and complex biological samples: bone and keratin in bird beaks. J R Soc Interface. 2012; doi: 10.1098/rsif.2011.0667

23. Soons J, Herrel A, Genbrugge A, Adriaens D, Aerts P, Dirckx J. Multi-layered bird beaks: a finite-element approach towards the role of keratin in stress dissipation. J R Soc Interface. 2012; 9(73), 1787-96. doi: 10.1098/rsif.2011.0910 PMID: 22337628

24. Soons J, Herrel A, Genbrugge A, Aerts P, Podos J, Adriaens D, et al. Mechanical stress, fracture risk and beak evolution in Darwin's ground finches (Geospiza). Philos T Roy Soc B. 2010; 365:1093-8. doi: 10.1098/rstb.2009.0280 PMID: 20194171

25. Herrel A, Soons J, Aerts P, Dirckx J, Boone M, Jacobs P, et al. Adaptation and function of the bills of Darwin's finches: divergence by feeding type and sex. Emu. 2010; 110(1), 39-47.

26. Dumont E R, Grosse I R, Slater G J. Requierements for comparing the performance of finite element models of biological structures. J Theor Biol. 2009; 256:96-103. doi: 10.1016/j.jtbi.2008.08.017 PMID: 18834892
27. Genbrugge A, Herrel A, Boone M, Van Hoorebeke L, Podos J, Dirckx J, et al. The head of the finch: the anatomy of the feeding system in two species of finches (Geospiza fortis and Padda oryzivora). J Anat. 2011; 219,6: 676-95. doi: 10.1111/j.1469-7580.2011.01437.x PMID: 21999913

28. Loeb GE, Gans C. Electromyography for experimentalists. University of Chicago Press, Chicago; 1986.

29. Mendez J, Keys A. Density and composition of mammalian muscle. Metabolism. 1960; 9, 184-8

30. Nigg B, Herzog W. Biomechanics of the Musculo-skeletal System. John Wiley & Sons Ltd, University of Calgary, 3rd edition; 2007. PMID: 18001971

31. Nuijens FW, Zweers GA. Characters discriminating two seed husking mechanisms in finches (Fringillidae: Carduelinae) and estrildids (Passeridae: Estrildinae). J. Morp. 1997; 232(1):1-33

32. Si H. TetGen: A Quality Tetrahedral Mesh Generator and three-dimensional Delaunay Triangulator v1.4. Control at <http://tetgen.berlios.de>; 2006.

33. Maas S, Weiss JA. FEBio: Finite Elements for Biomechanics. User's Manual, Version 1.0. Online publication; 2008.

34. Carter DR, Hayes WC. Bone compressive strength: the influence of density and strain rate. Science. 1976; 194:1174-6 PMID: 996549

35. Rice JC, Cowin SC, Bowman JA. On the dependence of the elasticity and strength of cancellous bone on apparent density. J biomech. 1988; 21.2 155-68. PMID: 3350829

36. Ford CM, Keaveny TM. The dependence of shear failure properties of trabecular bone on apparent density and trabecular orientation. J biomech. 1996; 29.10: 1309-17. PMID: 8884476

37. Fyhrie DP, Vashishth D. Bone stiffness predicts strength similarly for human vertebral cancellous bone in compression and for cortical bone in tension. Bone. 2000; 26:169-73. PMID: 10678412

38. Farren L, Shayler S, Ennos AR. The fracture properties and mechanical design of human fingernails. J Exp Biol. 2004; 207.5: 735-41.

39. Chou SF, Overfelt RA, Miller ME. Anisotropic mechanical behavior of keratin tissue from quill shells of North American porcupine (Erethizon dorsatum). Mater Sci Eng A-Struct. 2012; 557:36-44.

40. Li BW, Zhao HP, Feng XQ, Guo WW, Shan SC. Experimental study on the mechanical properties of the horn sheaths from cattle. J exp biol. 2010; 213(3):479-86. doi: 10.1242/jeb.035428 PMID: 20086133

41. Genbrugge A, Adriaens D, Kegel B, Brabant L, Hoorebeke L, Podos J, et al. Structural tissue organization in the beak of Java and Darwin's finches. J anat. 2012; 221, 5: 383-93. doi: 10.1111/j.1469-7580.2012.01561.x PMID: 22938039

42. Yamada H. Strength of Biological Materials. Williams & Wilkins: Baltimore; 1970.

43. Vogel S. Comparative Biomechanics: Life's Physical World. Princeton University Press: Princeton; 2003.

44. Currey JD. Bones: Structure and Mechanics. Princeton University Press: Princeton; 2006.

45. Panagiotopoulou O, Kupczik K, Cobb SN. The mechanical function of the periodontal ligament in the macaque mandible: a validation and sensitivity study using finite element analysis. J. Anat. 2011; 218,75-86. doi: 10.1111/j.1469-7580.2010.01257.x PMID: 20584094

46. Berthaume MA, Dechow PC, Iriarte-Diaz J, Ross CF, Strait DS, Wang Q, et al. Probabilistic finite element analysis of a craniofacial finite element model. j theor boil. 2012; 300: 242-53. doi: 10.1016/j.jtbi. 2012.01.031 PMID: 22306513

47. Kupczik K, Dobson CA, Fagan MJ, Crompton RH, Oxnard CE, O'Higgins P. Assessing mechanical function of the zygomatic region in macaques: validation and sensitivity testing of finite element models. J Anat. 2007; 210,1: 41-53. PMID: 17229282

48. Fitton LC, Proa M, Rowland C, Toro-ibacache V, O'Higgins P. The impact of simplifications on the performance of a finite element model of a Macaca fascicularis cranium. The Anat Rec. 2015; 298,1: 107-21. doi: 10.1002/ar.23075 PMID: 25339306

49. Grant P R, Grant B R. How and why Species Multiply: The radiation of Darwin's Finches. Princeton University Press: Princeton & Oxford; 2008.

50. Genbrugge A, Heyde AS, Adriaens D, Boone M, Van Hoorebeke L, Dirckx J, et al. Ontogeny of the cranial skeleton in a Darwin's finch (Geospiza fortis). J anat. 2011; 219, 2: 115-31. doi: 10.1111/j.1469-7580.2011.01388.x PMID: 21599660