Introduction
Although amphibian enigmatic declines had been identified by herpetologists as early as the 1970s, they were only recognized two decades later as a global phenomenon that in some cases could not be explained by environmental changes or other expected anthropogenic factors alone (1–3). The discovery of the amphibian-killing fungus Batrachochytrium dendrobatidis [hereafter Bd (4, 5)] was a turning point in understanding why many amphibian species have been in steep decline. The emergence of Bd, which causes the lethal disease, amphibian chytridiomycosis, has been associated with amphibian population declines of more than 500 species, including the presumed extinction of at least 90 species (6). Evidence suggests that Bd recently spread across the globe from an endemic focus, with East Asia as the most likely source from where it expanded to other continents during the past century (7–10). Its global spread has been mainly facilitated by the international trade of amphibians, particularly the North American bullfrog (Lithobates catesbeianus), the most intensively farmed frog worldwide (11–15).
Novel genomic techniques, including whole-genome sequencing and multilocus sequence typing have shown the existence of at least five major lineages of Bd: BdGPL, BdCAPE, BdASIA-1 (including BdCH), BdASIA-2/BdBRAZIL and BdASIA-3 (9, 10, 16). Of these, the global panzootic lineage (BdGPL) is the most widespread variant of Bd and is responsible for most known cases of amphibian population declines due to chytridiomycosis (9). Although BdGPL is highly virulent, its impacts are context-dependent (17), and under some conditions other lineages may be responsible for lethal disease and population declines (18). Additionally, multiple introductions of Bd have led to different lineages of Bd coming into contact, resulting in the formation of interlineage recombinants (e.g., through the co-infection of amphibians) that may have higher pathogenicity or transmissibility (9, 10, 14, 16, 19, 20). For instance, interlineage recombinants have been reported for BdGPL with BdASIA-2/BdBRAZIL, and for BdGPL with BdCAPE (9, 10, 14, 19).
Despite South America being the region with the greatest loss of biodiversity due to Bd (6), only low numbers of Bd isolates from this region have been genetically characterized (14, 21–26). This limits our capacity to adequately understand the epidemiological processes that have led to impacts of chytridiomycosis on South American native amphibians. Our focus here, the Chilean giant frog (Calyptocephalella gayi), is endemic to Chile and is considered as a living fossil since its family represents an old neobatrachian clade that diverged during the Cretaceous around 100–120 Mya (27). With females reaching up to 2 kg, this is the second largest anuran species worldwide, which has led this species to be of economic interest as a food source (28). This highly aquatic species is considered Vulnerable by the IUCN Red List, and is threatened by overconsumption, habitat loss due to agriculture, and invasive species including several introduced fish and the African clawed frog [Xenopus laevis (29)]. Chytridiomycosis has been suspected to be a contributing factor in its steep decline (30), but to date, no evidence has been found linking Bd to lethal effects in C. gayi. Based on histology and whole-genome sequencing, the aim of this study is to describe the lesions of chytridiomycosis in C. gayi and the genetic characterization of Bd isolates obtained from a chytridiomycosis outbreak that occurred in a C. gayi captive breeding program in Chile. In addition, the genomics of Bd isolates from captive C. gayi together with previously obtained Chilean Bd isolates are compared with a global panel of Bd.
Materials and Methods
Chilean Giant Frog Captive Breeding Center
The Chilean giant frog captive breeding center (Resolution N°2358/2013 by the Chilean Agriculture and Livestock Service) in Santiago has been functioning since 2013, with its objective to generate reproductive knowledge and to support the conservation of C. gayi. The center was built in an area of 70 m2 and was composed of 10 large and 10 small tanks for tadpoles (100 and 30 L each, respectively), 10 small tanks for recently metamorphosed frogs (30 L) and 20 medium size tanks for adult frogs (50 L). Water used in tanks was from the mains supply but had been left to stand for 2 days to allow chlorine evaporation. Around 50% of water in the tanks was changed twice a week. Tadpoles were fed on spirulina algae supplemented with lettuce given ad libitum, while postmetamorphs and adults were fed twice a day with dried amphipod crustaceans (Orchistoidea spp.) supplemented with chicken protein, vitamins, and minerals. By August 2016, the captive breeding program comprised ~400 1-year-old tadpoles, six juveniles, and 86 breeding adults (43 females and 43 males).
Mortality and Pathological Analyses
In September 2016, the program received new individuals from a separate C. gayi captive breeding program that had been terminated. The newly incoming individuals consisted of 800 2-year-old tadpoles and 18 breeding adults (9 females, 9 males), which were maintained in separate tanks from the resident animals. By early November 2016, 40 of the new tadpoles completed metamorphosis. From December 2016 to January 2017, 75 of the new individuals died: 37 tadpoles, 35 postmetamorphs and three reproductive adults. All dead postmetamorphs had not consumed any food after metamorphosis, and some of the tadpoles were observed to be lethargic prior to death and exhibiting partial depigmentation of the mouthparts, a finding consistent with amphibian chytridiomycosis (31). Freshly dead animals (4 tadpoles and 29 postmetamorphs) were transported in refrigerated conditions to the laboratory for postmortem examination and Bd isolation. Necropsies of the 29 postmetamorphs were performed according to standard protocol (32). Tissue sections of postmetamorphs were collected in neutral-buffered 10% formaldehyde from any organ displaying gross lesions and from lung, liver, spleen, kidney, skeletal muscle, heart, skin, stomach, small and large intestines. For histopathological analyses, tissues were embedded in paraffin wax, sectioned (4–5 μm) and stained with hematoxylin and eosin (H&E).
<italic>Bd</italic> Sampling and qPCR Assay
Non-invasive skin swabs (MW100, Medical & Wire Equipment Co.) were obtained from postmetamorphic amphibians (n = 29) by firmly running the swabs five times each over the ventral abdomen and pelvis, each ventral hind limb (femur and tibia) and the plantar surface of each hind foot (33). Also, from dead tadpoles (n = 4), samples were obtained from the oral disc by rotating the swab 10 times around the mouth opening. Swabs were kept in a cool box until freezing at −80°C once back at the laboratory. Briefly, DNA extraction from skin and oral swabs and subsequent detection of Bd DNA using a specific real time qPCR assay was done following Soto-Azat et al. (34). For each sample, diagnostic assays were performed in duplicate, and standards of known zoospore concentration (obtained from a previous Bd culture) were included within each PCR plate as positive controls. We assumed that a Bd-positive swab was indicative of Bd infection. By including known concentrations of Bd DNA in serial diluted positive control wells on each PCR plate, we were able to quantify infection intensity, which we defined as the number of zoospore equivalents/swab (ZE). To quantify and correct the infection intensity per swab, each genomic value obtained from the qPCR assay was multiplied by 120 to account for sample dilution (35).
<italic>Bd</italic> Isolation
Freshly dead tadpoles with suspected Bd infection (n = 4) were used for Bd isolation following Longcore et al. (5) and Fisher et al. (36). Subsequent confirmation of Bd infection status and load by qPCR, served to guide Bd isolation efforts. Within 8 h. after death, the whole mouthparts of Bd-positive dead tadpoles were removed and sectioned into small pieces and deposited in a fungal growth TGhL medium (8 g. tryptone, 2 g. gelatin hydrolysate, 4 g. lactose, 10 g. agar). Cultured sections were first cleaned using an agar plate with antibiotics (200 mg/L penicillin-G and 400 mg/L streptomycin sulfate), and then placed singly into TGhL agar plate with antibiotics incubated at 15–20°C. Because zoospore release may occur immediately, especially from tadpole mouthparts, cultures were examined with an inverted microscope for the presence of active zoospores every day for up to 1 week. Once growth of zoospores and/or zoosporangia was observed, part of the agar was transferred to a new TGhL agar plate without antibiotics and incubated at 15–20°C up to 1 week. Isolates were then passaged no more than three times in order to lessen the chance of genomic change due to prolonged laboratory culture (37).
DNA Extraction, Sequence Library Preparation and Phylogenomic Analyses
We performed DNA extraction using the MasterPureTM Yeast DNA Purification Kit (Epicentre, Wisconsin, USA) from all obtained purified Bd cultures. DNA extractions were first quantified using a Tapestation 2200 (Agilent Technologies, California, USA) and Qubit 2.0 fluorimeter (Thermo Fisher Scientific, Massachusetts, USA), and then sequenced using an Illumina HiSeq 2000 (Illumina, California, USA). Subsequently, TruSeq Nano 350 gel-free sequencing libraries were prepared for 125 + 125 bp paired-end sequencing using Illumina HiSeq high output v4 chemistry, and sequencing reads cleaned of adapter sequences and quality trimmed using cutadapt v1.10 (38). Reads of the JEL423 reference genome (GenBank assembly accession: GCA_000149865.1) were mapped using Burrows-Wheeler Aligner v0.7.8 (39).
We processed the resulting sequence alignment map (SAM) files using SAMtools v1.3.1 with the “fixmate” and “sort” programs to read the files for variants discovery. We performed the variant detection in a two-step process using Freebayes version dbb6160 (40). First, sorted SAM files for each isolate in the phylogeny were independently called to find variant positions and merged into a single variant call format (VCF) file. Second, each of the samples was re-called using the positions in VCF file, to produce a squared-off call set (a genotype call was made at every locus for each isolate, including for missing data). All VCF files were processed using vcflib (41) to break complex variants into allelic primitives and vt (42) to normalize short insertion and deletion sequences. VCF files were then quality filtered using bcftools v1.3.1 (43) to accept only variants with sufficient supporting evidence. Potentially polymorphic sites were filtered using the settings in the bcbio.variation.recall squaring-off pipeline (44), with sites not passing these filters set to homozygous reference (there was not enough evidence to call a variant at that position). Then, the processed VCF files were merged into a single multisample VCF and extracted a FASTA file of the SNP variant calls. Phylogenomic analyses were conducted using RAxML v8.2.9 with GTRCAT model with 500 bootstrap runs. Weir and Cockerham's estimator was performed using a sliding-window comparison of FST of Chilean Bd isolates against a global panel of representative global diversity of Bd in vcftools. Single nucleotide polymorphisms (SNPs) that were in high linkage disequilibrium were pruned from the dataset using the SNPRelate package version 1.10.2 in R v3.4.0 (45). After pruning using a sliding-window based analysis and a linkage disequilibrium threshold of 0.125, 3,900 SNPs positions remained which were analyzed using SNPRelate and plotted with ggplot2 (46). Finally, the clustering of Chilean Bd isolates against a global panel of Bd was investigated using principal component analysis (PCA) with adegenet package (47) and plotted with ggplot2 (46).
Discussion
The now globalized Bd has caused the greatest loss of biodiversity known due to a single pathogen (6). The impacts of chytridiomycosis have likely been underestimated, as the affected amphibian species are often difficult to study, particularly in endangered cryptic species that occur in remote locations (48). In addition, not all Bd lineages have the same impact on infected amphibian populations and species, therefore a better understanding of the genetic diversity of Bd is critical to understanding the risk presented by this pathogen and to informing mitigation actions (16). In this study, we describe a mass mortality event due to chytridiomycosis in an endangered species of amphibian in a captive breeding program. We genetically characterized two Bd isolates from this outbreak, showing that they nested within the BdGPL clade and were highly related to Bd genotypes previously isolated from wild amphibians in Chile.
The presence of agnathia and brachygnathia associated with Bd infection in postmetamorphic amphibians has not been reported before. Tadpole oral malformations have been associated with low temperatures (49), water contamination (50), nutrition (51), or ecological factors (52, 53). Although, this malformation might have been due to an unknown environmental or other cause, it is likely associated with Bd infection of the oral discs of tadpoles (54–56). Absence or reduction in development of the lower jaw may have had a profound impact on the ability of postmetamorphic amphibians (and tadpoles) to acquire food, contributing to death along with chytridiomycosis. The presence of distended gall bladders and the absence of gastrointestinal content in all the animals examined suggests a lack of feeding. Although oral deformations in postmetamorphic amphibians have not been used before as an indication of chytridiomycosis, they might be still important and may indicate an unknown sequela of Bd infection that may have been overlooked previously. This study highlights the need to use accurate diagnostic techniques such as qPCR or histology to be able to complement these observations (57).
Captive breeding has increasingly been used as a tool for amphibian conservation, but for these initiatives to be successful, several aspects must be considered, including genetic management and biosecurity protocols (58, 59). In our case, newly admitted individuals of C. gayi came from a semi-open captive breeding program that was supplied with water from an agricultural canal, in which X. laevis had previously been recorded (34). Despite the implementation of quarantine, this was not enough to prevent the introduction of Bd to the captive breeding program, causing mortality in the newly admitted animals. This highlights the importance of implementing strict biosecurity protocols against Bd (and other pathogens), such as Bd testing prior to admittance, preventive antifungal treatment or disinfection of water and materials (59).
Susceptibility to Bd is also influenced by environmental factors, such as climate (20). Immune function in amphibians is closely dependent on environmental temperature (60). For instance, low temperature has been associated with a lower survival in Bd-exposed frogs under laboratory conditions (60, 61) and chytridiomycosis die-offs have been often associated with higher elevation, lower temperature and winter season (62–64). It is possible that in our case, stress associated with transportation might have induced immunosuppression, facilitating the development of chytridiomycosis. Ramsey et al. (65) demonstrated that natural resistance to Bd in X. laevis can be reversed with the implementation of sub-lethal immunosuppressive treatments, such as exposure to X-irradiation or norepinephrine injections.
The Bd lineage that we isolated from infected individuals in this chytridiomycosis outbreak, BdGPL, has been associated with catastrophic mass mortalities and population declines in multiple continents (22, 48, 66, 67). Over the last four decades, a severe population decline of two species of Darwin's frogs (Rhinoderma darwinii and R. rufum) populations has occurred across central and south Chile. Retrospective and cross-sectional Bd studies and population monitoring data suggest that chytridiomycosis has contributed to these extirpations (33, 48) and that Bd has been present in Chile at least since 1970 (33). Although characterization of Bd isolates infecting R. darwinii has not been achieved so far, the identification of BdGPL over a large area of Chile, provides support to the hypothesis that BdGPL is causing these declines (25).
In our phylogenomic whole-genome analysis, all the Chilean Bd isolates (including the two isolates from captive C. gayi) grouped together with a genotype isolated in 2009 from the United Kingdom (UKTvB). A similar phylogenetic relationship was observed when restricting the analysis to a subset of the genome spanning a heterozygosity loss event shared by all BdGPL isolates, but in this case, isolates from other European countries and a Canadian isolate also group with the Chilean isolates (25). While there has been no report of mortality in the wild caused by UKTvB, a challenge with this isolate in the Mallorcan midwife toad (Alytes muletensis) under laboratory conditions caused a 73% mortality rate (18). It is likely, therefore, that the BdGPL genotypes in Chile are virulent and have caused amphibian mortalities in nature (25). However, detecting mortalities in the wild is often difficult, particularly in cryptic species, for which better surveillance is needed. In addition, we also evaluated the presence of Frog virus 3 (FV3) ranavirus in the same individuals as a possible cause of mortality, however all tested negative using a qPCR assay. Previously, a mass mortality of adult C. gayi, allegedly due to a drought, was described in central Chile, however whether Bd was involved in this mortality event is unknown as no fresh carcasses were available for necropsy, Bd detection or other diagnostics to be performed (68).
The low number of segregating sites exclusive to the Chilean Bd isolates, compared to the total number of sites where the BdGPL isolates are polymorphic, suggests a single and recent introduction of Bd into Chile, possibly through the international movement of amphibians, other aquatic animals, fomites, or tourism (9, 14, 25, 66). Molecular characterization of further isolates from non-surveyed areas in Chile and neighboring countries (e.g., Argentina and Peru), along with the calibration of a genome-wide molecular clock, is required to confirm this hypothesis. The existence of an apparently unique and recently introduced lineage of Bd in Chile differs with the known history of this pathogen in Brazil, where both BdGPL and BdAsia2/BdBrazil coexist, with evidence of multiple interlineage recombination events between them (9, 24). Therefore, together with the potential introduction of novel Bd genotypes, interlineage recombination can potentially arise facilitated by the globalization of human and animal transport (16). This highlights the importance of biosecurity measures at the national and local level, to prevent the introduction and establishment of further pathogenic Bd lineages, as this pathogen has the capability to increase its genomic diversity through the exchange of haplotypes among lineages (25).
Our work describes for the first time a mass mortality event in the endangered giant Chilean frog from a captive breeding program. We also provide new data on the potential susceptibility of C. gayi to the impacts of chytridiomycosis, a species that has been declining fast across its distribution in Chile. The high mortality observed in C. gayi with postmetamorphs exhibiting agnathia or brachygnathia as a possible consequence of oral infection with Bd in tadpoles has not been described previously and is a condition that can be considered in the monitoring of amphibians maintained in captivity, such as farms, zoos and ex situ conservation programs. We described two new isolates of Bd in Chile, belonging to BdGPL and clustering in a single group with another three previously isolated Bd isolates from central and south Chile (9, 25), for which evidence as a cause of amphibian mortality and population declines is growing.