
==== Front
Nat Ecol Evol
Nat Ecol Evol
Nature Ecology & Evolution
2397-334X
Nature Publishing Group UK London

38383849
2337
10.1038/s41559-024-02337-4
Article
Genomics of post-bottleneck recovery in the northern elephant seal
http://orcid.org/0000-0002-7265-4180
Hoelzel A. Rus a.r.hoelzel@dur.ac.uk

1
http://orcid.org/0000-0002-6410-682X
Gkafas Georgios A. 2
Kang Hui 34
http://orcid.org/0000-0002-5842-6114
Sarigol Fatih 1
Le Boeuf Burney 5
http://orcid.org/0000-0002-0233-5782
Costa Daniel P. 5
Beltran Roxanne S. 5
http://orcid.org/0009-0009-1215-940X
Reiter Joanne 5
Robinson Patrick W. 5
http://orcid.org/0000-0002-6519-7671
McInerney Nancy 6
http://orcid.org/0000-0001-8594-7217
Seim Inge 37
http://orcid.org/0000-0003-1092-3186
Sun Shuai sunshuai@genomics.cn

8
http://orcid.org/0000-0001-7365-1590
Fan Guangyi fanguangyi@genomics.cn

8
http://orcid.org/0000-0003-4977-1722
Li Songhai lish@idsse.ac.cn

34
1 https://ror.org/01v29qb04 grid.8250.f 0000 0000 8700 0572 Biosciences, Durham University, Durham, UK
2 https://ror.org/04v4g9h31 grid.410558.d 0000 0001 0035 6670 Department of Ichthyology and Aquatic Environment, University of Thessaly, Volos, Greece
3 grid.9227.e 0000000119573309 Marine Mammal and Marine Bioacoustics Laboratory, Institute of Deep-sea Science and Engineering, Chinese Academy of Sciences, Sanya, China
4 grid.9227.e 0000000119573309 Innovation Research Center for Aquatic Mammals, Institute of Hydrobiology, Chinese Academy of Sciences, Wuhan, China
5 https://ror.org/05t99sp05 grid.468726.9 0000 0004 0486 2046 Ecology and Evolutionary Biology, University of California, Santa Cruz, CA USA
6 grid.1214.6 0000 0000 8716 3312 Center for Conservation Genomics, National Zoo and Conservation Biology Institute, Smithsonian Institution, Washington, DC USA
7 https://ror.org/036trcv74 grid.260474.3 0000 0001 0089 5711 Integrative Biology Laboratory, College of Life Sciences, Nanjing Normal University, Nanjing, China
8 https://ror.org/05gsxrt27 BGI Research, Qingdao, China
21 2 2024
21 2 2024
2024
8 4 686694
29 6 2023
19 1 2024
© The Author(s) 2024
https://creativecommons.org/licenses/by/4.0/ Open Access This article is licensed under a Creative Commons Attribution 4.0 International License, which permits use, sharing, adaptation, distribution and reproduction in any medium or format, as long as you give appropriate credit to the original author(s) and the source, provide a link to the Creative Commons license, and indicate if changes were made. The images or other third party material in this article are included in the article’s Creative Commons license, unless indicated otherwise in a credit line to the material. If material is not included in the article’s Creative Commons license and your intended use is not permitted by statutory regulation or exceeds the permitted use, you will need to obtain permission directly from the copyright holder. To view a copy of this license, visit http://creativecommons.org/licenses/by/4.0/.
Populations and species are threatened by human pressure, but their fate is variable. Some depleted populations, such as that of the northern elephant seal (Mirounga angustirostris), recover rapidly even when the surviving population was small. The northern elephant seal was hunted extensively and taken by collectors between the early 1800s and 1892, suffering an extreme population bottleneck as a consequence. Recovery was rapid and now there are over 200,000 individuals. We sequenced 260 modern and 8 historical northern elephant seal nuclear genomes to assess the impact of the population bottleneck on individual northern elephant seals and to better understand their recovery. Here we show that inbreeding, an increase in the frequency of alleles compromised by lost function, and allele frequency distortion, reduced the fitness of breeding males and females, as well as the performance of adult females on foraging migrations. We provide a detailed investigation of the impact of a severe bottleneck on fitness at the genomic level and report on the role of specific gene systems.

Genome analysis of modern and historical elephant seals reveals impacts of a severe bottleneck on the genomes and fitness of individual seals, and the implications for recovery.

Subject terms

Conservation biology
Population genetics
Molecular evolution
National Natural Science Foundation of ChinaNational Science FoundationOffice of Naval Researchissue-copyright-statement© Springer Nature Limited 2024
==== Body
pmcMain

Iconic species are vanishing. Our past and present activities have reduced populations to the point where they may go extinct by demographic processes alone1. When they survive, inbreeding and genetic drift may reduce the fitness of individuals and the survival potential of populations2. Nevertheless, some species survive and apparently thrive. After heavy exploitation led to a severe population bottleneck in 1892, reducing the population to ~20 individuals3, the northern elephant seal (Mirounga angustirostris; hereafter NES) recovered nearly exponentially to over 220,000 today4,5. Theoretically, a rapid demographic recovery may reduce the negative impact on genetic diversity, but this was not the case for the NES. Census data indicate a rapid recovery4,5, but years of genetic studies show profoundly reduced genetic diversity6. In one of the earliest indications of the impact of population bottlenecks on the health of a species, ref. 7 reported the lack of protein (allozyme) diversity in the surviving population of NES. Further studies on genetic diversity at a range of markers (allozymes, mitochondiral (mt)DNA, minisatellite DNA, microsatellite DNA and immune system genes) showed reduced diversity compared with southern elephant seals (Mirounga leonina, a sibling species not impacted by a similar population bottleneck3,6,8). The loss of diversity was also evident in comparisons of pre- and post-bottleneck NES DNA9,10. Small population size compounded by extreme polygyny8 probably contributed to the loss of genetic diversity. Of 150,388 species reviewed for the IUCN Red List, 42,108 are threatened and 25,615 are endangered or critically endangered (http://www.iucnredlist.org). Many of the endangered or threatened species have previously existed or currently exist as small populations. Since the survival potential of small populations can be influenced by genetic diversity2, we wonder whether the loss of genetic diversity for NES has had a measurable effect on their fitness, even though population growth has been robust during the 132 years (~22 generations) since the bottleneck. If so, they could still be vulnerable to some new environmental stress.

Here we sequenced 260 modern and 8 historical genomes, showing that inbreeding, loss of function and the distortion of allele frequencies have reduced the fitness of breeding males and females, as well as the performance of adults on foraging migrations. The loss of fitness associated with inbreeding is well documented11, but the importance of allele frequency distortion and the presence of loss of function (LOF) alleles in specific gene systems is less well understood and provides new inference about the general and lasting impact of population bottlenecks. Ecosystem function depends on biodiversity and the contribution of species within that system. Effective conservation management requires an understanding of the scope of impact from depleted populations on specific functions in ecosystem communities.

Results and Discussion

All of the modern NES samples investigated in our study were from the breeding colony at Año Nuevo, California, USA. They were chosen due to their inclusion in studies on reproductive success or diving profiles. The source of historical samples is given in Supplementary Table 1. We consider diversity across NES genomes, comparing them before and after the bottleneck, followed by inference from the genomic analyses about demographic history. We then consider fitness, first associated with reproductive success, then with diving performance.

Diversity

Most of our modern nuclear genome sequences had greater than 30X read depth (see Extended Data Fig. 1 for full range), while 8 historical genomes with degraded DNA had a broader range of coverage (14.5X ± 12.7X after mapping; Supplementary Table 1). Average heterozygosity per genome was 0.00142 ± 0.00092 (s.d.) before the bottleneck (N = 5) and 0.000176 ± 0.000013 in the modern population (N = 180 adults; Fig. 1). The size range and number of runs of homozygosity (ROH) fragments >100 kb and >1 Mb are shown for male and female modern samples in Extended Data Fig. 2. Coverage was limited for most of the historical samples, restricting the potential for accurately estimating ROH. However, a pairwise comparison for both heterozygosity per sliding 50 kb window and ROH greater than 100 kb is shown in Extended Data Fig. 3, comparing the genomes of two individuals, one from 1884 (Hist8) that sequenced well enough for these analyses and a randomly chosen modern sample from 2009. The individual from 1884 was an adult male. His diversity reflects the population of his parents, who were alive when the species was sufficiently abundant to support a commercial hunt. Each comparison shows a loss of variation after the bottleneck. For this historical genome, the proportion of genotypes homozygous for the LOF allele (among all LOF loci identified by snpEff, see Methods) was 0.0011, while the proportion of heterozygotes was 0.287. Among 180 modern adult genomes, 0.379 ± 0.016 (s.d.) of the identified loci were homozygous for the LOF allele, while 0.193 ± 0.024 were heterozygous.Fig. 1 Diversity and demographics.

a, Mean (±1 s.d.) nuclear heterozygosity for pre (N = 5), immediately post (1906–1923; N = 3) bottleneck (BN) and modern genomes (N = 180). b, Mitogenome DNA neighbour-joining tree comparing historical (blue) and modern (green) haplotypes for all historical and all modern adult samples. For the historical samples (Hist), numbers in brackets indicate dates; for the modern samples, numbers in brackets indicate the numbers of individuals with a given haplotype. Bootstrap values reflect proportional support (%) out of 1,000 replicates. c, PSMC plot for 14 modern samples showing historical demography (pale lines showing bootstrap replicates). Inset: effective population size (Ne) estimate based on linkage disequilibrium.

Source data

A mitogenome tree was constructed with the 5 pre-bottleneck samples, compared to 3 from shortly after the bottleneck (1906–1924) and 180 modern genomes (Fig. 1). Two major lineages remained after the bottleneck, and some mutations were gained in the modern population within those two lineages. The pre-bottleneck mitochondrial genomes had an average pairwise genetic distance (uncorrected percentage) between them of 0.00219 ± 0.00146 (s.d.). The two main modern haplotypes remaining (haplotypes 1 and 2; Fig. 1) differ by 0.00232, while diversity within each lineage was 0.0000048 in the haplotype 1 lineage and 0.0000081 in the haplotype 2 lineage. The tree indicates that two mtDNA lineages survived the bottleneck, consistent with earlier reports3.

Demographic history

The demographic history of the species can be estimated in deep time on the basis of coalescent analyses using genomes (with the pairwise sequential Markovian coalescent (PSMC), see Methods). The 14 modern genomes we chose at random (while ensuring 7 males and 7 females, all from the 1980s) all showed essentially the same pattern (Fig. 1). The effective population size (Ne) was ~40,000 during the last interglacial (Eemian; ~130–115 Ka) but fell to ~2,000–4,000 during the last glacial period, reaching a nadir around the last glacial maximum (~20 Ka; Fig. 1). Using the approximate effective to census population size (Ne/Nc)ratio published from a meta-analysis12 of wildlife species (~0.1), this would suggest a census population of ~20,000, perhaps explaining why the population was so quickly nearly eradicated during the nineteenth century. The population may have been closer to its current size during the Eemian. For an estimate of current Ne (post-bottleneck), we used a method13 based on linkage disequilibrium (SNeP; Fig. 1c), which indicated an Ne of 100, suggesting a very low Ne/Nc ratio (~0.00045), consistent with low post-bottleneck diversity and rapid regrowth. This result was replicated with an alternative method (GONE)14, which indicated the same current Ne (see Extended Data Fig. 4).

Reproductive fitness

The average number of successfully weaned pups per year over a female’s lifetime was known for 40 females in our dataset. This fitness estimate reflects the lifetime contribution to future generations, and numerous studies have considered correlations between heterozygosity and fitness to assess the impact of inbreeding or heterosis (see ref. 15). We tested for correlation with FROH (for ROH > 1Mb; r2 = 0.143, F = 6.35, P = 0.016) and with the number of affected alleles across 151 LOF loci (restricted to loss of start or new stop mutations; r2 = 0.114, F = 4.88, P = 0.033; Fig. 2). The regression was also significant for all 328 LOF loci found (r2 = 0.146, F = 6.49, P = 0.015; Extended Data Fig. 5). The correlation with FROH was marginally stronger when the total number of weaned pups over a lifetime was considered (N = 43; r2 = 0.19, F = 9.59, P = 0.0035; Extended Data Fig. 5), which does not control for variation in longevity, but longevity itself was not significantly correlated (N = 44; r2 = 0.024; F = 1.01; P = 0.32). There was also no significant correlation between FROH and the number of pups a female produced (fecundity; N = 44, r2 = 0.07, F = 3.16, P = 0.083). We considered possible environmental effects, although females are all from a very similar time range, born between 1981 and 1988, and first weaning between 1985 and 1991. A multiple regression with lifetime reproductive success (weaned per year) as the response variable, and female birth year and year of first weaning as explanatory variables, was not significant (adjusted r2 = −0.0096, F = 0.81, P = 0.45). When we included FROH as a third explanatory variable, FROH was significantly correlated with the response variable (P = 0.03), but none of the interaction factors between explanatory variables were significant (Supplementary Table 2).Fig. 2 Reproductive fitness, inbreeding and LOF.

Pups weaned annually over the female’s lifetime as a measure of lifetime success is compared by linear regression with a, proportional runs of homozygosity (FROH) and b, the number of affected alleles across 200 LOF loci. Regression details and significance are shown for each. The 95% confidence intervals (shown in grey shading) were calculated using the R package ggplot.

Source data

A significant correlation for pups weaned, but not for pups produced, suggests that the impairment may be more about the fitness of offspring affected by inbreeding16 rather than fecundity; however, more data on the genotypes and phenotypes of the relevant pups would be needed to test this further. The linear regression between the number of LOF loci and FROH (x axis) was positive and significant (F = 9.78, P = 0.0034, r2 = 0.3042; Extended Data Fig. 6), as expected. We also used an additional programme (PROVEAN; https://www.jcvi.org/research/provean) to identify missense mutations and confirm a positive correlation with FROH (F = 16.61, P = 0.0002; Extended Data Fig. 6). Three LOF loci were associated with oogenesis (MARF1)17, oocyte growth (KMT2B)18 and embryonic development (NBAS)19, but the profile among these loci was the same for all but three females, and there was no association with fitness. An association with diversity across the genome rather than specific loci would be consistent with the ‘general effect’ hypothesis of heterozygosity–fitness correlations20.

We next consider a possible impact on male reproductive success. In an earlier study, the paternal success of the NES alpha male M12 was low compared with that expected by his frequency of observed copulations21. That study had investigated 10 NES harems, all from the same beach and same year, and 6 southern elephant seal (SES) harems. The average paternal success of alpha bulls was significantly lower than copulatory success for NES but not for SES. One NES male stood out as having especially low success (M12). Here we used genomic data to test the paternal success of alpha males at 4 of the same harems, including the harem held by M12. A total of 31 males and 77 pups were included in the paternity tests (including the 4 focal alphas) and there were 51 paternities detected. We found again that only M12 had significantly fewer paternities than expected based on observed copulations (Table 1). M12 also had the highest FROH, although all 4 were near the middle of the distribution for all adult NES (see Extended Data Fig. 7). We found a total of 5 loci with LOF (gained a stop codon), known to be associated with male fertility (associated with sperm production or function). These were LRGUK22, MNS1 (ref. 23), TUBB4B24, SRSF3 (ref. 25) and EZR26. We focused first on homozygotes for the LOF allele in case there is some function from co-dominance. M12 was homozygous for the non-functional version of 4 of these 5 loci, while the other alphas were homozygous for the non-functional version at 1 or 2 loci (Table 1). The negative relationship between paternal success and the frequency of the affected alleles was also seen for 27 non-alpha males (Extended Data Fig. 8 and Table 1). Males gain access to females in this highly polygynous species through competition among males, hence age and overall size are important factors associated with successes. For this reason, a general relationship with inbreeding or loss of function across the genome may not be expected. We tested this by looking for a correlation between non-alpha male paternal success and both FROH and relative LOF allele frequency (out of 151 LOF loci found that generated a new start or a stop codon). Neither were significant (FROH: r2 = 0.031, F = 0.79, P = 0.38; LOF: r2 = 0.057, F = 1.5, P = 0.23).Table 1 Alpha male reproductive success

Harema	BBS (M144)	TS (M160)	MBB (M12)	BMS (M4)	
% Copulations	31.5	25.9	90.5	42.7	
% Paternity	26.7	7.7	23.5	26.9	
N (genomes)	15	13	17	26	
χ2	0.13	2.0	17.6	1.36	
P	0.72	0.16	0.00003	0.24	
FROH (>1Mb)	0.1067	0.1093	0.1145	0.0967	
Homozygosity at LOF loci	1/5	2/5	4/5	1/5	
LOF allele frequency	0.6	0.6	0.8	0.6	
aFour focal northern elephant seal harems (column headings) with the alpha male shown in brackets. The percentage of observed copulations from earlier studies is compared with calculated percent paternities determined by genomic analysis. The significance of the difference was tested using two-sided Chi-square test. The proportional FROH (indicating inbreeding) and loss of function at five identified loci that are associated with male reproductive health are given. Bold text highlights the results for M12.

These data suggest that males, including M12 and some of the 27 non-alpha males, were impacted by specific LOF loci, reducing their fertility. Deleterious alleles can be retained after allele frequencies are distorted and purifying selection is weakened by strong genetic drift. The most successful of the non-alpha males (two males with five paternities each; Extended Data Table 1) were homozygous at two or three LOF loci. Both of these males and all alphas apart from M12 had functional alleles at EZR. M12 was homozygous for the LOF allele. This locus is dysfunctional in male humans with asthenozoospermia (reduced sperm motility)26. We compared the number homozygotes for the LOF allele at these five loci (average = 2.14) for all non-alpha males that had at least one offspring with those that had no offspring (average = 2.77). The difference was marginally significant (Mann–Whitney U = −1.72, P = 0.0427).

Diving performance

A critical aspect of elephant seal life history is their extensive deep-diving foraging excursions when they accumulate fat stores to facilitate fasting during the breeding season27. We acquired dive performance data from 92 females with complete dive profile data for up to three foraging trips per female lasting 70–220 d per trip. We generated a ‘performance’ metric as the product of the deepest dive (in metres, averaged over multiple trips) times the proportion of dives deeper than 516 m (median dive depth) times the relative dive duration (relative to the longest duration among individuals). Diving performance was not associated with metrics of genomic diversity (see Extended Data Fig. 9 for a correlation test with FROH). We then categorized each seal as having high or low performance (either side of the midpoint). These categories were not differentiated with respect to LOF allele frequency (t = 1.20, d.f. = 90, P = 0.232 for the LOF loci resulting in lost start or new stop codons and t = 0.126, P = 0.90 for all LOF loci). We then compared these two categories at the five variable loci with either non-synonymous mutations or upstream mutations (that may be associated with transcriptional regulation) associated with hypoxia found in our genomes and supported by publications about their relevant function. Two hypoxia loci were variable at non-synonymous sites within the gene: HIF3A (involved in hypoxia gene expression)28 and SETX (protects cells from DNA damage induced during transcription in hypoxia)29. Three other loci had segregating sites upstream of the genes: MB (myoglobin; oxygen provision)30, HIF1A (transcription regulator in response to hypoxia)30 and HYOU1 (cyto-protection during oxygen deprivation)31. For all five of these loci, minor allele frequency (MAF) was higher for those individuals that had lower performance scores (Χ2 = 9.789, P = 0.0017; Fig. 3). We compared 51 individuals born in the 1980s with 57 born between 2005 and 2015 to see whether there had been a change in MAF at these loci over time. All but one locus had diminished MAF over the intervening four or five generations (Table 2). This suggests that there may be purifying selection against these alleles (comparing MAF for all five loci combined, Χ2 = 2.76, P = 0.0483; see Table 2 for details for each locus), but tests involving a longer timeframe would help confirm this. We have no data on the dive performance of the females born in the 1980s, but the range is probably comparable. For our historical genomes, it was possible to genotype eight single-nucleotide polymorphisms (SNPs) among four of these loci (all but HYOU1) for three (4 SNPs|) or four (4 SNPs) seals. MAF ranged from 0.125 to 0.167 (average = 0.151) compared with an average of 0.222 among these loci for modern high-performance and 0.307 for low-performance divers. Thus, for our very limited sample, historical MAF is relatively low, as it is for high-performance modern seals, which suggests that it could be due to purifying selection, although drift remains a possibility. More data would be required to make a strong inference about this.Fig. 3 Fitness impact associated with deep diving.

Minor allele frequencies at six SNP sites across five loci (x axis) associated with hypoxia (top) (see Table 2 for genotype frequencies) are compared with minor allele frequencies across six SNPs in five MHC loci (bottom) (Supplementary Table 3) for high and low dive performance.

Source data

Table 2 Genotype data for hypoxia loci in seals showing high and low dive performance

Locus	SNP-ref	Genotype	High P	Low P	Position	MAF-1980s	MAF- ~2010	
HIF1A	96892683	CC	16	11	Upstream variant	0.44	0.438	
CA	23	21	
AA	7	14	
96894173	CC	12	9	Upstream variant	0.45	0.5	
CT	22	20	
TT	12	17	
HIF3A	45676640	AA	38	37	Non-synonymous	0.147	0.052	
AG	7	9	
GG	1	0	
SETX	30169236	CC	42	37	Non-synonymous	0.118	0.052	
CT	4	9	
TT	0	0	
HYOU1	30169060	GG	28	21	Upstream variant	0.275	0.254	
GT	17	21	
TT	1	4	
MB	30169078	CC	37	28	Upstream variant	0.225	0.149	
GC	7	17	
GG	2	1	
Within each locus, the SNP reference and position relative to the locus is shown. Also shown are changes in MAF comparing females born in the 1980s (N = 51) with those born between 2005 and 2015 (N = 57).

As a control to test the hypoxia genes results, we used loci of known function (immune system) but for which we expected no correlation with diving performance (upstream SNPs at major histocompatibility complex (MHC) Class 1 B-alpha, CD74, DMA and DPA1, and a non-synonymous change in DQB1). Five MHC loci were found to have either a variable non-synonymous mutation, or a variable mutation in an upstream region. As expected, no significant correlation was found (Χ2 = 1.35, P = 0.24; Fig. 3). Stochastic distortion of allele frequencies during a bottleneck can increase the frequency of less-fit minor alleles (although none of these showed LOF). These alleles could be retained in the post-bottleneck population if dominance or co-dominance is protective.

Conclusions

Our results show that the fitness of post-bottleneck northern elephant seals is impacted by stochastic effects and reduced diversity, even though recovery has been rapid, rebounding to a population size comparable to historical maxima. Population bottlenecks are known to distort allele frequency distributions, and distortions have been used to detect bottlenecks32. If alleles with stochastically increased frequency are deleterious, they can be maintained in small populations where purifying selection is relatively weak. This can lead to ‘mutational meltdown’ and extinction in asexual species (Müller’s ratchet33), but potentially in small populations of sexual species as well34. Genomic studies of endangered species have shown the accumulation of LOF loci at levels comparable to what we have found in the post-bottleneck northern elephant seal. For example, in the pygmy hog (Porcula salvania), of which there are only a few hundred left, there are substantially more frameshift, stop-gained and missense mutations than in related species35.

For the northern elephant seal, we found three categories of post-bottleneck impact. There was a reduction in diversity and an increase in ROH comparing pre- and post-bottleneck genomes, and the loss of diversity was correlated with lower female lifetime reproductive success. Having details on lifetime reproductive success in long-lived animals (as we have here) is rare, but examples of correlations between diversity and fitness proxies are common36. The frequency of LOF alleles was also negatively correlated with female lifetime reproductive success. At specific loci associated with reproductive health, LOF allele frequency was correlated with male reproductive success. However, the lack of a genome-wide heterozygosity–fitness correlation and instead a difference associated with MAF at relevant loci for dive performance was especially striking. There was variation within (non-synonymous change) or upstream of five loci associated with hypoxia, not identified as generating LOF. However, those individuals with higher MAF showed lower dive performance, suggesting an impact from post-bottleneck allele frequency distortion. Historical genomes showed lower average MAF than either high- or low-performance divers, and there was an indication of selection against higher MAF over time. We propose that together, these impacts leave the species vulnerable to environmental stresses (such as climate-induced resource bottlenecks)37 that an uncompromised population may be able to overcome. In conclusion, our data show that despite rapid recovery and apparent stability, the northern elephant seal has reduced fitness impacting their reproductive output and ability to forage efficiently. Important aspects of impact included the stochastic distortion of allele frequencies and the retention of LOF alleles at critical loci.

Methods

Field observation and sampling

Field work was conducted at Año Nuevo State Reserve in California, USA. Details of elephant seal harem observation and tissue sample collection are given in ref. 3. Harem observations were for 6–7 hours per day at all harems. A harem was defined as a group of females associated with a single alpha male. The alpha was the highest-ranking male in the dominance hierarchy associated with the harem. Copulations provide a useful measure of reproductive success38. Observations at night with a photomultiplier video camera revealed the same rate and pattern of copulations by individuals as during the day38. Ronguers were used to take tissue samples from the outer edge of the hind flippers. The samples were stored in 20% dimethylsulfoxide saturated with NaCl39. Samples for estimating reproductive success analyses were collected in 1990 and 1991 and those individuals were tagged and tracked over time (between 1981 and 2005), permitting lifetime reproductive success estimates for a subset of females. Details are provided in ref. 40. Paternity testing was possible for females attending a given harem that were sampled the following year (a return rate of 40–60% of tagged females).

Data on diving performance were collected between 2004 and 2018. Details on the tags and data collected are provided in ref. 41. Tissue samples were again collected from hind flippers and stored in salt/dimethylsulfoxide. All subjects were marked and tagged and recognizable as individuals42. Satellite platform transmitter terminals (model ST-6, Telonics) were affixed on the head with marine epoxy, with the antenna angled so that it would be exposed when the seal was at the surface42. The transmitter interacted with the ARGOS satellite system to generate locations, which were filtered using standard methods43. Data on position, date and time permitted a record of the distance, track and duration of time spent on the foraging trip. A separate tag was affixed to record the depth profiles (a time-temperature-depth-recorder; the MK7 by Wildlife Computers). Data recorded included the age of the seal, the departure and return dates, the mean depth, the maximum depth, the duration and the mass on departure and return. A single foraging trip was recorded for 8 seals, two trips for 68 seals and three trips for 16 seals. We calculated correlations between genomic diversity and a composite measure of diving behaviour. The composite measure was termed dive ‘performance’ and was generated as the average maximum depth (among all trips recorded) times the proportion of dives deeper than the median depth (516 m) times the relative mean duration (compared to the duration of all other dives from the set of 92 individuals). This metric was devised for this study to account for the various aspects of dive endurance recorded.

Genome sequencing and SNP detection

A total of 260 modern samples were subjected to sequencing. Total genomic DNA was initially extracted using standard phenol/ chloroform extraction and subsequently stored in TE buffer. Short-read libraries for whole-genome sequencing were then constructed as follows. The extracted DNA was sheared into 50–800 bp fragments. Fragments ranging from 150 bp to 500 bp were selected and treated with T4 DNA polymerase (Enzymatics, P7080L) to obtain blunt ends. T-tailed adapters were then ligated to repair these blunt ends. PCR amplification was performed and AMPure XP beads (Agencourt, A63881) were used to purify the PCR products. These short-read libraries were then sequenced on DNBSEQ-T1 sequencers at BGI-Shenzhen, to generate paired-end 100 bp reads. The reference for assembling re-sequencing reads was Mirounga angustirostris from DNA Zoo44,45 (HI-C, chromosome level from https://www.dnazoo.org/assemblies/Mirounga_angustirostris, total length: 2,366,206,800 bp; scaffold N50: 139,676,048 bp; scaffold N90: 54,920,518 bp, GC content: 41.52%; Supplementary Table 4). The average number of cleaned reads was 1,026.11 million (range: 651.91–7,516.79 million). Percent Q30 was 94.41 on average (range 89.76–97.17; Supplementary Table 5) and average genome coverage was 99.59% (Supplementary Table 6). The average sequencing depth was 33.8X (see Supplementary Fig. 1 for full range). SOAPnuke (v.2.1.5)46 was employed for quality control of sequencing reads with parameter ‘-J -l 10 -q 0.1 -n 0.05’ and a config file containing the settings ‘trimBadHead=13,100 trimBadTail=13,100 trim=0,0,0,0 qualSys=2 seqType=0 outQualSys=2’. The clean reads were mapped to the reference genome using BWA47 and Realigner in Sentieon48, while SNP detection utilized GVCFtyper in Sentieon with default parameters to generate joint-calling of raw SNPs. Sentieon was also used for duplicate read identification and removal (see Supplementary Table 6). The hard filtering of SNPs was performed using VariantFiltration in GATK (v.3.8.1)49, with the criterion ‘QD < 2.0 || MQ < 40.0 || FS > 60.0 || ReadPosRankSum < −8.0 || MQRankSum < −12.5 || SOR > 3.0’ as recommended by GATK. Three replicate samples were removed from further analysis.

Historical samples were acquired on loan from the Smithsonian Museum of Natural History, the Harvard Museum of Comparative Zoology and the San Diego Natural History Museum. Samples were bone, tooth or dried dermal tissue. DNA extraction and library construction were done at the Smithsonian Institution’s Center for Conservation Genomics. Work was done in their contained ancient DNA facility using standard precautions to avoid cross-contamination. Extraction followed the protocol in ref. 50. Briefly, samples were powdered using a Dremel drill and the calcium extracted from bone and tooth with EDTA (0.5 M, pH 8.0). Samples were then extracted for DNA in extraction buffer (after ref. 50 Supplement SD3) overnight at 37 °C with shaking. Solutions were cleaned on Qiagen MinElute columns (28004). DNA quantity was checked on a Qubit 4 fluorometer. DNA libraries were constructed using a KAPA Hyper Prep kit (Roche) following manufacturer protocols. Samples were dual indexed with different pairs of P5 and P7 adapters. Libraries were quantified on Qubit and Tapestation (Agilent) and the results compared.

Historical samples were first sequenced on a MiSeq system (Illumina) to balance quantities on the basis of the number of reads. Reads were mapped against the DNA Zoo NES reference (see above) with Bowtie2 (v.2.3.4.1)51 using the ‘–very-sensitive’ setting. After balancing quantities, libraries were pooled for sequencing in two Novaseq 6000 (Illumina) S4 lanes (2 × 100 bp). All historical sample sequencing was done at the Baur Core facility at Harvard University. Total reads per sample are shown in Supplementary Table 1. Sequences were trimmed with BBDuk (sourceforge.net/projects/bbmap/) using these settings: ‘k = 25 mink = 15 edist = 1 ktrim = r rcomp = f t = 8 qtrim=rl trimq = 20 maq = 20’. Quality was on average Q36 and the peak fragment size was in the 50–60 bp range, calculated using FastQC (v.0.11.9)52. Using clumpify dedupe in BBtools53, duplicate reads were removed. Trimmed and filtered reads were mapped to the NES DNA Zoo reference genome using Bowtie2 with the ‘–very-sensitive’ setting. Mapped reads with a quality score below 10 were removed using SAMtools (v.1.12)54,55 ‘view’, and the filtered reads were sorted using SAMtools ‘sort’. Variant sites were called using the ‘mpileup’ and ‘call’ commands in BCFtools (v.1.12)55. Using SnpSift56 ‘filter’, variant calls with a quality score of 20 and depth below 5 were removed. VCF files were further filtered using BCFtools ‘view’ to include only the 17 identified chromosomal scaffolds.

Mitogenomes were generated for 180 modern adults and all 8 historical samples. The reads with low quality, duplications or adaptors were removed using SOAPnuke (v.2.1.5) (https://github.com/BGI-flexlab/SOAPnuke)46, leaving clean reads for the final mitogenome assembly. To normalize the samples, randomly resampling the sequences to 40,000,000 reads from the clean reads of each sample was performed using Seqtk (https://github.com/lh3/seqtk). The mitogenome assembly of the NES was carried out using NOVOPlasty v.3.7, which is a de novo seed-extend-based assembler for organelle genomes57. The mitochondrial genomes of the NES (Mirounga angustirostris, CM055130.1) and the SES (Mirounga leonine, NC_008422.1) were used as seed sequences for the mitochondria assembly and the assembly parameters were set as follows: PE mode, read length = 100, k-mer = 39, genome range = 12,000–22,000 and type = mito. The mitogenome tree was based on all historical (8) and all modern adult (180) samples and generated in PAUP (https://paup.phylosolutions.com/) by the neighbour-joining method using the Tamura Nei evolution rate correction and 1,000 bootstrap replicates. It was routed with the outgroup, Mirounga leonina.

Genome analysis

A sliding-window analysis using VCFtools54 ‘–window-pi’ measured heterozygous sites every 50,000 bp. The total numbers of heterozygous sites across the genome were measured using the command grep on VCF files and heterozygosity determined by dividing by the total number of sites passing the minimum quality and depth filters for each sample. ROHs were measured for runs greater than 100 kb and 1 Mb long using (for 100 kb) the ‘–bed infile.bed–fam infile.fam–bim infile.bim–homozyg–noweb–allow-extra-chr–homozyg-kb 100–out data’ command in PLINK58. Before assessing ROH in modern samples, the VCF file was filtered to remove singletons and minor allele frequencies lower than 0.01.

Demographic history was estimated using the coalescent method implemented in PSMC59. Aligned mapped reads (BAM files) for 14 samples (7 males and 7 females chosen randomly within sex and among samples from the 1980s) were converted to consensus sequence in FASTQ format using the samtools/bcftools pipeline. First, the samtools ‘mpileup’ command was used to produce the VCF file from the BAM files, and then through bcftools the consensus sequence was generated with the original consensus calling model. Following that, the vcfutils.pl script was used for the FASTQ conversion: bcftools mpileup -Q 30 -q 30 -f NES.fa sample.bam | bcftools call -c | vcfutils.pl vcf2fq -d 10 -D 100 | gzip > diploid_sample.fq.gz. Mapped reads were filtered for a minimum mapping quality (q) of 30 and a minimum base quality (Q) of 30. The minimum (d) and maximum (D) coverages were calculated to allow for vcf2fq, and were set to 10 and 100, respectively (-d value to a third of the average depth and -D value to twice). The FASTQ file was then converted to the input format for PSMC using the following: fq2psmcfa sample.fq.gz > sample.psmcfa. For the final PSMC command, we used 64 atomic time intervals and 28 (=1 + 25 + 1 + 1) free interval parameters: psmc -p ‘4 + 25*2 + 4 + 6’ -o sample.psmc sample.psmcfa. Finally, the PSMC plot was drawn using the command ‘psmc_plot.pl’, with the per-generation mutation rate ‘-u’ and the generation time in years ‘-g’ set to 2.2 × 10−09 and 2.2 × 10−06, respectively. One hundred bootstrap replicates were run for each sample and 15/100 samples chosen at random to be included in the plot.

Current effective population size was estimated using the linkage disequilibrium method implemented in SNeP13 and GONE14. SNeP estimates the historical effective population size on the basis of the relationship between r2, Ne and c (recombination rate) via linkage disequilibrium estimates using the standard PLINK input file format (.ped and .map files). The squared Pearson’s product-moment correlation coefficient between pairs of loci was used due to the unknown phase of the genotypes. The software uses the physical distance (δ) between two loci as a reference and translates it into linkage distance (d). We used the default values for minimum distance in bp between SNPs to be analysed ‘-mindist’ and maximum distance in bp between SNPs to be analysed ‘-maxdist’ to 50,000 and 4,000,000, respectively. To infer the recombination rate, we used the ‘-svedf’ flag58 as a recombination rate modifier and the default MAF < 0.05, as it has been shown that accounting for MAF results in unbiased r2 estimates irrespective of sample size60,61.

The software GONE14 was used to obtain an additional estimate of Ne based on the linkage disequilibrium method. The VCF file was converted to .ped and .map format using the PLINK software. The parameters for the analysis were set as follows: number of generations 2,000, number of bins 400, no MAF pruning was applied and the maximum value of recombination rate (c) was set to 0.05. The number of internal replicates was set to 40 for the programme to provide the geometric mean of the consensus estimate of the historical Ne out of these replicates.

Correlations between pups weaned per year (a measure of lifetime reproductive success) or diving performance and FROH were done using linear regression. A multiple regression was run using the software at https://stats.blue/Stats_Suite/multiple_linear_ regression_calculator.html. FROH was calculated as the proportion of the genome represented by ROH at the given filter length. LOF loci were identified using SnpEff and the Mirounga angustirostris annotation GFF file. There were no good alternatives to using the DNA Zoo NES reference, but it is a high-quality and annotated genome, and the results showed credible patterns. Analyses were split between using only those loci that generated a new stop or caused loss of a start codon and all putative LOF loci discovered (including those of minor effect). For graphing, we coded the homozygote of the LOF allele as ‘1’, the heterozygote as ‘0.5’ and unaffected genotypes as ‘0’, and totalled across all identified loci. SnpEff was also used to identify variation at specific loci of interest, such as those associated with reproduction and hypoxia. We used the following pipeline: java -jar snpEff.jar build -gff3 -v NES.reference based on the annotation file GFF. Following that, we used the command ‘ann’ to annotate the variant file. We used the snpSift toolbox as implemented in snpEff to filter the annotated VCF file on the basis of the effect of the genetic variants, for example, LOF or non-synonymous substitutions (missenses). Although we focused on snpEff, which provides specificity on the nature of the change, we also used the Protein Variation Effect Analyzer (PROVEAN; https://github.com/MonashBioinformaticsPlatform/provean) as an alternative method to test the results. PROVEAN calculates the alignment score of a protein variant by comparing a query sequence with the existing database. To run PROVEAN, we configured it with BLAST v.2.2.31+, CD-HIT v.3.1.2 and the NCBI nr database containing 15,006,980 sequences from 15 August 2011. The tool was run with the default parameters, whereby missense variants with PROVEAN scores greater than or equal to −2.5 were classified as neutral, while those with scores less than −2.5 were considered deleterious on the basis of the default threshold.

To assess kinship, we used the relatedness2 function54 implemented in VCFtools. We included 77 pups born in the following year to females, most of whom had been in the focal harems, and determined the number of paternities achieved by each of the 31 males in the sample (including the 4 focal alpha males). Manhattan plot comparisons of high- and low-diving-performance females were generated by writing a trait file and generating the plot in R and TASSEL62. There were no clear outliers (data not shown). Data relevant to the LRS and diving analyses are provided in Supplementary Tables 7 and 8.

Reporting summary

Further information on research design is available in the Nature Portfolio Reporting Summary linked to this article.

Supplementary information

Reporting Summary

Supplementary Tables 1–8 Supplementary Tables 1–8.

Source data

Source Data Figs. 1–3 and Extended Data Figs. 1–9 Source data for all figures in one Excel file.

Extended data

Extended Data Fig. 1 Histogram showing the distribution of read depths for the 260 modern samples.

Source data

Extended Data Fig. 2 ROH number of segments compared to total bp covered.

a) fragments > 1MB, and b) fragments greater than 100KB.

Source data

Extended Data Fig. 3 Comparison of a historical (Hist8) and a modern (sample 1234) nuclear genomes.

a) The number of segments of ROH of different sizes greater than 100Kb. b) The number of heterozygous sites in 50,000 bp sliding windows. For each illustration, the box represented the inter-quartile range (50% of the data), the horizotal line through the box is the median, the whiskers include 75% of the data, and the remaining dots are outside that range.

Source data

Extended Data Fig. 4 Historical demographic profile estimated using the program GONE.

Source data

Extended Data Fig. 5 Fitness correlations with total LOF and all weaned.

a) Correlation between all LOF loci reported from snpEff (weak and strong effect) against weaner success per year. b) Correlation between FROH and the total pups weaned over a lifetime.

Source data

Extended Data Fig. 6 Correlations between LOF and FROH.

Linear regression correlations between FROH ( > 1 MB) and a) loss of function based on the analysis in snpEff, and b) missense mutations estimated in PROVEAN.

Source data

Extended Data Fig. 7 Relative proportion of runs of homozygosity greater than 1 MB in different range categories (in brackets) for all adult individuals, showing level found for male M12 (arrow).

Source data

Extended Data Fig. 8 Relationship between the number of paternities achieved and the frequencuy of loss of function alleles at loci associated with male health in non-alpha males.

Source data

Extended Data Fig. 9 Relationship between the proportion of runs of homozygosity greater than 1 MB and the diving performance metric (no significant correlation).

Source data

Extended Data Table 1 Paternal success of all males

Paternal success of all males

Extended data

is available for this paper at 10.1038/s41559-024-02337-4.

Supplementary information

The online version contains supplementary material available at 10.1038/s41559-024-02337-4.

Acknowledgements

We thank S. Edwards for hosting A.R.H. on sabbatical at OEB, Harvard University; R. Fleischer for hosting A.R.H. at the National Zoo, Smithsonian Institution; P. Byerly for great help at the National Zoo lab; A. Gusick and S. Robson at the Natural History Museum, Los Angeles County, M. McGowen at the Smithsonian Institution, and M. Omura and H. Hoekstra at the Museum of Comparative Zoology, Harvard University for all their help in obtaining historical samples; P. Unit for help from the San Diego Natural History Museum; S. Gaughran for helpful comments on the manuscript; C. Hartmann and C. Daly at the Baur Core, Harvard University for helpful bioinformatic advice; Z. Yang, Institute of Deep-sea Science and Engineering, Chinese Academy of Sciences for help in delivering and curating data; R. Condit for help and expertise with the maternity data; and the DNA Zoo team for establishing the Mirounga angustirostris reference genome. Approval for elephant seal handling and instrumentation for performance metrics was provided by the University of California Santa Cruz Institutional Animal Care and Use Committee and under National Marine Fisheries Service permits #786-1463, 87-143, 14636, 14535 and 19108. This work made use of the Hamilton HPC Service of Durham University. We acknowledge funding from the National Natural Science Foundation of China (Grant numbers 42225604 and 41422604) to S.L.; ‘One Belt and One Road’ Science and Technology Cooperation Special Program of the International Partnership Program of the Chinese Academy of Sciences (183446KYSB20200016) to S.L.; and the Specially-Appointed Professor Program of Jiangsu Province to I.S. Field work was supported by the National Science Foundation to BJL from 1970 to 2000, the National Geographic Society, and by the Office of Naval Research (N00014–18-1-2822 and N000014–13-1-0134 to D.P.C. and D. Crocker); the Strategic Environmental Research Development Program (RC20-2–1284 to D.P.C. and D. Crocker); and the Tagging of Pacific Predators Program including support from the Gordon and Betty Moore Foundation, the David and Lucile Packard Foundation, Alfred P Sloan Foundation and the Animal Telemetry Network.

Author contributions

A.R.H. conceptualized the project, conducted field work, lab work, sample sequencing and analysis, and wrote the original draft. S.L. and G.F. conceptualized the project, acquired funding and revised the manuscript. S.S. and H.K. performed sample sequencing and analysis. I.S. revised the manuscript. G.A.G. and F.S. conducted analysis and revised the manuscript. B.L.B., D.P.C. and J.R. conceptualized the project, and performed field work, analysis and manuscript revision. R.S.B. and P.W.R. collected and analysed data. N.M. performed lab work and sample sequencing.

Peer review

Peer review information

Nature Ecology & Evolution thanks the anonymous reviewers for their contribution to the peer review of this work.

Data availability

Sequences are deposited at the National Center for Biotechnology Information (NCBI) under Bioproject PRJNA1060307. Whole-genome data were uploaded as raw reads (Sequence Read Archive), including all historical samples and six high-coverage modern samples (Biosamples SAMN39224291–SAMN39224298; SAMN39305437, SAMN39305439, SAMN39305440, SAMN39305442, SAMN39305447 and SAMN39305448). Variant information for all modern sequences is included as a VCF file (European Variation Archive; PRJNA1060307). All raw nuclear genome data for the modern samples have also been deposited to the China National GeneBank Sequence Archive (CNSA) of the China National GeneBank DataBase (CNGBdb) with accession number CNP0005170. Source data are provided with this paper.

Competing interests

The authors declare no competing interests.

Publisher’s note Springer Nature remains neutral with regard to jurisdictional claims in published maps and institutional affiliations.

These authors contributed equally: Georgios A. Gkafas, Hui Kang, Fatih Sarigol.
==== Refs
References

1. Collen B Predicting how populations decline to extinction Phil. Trans. R. Soc. B. 2011 366 2577 2586 10.1098/rstb.2011.0015 21807738
2. Speilman D Brook BW Frankham R Most species are not driven to extinction before genetic factors impact them Proc. Natl Acad. Sci. USA 2004 101 15261 15264 10.1073/pnas.0403809101 15477597
3. Hoelzel AR Elephant seal genetic variation and the use of simulation models to investigate historical population bottlenecks J. Hered. 1993 84 443 449 10.1093/oxfordjournals.jhered.a111370 7505788
4. Abadia-Cardoso A Freimer NB Deiner K Garza JC Molecular population genetics of the northern elephant seal Mirounga angustirostris J. Hered. 2017 108 618 627 10.1093/jhered/esx053 28821186
5. Lowry MS Abundance, distribution and population growth of the northern elephant seal (Mirounga angustirostris) in the United States from 1991 to 2010 Aquat. Mamm. 2014 40 20 31 10.1578/AM.40.1.2014.20
6. Hoelzel AR Impact of a population bottleneck on genetic variation and the importance of life history; a case study of the northern elephant seal Biol. J. Linn. Soc. 1999 68 23 39 10.1111/j.1095-8312.1999.tb01156.x
7. Bonnell ML Selander RK Elephant seals: genetic variation and near extinction Science 1974 184 908 909 10.1126/science.184.4139.908 4825892
8. Le Boeuf BJ Male–male competition and reproductive success in elephant seals Am. Zool. 1974 14 163–176 10.1093/icb/14.1.163
9. Hoelzel AR Fleischer RC Campagna C Le Boeuf BJ Direct evidence for the impact of a population bottleneck on symmetry and genetic diversity in the northern elephant seal J. Evol. Biol. 2002 15 567 575 10.1046/j.1420-9101.2002.00419.x
10. Weber DS Stewart BS Garza JC Lehman N An empirical genetic assessment of the severity of the northern elephant seal population bottleneck Curr. Biol. 2000 10 1287 1290 10.1016/S0960-9822(00)00759-4 11069110
11. Keller LF Waller DM Inbreeding effects in wild populations Trends Ecol. Evol. 2002 17 230 241 10.1016/S0169-5347(02)02489-8
12. Frankham R Effective population size/adult population size ratios in wildlife: a review Genet. Res. 1995 66 95 107 10.1017/S0016672300034455
13. Barbato M Orozco-ter Wengel P Tapio M Bruford MW SNeP: a tool to estimate trends in recent effective population size trajectories using genome-wide SNP data Front. Genet. 2015 6 109 10.3389/fgene.2015.00109 25852748
14. Santiago E Recent demographic history inferred by high-resolution analysis of linkage disequilibrium Mol. Biol. Evol. 2020 37 3642 3653 10.1093/molbev/msaa169 32642779
15. Szulkin M Bierne N David P Heterozygosity–fitness correlations: a time for reappraisal Evolution 2010 64 1202 1217 20148954
16. Ralls K Ballou JD Templeton A Estimates of lethal equivalents and the cost of inbreeding in mammals Conserv. Biol. 1988 2 185 193 10.1111/j.1523-1739.1988.tb00169.x
17. Su T-Q MARF1 regulates essential oogenic processes in mice Science 2012 335 1496 1499 10.1126/science.1214680 22442484
18. Bilmez Y Talibova G Ozturk S Dynamic changes of histone methylation in mammalian oocytes and early embryos Histochem. Cell Biol. 2021 157 7 25 10.1007/s00418-021-02036-2 34599660
19. Anastasaki C Longman D Capper A Patton EE Ca´ceres JF Dhx34 and Nbas function in the NMD pathway and are required for embryonic development in zebrafish Nucleic Acids Res. 2011 39 3686 3694 10.1093/nar/gkq1319 21227923
20. Hansson B Westerberg L On the correlation between heterozygosity and fitness in natural populations Mol. Ecol. 2008 11 2467 2474 10.1046/j.1365-294X.2002.01644.x
21. Hoelzel AR Le Boeuf BJ Reiter J Campagna C Alpha male paternity in elephant seals Behav. Ecol. Sociobiol. 1999 46 298 306 10.1007/s002650050623
22. Liu Y LRGUK-1 is required for basal body and manchette function during spermatogenesis and male fertility PLoS Genet. 2015 11 e1005090 10.1371/journal.pgen.1005090 25781171
23. Leslie JS MNS1 variant associated with situs inversus and male infertility Eur. J. Hum. Genet. 2020 28 50 55 10.1038/s41431-019-0489-z 31534215
24. Feng M Tubulin TUBB4B is involved in spermatogonia proliferation and cell cycle processes Genes 2022 13 1082 10.3390/genes13061082 35741845
25. Feng S hnRNPH1 recruits PTBP2 and SRSF3 to modulate alternative splicing in germ cells Nat. Commun. 2022 13 3588 10.1038/s41467-022-31364-7 35739118
26. Salvolini E Involvement of sperm plasma membrane and cytoskeletal proteins in human male infertility Fertil. Steril. 2013 99 697 704 10.1016/j.fertnstert.2012.10.042 23174138
27. Le Boeuf, B. J. Elephant Seals: Pushing the Limits on Land and at Sea (Cambridge Univ. Press, 2021).
28. Ebersole JL Hypoxia-inducible transcription factors, HIF1A and HIF2A, increase in aging mucosal tissues Immunology 2018 154 452 464 10.1111/imm.12894 29338076
29. Ramachandran S Hypoxia-induced SETX links replication stress with the unfolded protein response Nat. Commun. 2021 12 3686 10.1038/s41467-021-24066-z 34140498
30. De Miranda MA Schlater AE Green TL Kanatous SB In the face of hypoxia: myoglobin increases in response to hypoxic conditions and lipid supplementation in cultured Weddell seal skeletal muscle cells J. Exp. Biol. 2012 215 806 813 10.1242/jeb.060681 22323203
31. Rao S Biological function of HYOU1 in tumors and other diseases OncoTargets Ther. 2021 14 1727 1735 10.2147/OTT.S297332
32. Luikart G Allendorf FW Cornuet J-M Sherwin WB Distortion of allele frequency distributions provides a test for recent population bottlenecks J. Hered. 1998 89 238 247 10.1093/jhered/89.3.238 9656466
33. Muller HJ The relation of recombination to mutational advance Mutat. Res. 1964 1 2 9 10.1016/0027-5107(64)90047-8
34. Lynch M Conery J Burger R Mutational meltdowns in sexual populations Evolution 1995 46 1067 1080 10.2307/2410432
35. Liu L Genetic consequences of long-term small effective population size in the critically endangered pygmy hog Evol. Appl. 2020 14 710 720 10.1111/eva.13150 33767746
36. Hoffman JI High-throughput sequencing reveals inbreeding depression in a natural population Proc. Natl Acad. Sci. USA 2014 111 3775 3780 10.1073/pnas.1318945111 24586051
37. Maron M Climate-induced resource bottlenecks exacerbate species vulnerability: a review Divers. Distrib. 2015 21 731 743 10.1111/ddi.12339
38. Le Boeuf BJ Sexual behavior in the northern elephant seal, Mirounga angustirostris Behaviour 1972 41 1 26 10.1163/156853972X00167 5062032
39. Hoelzel AR Dover GA Molecular techniques for examining genetic variation and stock identity in cetacean species IWC Spec. Issue 1989 11 81 120
40. Le Boeuf BJ Condit R Reiter J Lifetime reproductive success of northern elephant seals (Mirounga angustirostris) Can. J. Zool. 2019 97 1203 1217 10.1139/cjz-2019-0104
41. Boedlert GW Autonomous pinniped environmental samplers: using instrumented animals as oceanographic data collectors J. Atmos. Ocean. Technol. 2001 18 1882 1893 10.1175/1520-0426(2001)018<1882:APESUI>2.0.CO;2
42. Le Boeuf BJ Foraging ecology of northern elephant seals Ecol. Monogr. 2000 70 353 382 10.1890/0012-9615(2000)070[0353:FEONES]2.0.CO;2
43. Robinson PW Foraging behavior and success of a mesopelagic predator in the northeast Pacific Ocean: insights from a data-rich species, the northern elephant seal PLoS ONE 2012 7 e36728 10.1371/journal.pone.0036728 22615801
44. Dudchenko O De novo assembly of the Aedes aegypti genome using Hi-C yields chromosome-length scaffolds Science 2017 356 92 95 10.1126/science.aal3327 28336562
45. Dudchenko, O. et al. The Juicebox Assembly Tools module facilitates de novo assembly of mammalian genomes with chromosome-length scaffolds for under $1000. Preprint at bioRxiv10.1101/254797 (2018).
46. Chen Y SOAPnuke: a MapReduce acceleration-supported software for integrated quality control and preprocessing of high-throughput sequencing data GigaScience 2018 7 1 6 10.1093/gigascience/gix120 29659813
47. Li H Durban R Fast and accurate short read alignment with Burrows–Wheeler transform Bioinformatics 2009 25 1754 1760 10.1093/bioinformatics/btp324 19451168
48. Freed, D., Aldana, R., Weber, J. A. & Edwards, J. S. The Sentieon Genomics Tools – a fast and accurate solution to variant calling from next-generation sequencing data. Preprint at bioRxiv10.1101/115717 (2017).
49. McKenna A The genome analysis toolkit: MapReduce framework for analyzing next-generation DNA sequencing data Genome Res. 2010 20 1297 1303 10.1101/gr.107524.110 20644199
50. McDonough MM Parker LD McInerney NR Campagna MG Maldonando JE Performance of commonly requested destructive museum samples for mammalian genomics studies J. Mammal. 2018 99 789 802 10.1093/jmammal/gyy080
51. Langmean B Saltzberg SL Fast gapped-read alignment with Bowtie 2 Nat. Methods 2012 9 357 359 10.1038/nmeth.1923 22388286
52. Andrews, S. FastQC: A Quality Control Tool for High Throughput Sequence Data (Babraham Bioinformatics, 2010).
53. Bushnell B Rood J Singer E BBMerge – accurate paired shotgun read merging via overlap PLoS ONE 2017 12 e0185056 10.1371/journal.pone.0185056 29073143
54. Li H A statistical framework for SNP calling, mutation discovery, association mapping and population genetical parameter estimation from sequencing data Bioinformatics 2011 27 2987 2993 10.1093/bioinformatics/btr509 21903627
55. Danacek P Twelve years of SAMtools and BCFtools GigaScience 2021 10 1 4
56. Cingolani P Using Drosophila melanogaster as a model for genotoxic chemical mutational studies with a new program, SnpSift Front. Genet. 2012 3 35 10.3389/fgene.2012.00035 22435069
57. Dierckxsens, N., Mardulyn, P. & Smits, G. NOVOPlasty: de novo assembly of organelle genomes from whole genome data. Nucleic Acids Res. 45, e18 (2017).
58. Purcell S PLINK: a tool set for whole-genome association and population-based linkage analyses A. J. Hum. Genet. 2007 81 559 575 10.1086/519795
59. Li H Durbin R Inference of human population history from individual whole-genome sequences Nature 2011 475 493 496 10.1038/nature10231 21753753
60. Sved JA Feldman MW Correlation and probability methods for one and two loci Theor. Popul. Biol. 1973 4 129 132 10.1016/0040-5809(73)90008-7 4726005
61. Manichaikul A Robust relationship inference in genome-wide association studies Bioinformatics 2010 26 2867 2873 10.1093/bioinformatics/btq559 20926424
62. Bradbury PJ TASSEL: software for association mapping of complex traits in diverse samples Bioinformatics 2007 23 2633 2635 10.1093/bioinformatics/btm308 17586829
