Selective constraint acting on TLR2 and TLR4 genes of Japanese Rana frogs

Introduction

Toll-like receptors (TLRs) are a type of pattern recognition receptor that recognize pathogen-associated molecular patterns (PAMPs) such as bacterial cell walls and nucleic acids (Medzhitov, 2001). The signalling of TLRs triggers the synthesis and release of pro-inflammatory cytokines, and thus TLRs have an important role in innate immunity and activation of adaptive immunity. TLRs are type 1 membrane glycoproteins comprised of extracellular and cytoplasmic domains; the extracellular domain is also considered as a ‘pathogen-recognition domain’ with a variable number of leucine rich repeats (LRR, protein motifs in the ectodomain inferred to be important for recognising molecules), while the cytoplasmic signalling domain is a conserved toll/IL-1 domain (Mikami et al., 2012). The TLR repertoire differs between vertebrate groups ranging from 10 loci in humans and great apes to 21 in urodele amphibians (Takeda, Kaisho & Akira, 2003; Quach et al., 2013; Babik et al., 2014). Specific TLR loci are generally considered to recognize different groups of PAMPs, for example, TLR1 and TLR6 associate with TLR2 to recognize acylated lipopeptides, TLR4 recognize lipopolysaccharides in gram-negative bacteria, and TLR5 recognize bacterial flagellin (Poltorak et al., 1998; Voogdt et al., 2016). In addition, TLR2 and TLR4 can recognize PAMPs derived from different fungi species (reviewed by Roeder et al., 2004; Luther & Ebel, 2006).

In amphibians, TLRs have been described in Xenopus frogs (20 loci from 14 families, Ishii et al., 2007) and newts (21 loci from 13 families, Babik et al., 2014). TLR families characterized in both frogs and newts include TLR01, TLR02, TLR03, TLR05, TLR07, TLR08, TLR09, TLR12, TLR13, TLR14, TLR21, and TLR22. However, the prevalence of TLR4 across amphibian taxa is uncertain: no TLR4 orthologs were found in newts, but putative TLR4 were identified in Xenopus (Ishii et al., 2007; Babik et al., 2014). In addition, TLR4 was one of the 11 TLR genes that were isolated from transcriptome data of Bombina maxima frogs (Zhao et al., 2014).

Most TLR genes are functionally constrained, and typically have high sequence conservation and slow evolutionary rates to maintain a functional role of recognizing conserved PAMPs (Roach et al., 2005), although extracellular LRR domains have higher evolutionary rates compared to intracellular domains (Mikami et al., 2012). An overall signature of purifying selection was identified in TLRs of urodele amphibians, although a few individual codons were found to be evolving under positive selection (Babik et al., 2014), as has been observed in other vertebrates (Wlasiuk & Nachman, 2010; Shang et al., 2018). The identification of TLR-disease associations in vertebrates (Tschirren et al., 2013; Noreen & Arshad, 2015) also supports that episodic selection can occur in TLR genes as a response to changes in pathogen diversity. In addition, human TLR4 displayed significantly negative Tajima’s D values in nonsynonymous variants (Smirnova et al., 2001), and subsequent evidence supporting selection for rare TLR4 variants was found (Smirnova et al., 2003). Therefore, while TLRs may be under functional constraints, there is evidence across vertebrates that positive selection may act on these genes in response to local pathogens.

Chytridiomycosis is a disease in amphibians caused by the fungal pathogen Batrachochytridium dendrobatidis (Bd). This disease has been linked to the decline of amphibian populations worldwide (Daszak, Cunningham & Hyatt, 2003; Longcore et al., 2007; Skerratt et al., 2007; Wake & Vredenburg, 2008). Despite Bd being prevalent in Korea and Japan (Goka et al., 2009; Bataille et al., 2013), within endemic East Asian frogs there is no evidence of Bd-related declines and no published reports of Bd susceptibility following experimental infection; this supports that such frogs could be Bd-tolerant. Additionally, genetic evidence for high Bd genetic diversity and endemism in this region indicates that Bd is endemic to Asia (Fisher, 2009; Bataille et al., 2013), suggesting a long co-evolutionary history between the Bd pathogen and Asian amphibians.

While adaptive immune genes of Japanese frogs have been studied in the context of diseases like chytridiomycosis (Lau et al., 2016; Lau et al., 2017), there are limited studies involving innate immunity genes including TLRs. Since TLR2 and TLR4 have been shown to play a role in innate immune responses to various fungi (Roeder et al., 2004; Luther & Ebel, 2006), TLRs therefore may be involved in Bd resistance (Richmond et al., 2009).

Basal expression of TLRs can provide an overview of the function of these genes, and have been examined in a few anuran species including Bombina and Xenopus frogs. In B. maxima, expression levels of TLR2 and TLR4, measured using quantitative RT-PCR, were variable between different adult tissues (Zhao et al., 2014). In X. laevis, ubiquitous expression of both TLR2 and TLR4 was detected in adults and tadpoles using standard PCR, but expression levels were not quantified (Ishii et al., 2007). This current study provides the opportunity to expand the knowledge about basal TLR expression in amphibians.

To better understand the gene complements of anuran innate immunity, here we characterize the genetic diversity and selection patterns of two candidate TLR genes in three Japanese Rana species using molecular cloning and sequence analyses. Given the potential immunological importance of TLRs against fungal pathogens, we hypothesized that TLRs would be subjected to purifying selection in species which show marked resistance to Bd. In addition, we also conducted preliminary investigation of TLR expression using published transcriptomic data, to determine whether basal expression is similar across different life stages and tissue types.

Methods

Results

We characterized a total of 27 TLR2 and 20 TLR4 alleles from the three study species, all of which were distinctive at the amino-acid level. Rana japonica had seven TLR2 and five TLR4 alleles, R. ornativentris had ten TLR2 and seven TLR4 alleles, and R. tagoi tagoi had ten TLR2 and eight TLR4 alleles (Table 1, Figs. S3 and S4, GenBank accession numbers MG999527–MG999573). All alleles clustered phylogenetically into species-specific clades (Fig. 1) and had similar domain structure to that of other frogs (Fig. S5). Allelic diversity among species was high (Table 1, Table S3), especially TLR2 in R. ornativentris and R. t. tagoi, where all individuals were heterozygous with two unique alleles.

Selection tests over the entire alignment indicated that TLR2 and TLR4 in the three focal species are under purifying selection (Z-value = 2.66 − 3.83, p < 0.01, Table 2) with the exception of TLR4 in R. t. tagoi (Z-value = 1.28, p = 0.102). In addition, Tajima’s D, normalized Fay and Wu’s H, and the McDonald-Kreitman test showed no significant support for selection (Table 1). However, omegaMap analyses identified six to nine positively selected sites (PSSs) in either TLR2 or TLR4 of each of the three focal species (Posterior probability >99%). Of these, two (22.2%) to six (75.0%) PSSs were located in predicted leucine rich repeat domains (Table 2, Figs. S3 and S4). There were no PSSs common in all three species studied, but two PSSs of TLR4 (sites D56 and S65) were shared across two species (Table 2, Fig .S4).

Pairwise comparisons of d_N∕d_S ratios between Rana species were low and ranged from 0.188 to 0.398 (Table 3). These ratios were comparable to that of MHC class 1 α3 and MHC class II β2 domains as well as other transcribed genes within the species (mean d_N∕d_S of over 3,000 genes = 0.380, Fig. 2). In contrast, these ratios were lower than that of MHC class I α1 and α2 domains and class II β1 domain.

When we checked published transcriptome data for TLR2 and TLR4 expression, preliminary examination of adult tissues indicated no significant differences between tissue types (blood, spleen, skin) at FDR cut-off of 0.0001. Nevertheless, overall expression of TLR2 and TLR4 was low and there were a few consistent trends seen across all three species. This included higher expression of TLR2 in adult skin relative to blood, and lower expression of TLR4 in adult skin relative to either blood or spleen (Table 4). Across life stages, expression of TLR2 and TLR4 was seemingly low in tadpoles relative to adult samples (FDR > 0.0001, n.s., Table 4).

Discussion

The preliminary characterization of TLR2 and TLR4 in the three focal species here provides a platform for future population genetics studies across the species’ distributions, uncovering the full TLR diversity with more targeted PCR and sequencing approaches. Overall allelic diversity of TLR2 and TLR4 appeared to be high, whereby only a few alleles were shared between more than one individual (Table S3). Commonly shared alleles were more apparent in R. japonica, likely due to a population bottleneck history in the source population of Etajima (Lau et al., 2016). The presence of the TLR4 gene in amphibians was previously unclear, described as ‘putative’ in Xenopus frogs (Ishii et al., 2007) and ‘predicted’ from genomic data in Lithobates catesbeianus and Nanorana parkeri (GenBank accessions XP_018421367 and PIO23183; Sun et al., 2015; Hammond et al., 2017), and undetected in newts (Babik et al., 2014). The characterization of TLR4 genes in this study supports the existence of this gene family in anurans, whereby TLR4 alleles of the three Rana species were similar in phylogeny and domain structure to that of other frogs.

From selection tests, we found overall evidence of purifying selection and no support for sequence-wide positive selection. This agrees with data from other vertebrates, including newts, where TLRs are regarded as conserved with their evolution predominated by purifying selection (Roach et al., 2005; Babik et al., 2014). The d_N∕d_S ratios in TLR2 and TLR4 of the Rana species studied here were remarkably low compared to that of external domains of MHC class I (α1 and α2 domains) and class II (β1 domain), which are considered to be under balancing selection. However, the low d_N∕d_S of TLR2 and TLR4 was comparable to d_N∕d_S of MHC class 1 α3 and MHC class II β2 domains which are intracellular or not involved in peptide recognition, as well as those of over 3,000 transcribed genes within the species. These findings further support that TLR2 and TLR4 are under functional constraint.

Although we found that most of the TLR2 and TLR4 sequences of the Japanese Rana frogs were evolutionarily constrained, we identified evidence of adaptive evolution occurring at individual codon sites in our alignment, similar to other vertebrates studied (Wlasiuk & Nachman, 2010; Shang et al., 2018). When comparing with codon sites predicted to be important for binding of non-fungal ligands (Figs. S3 and S4), two PSSs identified in R.ornativentris TLR2 (Q284 and V299) corresponded to sites in human TLR2 predicted to be involved in ligand binding of lipopeptides (N294 and L312, Jin et al., 2007). In addition, one PSS each of R. tagoi tagoi (T128) and R. ornativentris (Q253) corresponded to human TLR4 sites predicted to be involved in secondary (N268) and phosphate (K388) binding, respectively, of bacterial lipopolysaccharides (Park et al., 2009). Positive selection at identical codon sites across different species, as observed in MHC adaptive immune genes (Lau et al., 2016; Lau et al., 2017), could be driven by a single selective force that is pathogen-related. However, in Rana TLRs there were no PSSs shared across all three species studied. As we did not examine the patterns of selection in species which are susceptible to Bd, any link is currently speculative; nevertheless, we cannot rule out the possibility of adaptive evolution, potentially driven by pathogens such as Bd, acting on TLR of the study species in recent evolutionary history.

Preliminary examination of TLR2 and TLR4 expression levels extracted from transcriptome data showed overall low expression. While the expression data is derived from single individuals that were housed in disease-free environments, it appears that skin of healthy frogs that are not immune-challenged express TLR2 more so than TLR4. However, expression of immune-related genes could be modulated following immune or stress challenges, and future studies should monitor immune gene expression following experimental infection with pathogens like Bd. The adult tissue-specific differences in TLR2 and TLR4 expression from R. ornativentris in this study were distinct from B. maxima (Zhao et al., 2014), but sample size should be increased for both species before further inferences can be made. A previous study in Xenopus frogs detected ubiquitous expression of both TLR2 and TLR4 in adults and tadpoles using PCR, but expression levels were not quantified (Ishii et al., 2007). Although we found low TLR expression in tadpoles in this study, further conclusions cannot be made due to limited sampling and overall low TLR expression across the samples. Future quantitative studies can investigate expression level changes of TLRs during development from tadpole to adults, as well as that of other innate and adaptive immune genes extracted from the transcriptome data set (Lau et al., 2017).

Conclusion

In this study, we characterized TLR2 and TLR4 genes from three Japanese Rana species. We provide strong evidence of purifying selection acting across TLR2 and TLR4, and evidence of a few specific codon sites under positive selection. Further research is necessary to determine if the positive selection we detected is due to pathogen-driven selection. Since immunity to infectious diseases is usually polygenetic, our study adds to the growing body of literature related to genes that potentially impact resistance to Bd and other pathogens in amphibians.

Supplemental Information