Skip to main content

An annual cycle of gene regulation in the red-legged salamander mental gland: from hypertrophy to expression of rapidly evolving pheromones



Cell differentiation is mediated by synchronized waves of coordinated expression for hundreds to thousands of genes, and must be regulated to produce complex tissues and phenotypes. For many animal species, sexual selection has driven the development of elaborate male ornaments, requiring sex-specific differentiation pathways. One such male ornament is the pheromone-producing mental gland of the red-legged salamander (Plethodon shermani). Mental gland development follows an annual cycle of extreme hypertrophy, production of pheromones for the ~ 2 month mating season, and then complete resorption before repeating the process in the following year. At the peak of the mating season, the transcriptional and translational machinery of the mental gland are almost exclusively redirected to the synthesis of rapidly evolving pheromones. Of these pheromones, Plethodontid Modulating Factor (PMF) has experienced an unusual history: following gene duplication, the protein coding sequence diversified from positive sexual selection while the untranslated regions have been conserved by purifying selection. The molecular underpinnings that bridge the processes of gland hypertrophy, pheromone synthesis, and conservation of the untranslated regions remain to be determined.


Using Illumina sequencing, we prepared a de novo transcriptome of the mental gland at six stages of development. Differential expression analysis and immunohistochemistry revealed that the mental gland initially adopts a highly proliferative, almost tumor-like phenotype, followed by a rapid increase in pheromone mRNA and protein. One likely player in this transition is Cold Inducible RNA Binding Protein (CIRBP), which selectively and cooperatively binds the highly conserved PMF 3′ UTR. CIRBP, along with other proteins associated with stress response, have seemingly been co-opted to aid in mental gland development by helping to regulate pheromone synthesis.


The P. shermani mental gland utilizes a complex system of transcriptional and post-transcriptional gene regulation to facilitate its hypertrophication and pheromone synthesis. The data support the evolutionary interplay of coding and noncoding segments in rapid gene evolution, and necessitate the study of co-evolution between pheromone gene products and their transcriptional/translational regulators. Additionally, the mental gland could be a powerful emerging model of regulated tissue proliferation and subsequent resorption within the dermis and share molecular links to skin cancer biology.


The processes of cellular differentiation and tissue remodeling are ubiquitous for multicellular organisms. The transition of cells from totipotency to terminal differentiation are the result of highly coordinated gene networks and changes in expression patterns [1]. Cell differentiation can be induced by a range of signals, including cytokines, hormones, cell contact, external stressors, and extracellular matrix composition [2,3,4,5,6]. Generally, gene expression is regulated through activation of specific transcription factors and/or post-transcriptional mechanisms (e.g. microRNAs and RNA-binding proteins) which permit the regulated expression of particular proteins [7, 8]. An extreme instance of targeted gene expression and differentiation is in the male ornaments of many animal species, such as ornate plumage in birds, vibrant coloration of fishes, and large antlers in deer [9, 10]. These elaborate ornaments are thought to arise from intense sexual selection imposed by female preferences and ultimately mate choice [11]. Such sexual selection can also influence molecular traits and has resulted in many reproductive proteins that evolve at extraordinary rates [12, 13]. For example, in the marine gastropod abalone, males produce billions of sperm that overexpress the rapidly evolving protein lysin (up to ~ 1 g lysin per male abalone) that facilitates species-specific fertilization [14,15,16]; in fruit flies, male accessory glands synthesize complex mixtures of rapidly evolving seminal fluid proteins that restrict female re-mating by reducing her viability and survival [17, 18]; similarly, in Pieris butterflies, males produce enormous spermatophores (~ 13% of their body mass) that are encased in a nearly indestructible protein shell that slows spermatophore clearance and prevents female re-mating [19]. In these examples, however, the molecular mechanisms underlying the regulated expression of these unusual reproductive proteins are not fully understood.

Protein sex pheromones are another example of rapidly evolving reproductive proteins likely shaped by sexual selection through interaction with receptors in the other sex [13, 20]. For more than 60 million years, male plethodontid salamanders have utilized a system of nonvolatile protein courtship pheromones to regulate female reproductive behavior [21]. In the red-legged salamander (Plethodon shermani), during a courtship behavior known as tail-straddling walk, male salamanders privately deliver pheromones by “slapping” a large pad-like gland on his chin (the mental gland) to the female’s nares [22]. Pheromones diffuse into the female nasal cavity, bind to receptors on neurons in the vomeronasal organ, activate regions of the brain involved in pheromone response, and regulate female mating behavior [23,24,25,26,27]. Chemical analysis of the pheromone extract revealed two major components: Plethodontid Receptivity Factor (PRF) and Plethodontid Modulating Factor (PMF). PRF is a 22-kDa protein with sequence similarity to IL-6 cytokines [25], while PMF is a 7-kDa protein related to the highly diverse three-finger protein (TFP) superfamily that includes snake venom neuro- and cytotoxins, the complement receptor CD59, the human Ly6 antigen, the urokinase receptor uPAR, and the amphibian regeneration factor Prod1 [28]. When experimentally applied to female salamanders, both PRF and PMF altered the length of courtship time [25, 29,30,31]. Analysis by high performance liquid chromatography (HPLC) and mass spectrometry (MS) revealed multiple isoforms of both PRF and PMF; however, compared to 3 highly conserved PRF isoforms (> 95% identity), individual male P. shermani expressed more than 30 diverse PMF isoforms (~ 30% identity) [32]. The ratios of different PRF and PMF isoforms are quite variable between male salamanders [33], and the source of isoform sequence diversity is primarily from gene duplication within the diploid genome [32]. Examination of PRF and PMF sequences from 28 plethodontid species revealed that both genes have repeatedly experienced positive selection [34, 35]. Sampling from these many species by RT-PCR was facilitated by the unique quality that both PRF and PMF have unusually conserved, AU-rich untranslated regions (UTRs). The contrast is most striking for PMF: compared to the ~ 30% amino acid identity between isoforms, the average conservation for both the 5′ and 3′ UTRs is ~ 98%. We hypothesized that the coding regions of the many PMF gene copies have been under positive selection to expand the functional breadth of PMF as a pheromone, while purifying selection on the UTRs permitted coordinated, synchronized expression of the many PMF isoforms. The mechanism(s) by which these UTRs mediate such expression remains unknown, but we postulated that RNA binding proteins were likely involved [32].

As with many elaborate male ornaments, the mental gland of male P. shermani is seasonally regulated. During the non-breeding season, it is mostly absent from male salamanders; however, presumably in response to elevated plasma androgens [36, 37], the gland hypertrophies over ~ 2 months and develops into a large pad-like structure solely dedicated to the production of protein pheromones (Fig. 1). Once the gland has fully developed, PRF and PMF represent ~ 85% of the secreted protein [38]. Similarly, cDNA sequence analysis revealed that ~ 70% of the total mRNA coded for pheromones [39]. Following the end of the courtship season, the gland resorbs and a new one forms each subsequent year. It is noteworthy that surgical removal of the mental gland is followed by rapid wound healing and prevents gland regrowth in subsequent years (Lynne D. Houck and R.C. Feldhoff, personal communication), suggesting the existence of androgen-sensitive precursor cells embedded in the dermis. Given that the transcriptional and translational machinery of fully developed mental glands are directed almost exclusively towards pheromone production, there likely exist some earlier developmental phase characterized by greater mitosis and/or general growth to form the glandular structure. We hypothesized that the unusually conserved pheromone UTRs may be critical for regulating the transition from gland hypertrophication to pheromone synthesis. In this study, we applied transcriptome sequencing to characterize the P. shermani mental gland developmental profile at ~ 3 week intervals from an early precursor stage to the active pheromone producing gland. We identified an RNA binding protein that binds to the highly conserved pheromone UTRs and likely contributes to the synchronized transition from gland hypertrophication to the overexpression of rapidly evolving pheromones during the short annual breeding season.

Fig. 1

Mental gland hypertrophication. Comparison of male P. shermani from (a) late May (non-breeding condition) and (b) mid-August (breeding condition), with the mental gland being the large pad-like structure on the male’s lower jaw in panel b


Qualitative observation of mental gland development

In a previous study by Woodley [36], mental gland development was examined from mid-June to late August: during this period, plasma testosterone increased > 3 fold and mental gland volume significantly increased. Based on this work, we collected male P. shermani approximately every 3 weeks from late May through mid-September to span all phases of gland development. Our qualitative observations were similar to those of Woodley [36]. In late May, male salamanders had visibly different skin pigmentation near the mentum compared to females, including 2 of 15 males with extremely faint “outlines” of a mental gland. By mid June, 8 of 10 collected males had visible outlines and/or thin mental glands; and in mid July, all males collected had visible and protruding mental glands. At this mid-July point, pheromone was collected from 8 mental glands by incubation in an Amphibian Ringers buffer containing acetylcholine. Subsequent analysis by reverse-phase high performance liquid chromatography (RP-HPLC) revealed normal proportions of PRF and PMF, but with protein concentrations at ~ 33–50% of levels normally observed in mid-August. Glands from the last three time points (early August, late August, and mid September) were visibly well-developed and had normal levels of pheromone. In a second time series we collected animals from mid-June and early August, and observed similar progression of gland development, suggesting that mental gland development is tightly seasonally regulated in our study population. Histological analysis by hematoxylin/eosin staining showed that the glandular tissue was immediately under the ~ 2–3 cell layer thick epidermis. The mental gland was structured as cylindrically shaped bundles with nuclei localizing exclusively near the periphery. Strong eosin staining was observed in the center of these bundles during early August, contrasting with mid June where only nuclear staining was visible for the smaller, more tightly packed bundles (Fig. 2a). Higher resolution images were obtained using confocal fluorescence microscopy: for early August, actin staining was strongest along the periphery (near the nuclei), with light, diffuse staining and a few small fibers visible in the eosin-stained space, whereas actin was only found adjacent to DAPI-stained nuclei in mid June (Fig. 2b, c). Lectin staining suggested condensation and/or degradation of much of the ECM as the gland expanded (Fig. 2b). Immunohistochemical labelling with anti-PRF produced a strong, punctate pattern throughout the eosin-positive space, reflective of secretory vesicles observed in scanning electroscopy microscopy [40]. There was minimal PRF staining for mid June, both in intensity and volume (Fig. 3). These data suggest that the mental gland initially forms as a tightly packed mass of cells with little cytoplasm, and upon induction of pheromone synthesis, cells swell with large volumes of pheromone, adopt a columnar shape, and the ECM condenses and/or degrades in order to support the enlarged cells.

Fig. 2

Mental gland histology. Comparison of mental glands from male P. shermani from two stages in mental gland development (mid June and early August) using (a) hematoxylin and eosin staining, (b, c) fluorescent confocal microscopy with dyes labelling the nucleus (blue), actin cytoskeleton (red), and ECM (green) at 100X (b) and 400X (c) magnifications

Fig. 3

Pheromone immunohistochemistry. Comparison of pheromone expression and localization for mental glands at two stages of development (mid June and early August) by immunohistochemistry (using anti-PRF; red), with fluorescent dyes labelling the nucleus (blue) and ECM (green)

Mental gland transcriptome and changes in gene expression

A transcriptome of mental gland development was constructed from the six time points spanning late May through mid-September. It should be noted that “mental gland” is used here to describe the tissue where the mental gland is normally found; for the earlier time points, this is largely comprised of skin and connective tissue. Interestingly, bands that correspond to PMF and PRF were visible in an agarose gel of full-length cDNA for all six time points, even in the early points when there was little-to-no pheromone protein (Additional file 1: Figure S1). The intense ~ 850 bp band in the late May and mid June samples was later identified to encode a 15.3 kDa secreted protein with no significant BLAST hits in Genbank. A de novo transcriptome from Illumina reads was constructed using Trinity (Additional file 2: Table S1), resulting in ~ 55,000 groups of reconstructed transcripts – herein referred to as genes for simplicity – which were annotated using several publicly available bioinformatics tools (results included in Additional file 3). Given the unique biology of the mental gland and the evolutionary distance of P. shermani from species with well curated genomes, assigning specific functions to many of our newly identified genes was highly speculative. As such, we chose to focus our analyses principally on genes with evolutionarily deep functional conservation, very high amino acid similarity to annotated homologs in other taxa, or functions described from previous studies in plethodontid salamanders. Nearly all assembled genes were detectable at all six points (at varying levels), and cluster analysis performed to identify genes with correlated expression patterns. Using a hierarchical approach, five major groups of genes were identified with different expression profiles and numbered 1 to 5 in order of decreasing gene count (Fig. 4). Cluster 1 (which contained 87.7% of genes) showed maximum expression in May/June, and included many housekeeping genes (e.g. β-actin, GAPDH, PCNA, ferritin, ribosomal proteins). In relative terms, cluster 2 had the most stable expression patterns (~ 1-4X fold difference between time points), and included a range of genes from different biological pathways, including ribosomal proteins, lysosomal proteases, signal peptidase complex members, and lipid biogenesis enzymes. Cluster 3 included genes almost exclusively found in the earliest time point (with slightly elevated expression in the late August and mid September time points, possibly suggesting a cyclical response as the gland begins to resorb). Some of the most highly expressed genes in cluster 3 included ribosomal proteins (S6, S15, S17, S23, L14, L24, L32, L37a) and histone proteins (H1E, H2A, H3). Clusters 4 and 5 together include the genes most highly expressed in the later phases of gland maturation (0.8% of all genes). As expected, the majority of transcripts coded for pheromone, including PRF and PMF, but also (in lower abundance) many putative pheromones that were identified in other plethodontid species (natriuretic peptide, vasoactive intestinal peptide, sodefrin precursor-like factor, cysteine rich secretory protein) [39, 41, 42]. Included in these sequences was a predicted protein related to the tissue inhibitor of metalloproteinase (TIMP) family that included an extraordinarily long 3′ UTR (~ 3700 nt); by mass spectrometry, this sequence matched a protein previously termed C3 that constitutes ~ 10% of the pheromone extract [33]. The function of this protein is still unknown, but given this new information as to its likely homology, we now refer to it as Plethodontid TIMP-like Protein (PTP) [20, 41]. Multiple other protease inhibitor-like proteins were identified, included cystatin C and multiple Kazal-type inhibitors. Clusters 4 also included several retrotransposon and reverse transcriptase-like sequences. The biological importance of these sequences is unclear, but may explain the presence of intron-less processed PMF pseudogenes in the P. shermani genome [32]. Three other proteins of interest in cluster 4 included acetylcholinesterase (AChE), nuclear protein 1 (NP1), and vascular endothelial growth factor (VEGF). Incubation with acetylcholine is the standard methodology by which pheromones are extracted from the mental gland [25, 43]; under normal biological conditions, co-secretion of AChE would allow for tightly controlled pheromone release [44]. VEGF, an angiogenic factor, may be necessary for mental gland maintenance: as the mental gland enlarges, the cells closest to the dorsal surface will be distant from dermal capillaries and may not receive sufficient nutrition without recruitment of new blood vessels. Related, NP1 classically functions in chromatin remodeling as part of stress responses (such as nutrient starvation) to prevent apoptosis [45, 46]. In mental gland development, NP1 is one of the few genes to steadily increase in mRNA abundance over the time course, and may play an important role in ensuring that the gland persists throughout the courtship season after it has transitioned to pheromone synthesis. As a quality control, qRT-PCR was performed for 16 genes of interest (Additional file 1: Figure S2), and similar expression patterns were observed between RNASeq- and qPCR-based estimates, with no significant biases from gene or time point.

Fig. 4

Cluster analysis of mental gland gene expression. (a) Cartoon of mental gland development for each of the six time points from which RNA was extracted for transcriptome sequencing. Mental glands were not observed in May, partially developed in June, and fully developed by early August. (b) Cladogram representing the ~ 55,000 “genes” organized into 6 clusters. Shades of grey represent log fold changes between time points. (c) Line graphs of cluster means vs time. (d) The six most abundant genes in each cluster (at the time point with the highest expression levels per cluster)

Cold inducible RNA binding protein (CIRBP) binds the PMF 3′ UTR

A perplexing feature of the mental gland developmental transcriptome was the high levels of pheromone mRNA at all time points: from May to August, PMF mRNA levels increased > 100-fold and was still the third and first most abundant transcript in May and June, respectively. However, even at proportionally high levels, pheromone protein was not detected in these early time points. The unusual evolutionary history of PMF – rapidly evolving coding sequence with conserved UTRs – led us to hypothesize that the UTRs may recognize RNA binding proteins (RNA-BPs) to coordinate the expression of the many diverse PMF isoforms [32]. While no RNA-BPs were detected as differentially expressed with a 5% false discovery rate, we found one highly abundant RNA-BP by manual inspection: Cold Inducible RNA-BP (CIRBP), which represented ~ 0.13% of transcripts in mid June. Analysis by qRT-PCR confirmed differential mRNA expression over the six time points, with maximum expression at mid June (Additional file 1: Figure S2). To test for biological activity in vitro, recombinant CIRBP fused to the enhanced cyan fluorescent protein (rCIRBP/ECFP) was expressed in E. coli, and electrophoretic mobility shift assay (EMSA) experiments were performed with rCIRBP/ECFP and in vitro transcribed RNA. rCIRBP/ECFP was titrated against a nearly full length PMF 3′ UTR (nucleotides 26–667), and a very clear shift was observed in both the RNA and protein bands (Fig. 5). Interestingly, there was visible RNA smearing at lower concentrations of rCIRBP/ECFP, suggesting possible dissociation of the RNA-protein complex during electrophoresis. Simultaneously, the altered position of the RNA band in the presence of greater rCIRBP/ECFP suggested a non-equimolar stoichiometry. To narrow the potential binding sequences for CIRBP, four overlapping ~ 250 nt segments of the PMF 3′ UTR were prepared, along with a similar sized segment of the keratin 3′ UTR as a negative control. When these five different RNAs were analysed by EMSA, all showed visible gel shifts in the presence of increasing rCIRBP/ECFP, yet the overlap between RNA/protein was most intense in the PMF 3′ UTR 26–288 and 99–368 fragments (Additional file 1: Figure S3). Variability in band number and intensity in zero-protein control lanes suggested different degrees of RNA secondary structure between the different sequences, which may have had an impact on CIRBP binding. Specificity of CIRBP towards the PMF 3′ UTR was assessed using a competition assay where fluorescently tagged versions of the different RNA molecules were observed in the presence and absence of excess unlabelled RNAs. With PMF 3′ UTR 99–368, addition of a 100-fold excess of unlabeled PMF 3’UTR 99–368 eliminated the gel shift, while a 100-fold excess of unlabeled keratin 3′ UTR only reduced the gel shift to a smear (Fig. 6a). These data suggest that CIRBP has relatively greater affinity for the PMF 3′ UTR, with some lower affinity for other RNA molecules. In a similar competition assay using fluorescent PMF 3′ UTR 26–667, 100X unlabeled RNA was added for three different lengths of the PMF 3′ UTR (26–288, 26–565, 26–667) and keratin 3′ UTR. Only the full length PMF 3′ UTR 26–667 was able to fully eliminate the observed gel shift, suggesting a positive correlation between RNA length and strength of CIRBP binding (Fig. 6b).

Fig. 5

CIRBP – PMF 3′ UTR interaction. EMSA using a constant amount of PMF 3′ UTR in vitro RNA (200 ng) with increasing concentrations of rCIRBP/ECFP (μM; RBP), and a protein-only control. Protein fluorescence was detected by ECFP (green), and RNA was stained using Sybr Green II (red)

Fig. 6

Competition EMSA with CIRBP. (A) EMSA between TAMRA-labeled PMF 3′ UTR 99–368 RNA (30 ng; red) and rCIRBP/ECFP (1.5 μg; RBP; green), with competition using 100X (3 μg) of a specific competitor (SC; unlabeled PMF 3′ UTR 99–368) or a non-specific competitor (NC; unlabeled Keratin 3′ UTR 205–441). (B) EMSA with TAMRA-labeled PMF 3′ UTR 26–668 RNA (30 ng; red) and rCIRBP/ECFP (1.5 μg; RBP; green), with competition using 100X (3 μg) of four different unlabeled competitors: SC1 = PMF 3′ UTR 26–288, SC2 = PMF 3′ UTR 26–565, SC3 = PMF 3′ 26–668, or NC = Keratin 3′ UTR 205–441

Dynamics of CIRBP-PMF 3′ UTR interactions

CIRBP contains two structural domains: a N-terminal RNA recognition motif (RRM) and a C-terminal glycine-rich, low complexity domain (LCD). Studies on the human homolog of CIRBP suggested that both domains bind RNA, with the RRM and LCD having specific and non-specific interactions, respectively [47]. For P. shermani, each domain was expressed as a separate fusion protein to enhanced cyan fluorescent protein (rCIRBP-RRM/ECFP and rCIRBP-LCD/ECFP), and neither domain in isolation induced strong gel shifts with fluorescent PMF 3′ UTR (99–368) (Additional file 1: Figure S4). A small amount of smearing that occurred in the highest concentrations of rCIRBP-LCD/ECFP, suggested a weak interaction. As similar smearing was observed at lower concentrations with rCIRBP/ECFP, EMSAs with rCIRBP/ECFP and PMF 3′ UTR 99–368 were repeated with and without formaldehyde treatment to crosslink protein-RNA complexes. Crosslinking successfully reduced the amount of visible RNA smearing in the gel (Additional file 1: Figure S5), suggesting that under sufficiently low stoichiometry, rCIRBP/ECFP (and likely rCIRBP-LCD/ECFP) forms unstable complexes that readily dissociate during electrophoresis.

When using fluorescently labeled RNAs, interaction with CIRBP caused fluorescence quenching of bound RNA – likely due to shielding of the fluorophores attached to the uracil bases (Fig. 6). This fluorescence quenching was exploited to estimate binding affinities of all three CIRBP constructs to PMF 3′ UTR 99–368 in solution (Fig. 7). While there was no detectable fluorescence quenching with rCIRBP-RRM/ECFP, both rCIRBP/ECFP and rCIRBP-LCD/ECFP yielded sigmoidal curves characteristic of cooperative binding. Fitting of this data to the Hill equation by nonlinear regression yielded both significant equilibrium constants (KD) and Hill coefficients (n), but with rCIRBP-LCD/ECFP having both lower affinity (higher KD) and weaker cooperativity (lower n) compared to rCIRBP/ECFP. Therefore, there likely exists synergism between the two domains to promote binding to the PMF 3′ UTRs.

Fig. 7

CIRBP titration curve. TAMRA-labeled PMF 3′ UTR 99–368 RNA was titrated with increasing concentrations of rCIRBP/ECFP (black), rCIRBP-RRM/ECFP (blue), and rCIRBP-LCD/ECFP (red), and binding measured by fluorescence quenching. Data were fit to the Hill equation by nonlinear modeling to obtain measures of binding affinity (KD) and cooperativity (n). No significant change was detected for rCIRBP-RRM/ECFP, denoted by a dashed line

Recent models of other RNA-BPs with LCDs suggested that, upon binding to a proper catalyst, unstructured LCDs adopt regular β-sheet structure that drive aggregation and formation of stress granules, processing bodies, or other macromolecular RNA-protein complexes [48]. To test if the PMF 3′ UTR may be acting as such a catalyst, rCIRBP/ECFP was analysed by circular dichroism (CD) and titrated with increasing amounts of PMF 3′ UTR 26–667 (Fig. 8). There was a detectable increase in the CD absorbance, particularly near ~ 215 nm where β-sheet can be measured. While the percentage change is relatively small, the majority of CD signal likely originates from ECFP (a highly structured β-barrel that comprises ~ 60% of the fusion protein), and since CD only reports on the average secondary structure content, the “induced” β-sheet in CIRBP may only be occurring in a small proportion of the available molecules. These data support that binding of CIRBP to the PMF 3′ UTR promotes a conformational change and increased secondary structure, likely in the LCD shifting from random coil to β-sheet.

Fig. 8

PMF 3′ UTR induces secondary structure changes in CIRBP. CD spectra of rCIRBP/ECFP with increasing concentrations of PMF 3′ UTR 26–668 RNA, with changes in CD suggesting higher levels of secondary structure (likely β-sheet)

CIRBP expression and correlation to PMF in vivo

Polyclonal antibodies to rCIRBP/ECFP were prepared and affinity purified specifically against CIRBP-RRM. Immunohistochemical staining using anti-CIRBP-RRM revealed that CIRBP protein was localized to the cytoplasm, and more abundant mid June compared to early August. The staining was uniformly distributed, uncharacteristic of cytologically visible RNA granules such as stress granules, Cajal bodies, or nuclear speckles (Fig. 9). CIRBP protein levels for individual mental glands were estimated by western blot analysis. However, multiple bands were observed near and below the approximate 17 kDa expected molecular weight (Fig. 10a-d). We postulated that these bands were CIRBP degradation products generated by a protease sequentially processing the C-terminal LCD. In support of this hypothesis, when proteins from immuno-pulldown were analysed by mass spectrometry, peptides were identified that fully span the N-terminal RRM (Fig. 10b). Addition of either a broad protease inhibitor cocktail or 0.1 mM iodoacetamide (to specifically inhibit cysteine proteases) limited the extent of degradation, but did not completely ablate it (Fig. 10a). The same samples were examined by western blot analyses over multiple days, and even without a protease inhibitor, degradation was incomplete (data not shown). Thus, some of this visible degradation may be a result of natural CIRBP turnover and perhaps play a biological role. Interestingly, neither inhibitor altered the number of observed bands, but rather their relative intensities, such that a single protease may be contributing to both natural and experimentally-imposed degradation. While we were unable to determine the specific protease involved, multiple lines of evidence suggested that it may be cathepsin S: (1) iodoacetamide performed as well or better than the diverse mixture of protease inhibitors, implicating a cysteine protease; (2) most of the identified mental gland proteases in the transcriptome were lysosomal, with cathepsin S having the highest expression; (3) in contrast to most lysosomal proteases, cathepsin S is active only at near-neutral pH, and the lysis buffer was at pH 8; (4) the CIRBP LCD is enriched for the cathepsin S target sequences (aliphatic or aromatic residues followed by Gly [49], specifically YG in the CIRBP-LCD), and cleavage at these sites would generate proteins of similar molecular weight to those observed by western blot (Fig. 10c).

Fig. 9

CIRBP immunohistochemistry. Comparison of CIRBP expression and localization for mental glands at two stages of development (mid June and early August) by immunohistochemistry (using anti-CIRBP-RRM; red), with fluorescent dyes labeling the nucleus (blue) and ECM (green)

Fig. 10

CIRBP protein analysis. (a) Western blot analysis comparing CIRBP degradation with and without different protease inhibitors. Mental glands were dissected into approximate halves, and incubated with RIPA extract containing no protease inhibitors, a commercially available protease inhibitor cocktail (PIC), or 100 μM iodoacetamide (IAA). (b) SDS-PAGE with SYPRO Ruby stain of immunopulldown products using either a rabbit IgG mixture vs anti-CIRBP-RRM for two time points in gland development. (c) Estimated molecular weights (in kDa) for different CIRBP degradation products compared with masses predicted by cleavage of C-terminal Cathepsin S sites (YG). (d) Western blot analysis demonstrating the diversity of CIRBP abundance and degradation state for individual mental glands from mid June (5 μg total protein) vs early August (30 μg total protein) when extracted in the presence of IAA

To ascertain a potential role for CIRBP in regulating mental gland development and PMF synthesis, several variables were correlated from individual mental glands collected at two time points (mid June and early August). Pheromone was extracted by incubation in acetylcholine and protein concentrations measured. Because PMF consistently comprises ~ 50% of the total pheromone [33], we used pheromone concentration as a proxy for PMF protein expression (PPMF). Following pheromone extraction, mental glands were stored in RNAlater, and later dissected into two approximately equal halves, permitting independent isolation of total RNA and cellular protein. Using standardized amounts of cellular protein, CIRBP expression was measured by western blot and densitometry for both total CIRBP (PTotal-CIRBP) and intact CIRBP (PIntact-CIRBP) (Fig. 10d). Using total RNA, expression levels were measured for three genes by qRT-PCR: PMF (RPMF), CIRBP (RCIRBP), and cathepsin S (RCath). Using these 6 variables plus a time covariate (mid June vs early August), a series of nested MANOVAs were performed to identify potentially meaningful correlations (Additional file 2: Table S3). Time was a significant covariate for all variables, and the only significant factor for all three measured mRNA levels. Total CIRBP protein (including degradation products) increased proportionally to CIRBP mRNA levels, but with a higher ratio in June compared to August. Total CIRBP was the best predictor for intact CIRBP levels, however, there were significant interaction terms between time/RCath and time/RCath/RCIRBP/PTotal-CIRBP such that only in June, intact CIRBP levels decreased as Cathepsin S mRNA levels increase (both independently and in proportion to total CIRBP levels). Finally, while neither CIRBP protein variable affected PMF mRNA levels, pheromone protein levels were negatively correlated with intact CIRBP (reflected by the time/PIntact-CIRBP interaction term, since no pheromone protein detected in June). While some of these correlations may be a consequence of natural gland progression in response to other unmeasured variables, these data provide additional evidence to support that Cathepsin S regulates steady-state levels of CIRBP. As intact CIRBP was negatively correlated with PMF protein levels but not mRNA levels, we hypothesize that CIRBP may be acting as a translational repressor.


The process of mental gland development is an exciting example of sex-specific differentiation and organogenesis that yields an elaborate male ornament which expresses numerous rapidly evolving reproductive proteins. The last common ancestor of all plethodontid salamanders likely had a large pad-shaped mental gland similar to P. shermani [50]. Over the last ~ 66 million years [51], it has adopted multiple morphological states (large pad, small pad, bifurcated, fan-shaped, and anterior protrusion) that accommodate varying courtship behaviors to deliver pheromone, including slapping the female’s snout or scratching her dorsum [52]. Coupled with this morphological evolution is accelerated molecular evolution where numerous gene families have been co-opted for pheromone activity and subjected to rapid positive sexual selection [20]. During the non-breeding season, the mentum of male P. shermani appears as normal skin and is visually indistinguishable from that of a female. The presented molecular and histological data support our original hypothesis that mental glands begin in a highly mitogenic state, and then transition into a veritable pheromone factory. Cells first rapidly divide and proliferate, evidenced by elevated PCNA levels at the earliest time points in gland development (Additional file 1: Figure S2). As these cells divide, many housekeeping genes (representing over 80% of the expressed genes) are activated and construct the basic architecture of the gland (Fig. 4). The overall glandular structure seems to consist of many dozens to hundreds of cells surrounding an open lumen with a channel leading to the epidermis (Fig. 2).

While the mental gland is an organ unique to plethodontid salamanders, the initial phase of rapid mitosis likely shares the same molecular players as other highly mitotic systems, such as invasive skin cancers and/or early stages of tumourogenesis. At early stages of development, the mental gland cells consist volumetrically of little more than a nucleus, characteristic of rapid mitosis and mimic a tumour-like structure [53,54,55]. While describing the early mental gland as “tumour-like” may be allegorical, there is a reasonable degree of analogy: the mental gland seemingly “invades” the dermis, compresses and/or degrades the ECM, rapidly divides to form a larger structure, and eventually releases angiogenic factors to increase blood supply to extend its viability. While ECM composition was only qualitatively examined through lectin staining (Figs. 2, 3, 9), it is noteworthy that several different protease inhibitor families (TIMPs, cystatins, Kazal-type serine protease inhibitors) are overexpressed at later stages of development. Many of these proteins have signal peptides, and may be secreted to inhibit proteases that are initially required for mental gland expansion, but must be later inactivated to prevent cell degradation. Notably, Plethodontid TIMP-like Protein is part of the pheromone mixture when extracted with acetylcholine [33], and is presumably released when male salamanders slap to deliver pheromone. It is unclear if the other proteins are released in the same manner, or secreted in response to other stimuli. Alternatively, these protease inhibitors may be necessary to protect synthesized pheromones from premature degradation. What makes this tumour-like phenotype particularly interesting is that, unlike normal cancerous tissue, the mental gland completely resorbs into a non-precursor state following the courtship season. Detailed molecular characterization of both the proliferation and resorption processes could have powerful implications towards cancer biology and identifying potential factors with which to reprogram tumour cells.

For the gland to transition towards pheromone synthesis, extensive changes in gene expression are required, including overexpression of both PRF and PMF mRNA. While the proportions changed, nearly all genes in the transcriptome were detectable at all 6 time points. On average, housekeeping genes in Cluster 1 were still represented at ~ 33% of maximum levels after the upregulation of pheromone synthesis. Despite minimum pheromone translation at early time points (Fig. 3 and Additional file 2: Table S3), PRF and PMF together comprised ~ 1.5% and ~ 10% of all mRNA in May and June, respectively. One possibility may be that androgen stimulation uniformly activates transcription of all mental gland-associated genes. Subsequently, additional variables such as chromatin remodelling and mRNA stability could regulate steady state levels of different transcripts to control the timing and expression of key proteins. Models with this type of gene regulation already exist within the scope of reproductive biology: gametogenesis and the post-fertilization maternal-to-zygotic transition. Chromatin condensation occurs throughout both spermatogenesis and oogenesis, with complete transcriptional silencing in mature spermatids and oocytes [56,57,58,59]. Consequently, mRNAs are commonly regulated by cytoplasmic polyadenylation to control poly(A) tail length, recruitment of the poly(A) binding protein (PABP), and formation of the translation initiation complex [60,61,62]. However, additional RNA-BPs provide further regulation. One example is the Deleted-in-Azoospermia (DAZ) protein and its autosomal homolog, DAZ-like (DAZL) protein. Both DAZ and DAZL recognize target sequences within the 3′ UTR of select mRNAs and can recruit PABP through protein-protein interactions, allowing formation of the initiation complex independently of a poly(A) tail [63,64,65]. However, the stoichiometry and relative spacing of DAZ/DAZL molecules on a 3′ UTR can recruit repressor proteins (DAZAP, PUM2) whose activity is dependent on phosphorylation state [66,67,68]. Genes such as the meiosis-associated sycp3 are regulated by all of the aforementioned regulatory layers, but the net result is a process of ordered and synchronized gene expression independent of changes in transcription [64, 69, 70]. In this context, NP1 and CIRBP may be playing critical roles in chromatin remodelling and controlling mRNA stability, respectively. While we provided direct evidence for CIRBP selectively recognizing the PMF 3′ UTR, further biochemical experiments such as crosslinking mass spectrometry and RNA protein assays will be necessary to fully characterize its molecular function. As genomic technologies expand to salamanders, ATAC-Seq, ChIP-Seq, and other assays of DNA structure will facilitate understanding the roles of nuclear proteins such as NP1.

In mammalian systems, CIRBP plays many diverse, yet highly integrated roles. It is classically recognized for stress granule formation in response to cellular stressors, including heat shock, UV irradiation, hypoxia, oxidative stress, osmotic shock, or arsenic exposure [47, 71]. CIRBP is normally stored in the nucleus, and only translocates to the cytosol upon stress induction where it binds target mRNAs and associates into stress granules [47]. As “cold inducible” RNA binding protein, CIRBP received its name because of its role in testes, which operate ~ 3 °C below normal human body temperature [72]. It was later identified as a major cold shock protein in Arabidopsis, and the Arabidopsis homolog could functionally replace non-homologous cold shock proteins in E. coli [73]. In humans, the temperature-dependent expression is a product of alternative promoter usage, producing a transcript with an elongated 5′ UTR that contains an internal ribosome entry site [74]. CIRBP has recently been identified as a candidate gene for temperature-dependent sex determination in the common snapping turtle [75]. Echoing the linkage between the mental gland and cancerous tumours, CIRBP is also part of the human telomerase complex and plays a crucial role in telomere maintenance [76]. Because of the seasonal nature of the mental gland and its presence in a cold-blooded animal, temperature-dependent expression in P. shermani is a compelling hypothesis that will require further investigation. In one study [48], the low complexity domains of multiple RNA-BPs were found to form hydrogels composed of β-sheet cross strands similar to amyloid fibers. Repeats of the tripeptide Gly/Ser-Tyr-Gly/Ser were determined to be essential for formation, and likely involved π-π overlap between tyrosine residues. Using mutagenesis to compare the effect of 0, 12, 18, 22, or 27 tripeptide repeats, hydrogel size and stability were positively correlated with the number of repeats; in contrast, P. shermani CIRBP only has 7 repeats, such that it may form smaller aggregates that are not visible microscopically. Given that CIRBP-RNA interactions induced a conformational change in CIRBP secondary structure and recruited additional CIRBP molecules through cooperativity, it is plausible that CIRBP may coat the entire PMF 3′ UTR (possibly the whole PMF mRNA), and interfere with ribosome binding and translation. As PMF mRNA abundance changes dramatically over gland development, this binding may also facilitate mRNA degradation. Stress granules are generally thought to protect mRNA molecules from degradation [77, 78], but their close proximity and possible association with processing bodies (P-bodies) has led to the hypothesis that there may be mRNA exchange between these macromolecular complexes, leading to mRNA degradation under the proper cellular conditions [79, 80].

Similar to reports on the human homolog [47], we determined that both the RRM and LCD were required for stable, highly cooperative binding of CIRBP to the PMF 3′ UTR. Structural comparison of different RRM-RNA complexes from other RNA-BPs has revealed that the RRM is a highly plastic structure with respect to nucleotide binding. The αβ sandwich structure contains two conserved regions of 7 and 6 amino acids (termed RNP1 and RNP2) which bind single stranded dinucleotides on either DNA or RNA through a range of interactions, including π-π overlap between aromatic residues and the nitrogen bases, hydrophobic interactions between additional aromatics and the sugar moieties, and/or salt bridges with the phosphate between the two nucleotides. Despite these conserved structural elements, sequence specificity is dictated by less conserved residues forming additional interactions [81]. RNA secondary structure also plays a critical role, with some RNA-BPs requiring dinucleotides to be part of stem loops or other types of internal loops [82,83,84]. When the dinucleotide frequency within the PMF 3′ UTR was examined relative to predicted values based on nucleotide abundance, five of the sixteen combinations were observed at higher than expected frequencies: UC (+ 26%), UG (+ 58%), CA (+ 7%), CU (+ 29%), and GA (+ 8%). When these distributions were compared by χ2 test, UG contributed much more to the test statistic relative to the other dinucleotides with higher than expected frequencies (χ2UG = 12.02, compared to χ2CU = 3.30). Future experiments examining both CIRBP-RRM specificity as well as PMF 3′ UTR secondary structure will be required to determine if these altered frequencies hold functional significance.

In addition to CIRBP, many other cold shock proteins were identified in the mental gland and had similar expression patterns (Figs. 4, Additional file 1: Fig. S2). All evidence suggests that CIRBP likely plays some critical role in facilitating the transition of the gland from growth and development into pheromone synthesis by interacting with PMF and possibly other mRNAs (likely relating to translational repression and/or mRNA degradation). While CIRBP was originally characterized based on its roles in temperature-dependent expression, more recent studies have demonstrated its roles in systemic and general stress response [47, 85]. In addition to facilitating proliferation and development, the role of CIRBP as a regulator of pheromone synthesis is also interesting within the scope of gene co-option. A common theme among plethodontid pheromone genes – as well snake toxins and the products of other exocrine tissues – is gene duplication followed by rapid evolution which drives neofunctionalization [32, 86, 87]. Through the same basic processes, it is plausible that novel regulators could be recruited and co-opted to provide tight control over gene expression of exocrine products in these systems. With many of these gene products having evolutionary histories of positive selection, it would likely be of value to explore both the products and their regulators in a co-evolutionary framework. At least in the case of PMF, CIRBP may have been the causative agent driving purifying selection on the UTRs while sexual selection from the female receptors promoted positive selection on the coding regions. The putative role of CIRBP directing the selective forces on the noncoding regions of PMF might be relevant to future studies of the regulatory elements of other rapidly evolving genes.


For more than 60 million years, plethodontid salamanders have utilized a rapidly evolving system of non-volatile protein courtship pheromones that regulate female mating behaviour and facilitate reproduction. In the red-legged salamander (P. shermani), these pheromones are synthesized in a submandibular mental gland that annually hypertrophies into a large, pad-like structure whose transcriptional and translational machinery are almost exclusively programmed for pheromone synthesis. We now report that this highly effective system results from a tightly coordinated gene expression cascade, allowing for annual organogenesis followed by rapid conversion to an efficient pheromone factory. Glandular cells initially have a phenotype characterized by rapid proliferation and ECM dissolution, followed by a tremendous increase in pheromone mRNA levels. A key regulator in this process is Cold Inducible RNA Binding Protein: a stress-responsive RNA binding protein used by both animals and plants to store select mRNAs in stress granules and promote cell survival. For at least one pheromone (Plethodontid Modulating Factor), CIRBP selectively binds the 3′ UTR and cooperatively recruits additional molecules, with protein-protein and protein-RNA interactions likely stabilized through induced intermolecular β-sheets. This interaction may inhibit translation of PMF mRNA and/or promote its degradation through association with P-bodies. The net result is suppression of pheromone translation until the gland is sufficiently large to support the storage of 10s to 100s of micrograms of pheromone, which can be used to increase receptivity in almost any female in the breeding population. CIRBP may be one player that has exerted purifying selection on the PMF UTRs, creating a set of genes with highly conserved noncoding segments and rapidly evolving coding regions. The mechanisms behind this exciting dichotomy may serve as a model for future studies of gene regulation on rapidly evolving proteins.


Animal collection, gland removal, and pheromone extraction

P. shermani males were collected during their breeding season from a single site in Macon Co., North Carolina, USA (35°10′48″ N, 83°33′38″ W). Males were anesthetized in a mixture of 7% (v/v) diethyl ether in water. For analysis of total RNA, pheromone extract, or cellular proteins, mental glands were surgically removed using iridectomy scissors. For RNA analysis, glands were incubated overnight in RNAlater (Ambion, Austin, TX) at 4 °C before long term storage at − 20 °C. Pheromone was extracted from mental glands based on the protocol in Chouinard et al. [33]: briefly, mental glands were individually incubated in 0.2 mL acetylcholine chloride (0.8 mM in Amphibian Ringer’s solution) for 30 min, centrifuged at 14,000 x g for 10 min, the supernatant collected, and centrifugation repeated before storage at − 80 °C. Following extraction, mental glands were stored in RNAlater to allow preservation of RNA and cellular proteins. Pilot experiments confirmed that 30 min in the acetylcholine solution was sufficient to extract > 90% of the total pheromone with no detectable RNA degradation (Wilburn and Feldhoff, unpublished data). Following mental gland removal, each male was placed in a clean plastic box lined with a damp paper towel and his chin rested on a small piece of gauze containing antibiotic ointment. Males were allowed to heal before use in additional experiments outside of this study or returned to the field site. For histological analyses, salamanders were anesthetized in 7% (v/v) diethyl ether in water and sacrificed by rapid decapitation.

cDNA preparation and transcriptome sequencing

Mental glands were collected from male P. shermani at six time points approximately every 3 weeks during 2010 (5/29, 6/19, 7/10, 8/1, 8/21, and 9/11). This range preceded and spanned the principal August mating season. Five glands were collected at each time point and immediately stored in RNAlater. Total RNA was extracted from individual glands using the RNeasy kit (Qiagen, Valencia, CA) following the manufacturer’s instructions. RNA concentrations were estimated by 260 nm absorbance. For each time point, standardized amounts of RNA were pooled from all five glands, and double-stranded cDNA prepared by oligo-dT priming using the SMARTer cDNA synthesis kit (Clontech, Palo Alto, CA). cDNA was supplied to Otogenetics Corporation (Norcross, GA) for library preparation and sequenced using the Illumina HiSeq 2000 platform (> 20 million reads per time point, 100-bp paired end reads). Data was deposited in the NCBI Sequence Read Archive (BioProject ID PRJNA427596).

Transcriptome bioinformatics analysis

Illumina reads from all six time points were pooled and assembled into a single transcriptome using Trinity (r2012-10-05) [88]. Initial assemblies with default settings resulted in over-compaction of deBruijn graphs for PMF, limiting both isoform detection and full-length mRNA re-construction. Butterfly parameters were then optimized, and the final assembly included the additional settings --min_kmer_cov 2 --bfly_opts “-path_reinforcement_distance= 25 -min_per_id_same_path= 98”. Reads from each time point were re-aligned to the full transcriptome using RSEM (v1.2.5) [89], with differential expression analysis conducted using EBSeq [90]. Expression differences were compared between adjacent time points (5/29 to 6/19, 6/19 to 7/10, etc.). Based on visual observations and analysis of pheromone extract from additional glands collected at each time point, the separate phases of gland development were best characterized by the 6/19 time point (growth/development) and the 8/1 time point (pheromone production); therefore, an additional comparison with EBSeq was performed between the 6/19 and 8/1 time points. For putative gene annotation, the Trinotate package (r2013-02-25) was used to determine (1) putative open reading frames (using TransDecoder, [88]), (2) protein BLAST (blastp) against both the SwissProt and TrEMBL databases [91,92,93,94], (3) orthologous group identification with eggNOG [95], (4) gene ontology assignment [96], (5) signal peptide prediction with SignalP [97], (6) protein family assignment with Pfam [98], and (7) transmembrane domain prediction with TmHMM [99]. As most of these databases searches relied on proper open reading frame assignment by TransDecoder, an additional BLAST search was performed with blastx using assembled nucleotide sequences against the full TrEMBL database [92, 94]. There were multiple cases of the TransDecoder proteins having no blastp hits in SwissProt or TrEMBL, yet the nucleotide sequence produced strong blastx hits (e-value < 0.001), suggesting that TransDecoder identified the wrong open reading frame. In these cases, an alternative open reading frame was selected based on the longest amino acid sequence that contained the aligned region of the blastx hit. Hierarchical clustering of RNASeq expression profiles was performed using the MBCluster package in R with a negative binomial model and FPKM values reported by RSEM. The number of clusters was varied, and based on manual inspection of both branch lengths as well as biology of the genes being grouped together, five clusters was found to be optimal.

qRT-PCR analysis of differentially expressed genes

For select genes (see Additional file 2: Table S4 for primers), qRT-PCR analysis was individually performed on RNA samples isolated from the mental gland that were pooled to construct the transcriptome (six time points each with 5 glands, 30 glands total). Total RNA was diluted to ~ 5–20 ng/uL, and accurate concentrations were determined using Quant-iT RiboGreen RNA assay kit (Invitrogen, Carlsbad, CA). qRT-PCR reactions were performed in triplicate using 20 μL reactions with the Power SYBR Green RNA-to-CT 1-Step kit (Ambion) containing 1 μL diluted total RNA. Expression levels were calculated by the pcrfit function in the R package qpcR using the cm3 mechanistic model [100]. Based on the gross morphological changes of the mental gland, it was expected that few genes would be stably expressed across mental gland development; therefore, RNA input was used to normalize qPCR measurements, with the literature supporting that input is often a more robust reference than housekeeping genes [101, 102]. Expression levels for each gene were normalized based on the time point with the highest expression, and a single linear mixed effect model was fit by maximum likelihood using the lme function in the R package lmer. The model included fixed effects for gene and time, with male as a random effect; all variables, including the interaction between gene and time, were significant at p < 0.001.

Expression of recombinant CIRBP

To conduct in vitro protein-RNA binding assays, recombinant CIRBP was prepared using an E. coli expression system. Preliminary experiments revealed that CIRBP was insoluble at > 0.2 mg/mL in all tested buffers (data not shown); however, solubility was dramatically improved when recombinant CIRBP was expressed as a fusion protein with enhanced cyan fluorescent protein (rCIRBP/ECFP). Simultaneously, fusion proteins were prepared with only the RRM-containing N-terminus (residues 1–84; rCIRBP-RRM/ECFP) and the glycine-rich C-terminus (residues 85–165; rCIRBP-LCD/ECFP). All constructs included ECFP on the N-terminus, a short hydrophilic linker (SGAAAAGGSDP), and the CIRBP element at the C-terminus. The CIRBP coding regions were amplified from mental gland cDNA using the Accuprime High Fidelity Taq polymerase system (Invitrogen) (see Additional file 2: Table S4 for primers). ECFP was amplified from a pcDNA3.1-based vector (supplied by Dr. Ronald Gregg, University of Louisville), and fusion genes were prepared by modified assembly PCR [103]. Fusion PCR products were purified using the QIAquick PCR cleanup system (Qiagen), and cloned the pET45b vector (Novagen, San Diego, CA) following restriction digest with KpnI and HindIII (New England Biolabs, Ipswisch, MA), gel purification of cleavage products (GFX purification system, GE Healthcare, Piscataway, NJ), ligation using T4 DNA Ligase (New England Biolabs), and transformation into T7 Express lysY/Iq chemically competent E. coli (New England Biolabs). Sequence of transformed DNA was validated by colony PCR and Sanger sequencing by the University of Louisville DNACore facility. For protein expression, clones were cultured overnight at 28 °C in LB media with 100 μg/mL ampicillin (LB/Amp), diluted in 1 L LB/Amp with a 600 nm optical density (OD600) equal to ~ 0.05, incubated with shaking until the OD600 equalled ~ 0.7, IPTG added to a final concentration of 0.1 mM, and incubated overnight. E. coli were harvested by centrifugation, resuspended in 50 mM NaCl/0.1% Triton X-100/2 mM EDTA/50 mM Tris, pH 8, and lysed by sonication followed by lysozyme treatment (final concentration 1 mg/mL) for 1 h. Insoluble material was removed by centrifugation, and proteins concentrated by ammonium sulfate precipitation (final concentration 70%). Ammonium sulfate pellets were resolubilized in Ni-NTA binding buffer (0.5 M NaCl/5 mM Imidazole/20 mM Tris, pH 8), and passed through a 15 mL Ni-NTA column (Thermo-Pierce, Rockford, IL) at 1 mL/min. The column was then washed and eluted with increasing concentrations of imidazole (all in 0.5 M NaCl/20 mM Tris, pH 8): 10 column volumes (CVs) at 5 mM, 2 CVs at 20 mM, 1 CV at 40 mM, 1 CV at 60 mM, and finally 3 CVs at 200 mM (collected in 5 mL fractions). Highly fluorescent fractions were pooled, concentrated, and buffer exchanged to 1X EMSA Buffer (100 mM KCl/2 mM EDTA/0.05% Tween-20/20 mM HEPES, pH 7.4) by ultrafiltration (YM-10 Centiprep, Millipore). All fusion products were standardized to equimolar concentrations (rCIRBP/ECFP: 1.4 mg/mL, rCIRBP-N/ECFP: 1.1 mg/mL, rCIRBP-C/ECFP: 1.0 mg/mL).

Electrophoretic mobility shift assays

CIRBP-RNA interactions were characterized by electrophoretic mobility shift assays (EMSAs). RNA based on the PMF Class I 3′ UTR was prepared by in vitro transcription using the T7 High Yield RNA Synthesis Kit (New England Biolabs) with purified PCR products amplified from P. shermani 8/1 cDNA using primers that included engineered T7 promoters (Additional file 2: Table S4). Synthesized RNA was subsequently treated with TURBO DNase I (Ambion), and purified using the RNeasy Kit (Qiagen). For fluorescently labeled RNA, in vitro transcription reactions were adjusted to include 7.5 mM UTP and 2.5 mM aminoallyl-UTP (Ambion), treated with TURBO DNase I (Ambion), purified using the RNeasy kit (Qiagen), and adjusted to 0.4 mg/mL with 2.5 mg/mL TAMRA-carboxylic acid (Invitrogen) in 0.1 M MES (pH 5). One-fourth volume of EDAC (0.1 mM) was added to the reaction, incubated for 2 h in the dark, and then labelled RNA was purified using the RNeasy Kit (Qiagen). All EMSAs were performed as 15 μL reactions in 1X EMSA Buffer with 3% glycerol using either 300 ng unlabeled or 45 ng labeled RNA. rCIRBP/ECFP (or the domain truncations), competitors, and other reagents were added at the listed concentrations, and incubate at room temperature for 20 min prior to gel loading. For crosslinking assays, formaldehyde was added to the protein-RNA mixture at a final concentration of 1%, incubated for 5 min, and the reaction quenched by addition of glycine to a final concentration of ~ 330 mM prior to gel loading. RNA was separated using agarose gels of 2% (for PMF 3′ UTR 26–667) or 3% (for ~ 250 nt products), and electrophoresis performed at 80 V for 3 h. For unlabeled RNA, gels were stained with Sybr Green II RNA Stain (Invitrogen) for 30 min. Gels were imaged using a Typhoon 9400 fluorescent bed scanner (GE Life Sciences, Piscataway, NJ) with the appropriate laser and filter settings. For fluorescence quenching assay, 0.5 μg of TAMRA-labeled RNA was incubated with increasing concentrations of different rCIRBP constructs in a 200 μL reaction, and fluorescence measured using the Synergy2 plate reader (Biotek, Winooski, VT). Data were fit to the Hill equation (\( \theta =\frac{1}{{\left(\frac{K_D}{\left[L\right]}\right)}^n+1} \)) using the nlsLM function in the R package minpack.lm.

Antibody preparation

PRF antisera was prepared by immunizing two rabbits with PRF fractions collected by anion-exchange chromatography (performed by RCF at U. Louisville). The antigen was > 90% enriched for PRF, but contained ~ 5% PTP such that there was significant immunoreactivity for both proteins. CIRBP antisera was prepared by immunizing two rabbits with purified rCIRBP/ECFP (~ 98%) by Ni-NTA chromatography (performed by the Proteintech Group, Chicago, IL). Highly purified antigens (> 99% by reverse phase HPLC) were coupled to agarose matrices for affinity purification of antibodies: briefly, 6% crosslinked agarose beads (CL6B; Sigma-Aldrich, St. Louis, MO) were activated with carbonyldiimidazole (CDI; Sigma-Aldrich) in dry acetone for 20 min, quickly rinsed under vacuum with cold distilled water, and then incubated in antigen solutions (~ 2 mg/mL in 100 mM KCl/0.05% Tween-20/100 mM NaCO3, pH 11) overnight with gentle agitation. Resins were packed into ~ 1-3 mL columns (depending on the starting amount of antigen), remaining active sites blocked with 50 mM Tris, pH 10, and equilibrated in 1X phosphate buffered saline (PBS). Antibodies were purified by incubating antisera with the resin overnight at 4 °C, washing with 10 CVs of 1X PBS to collect unbound antibodies, stringently washed with 5 CVs of 0.5 M NaCl/0.05% Tween-20/20 mM Tris, pH 7.5 (TTBS), and eluted with 6 CVs 0.1 M citric acid/0.1% Tween-20, pH 3 collected in 0.2 mL fractions. Elution fractions were neutralized with 1 M Na2HPO4, and antibody-containing fractions determined by absorbance measurements at 280 nm. Fractions were then pooled, concentrated, and buffer exchanged to 1X PBS by ultrafiltration (YM-30 Centriprep; Millipore, Billerica, MA). Four antigen columns were prepared: PRF, PTP, rPRF/ECFP, and rCIRBP-RRM/ECFP. PRF and PTP antibodies were serially purified from the same antisera by first passing whole antisera through the PRF column for purification of anti-PRF, and the flowthrough (i.e., unbound antibodies) then passed over the PTP column for purification of anti-PTP. Both solutions of affinity-purified antibodies were then adsorbed against the opposing column to ensure removal of non-specific or low-affinity antibodies. Anti-CIRBP-RRM was purified by passing rCIRBP/ECFP antisera over the rPRF/ECFP column to remove anti-ECFP, the flowthrough next incubated with the rCIRBP-N/ECFP column, anti-CIRBP-RRM purified, and the two-column process repeated to ensure removal of non-specific or low-affinity antibodies. Antibody specificity was validated by western blot.

Histological analysis

For histological analysis of mental glands, male salamanders were anesthetized in 7% ether, sacrificed by rapid decapitation, the lower jaw removed, immediately placed in 4% paraformaldehyde (in 150 mM NaCl/100 mM Na2HPO4, pH 7.4), and incubated overnight with gentle agitation. Fixed jaws were decalcified by incubation in 10% EDTA (pH 7.4, DEPC-treated) for 48 h, cryoprotected in 30% sucrose for 48 h, embedded in Optimum Cutting Temperature (OCT) media (Sakura-Finetek, Torrance, CA), and stored at − 80 °C prior to cryosectioning. Lower jaws were sectioned coronally at a thickness of 16 μm (for immunohistochemistry) or 40 μm (for structural morphology), and thaw mounted onto superfrost plus slides pre-coated with polylysine. Slides were stored at − 20 °C until analyzed. Hematoxylin/eosin staining was performed using standard protocols [104]. For immunohistochemical labeling, slides were first equilibrated to room temperature for 30 min and washed five times with 1X PBS for 5 min each. Antigen retrieval was conducted by incubating the slides in 10 mM citric acid (pH 6)/0.05% Tween-20 at 70 °C for 30 min, cooled to room temperature for 20 min, and washed twice with PBS with 0.05% Tween-20 (PBST) for 5 min each. Blocking was performed for 30 min with 1X PBS/0.1% BSA/0.5% Tween-20, and incubated overnight with primary antibody (anti-PRF at 7.5 μg/mL or anti-CIRBP-RRM at 3.5 μg/mL) in PBST with 0.1% BSA. Slides were then washed five times with PBST for 5 min each, incubated with secondary antibody (Alexa Fluor 633 goat anti-rabbit IgG (Invitrogen), 1 μg/mL in PBST) for 30 min, washed five times with PBST for 5 min each, and finally counterstained with Alexa Fluor 488-wheat germ agglutinin (Invitrogen; 10 μg/mL in PBST) and DAPI (Invitrogen; 3.6 μM). For gland morphology, thick sections (40 μm) were stained with Alexa Fluor 488-wheat germ agglutinin, Alexa Fluor 555-Phalloidin (Invitrogen; 0.005 U/μL), and DAPI. Slides were coverslipped with Prolong Gold Antifade reagent (Invitrogen), cured overnight in the dark, and stored at 4 °C. Imaging was accomplished using an Olympus Fluoview FV1000.

Circular dichroism

To measure secondary structure changes in rCIRBP/ECFP with different concentrations of RNA, far-UV CD spectra (185–255 nm) were acquired by averaging 10 scans across a 0.1-cm path at 0.2 nm intervals using a Jasco J-810 Spectropolarimeter. Stock rCIRBP/ECFP (1.4 mg/mL in 1X EMSA) buffer was diluted 10-fold in water (to reduce chloride background). Following initial measurements with no RNA, PMF 3′ UTR 26–667 (at 70 ng/μL in 0.1X EMSA Buffer) was titrated into the solution, incubated for 5 min, and CD spectra recorded. Spectra were adjusted for slight changes in volume/concentration following serial dilution throughout the experiment. At the concentrations used, PMF 3′ UTR 26–667 produced no measurable CD signal over buffer.

Quantification and modelling of CIRBP expression

For two separate time points in 2013 (6/13 and 8/3), 15 male P. shermani were collected, mental glands removed, pheromone extracted, and mental glands stored in RNAlater (Ambion) as previously described. One gland for 8/3 was inadvertently destroyed, such that n = 29. Pheromone concentrations were accurately determined by BCA Protein Assay (Thermo-Pierce), and normal proportions of pheromone components (PMF:PRF:PTP ~ 5:3:1) were validated by RP-HPLC. Stored mental glands were later dissected into approximately equal halves to extract both RNA (RNeasy Kit, Qiagen) and cellular protein (homogenization in RIPA Buffer, supplemented with 0.1 mM iodoacetamide). Pilot studies confirmed that the two halves were equivalent for levels of target proteins and mRNA. CIRBP protein levels were measured by western blot. Preliminary experiments suggested that glands from 6/13 had ~6X higher levels of CIRBP, such that, in order to maintain similar blot intensities for quantification, 5 μg was loaded per lane for 6/13 glands and 30 μg per lane for 8/3 glands. Briefly, proteins were separated using 15% Tris-Tricine gels with 4% stacking gels [105] at 50 V for 15 min followed by 100 V for 90 min, transferred to PVDF membranes by electroblotting in 15% methanol/0.025% SDS/10 mM CAPS, pH 11 at 55 V for 75 min. The membrane was blocked with 0.5 M NaCl/0.5% Tween-20/20 mM Tris, pH 7.5 for 1 h, incubated with anti-CIRBP-RRM (1.4 μg/mL) for 1 h, then incubated with alkaline phosphatase-conjugated goat anti-rabbit IgG (1 μg/mL; Sigma-Aldrich), and developed using BCIP/NBT. All membrane washes and antibody dilutions were performed using TTBS (0.5 M NaCl/0.05% Tween-20/20 mM Tris, pH 7.5). All blots were processed and developed simultaneously in order to minimize run-to-run variation, and a reference of 10 ng rCIRBP/ECFP was included on each membrane for normalization. Densitometry analysis was performed using ImageJ. Extracted RNA was used for qRT-PCR analysis for PMF, CIRBP, and Cathepsin S based on the previous methods (see Additional file 2: Table S4 for primers). Correlation between variables was determined by multivariate ANOVA (MANOVA).

Immunopulldown and mass spectrometry of CIRBP

To validate the identity of different bands visualized during western blot, immuno-pulldown experiments were performed. Anti-CIRBP-RRM was used for CIRBP pulldown, and as a negative control, an IgG-enriched fraction from pre-immunization serum of CIRBP immunized rabbits was prepared by ammonium sulfate precipitation and DE52 chromatography. Using Protein G coupled Dynabeads (Invitrogen), 10 μg of antibody was adsorbed to the beads, stringently washed with PBST, incubated with pooled RIPA extract from 6/13 (~ 180 μg) or 8/3 (~ 840 μg) for 30 min, washed with PBST, the beads transferred to a clean 1.7 mL tube, and incubated with 1X gel loading buffer (1% SDS/12% glycerol/0.005% bromophenol blue/50 mM Tris, pH 6.8) at 65 °C for 30 min. The complete samples were loaded into a 15% Tris-Tricine gel with 4% stacking gel, electrophoresed at 50 V for 15 min followed by 100 V for 90 min, and stained using SYPRO Ruby fluorescent gel stain (Invitrogen). The gel was imaged using the Typhoon 9400 fluorescent bed scanner (GE Life Sciences). All bands were then individually excised, proteins reduced and alkylated with dithiothreitol/iodoacetamide, treated overnight with trypsin, peptides isolated by acetonitrile extraction, and supplied to the University of Louisville Biomolecular Mass Spectrometry Core Laboratory for analysis by LC/MS-MS.





Circular dichroism


Cold inducible RNA binding protein






Electrophoretic mobility shift assay


High performance liquid chromatography


Low complexity domain


Mass spectrometry


Nuclear protein 1


Optical density at 600 nm


Poly(A) binding protein


Processing body


Phosphate buffered saline


Plethodontid modulating factor


Plethodontid receptivity factor


Plethodontid TIMP-like factor


Recombinant CIRBP fused to enhanced cyan fluorescent protein


RNA binding protein


Reverse phase HPLC


RNA recognition motif


Three-finger protein


Tissue inhibitor of metalloproteinase


untranslated region


Vascular endothelial growth factor


  1. 1.

    Tamayo P, Slonim D, Mesirov J, Zhu Q, Kitareewan S, Dmitrovsky E, Lander ES, Golub TR. Interpreting patterns of gene expression with self-organizing maps: methods and application to hematopoietic differentiation. Proc Natl Acad Sci. 1999;96(6):2907–12.

    Article  CAS  Google Scholar 

  2. 2.

    Reilly GC, Engler AJ. Intrinsic extracellular matrix properties regulate stem cell differentiation. J Biomech. 2010;43(1):55–62.

    Article  Google Scholar 

  3. 3.

    Sekiguchi H, Ii M, Jujo K, Thorne T, Ito A, Klyachko E, Hamada H, Kessler J, Tabata Y, Kawana M, et al. Estradiol promotes neural stem cell differentiation into endothelial lineage and angiogenesis in injured peripheral nerve. Angiogenesis. 2013;16(1):45–58.

    Article  CAS  Google Scholar 

  4. 4.

    Chadwick K, Wang L, Li L, Menendez P, Murdoch B, Rouleau A, Bhatia M. Cytokines and BMP-4 promote hematopoietic differentiation of human embryonic stem cells. Blood. 2003;102(3):906–15.

    Article  CAS  Google Scholar 

  5. 5.

    Johe KK, Hazel TG, Muller T, Dugich-Djordjevic MM, McKay RD. Single factors direct the differentiation of stem cells from the fetal and adult central nervous system. Genes Dev. 1996;10(24):3129–40.

    Article  CAS  Google Scholar 

  6. 6.

    Spradling A, Drummond-Barbosa D, Kai T. Stem cells find their niche. Nature. 2001;414(6859):98–104.

    Article  CAS  Google Scholar 

  7. 7.

    Glisovic T, Bachorik JL, Yong J, Dreyfuss G. RNA-binding proteins and post-transcriptional gene regulation. FEBS Lett. 2008;582(14):1977–86.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  8. 8.

    Krol J, Loedige I, Filipowicz W. The widespread regulation of microRNA biogenesis, function and decay. Nat Rev Genet. 2010;11(9):597–610.

    Article  CAS  Google Scholar 

  9. 9.

    Arnold SJ, Houck LD. Can the fisher-Lande process account for birds of paradise and other sexual radiations? Am Nat. 2016;187(6):717–35.

    Article  Google Scholar 

  10. 10.

    Zuk M, Thornhill R, Ligon JD, Johnson K, Austad S, Ligon SH, Thornhill NW, Costin C. The role of male ornaments and courtship behavior in female mate choice of red jungle fowl. Am Nat. 1990;136(4):459–73.

    Article  Google Scholar 

  11. 11.

    Andersson M. Sexual selection. Princeton, NJ: Princeton Univ. Press; 1994.

    Google Scholar 

  12. 12.

    Swanson WJ, Vacquier VD. The rapid evolution of reproductive proteins. Nature Review Genetics. 2002;3(2):137–44.

    Article  CAS  Google Scholar 

  13. 13.

    Wilburn DB, Swanson WJ. From molecules to mating: rapid evolution and biochemical studies of reproductive proteins. J Proteomics. 2016;135:12–25.

  14. 14.

    Lee Y-H, Ota T, Vacquier VD. Positive selection is a general phenomenon in the evolution of abalone sperm lysin. Mol Biol Evol. 1995;12(2):231–8.

    CAS  PubMed  Google Scholar 

  15. 15.

    Lewis CA, Talbot CF, Vacquier VD. A protein from abalone sperm dissolves the egg vitelline layer by a nonenzymatic mechanism. Dev Biol. 1982;92(1):227–39.

    Article  CAS  Google Scholar 

  16. 16.

    Wilburn DB, Tuttle LM. Klevit RE. Solution structure of sperm lysin yields novel insights into molecular dynamics of rapid protein evolution. Proceedings of the National Academy of Science USA: Swanson WJ; 2018;115(6):1310-15.

  17. 17.

    Wolfner M. The gifts that keep on giving: physiological functions and evolutionary dynamics of male seminal proteins in drosophila. Heredity. 2002;88(2):85–93.

    Article  CAS  Google Scholar 

  18. 18.

    Wigby S, Sirot LK, Linklater JR, Buehner N, Calboli FCF, Bretman A, Wolfner MF, Chapman T. Seminal fluid protein allocation and male reproductive success. Current biology : CB. 2009;19(9):751–7.

    Article  CAS  Google Scholar 

  19. 19.

    Meslin C, Cherwin TS, Plakke MS, Small BS, Goetz BJ, Morehouse NI, Clark NL. Structural complexity and molecular heterogeneity of a butterfly ejaculate reflect a complex history of selection. Proc Natl Acad Sci. 2017;114(27):E5406–13.

    Article  CAS  Google Scholar 

  20. 20.

    Wilburn DB, Arnold SJ, Houck LD, Feldhoff PW, Feldhoff RC. Gene duplication, co-option, structural evolution, and phenotypic tango in the courthsip pheromones of plethodontid salamanders. Herpetologica. 2017;73(3):206–19.

    Article  Google Scholar 

  21. 21.

    Houck LD, Arnold SJ. Courtship and mating behavior. In: Sever DM, editor. Phylogeny and Reproductive Biology of Urodela (Amphibia). Enfield: New Hampshire: Science Publishers; 2003. p. 383–424.

    Google Scholar 

  22. 22.

    Houck LD, Bell AM, Reagan-Wallin NL, Feldhoff RC. Effects of experimental delivery of male courtship pheromones on the timing of courtship in a terrestrial salamander, Plethodon jordani (Caudata: Plethodontidae). Copeia. 1998;1998:214–9.

    Article  Google Scholar 

  23. 23.

    Wirsig-Wiechmann CR, Houck LD, Feldhoff PW, Feldhoff RC. Pheromonal activation of vomeronasal neurons in plethodontid salamanders. Brain Res. 2002;952:335–44.

    Article  CAS  Google Scholar 

  24. 24.

    Laberge F, Roth G. Connectivity and cytoarchitecture of the ventral telencephalon in the salamander Plethodon shermani. J Comp Neurol. 2005;482:176–200.

    Article  Google Scholar 

  25. 25.

    Rollmann SM, Houck LD, Feldhoff RC. Proteinaceous pheromone affecting female receptivity in a terrestrial salamander. Science. 1999;285:1907–9.

    Article  CAS  Google Scholar 

  26. 26.

    Fontana MF, Houck LD, Staub NL. In situ localization of plethodontid courtship pheromone mRNA in foramlin-fixed tissue. Gen Comp Endocrinol. 2006;150:480–5.

    Article  CAS  Google Scholar 

  27. 27.

    Wirsig-Wiechmann CR, Houck LD, Wood JM, Feldhoff PW, Feldhoff RC. Male pheromone protein components activate female vomeronasal neurons in the salamander Plethodon shermani. BMC Neurosci. 2006;7(1):26.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  28. 28.

    Palmer CA, Hollis DM, Watts RA, Houck LD, McCall MA, Gregg RG, Feldhoff PW, Feldhoff RC, Arnold SJ. Plethodontid modulating factor, a hypervariable salamander courtship pheromone in the three-finger protein superfamily. FEBS J. 2007;274:2300–10.

    Article  CAS  Google Scholar 

  29. 29.

    Wilburn DB, Eddy SL, Chouinard AJ, Arnold SJ, Feldhoff RC, Houck LD. Pheromone isoform composition differentially affects female behaviour in the red-legged salamander, Plethodon shermani. Anim Behav. 2015;100:1–7.

    Article  Google Scholar 

  30. 30.

    Houck LD, Palmer CA, Watts RA, Arnold SJ, Feldhoff PW, Feldhoff RC. A new vertebrate courtship pheromone that affects female receptivity in a terrestrial salamander. Anim Behav. 2007;73:315–20.

    Article  Google Scholar 

  31. 31.

    Houck LD, Watts RA, Arnold SJ, Bowen KE, Kiemnec KM, Godwin HA, Feldhoff PW, Feldhoff RC. A recombinant courtship pheromone affects sexual receptivity in a plethodontid salamander. Chem Senses. 2008;33(7):623–31.

    Article  CAS  Google Scholar 

  32. 32.

    Wilburn DB, Bowen KE, Gregg RG, Cai J, Feldhoff PW, Houck LD, Feldhoff RC. Proteomic and UTR analyses of a rapidly evolving hypervaraible family of vertebrate pheromones. Evolution. 2012;66(7):2227–39.

    Article  CAS  Google Scholar 

  33. 33.

    Chouinard AJ, Wilburn DB, Houck LD, Feldhoff RC: Individual variation in pheromone isoform ratios of the red-legged salamander, Plethodon shermani. In: Chemical Signals in Vertebrates 12. Edited by East ML, Dehnhard M. New York: Springer; 2013: 99–115.

  34. 34.

    Palmer CA, Watts RA, Gregg RG, McCall MA, Houck LD, Highton R, Arnold SJ. Lineage-specific differences in evolutionary mode in a salamander courtship pheromone. Mol Biol Evol. 2005;22(11):2243–56.

    Article  CAS  Google Scholar 

  35. 35.

    Palmer CA, Picard AL, Watts RA, Houck LD, Arnold SJ. Rapid evolution of plethodontid modulating factor (PMF), a hypervariable salamander courtship pheromone, is driven by positive selection. J Mol Evol. 2010;70:427–40.

    Article  CAS  Google Scholar 

  36. 36.

    Woodley SK. Plasma androgen levels, spermatogenesis, and secondary sexual characteristics in two species of plethodontid salamanders with dissociated reproductive patterns. Gen Comp Endocrinol. 1994;96(2):206–14.

    Article  CAS  Google Scholar 

  37. 37.

    Sever DM. Induction of secondary sex characters in Eurycea quadridigitata. Copeia. 1976;1976(4):830–3.

    Article  Google Scholar 

  38. 38.

    Feldhoff RC, Rollmann SM, Houck LD: Chemical analyses of courtship pheromones in a Plethodontid salamander. In: Advances in Chemical Signals in Vertebrates. Edited by Johnston RE, Műller-Schwarze D, Sorensen P. New York: Kluwer Academic/Plenum; 1999: 117–125.

  39. 39.

    Kiemnec-Tyburczy KM, Watts RA, Gregg RG, von Borstel D, Arnold SJ. Evolutionary shifts in courtship pheromone composition revealed by EST analysis of plethodontid salamander mental glands. Gene. 2009;432:75–81.

    Article  CAS  Google Scholar 

  40. 40.

    Rupp AE, Sever DM. Histology of mental and caudal courtship glands in three genera of plethodontid salamanders (Amphibia: Plethodontidae). Acta Zool. 2018;99(1):20–31.

    Article  Google Scholar 

  41. 41.

    Wilburn DB, Bowen KE, Feldhoff PW, Feldhoff RC. Proteomic analyses of courtship pheromones in the redback salamander, Plethodon cinereus. J Chem Ecol. 2014;40(8):928–39.

    Article  CAS  Google Scholar 

  42. 42.

    Doty KA, Wilburn DB, Bowen KE, Feldhoff PW, Feldhoff RC. Co-option and evolution of non-olfactory proteinaceous pheromones in a terrestrial lungless salamander. J Proteome. 2016;135:101–11.

    Article  CAS  Google Scholar 

  43. 43.

    Rollmann SM, Houck LD, Feldhoff RC. Population variation in salamander courtship pheromones. J Chem Ecol. 2000;26:2713–22.

    Article  CAS  Google Scholar 

  44. 44.

    Rubin LL. Increases in muscle Ca2+ mediate changes in acetylcholinesterase and acetylcholine receptors caused by muscle contraction. Proc Natl Acad Sci. 1985;82(20):7121–5.

    Article  CAS  Google Scholar 

  45. 45.

    Chowdhury U, Samant R, Fodstad O, Shevde L. Emerging role of nuclear protein 1 (NUPR1) in cancer biology. Cancer Metastasis Rev. 2009;28(1–2):225–32.

    Article  CAS  Google Scholar 

  46. 46.

    Cano CE, Hamidi T, Sandi MJ, Iovanna JL. Nupr1: the Swiss-knife of cancer. J Cell Physiol. 2011;226(6):1439–43.

    Article  CAS  Google Scholar 

  47. 47.

    De Leeuw F, Zhang T, Wauquier C, Huez G, Kruys V, Gueydan C. The cold-inducible RNA-binding protein migrates from the nucleus to cytoplasmic stress granules by a methylation-dependent mechanism and acts as a translational repressor. Exp Cell Res. 2007;313(20):4130–44.

    Article  CAS  Google Scholar 

  48. 48.

    Kato M, Han Tina W, Xie S, Shi K, Du X, Wu Leeju C, Mirzaei H, Goldsmith Elizabeth J, Longgood J, Pei J, et al. Cell-free formation of RNA granules: low complexity sequence domains form dynamic fibers within hydrogels. Cell. 2012;149(4):753–67.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  49. 49.

    Biniossek ML, Nägler DK, Becker-Pauly C, Schilling O. Proteomic identification of protease cleavage sites characterizes prime and non-prime specificity of cysteine cathepsins B, L, and S. J Proteome Res. 2011;10(12):5363–73.

    Article  CAS  Google Scholar 

  50. 50.

    Sever DM, Siegel DS, Taylor MS, Beachy CK. Phylogeny of mental glands, revisited. Copeia. 2016;104(1):83–93.

    Article  PubMed  PubMed Central  Google Scholar 

  51. 51.

    Shen X-X, Liang D, Chen M-Y, Mao R-L, Wake DB, Zhang P. Enlarged multilocus dataset provides surprisingly younger time of origin for the Plethodontidae, the largest family of salamanders. Syst Biol. 2015;65(1):66–81.

    Article  Google Scholar 

  52. 52.

    Houck LD, Sever DM. The role of the skin in reproduction and behavior. In: Heatwole H, Barthalamus G, editors. Amphibian Biology, vol. 1. Chipping Norton: Australia: Surrey Beatty and Sons; 1994. p. 351–81.

    Google Scholar 

  53. 53.

    Sheikh MS, Rochefort H, Garcia M. Overexpression of p21WAF1/CIP1 induces growth arrest, giant cell formation and apoptosis in human breast carcinoma cell lines. Oncogene. 1995;11(9):1899–905.

    CAS  PubMed  Google Scholar 

  54. 54.

    Kellogg DR. Wee1-dependent mechanisms required for coordination of cell growth and cell division. Journal of Cell Science. 2003;116(24):4883–90.

    Article  CAS  Google Scholar 

  55. 55.

    Yokota J, Akiyama T, Fung YK, Benedict WF, Namba Y, Hanaoka M, Wada M, Terasaki T, Shimosato Y, Sugimura T. Altered expression of the retinoblastoma (RB) gene in small-cell carcinoma of the lung. Oncogene. 1988;3(4):471–5.

    CAS  PubMed  Google Scholar 

  56. 56.

    Thongkukiatkul A, Jungudomjaroen S, Ratanapahira C. Spermatogenesis and chromatin condensation in male germ cells of sea cucumber Holothuria leucospilota (Clark, 1920). Tissue Cell. 2008;40(3):167–75.

    Article  CAS  Google Scholar 

  57. 57.

    Chapman J, Michael S. Proposed mechanism for sperm chromatin condensation/decondensation in the male rat. Reprod Biol Endocrinol. 2003;1(1):20.

    Article  PubMed  PubMed Central  Google Scholar 

  58. 58.

    Swain JE, Ding J, Brautigan DL, Villa-Moruzzi E, Smith GD. Proper chromatin condensation and maintenance of histone H3 phosphorylation during mouse oocyte meiosis requires protein phosphatase activity. Biol Reprod. 2007;76(4):628–38.

    Article  CAS  Google Scholar 

  59. 59.

    Tan J-H, Wang H-L, Sun X-S, Liu Y, Sui H-S, Zhang J. Chromatin configurations in the germinal vesicle of mammalian oocytes. Mol Hum Reprod. 2009;15(1):1–9.

    Article  CAS  Google Scholar 

  60. 60.

    Tay J, Richter JD. Germ cell differentiation and synaptonemal complex formation are disrupted in CPEB knockout mice. Dev Cell. 2001;1(2):201–13.

    Article  CAS  Google Scholar 

  61. 61.

    Richter JD. CPEB: a life in translation. Trends Biochem Sci. 2007;32(6):279–85.

    Article  CAS  Google Scholar 

  62. 62.

    Belloc E, Pique M, Mendez R. Sequential waves of polyadenylation and deadenylation define a translation circuit that drives meiotic progression. Biochem Soc Trans. 2008;36:665–70.

    Article  CAS  Google Scholar 

  63. 63.

    Collier B, Gorgoni B, Loveridge C, Cooke HJ, Gray NK. The DAZL family proteins are PABP-binding proteins that regulate translation in germ cells. EMBO J. 2005;24(14):2656–66.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  64. 64.

    Brook M, Smith JWS, Gray NK. The DAZL and PABP families: RNA-binding proteins with interrelated roles in translational control in oocytes. Reproduction. 2009;137(4):595–617.

    Article  CAS  Google Scholar 

  65. 65.

    Hasegawa E, Karashima T, Sumiyoshi E, Yamamoto M. C. elegans CPB-3 interacts with DAZ-1 and functions in multiple steps of germline development. Dev Biol. 2006;295(2):689–99.

    Article  CAS  Google Scholar 

  66. 66.

    Moore FL, Jaruzelska J, Fox MS, Urano J, Firpo MT, Turek PJ, Dorfman DM, Pera RAR. Human Pumilio-2 is expressed in embryonic stem cells and germ cells and interacts with DAZ (deleted in AZoospermia) and DAZ-like proteins. Proc Natl Acad Sci U S A. 2003;100(2):538–43.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  67. 67.

    Tsui S, Dai T, Roettger S, Schempp W, Salido EC, Yen PH. Identification of two novel proteins that interact with germ-cell-specific RNA-binding proteins DAZ and DAZL1. Genomics. 2000;65(3):266–73.

    Article  CAS  Google Scholar 

  68. 68.

    Morton S, Yang H-T, Moleleki N, Campbell DG, Cohen P, Rousseau S. Phosphorylation of the ARE-binding protein DAZAP1 by ERK2 induces its dissociation from DAZ. Biochem J. 2006;399(2):265–73.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  69. 69.

    Yuan L, Liu J-G, Zhao J, Brundell E, Daneholt B, Höög C. The murine SCP3 gene is required for synaptonemal complex assembly, chromosome synapsis, and male fertility. Mol Cell. 2000;5(1):73–83.

    Article  CAS  Google Scholar 

  70. 70.

    Reynolds N, Collier B, Bingham V, Gray NK, Cooke HJ. Translation of the synaptonemal complex component Sycp3 is enhanced in vivo by the germ cell specific regulator Dazl. RNA. 2007;13(7):974–81.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  71. 71.

    Li S, Zhang Z, Xue J, Liu A, Zhang H: Cold-inducible RNA binding protein inhibits H2O2-induced apoptosis in rat cortical neurons. Brain Res 2012, 1441(0):47–52.

  72. 72.

    Nishiyama H, Danno S, Kaneko Y, Itoh K, Yokoi H, Fukumoto M, Okuno H, Millan JL, Matsuda T, Yoshida O, et al. Decreased expression of cold-inducible RNA-binding protein (CIRP) in male germ cells at elevated temperature. Am J Pathol. 1998;152(1):289–96.

    CAS  PubMed  PubMed Central  Google Scholar 

  73. 73.

    Kim JS, Park SJ, Kwak KJ, Kim YO, Kim JY, Song J, Jang B, Jung C-H, Kang H. Cold shock domain proteins and glycine-rich RNA-binding proteins from Arabidopsis thaliana can promote the cold adaptation process in Escherichia coli. Nucleic Acids Res. 2007;35(2):506–16.

    Article  CAS  Google Scholar 

  74. 74.

    Al-Fageeh MB, Smales CM. Cold-inducible RNA binding protein (CIRP) expression is modulated by alternative mRNAs. RNA. 2009;15(6):1164–76.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  75. 75.

    Schroeder AL, Metzger KJ, Miller A, Rhen T. A novel candidate gene for temperature-dependent sex determination in the common snapping turtle. Genetics. 2016;203(1):557–71.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  76. 76.

    Zhang Y, Wu Y, Mao P, Li F, Han X, Zhang Y, Jiang S, Chen Y, Huang J, Liu D, et al. Cold-inducible RNA-binding protein CIRP/hnRNP A18 regulates telomerase activity in a temperature-dependent manner. Nucleic Acids Res. 2016;44(2):761–75.

    Article  CAS  Google Scholar 

  77. 77.

    Hilgers V, Teixeira D, Parker R. Translation-independent inhibition of mRNA deadenylation during stress in Saccharomyces cerevisiae. RNA. 2006;12(10):1835–45.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  78. 78.

    Gowrishankar G, Winzen R, Dittrich-Breiholz O, Redich N, Kracht M, Holtmann H. Inhibition of mRNA deadenylation and degradation by different types of cell stress. Biol Chem. 2006;387(3):323–7.

    Article  CAS  Google Scholar 

  79. 79.

    Buchan JR, Parker R. Eukaryotic stress granules: the ins and outs of translation. Mol Cell. 2009;36(6):932–41.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  80. 80.

    Parker R, Sheth U. P bodies and the control of mRNA translation and degradation. Mol Cell. 2007;25(5):635–46.

    Article  CAS  Google Scholar 

  81. 81.

    Maris C, Dominguez C, Allain FHT. The RNA recognition motif, a plastic RNA-binding platform to regulate post-transcriptional gene expression. FEBS J. 2005;272(9):2118–31.

    Article  CAS  Google Scholar 

  82. 82.

    Oubridge C, Ito N, Evans PR, Teo CH, Nagai K. Crystal structure at 1.92 a resolution of the RNA-binding domain of the U1A spliceosomal protein complexed with an RNA hairpin. Nature. 1994;372(6505):432–8.

    Article  CAS  Google Scholar 

  83. 83.

    Allain FHT, Gubser CC, Howe PWA, Nagai K, Neuhaus D, Varani G. Specificity of ribonucleoprotein interaction determined by RNA folding during complex formation. Nature. 1996;380(6575):646–50.

    Article  CAS  Google Scholar 

  84. 84.

    Varani L, Gunderson SI, Mattaj IW, Kay LE, Neuhaus D, Varani G. The NMR structure of the 38 kDa U1A protein - PIE RNA complex reveals the basis of cooperativity in regulation of polyadenylation by human U1A protein. Nat Struct Mol Biol. 2000;7(4):329–35.

    Article  CAS  Google Scholar 

  85. 85.

    Qiang X, Yang W-L, Wu R, Zhou M, Jacob A, Dong W, Kuncewitch M, Ji Y, Yang H, Wang H, et al. Cold-inducible RNA-binding protein (CIRP) triggers inflammatory responses in hemorrhagic shock and sepsis. Nat Med. 2013;19(11):1489–95.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  86. 86.

    Watts RA, Palmer CA, Feldhoff RC, Feldhoff PW, Houck LD, Jones AG, Pfrender ME, Rollmann SM, Arnold SJ. Stabilizing selection on behavior and morphology masks positive selection on the signal in a salamander pheromone signaling complex. Mol Biol Evol. 2004;21(5):1032–41.

    Article  CAS  Google Scholar 

  87. 87.

    Fry BG. From genome to “venome”: molecular origin and evolution of the snake venom proteome inferred from phylogenetic analysis of toxin sequences and related body proteins. Genome Res. 2005;15(3):403–20.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  88. 88.

    Grabherr MG, Haas BJ, Yassour M, Levin JZ, Thompson DA, Amit I, Adiconis X, Fan L, Raychowdhury R, Zeng Q, et al. Full-length transcriptome assembly from RNA-Seq data without a reference genome. Nat Biotech. 2011;29(7):644–52.

    Article  CAS  Google Scholar 

  89. 89.

    Li B, Dewey C. RSEM: accurate transcript quantification from RNA-Seq data with or without a reference genome. BMC Bioinformatics. 2011;12(1):323.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  90. 90.

    Leng N, Dawson JA, Thomson JA, Ruotti V, Rissman AI, Smits BMG, Haag JD, Gould MN, Stewart RM, Kendziorski C. EBSeq: an empirical Bayes hierarchical model for inference in RNA-seq experiments. Bioinformatics. 2013;29(8):1035–43.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  91. 91.

    Altschul SF, Gish W, Miller W, Myers EW, Lipman DJ. Basic local alignment search tool. J Mol Biol. 1990;215(3):403–10.

    Article  CAS  Google Scholar 

  92. 92.

    Camacho C, Coulouris G, Avagyan V, Ma N, Papadopoulos J, Bealer K, Madden T. BLAST+: architecture and applications. BMC Bioinformatics. 2009;10(1):421.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  93. 93.

    Bairoch A, Apweiler R. The SWISS-PROT protein sequence data Bank and its new supplement TREMBL. Nucleic Acids Res. 1996;24(1):21–5.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  94. 94.

    Consortium TU. Ongoing and future developments at the universal protein resource. Nucleic Acids Res. 2011;39(suppl 1):D214–9.

    Article  CAS  Google Scholar 

  95. 95.

    Powell S, Szklarczyk D, Trachana K, Roth A, Kuhn M, Muller J, Arnold R, Rattei T, Letunic I, Doerks T, et al. eggNOG v3.0: orthologous groups covering 1133 organisms at 41 different taxonomic ranges. Nucleic Acids Res. 2012;40(D1):D284–9.

    Article  CAS  Google Scholar 

  96. 96.

    Ashburner M, Ball CA, Blake JA, Botstein D, Butler H, Cherry JM, Davis AP, Dolinski K, Dwight SS, Eppig JT, et al. Gene ontology: tool for the unification of biology. Nat Genet. 2000;25:25–9.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  97. 97.

    Petersen TN, Brunak S, von Heijne G, Nielsen H. SignalP 4.0: discriminating signal peptides from transmembrane regions. Nat Meth. 2011;8(10):785–6.

    Article  CAS  Google Scholar 

  98. 98.

    Punta M, Coggill PC, Eberhardt RY, Mistry J, Tate J, Boursnell C, Pang N, Forslund K, Ceric G, Clements J, et al. The Pfam protein families database. Nucleic Acids Res. 2012;40(D1):D290–301.

    Article  CAS  Google Scholar 

  99. 99.

    Krogh A, Larsson B, von Heijne G, Sonnhammer ELL. Predicting transmembrane protein topology with a hidden markov model: application to complete genomes. J Mol Biol. 2001;305(3):567–80.

    Article  CAS  PubMed  Google Scholar 

  100. 100.

    Carr AC, Moore SD. Robust quantification of polymerase chain reactions using global fitting. PLoS One. 2012;7(5):e37640.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  101. 101.

    Bustin S. Absolute quantification of mRNA using real-time reverse transcription polymerase chain reaction assays. J Mol Endocrinol. 2000;25(2):169–93.

    Article  CAS  Google Scholar 

  102. 102.

    Bustin S. Quantification of mRNA using real-time reverse transcription PCR (RT-PCR): trends and problems. J Mol Endocrinol. 2002;29(1):23–39.

    Article  CAS  Google Scholar 

  103. 103.

    Stemmer WPC, Crameri A, Ha KD, Brennan TM, Heyneker HL. Single-step assembly of a gene and entire plasmid from large numbers of oligodeoxyribonucleotides. Gene. 1995;164(1):49–53.

    Article  CAS  Google Scholar 

  104. 104.

    Lillie R. Histopathologic technic and practical Histochemistry. 3rd ed. New York: McGraw-Hill Book Co; 1985.

    Google Scholar 

  105. 105.

    Schägger H, von Jagow G. Tricine-sodium dodecyl sulfate-polyacrylamide gel electrophoresis for the separation of proteins in the range from 1 to 100 kDa. Anal Biochem. 1987;166:368–79.

    Article  Google Scholar 

Download references


We thank Pamela W. Feldhoff, Kari A. Doty, Kathleen Bowen, Lynne D. Houck, Stevan Arnold, Sarah Eddy, Adam Chouinard, and Sean Darrow for assistance in the field, and the Highlands Biological Station for providing laboratory space and housing. We also thank the University of Louisville Biomolecular Mass Spectrometry Core Laboratory for their continued support.


Funding was provided in part by National Science Foundation (Collaborative) grants IOS-1146899 (RCF), a NSF Graduate Research Fellowship to DBW, and two HBS Grant-in-Aid awards to DBW. The funding agencies did not contribute to the design of the study, participate in the collection, analysis, and interpretation of the data, or aid in the writing of the manuscript.

Availability of data and materials

Sequence data generated by the project was deposited in the NCBI Sequence Read Archive (BioProject ID PRJNA427596). Further data supporting the conclusions of this article are included within the article and its additional files. Bacterial strains and antibodies are available from the authors upon request.

Author information




The study was conceived and designed by DBW and RCF. All experiments and data analysis were performed by DBW. The manuscript was drafted by DBW and edited by RCF. Both authors read and approved the final manuscript.

Corresponding author

Correspondence to Damien B. Wilburn.

Ethics declarations

Ethics approval and consent to participate

Salamanders were collected under permits issued by the North Carolina Wildlife Resource Commission (to R.C. Feldhoff and D.B. Wilburn), and all animal protocols were approved by University of Louisville IACUC (#12041 to Pamela W. Feldhoff) and Highlands Biological Station IACUC (to D.B. Wilburn).

Consent for publication

Not applicable.

Competing interests

The authors declare that they have no competing interests.

Publisher’s Note

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

Additional files

Additional file 1:

Figure S1. P. shermani mental gland cDNA. Figure S2. qRT-PCR analysis of select mental gland genes. Figure S3. CIRBP affinity for different RNAs. Figure S4. EMSA with CIRBP domains. Figure S5. CIRBP – RNA interactions stabilized by formaldehyde crosslinking. (PDF 1358 kb)

Additional file 2:

Table S1. Summary of P. shermani mental gland transcriptome. Table S2. Select list of differentially expressed genes. Table S3. MANOVA results for CIRBP expression and regulation. Table S4. Primers used for qRT-PCR, CIRBP expression, and in vitro transcription. (PDF 168 kb)

Additional file 3:

Assembled gene data. Spreadsheet with Trinity assembled genes, associated annotations, and gene expression estimates over the time series. (XLSX 57561 kb)

Rights and permissions

Open Access This article is distributed under the terms of the Creative Commons Attribution 4.0 International License (, which permits unrestricted use, distribution, and reproduction in any medium, provided 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 Creative Commons Public Domain Dedication waiver ( applies to the data made available in this article, unless otherwise stated.

Reprints and Permissions

About this article

Verify currency and authenticity via CrossMark

Cite this article

Wilburn, D.B., Feldhoff, R.C. An annual cycle of gene regulation in the red-legged salamander mental gland: from hypertrophy to expression of rapidly evolving pheromones. BMC Dev Biol 19, 10 (2019).

Download citation


  • Salamander
  • Pheromone
  • Mental gland
  • Organogenesis
  • Molecular evolution
  • Untranslated region
  • RNA binding proteins