Changing diets over time: knock-on effects of marine megafauna overexploitation on their competitors in the southwestern Atlantic Ocean

Abstract. This study compares the δ15N values and the trophic position of two seabird species throughout the late Holocene in three regions in the southwestern Atlantic Ocean to assess the hypothesis that the decimation of megafauna led to changes in the trophic position of mesopredators. Modern and ancient mollusk shells were also analyzed to account for changes in the isotopic baseline through time. Results revealed that modern Magellanic penguins have higher δ15N values than their ancient conspecifics in the three regions, after controlling for changes in the isotopic baseline. This was also true for modern Imperial shags compared with ancient unidentified cormorants/shags from the two areas where ancient specimens were recovered (southern Patagonia and the Beagle Channel). Such temporal variability might be caused by three non–mutually exclusive processes: decreased availability of pelagic squat lobster resulting from decreasing primary productivity through the late Holocene, increased availability of small fishes resulting from the sequential depletion of other piscivores (South American fur seal and sea lion and Argentine hake) since the late eighteenth century, and modification of the migratory patterns of Magellanic penguins. Although disentangling the relative contribution of all those processes is impossible at this time, the results reported here demonstrate that the ecology of Magellanic penguins and Imperial shags has undergone major changes since the late Holocene.


Introduction
Dwindling numbers of large marine predators are one of the most pervasive signatures of ecosystem overfishing (Jackson et al. 2001), and the removal of top predators often spreads through food webs as trophic cascades (Frank et al. 2005;Mumby et al. 2006) and competitor release (Laws 1977;Aalto and Baskett 2013;Surma et al. 2014). The southwestern Atlantic Ocean is no exception, and the sequential development of industrial whaling, sealing, and fishing largely reduced the populations of many marine mammals and predatory fishes between the eighteenth and twentieth centuries (Vales et al. 2020 and references therein). In parallel, most coastal predators in the southwestern Atlantic Ocean shifted their diets following the decimation of their own populations, which resulted in a significant increase in their trophic positions and a reduction in the degree of individual trophic specialization (Drago et al. 2009(Drago et al. , 2017Zenteno et al. 2015;Vales et al. 2017;Bas et al. 2019Bas et al. , 2020b. In contrast to marine mammals and large fishes, seabirds have not been intensely exploited by humans in the southwestern Atlantic Ocean, although they were consumed regularly by hunter-fisher-gatherer people inhabiting the region since the middle Holocene (Tivoli and Zangrando 2011;Borella and Cruz 2012;Zangrando and Tivoli 2015), and European sailors and settlers hunted them for oil and collected their eggs (Armstrong 1994;Cruz et al. 2010;Grosso 2016 and references therein).
Magellanic penguins (Spheniscus magellanicus) and Imperial shags (Leucocarbo atriceps) are currently the most abundant coastal seabirds off Patagonia and nest all the way from latitude 42°S to latitude 54°S (Frere et al. 2005;Schiavini et al. 2005). Magellanic penguins have dramatically increased both their population size and geographic range during the twentieth century (Boersma et al. 1990;Bouzat et al. 2009;Raya Rey et al. 2014). Information about the Imperial shags is scarcer and does not reveal any consistent trend in the region, although most colonies increased since the 1990s (Frere et al. 2005;Raya Rey et al. 2014).
Magellanic penguins foraging off Patagonia feed mainly on small pelagic fish and juvenile Argentine hake, as well as squid and crustaceans (Frere et al. 1996;Scolaro et al. 1999;Schiavini et al. 2005;Scioscia et al. 2014). This pattern is reversed in Tierra del Fuego, where pelagic crustaceans prevail in the diet of Magellanic penguins, although this has been a recent change (Dodino et al. 2020). In contrast, Imperial shags feed primarily on benthic fishes, although small pelagic fishes and cephalopods are also consumed (Gosztonyi and Kuba 1998;Punta et al. 2003;Harris et al. 2016).
Interestingly, the expansion of Magellanic penguins in Atlantic Patagonia has paralleled the decline of otariids and hake (Boersma et al. 1990;Koen-Alonso and Yodzis 2005;Vales et al. 2015), and all three are preying largely on small pelagic fishes and squid (Angelescu and Prensky 1987;Koen Alonso et al. 2000;Baylis et al. 2014;Vales et al. 2015). Predator decimation has resulted in the increase in biomass of small pelagic fishes and squid (Koen- Alonso and Yodzis 2005;Sánchez et al. 2012), thus suggesting that the expansion of Magellanic penguins might have resulted from a reduction of competition and increased food availability (Boersma et al. 1990).
Stable isotope analysis is a valuable technique to reconstruct historical changes in the diet of predators over long periods, because the stable isotope ratios in their tissues integrate those of their prey (Bearhop et al. 2004), although turnover rates vary across tissues and integrate dietary information at different temporal scales (Bearhop et al. 2004). Bone has a slow turnover rate, and stable isotope ratios in the bone organic matrix average the isotopic signatures of prey during several years (Tieszen et al. 1983;Hobson and Clark 1992) and hence offer a proxy for individual variability comparable to repeated measurements of other tissues . This is particularly useful in birds, as diets may vary largely during breeding and nonbreeding seasons (Hobson and Clark 1992;Silva et al. 2014).
This study aims to test the hypothesis that Magellanic penguins and Imperial shags currently have increased their trophic position compared with their conspecifics living during the middle and late Holocene. To do this, we compare the stable isotope ratios of nitrogen in the bone tissue of modern Magellanic penguins and Imperial shags with those of ancient conspecifics recovered from archaeological sites throughout Atlantic Patagonia and the Beagle Channel coasts. Furthermore, we compare resource partitioning between ancient King penguins (Aptenodytes patagonicus) and Magellanic penguins from the late Holocene of southern Patagonia with that of their modern conspecifics in the nearby Malvinas/Falkland Islands, as they are not sympatric anymore on the mainland.

Materials and Methods
Study Area and Sample Collection.-Both modern and archaeological samples were collected along the Argentine coast and grouped in three regions: northern Patagonia (Río Negro and Chubut Provinces), southern Patagonia (Santa Cruz Province and the Atlantic coast of Tierra del Fuego region), and the Beagle Channel (Tierra del Fuego Province). These three regions present different oceanographic features and distinct isotopic baselines . Additionally, with the purpose of comparing stable isotope ratios of ancient seabirds (this study) with those in their modern counterparts (Cherel et al. 2002;Weiss et al. 2009) in southern Patagonia, a fourth area was defined: the Malvinas/Falkland Islands (Fig. 1).
Bone samples of ancient seabirds (Magellanic penguins, King penguins, and unidentified cormorants/shags) and shells of ancient limpets (Nacella magellanica), Chilean mussels (Mytilus chilensis), and ribbed mussels (Aulacomya atra) were collected from archaeological sites dating back to the middle and late Holocene ( Fig. 1 and Supplementary Table 1, respectively). Mollusk shells were used to reconstruct the isotopic baseline for each region and period. The stable isotope ratios for a few ancient fish species were obtained from the literature (Supplementary Table 1).
In order to avoid pseudoreplication, bones from the neurocranium (fishes) or long limb bones (birds) with the same laterality were used for both ancient and modern specimens. The skeletal elements analyzed for each species varied across archaeological sites, because of uneven occurrence, but intra-individual variability in the stable isotope ratios of skeletal elements of fish and birds is usually much smaller than interindividual variability (Bas and Cardona 2018;Hyland et al. 2021).
Stable Isotope Analysis.-Modern samples were stored in a freezer at −20°C until analysis. Soft tissues were removed from the seabird bone and mollusk shells, rinsed with water, and allowed to dry at room temperature. Fishes were thawed at room temperature, boiled between 5 and 10 minutes, and dissected to remove the selected bones. Shells and bones were latter dried in a stove at 60°C for 24 hours.
Once dry, each sample was ground to fine powder using a mortar and pestle. Powdered shell samples were first demineralized by soaking in 1 N HCl until no more CO 2 was released  Supplementary Table S1) with an extension showing the species sampled there (black and white animals). Gray dots denote the modern sampling locations with an extension showing the species sampled there (grayscale animals). (Saporiti et al. 2014a,b), rinsed with distilled water for 24 hours, and dried again for 24 hours at 50°C, and then lipids were removed through sequential rinses with a 2:1 chloroform:methanol solution until the solution was transparent (Folch et al. 1957). Samples were then dried again for 24 hours at 50°C, and 0.5 mg of each sample was weighed into a tin cup. Lipids were removed from dry, powdered bone samples as described above, but they were demineralized with 0.5 N HCl (Newsome et al. 2006;Bas and Cardona 2018). Then, samples were dried again for 24 hours at 50°C, and 0.3 mg of each sample was weighed into a tin cup. Acidification has no significant effect on the δ 15 N value of mollusk shells (Carmichael et al. 2008) or bone collagen (Tuross et al. 1988). Tin cups were combusted at 900°C and analyzed in a continuous-flow isotope ratio mass spectrometer (Flash 1112 IRMS Delta C Series EA, Thermo Finnigan).
Abundance of stable isotopes is expressed using the δ notation, where the relative variations of stable isotope ratios are expressed as per mil (‰) deviations from a predefined reference scale: atmospheric nitrogen for δ 15 N. Stable isotopic reference materials of known 15 N/ 14 N ratios, as given by the International Atomic Energy Agency (Vienna, Austria), were used for calibration. Isotopic reference materials were employed to recalibrate the system once every 12 samples and were analyzed in order to compensate for any measurement drift over time. The raw data were recalculated taking into account a linear regression previously calculated for isotopic reference materials (Skrzypek 2013).
Statistical Analysis.-The stable isotope ratios of modern and ancient organisms cannot be compared directly, because the isotopic baseline may vary temporally (Casey and Post 2011). Nonetheless, the proteins that make up the organic matrix of mollusk shells are preserved and offer suitable material to reconstruct the changes in the isotopic baseline (Casey and Post 2011;Drago et al. 2017;Misarti et al. 2017). First, the offset between the average δ 15 N of ancient limpets and mussels and that of their modern conspecifics was calculated when differences were statistically significant, and that amount was later subtracted from the δ 15 N values of ancient vertebrate samples to allow comparison with modern values (called the "correction factor"; Table 1).
Second, after correction for any baseline shift according to molluscan stable isotope ratios for each area and period, ancient and modern values of δ 15 N of Magellanic penguins and cormorants/shags were compared (see "Results"). It should also be noted that the published δ 15 N values of modern grenadier and eelpout were obtained from demineralized and delipided bone (Zangrando et al. 2016) and muscle samples (Riccialdelli et al. 2017), respectively. According to Ankjaerø et al. (2012), δ 15 N values from muscle in fishes did not differ from those to bone collagen, and hence we used δ 15 N values from muscle of eelpout to compare it with the stable isotope ratios of bones from ancient conspecifics. In addition, published δ 15 N values of modern King and Magellanic penguins from the Malvinas/Falkland Islands were obtained from feathers (Weiss et al. 2009) and blood cell samples (Cherel et al. 2002), respectively. Therefore, the δ 15 N values of blood cells from Magellanic penguins were converted to those expected for feathers according to the offset between these two tissues of King penguins to allow comparison (Cherel et al. 2005a,b).
Third, the trophic position of each predator (TP p ) was calculated as: where δ 15 N p is the δ 15 N average values of each predator; δ 15 N m is the δ 15 N average value of mollusks; "3" corresponds to the trophic discrimination factor; and mussels and limpets were considered herbivores at TP = 2 (Caut et al. 2009). Then, ancient and modern values for the trophic position (TP) of Magellanic penguins and cormorants/shags were compared. All ancient and modern δ 15 N and TP values were compared independently using general linear models (GLM) as run in IBM SPSS Statistics (v. 23.0.0.2 for Mac), with two fixed factors (species and period) for invertebrates and one fixed factor (period) for the seabird species, unless otherwise stated. Then, Tukey's (HSD) post hoc tests were run to assess the temporal variation of the δ 15 N and TP values in each area. The Bonferroni correction was used to adjust α levels per test depending on compared periods per area, respectively.

Results
The δ 15 N values of all groups (species × period) were normally distributed and fulfilled the homoscedasticity requirement. The shells of ancient mollusks were always enriched in 15 N compared with those of their modern conspecifics, and differences in the δ 15 N values of mollusks of all periods were statistically significant (Table 1, Supplementary Fig. 1, Supplementary  Table 3), thus revealing the existence of a parallel drop in the δ 15 N baseline in the three regions during the past 2000 years. Accordingly, the stable isotope ratios from ancient penguins and unidentified cormorants/shags were transformed using the correction factors to allow comparison with those of modern conspecifics from the same region (Tables 1, 2).
The δ 15 N values and the mean trophic position of modern Magellanic penguins were significantly higher than those of their ancient conspecifics in the three regions (Table 2, Fig. 2). In addition, statistically significant differences were also observed between ancient samples from different periods in the Beagle Channel (Table 2, Fig. 2). Likewise, the δ 15 N values and the mean trophic position of modern Imperial shags were significantly higher than those of ancient unidentified cormorants/shags in southern Patagonia and the Beagle Channel (Table 2, Fig. 2). Ancient Magellanic and King penguins from southern Patagonia differed in their average δ 15 N and hence foraged at a different trophic position (Student's t-test; Table 3), but the difference was not as large as currently exists in the Malvinas/Falkland Islands.
We lack bone samples from ancient fishes from northern Patagonia, but the δ 15 N values of ancient Magellanic penguins were only 4‰ higher than those of contemporary ribbed mussels (Fig. 3A), thus suggesting a lower trophic TABLE 1. Results of general linear model (GLM) with two fixed factors (species and period) performed to assess the temporal variation of the δ 15 N values in shells and, when necessary, compensate for any isotopic baseline shift between the periods considered. N is sample size; δ 15 N (‰) is reported as mean ± SD. Correction factor (CF) was calculated by difference between mean isotope values of mollusks of modern and ancient samples. *Statistically significant differences ( p < 0.05) between ancient and modern samples. † Stable isotope data from Bas et al. (2020b) position than that of their modern conspecifics. Currently, the δ 15 N values of Magellanic penguins from northern Patagonia are similar to those of banded cusk eel, pink cusk eel, and cod icefish; 1‰ higher than those of hake; 3‰ higher than those of Argentine anchovy; and 6‰ higher than the average value of limpets and ribbed mussels (Fig. 3B). The δ 15 N values of ancient unidentified cormorants/shags from southern Patagonia were highly variable, but they were always in between those of predatory fishes such as hake, eelpout, snoek, and pink cusk eel and only 3‰-4‰ higher than the average of limpets and mussels (Fig. 4A-C). However, ancient Magellanic penguins from southern Patagonia were always extremely depleted in 15 N, and their δ 15 N values were so low that they did not differ from those of contemporary limpets in Teis XI, Cabo Vírgenes, and Margen Sur ( Fig. 4A-C), thus suggesting that they were not foraging locally most of that time. The same was true for ancient King penguins (Fig. 4B). Current Magellanic penguins and Imperial shags from southern Patagonia have similar δ 15 N values and higher values than those of fishes and mollusks: 1‰ above those of eelpout; 3‰ higher than those of hake, pink cusk eel, and Patagonian blenny; 6‰ above those of Patagonian grenadier; and 7‰ higher than the average of limpets and mussels (Fig. 4D).
Finally, the δ 15 N values of ancient Magellanic penguins from the Beagle Channel were only 3‰-5‰ higher than those of mollusks ( Fig. 5A-C), and ancient unidentified cormorants/shags had extremely low values of δ 15 N (Fig. 5B). Conversely, both species of modern seabirds have δ 15 N values similar to those of eelpout and cod icefish and 8‰ higher than the average of limpets and mussels (Fig. 5D).

Discussion
The protocol used here would allow obtaining unbiased δ 13 C and δ 15 N values for bone  collagen, as they are measured using the standard procedure for this type of sample (Newsome et al. 2006;Guiry et al. 2016;Bas and Cardona 2018;Guiry and Hunt 2020). DeNiro (1985) reported carbon to nitrogen (C:N) atomic ratios of bone collagen to range from 2.9 to 3.6, and this has been the standard requirement for decades both in ecology and archaeology. Most of the samples analyzed here satisfied this requirement (Supplementary Table 2), but recently Szpak (2020, 2021) have reported a much narrower acceptable range (3.0-3.3). Accordingly, many samples in this study with C:N atomic ratios ranging from 3.4 to 3.6, may still contain some traces of lipid or humic acid, in modern and ancient samples, respectively, and hence might yield slightly biased δ 13 C values. For this reason, we discuss here only their δ 15 N values, as neither lipids nor humic acid contain nitrogen, and hence collagen is the only source of nitrogen in acidified bone samples (Bas and Cardona et al. 2018;Bas et al. 2020a;Guiry and Hunt 2020;Szpak 2020, 2021). It should also be noted that the organic matrix of mollusk shells is a mixture of proteins and chitin, a polysaccharide containing nitrogen (Furuhashi et al. 2009). As a result, the C:N ratio of the organic matrix of mollusk shells including equal amounts of protein and chitin is close to 5.5 and hence differs from that expected for collagen.
The results reported here reveal major changes in the δ 15 N of Magellanic penguins and cormorants/shags in the southwestern Atlantic Ocean since the middle Holocene. Certainly, sample size for some species and archaeological sites is small, but differences between ancient and modern conspecifics are so huge and consistent across areas, particularly for Magellanic penguins, that we believe that our conclusions are robust. It should be noted that male and female modern Magellanic penguins do not differ in their average δ 15 N values (Scioscia et al. 2014;Silva et al. 2014;Barrionuevo et al. 2020;Rosciano et al. 2020;Dodino et al. 2021), and the comparison of δ 15 N in ancient and modern Magellanic penguins is therefore unlikely to be affected by Arithmetic mean and standard deviation (mean ± SD) are shown for each species. Key: squares, invertebrates (INV); diamonds, pelagic fishes (PF); triangles, benthic fishes (BF); circles, air-breathing predators (ABP). See Supplementary Table 1 for acronyms. biased sex ratios. Furthermore, the use of stable isotope ratios from ancient and modern mollusk shells allowed us to account for temporal changes in the isotopic baseline (Casey and Post 2011;Misarti et al. 2017;Bas et al. 2019) and compare the stable isotope ratios of ancient and modern seabirds, although selecting the most suitable isotopic baseline for migratory species is challenging. The use of compoundspecific isotopic analyses of individual amino acids offers an alternative approach to identify the relative contribution of changes in diet or baseline to the variability of δ 15 N in consumers (Lorrain et al. 2009) but requires a much larger amount of protein (3 mg), which limits the technique to bulky samples.
The analysis of mollusk shells confirmed a drop in the δ 15 N baseline of the three regions considered during the past 2000 years. Such a pattern has already been suggested by Saporiti et al. (2014a) with a more limited data set and was likely caused by a decrease in the intensity of vertical mixing and coastal upwelling (Somes et al. 2010). If so, the primary productivity of the coastal ecosystems of the southwestern Atlantic Ocean south to latitude 40°S is currently lower than during most of the late Holocene (Saporiti et al. 2014a). Furthermore, the coastal areas of the southwestern Atlantic Ocean currently support smaller populations of top predators due to sequential sealing and fishing since the late eighteenth century (Vales et al. 2020 and references therein). Both processes might have operated synergistically to modify the diet of Magellanic penguins and cormorants/shags and might contribute to the changes in δ 15 N reported here.
The optimal foraging theory predicts an increase in the trophic position of predators far below carrying capacity as they experience higher per capita food availability. This is because predators close to satiation select prey to maximize net energy intake, whereas hungry predators close to carrying capacity are less selective and capture prey according to their encounter rates (Schoener 1971;Pulliam 1974;Stephens and Krebs 1986). On the other hand, prey density in aquatic ecosystems decreases with body size (Blanchard et al. 2017), whereas the trophic position of prey is positively correlated with body size (Jennings 2005). From here follows that aquatic carnivores close to carrying capacity rely largely on small, highly abundant prey with a low trophic level. However, in a scenario in which the populations of predators decrease as a result of harvesting, individuals can shift to scarcer but more rewarding prey with a larger body size and a higher trophic position, due to a relaxation of intraspecific competition. This mechanism may explain the increase in the δ 15 N values of Magellanic penguins and cormorants/shags following competitor release as a result of the sequential exploitation of otariids and hake since the late eighteenth century.
Currently, the diet of Magellanic penguins is dominated by fish over most of its range (Scolaro et al. 1999;Clausen and Pütz 2002;Schiavini et al. 2005;Weiss et al. 2009;Scioscia et al. 2014), although the squat lobster (Munida gregaria) is also a major component of the diet of Magellanic penguins in the Malvinas/Falkland Islands (Thompson 1993;Pütz et al. 2001;Clausen and Pütz 2002) and very recently has become their staple diet in the Beagle Channel (Dodino et al. 2020). Squat lobsters have two morphs, and the specimens of the pelagic one are highly depleted in 15 N compared with both benthic and pelagic fishes (Ciancio et al. 2008;Drago et al. 2009;Dodino et al. 2020). Accordingly, the low δ 15 N values of ancient Magellanic penguins and unidentified cormorants/shags might have resulted from increased consumption of pelagic squat lobsters during the late Holocene compared with modern conspecifics.
Squat lobsters have a much lower energy density than small pelagic fishes (Thompson 1993;Ciancio et al. 2007) and hence the former are expected to be consumed by Magellanic penguins only when the latter are scarce. This prediction fits the expectations of the competitive release hypothesis, because large populations of piscivorous fur seals, sea lions, and hake might have reduced dramatically the availability of small fishes for seabirds during the late Holocene. However, this is not the only possible explanation for the pattern reported here because of the changes in primary productivity reported earlier. Since 2012, the abundance of pelagic squat lobsters has increased markedly in the Beagle Channel, as a result of increased primary productivity (Diez et al. 2016), which in turn has triggered a major dietary shift in Magellanic penguins (Dodino et al. 2020). It is worth noting that the samples of Magellanic penguins from the Beagle Channel analyzed here were collected in 2010, before the increase in the population of pelagic squat lobsters. The trophic position of the Magellanic penguins from the Beagle Channel reported here is similar to that reported by Dodino et al. (2020) for penguins sampled in 2009, but much higher than those sampled from 2013 to 2017 and reported to rely largely on pelagic squat lobster (Dodino et al. 2020).
The abundance of pelagic squat lobsters is highly dependent on the intensity of vertical mixing and the level of primary productivity (Diez et al. 2016), and the high δ 15 N values of ancient limpets and mussels reported here suggest a more intense coastal upwelling and vertical mixing during the second half of the Holocene (Somes et al. 2010;Saporiti et al. 2014a). If so, the abundance of pelagic squat lobsters was perhaps much higher in the past, independent of the abundance of small fishes, and possibly due to more favorable environmental conditions. In this scenario, the lower trophic position of ancient Magellanic penguins and unidentified cormorants/shags might be just the result of a higher availability of pelagic squat lobsters resulting from a bottom-up processes.
The migratory behavior of Magellanic penguins (Stokes et al. 1998;Pütz et al. 2007;Barrionuevo et al. 2020;Dodino et al. 2021) and the low turnover of bones (Tieszen et al. 1983;Hobson and Clark 1992) add additional complexity to interpreting the temporal changes of δ 15 N values in these species. Currently, Magellanic penguins nesting in northern Patagonia overwinter in southern Brazil (Stokes et al. 1998), whereas those nesting in southern Patagonia and the Beagle Channel overwinter at much higher latitudes, with just a small fraction of the penguins from southern Patagonia reaching southern Brazil in winter (Pütz et al. 2007;Barrionuevo et al. 2020;Dodino et al. 2021). The baseline of δ 15 N values shifts latitudinally along the Atlantic coast of South America, from Brazil to Tierra del Fuego Vélez-Rubio et al. 2018), and hence the δ 15 N values of penguin bone collagen are influenced not only by diet but also by the foraging grounds exploited. This might also explain the differences observed in the δ 15 N values of ancient and modern Magellanic penguins, if migratory patterns changed throughout time. This is certainly the case for the ancient Magellanic penguins from southern Patagonia, as they were extremely depleted in 15 N compared with coastal contemporary species from the same area ( Fig. 4A-C), which resulted in an unrealistically low trophic position. They were likely vagrants that had foraged elsewhere for most of the time recorded in their bones and were captured by hunter-fisher-gatherer people when they approached the coast. Interestingly, the δ 15 N values of ancient Magellanic penguins from Margen Sur (885-838 cal yr BP) were only 2‰ above those of almost contemporary King penguins from Cabo Virgenes 20 (1131 cal yr BP), whereas the offset between the two species in the Malvinas/Falkland Islands is currently 4.4‰ (Fig. 4B, Table 3). Modern myctophids, which are the main prey for King penguins (Cherel et al. 1993(Cherel et al. , 1996Moore et al. 1998), are highly depleted in 15 N compared with small pelagic fishes (Ciancio et al. 2008), thus suggesting that the low δ 15 N values of ancient Magellanic penguins from southern Patagonia and their unrealistically low trophic positions might be because of an intense use of offshore foraging grounds. This is not true for modern Magellanic penguins from southern Patagonia, whose δ 15 N values fit well with those of coastal fish species (Fig. 4D).
In contrast to Magellanic penguins, the cormorants and shags inhabiting the southwestern Atlantic Ocean are sedentary, and temporal changes in their δ 15 N values are therefore likely the result of dietary shifts. Furthermore, local mollusks offer a reliable proxy to set the baseline for the calculation of the trophic positions of cormorants and shags. Interestingly, cormorants/shags from southern Patagonia and the Beagle Channel exhibited the same temporal pattern as Magellanic penguins, with a recent increase in both δ 15 N values and trophic position. It should be noted, however, that the ancient unidentified cormorants/shags from the Beagle Channel (Lanshuaia II archaeological site) were extremely depleted in 15 N compared with their contemporary marine potential prey. This divergence was interpreted by Bas et al. (2019) as an evidence of consumption of catadromous galaxid fishes, an abundant resource at that time in estuarine ecosystems of Tierra del Fuego, which are highly depleted in 15 N compared with marine fishes. Differences in the feeding habits between Magellanic penguins and unidentified cormorants/shags during the Holocene were already reported in the Beagle Channel by Kochi et al. (2018), although the δ 15 N values reported for unidentified cormorants/shags were not so low.
Although decreasing primary productivity and changes in the migratory patterns of Magellanic penguins may explain the changes in the δ 15 N values reported here, only competitor release simultaneously explains the changes in the distribution of Magellanic penguins observed during the twentieth century and the changes in the δ 15 N values of both Magellanic penguins and Imperial shags. Historical data suggest that Magellanic penguins in the southwestern Atlantic Ocean have been moving northward since the early twentieth century, when they started nesting along the northern continental coast of Patagonia (Boersma et al. 1990;Boersma 2008;Bouzat et al. 2009;Cruz et al. 2010). Genetic data also supported the hypothesis that Magellanic penguins have experienced a recent expansion during a favorable period (Bouzat et al. 2009). This expansion suggests that ecological conditions in this region during the second half of the twentieth century were better for this species than at the beginning of the last century (Boersma et al. 1990), which is not explained by declining productivity. However, competitor release resulting from the overharvesting of South American fur seals and sea lions and Argentine hake between the late eighteenth and the twentieth centuries (Lloris et al. 2005;Romero et al. 2017) might have contributed to improving the habitat conditions for Magellanic penguins and simultaneously increased the consumption of fish revealed by increased δ 15 N values. This is because the populations of those competitors, feeding on small pelagic fishes and squid (Scolaro et al. 1999;Punta et al. 2003;Schiavini et al. 2005;Baylis et al. 2014;Scioscia et al. 2014), are currently well below original numbers, despite the legal protection of pinnipeds Romero et al. 2017).

Conclusion
The trophic ecology of modern Magellanic penguins and Imperial shags in the southwestern Atlantic Ocean differs from that of conspecifics living in the same region during the late Holocene, as revealed by a recent increase in δ 15 N values. Three, non-mutually exclusive processes, namely competitor release, reduced primary productivity, and changes in migratory patterns between isotopically dissimilar regions, may explain that increase. Although disentangling the relative contribution of all those processes is not possible at this time, the results reported here demonstrate that the ecology of Magellanic penguins and Imperial shags has undergone major changes since the late Holocene, and competitor release remains as a plausible hypothesis.

Acknowledgments
We are very grateful to L. Silva for her assistance with sample processing and to M. Álvarez and J. Gómez Otero who kindly provided us with archaeological samples from Teis XI and Playa Las Lisas 2, respectively. P. Rubio helped us with isotopic analyses at Centres Científics i Tecnològics de la Universitat de Barcelona (Barcelona, Spain). I. Briz i Godino is member of the "María Zambrano" program at the University of Barcelona. This study was funded by project no. 309765 from Fundació Bosch i Gimpera. All biological samples included in this paper were obtained, transported, and analyzed following the legal terms and conditions of the Argentine government.

Data Availability Statement
The data used to support the findings of this article are available in the article and in its online Supplementary Material, which is available on Dryad at: https://doi.org/10.5061/ dryad.dbrv15f3k.