#PAGE_PARAMS# #ADS_HEAD_SCRIPTS# #MICRODATA#

Two transcriptionally distinct pathways drive female development in a reptile with both genetic and temperature dependent sex determination


Authors: Sarah L. Whiteley aff001;  Clare E. Holleley aff002;  Susan Wagner aff001;  James Blackburn aff003;  Ira W. Deveson aff003;  Jennifer A. Marshall Graves aff001;  Arthur Georges aff001
Authors place of work: Institute for Applied Ecology, University of Canberra, Canberra, Australia aff001;  Australian National Wildlife Collection CSIRO National Research Collections Australia, Canberra, Australia aff002;  Garvan Institute of Medical Research, Sydney, Australia aff003;  St. Vincent’s Clinical School, UNSW, Sydney, Australia aff004;  Latrobe University, Melbourne, Australia aff005
Published in the journal: Two transcriptionally distinct pathways drive female development in a reptile with both genetic and temperature dependent sex determination. PLoS Genet 17(4): e1009465. doi:10.1371/journal.pgen.1009465
Category: Research Article
doi: https://doi.org/10.1371/journal.pgen.1009465

Summary

How temperature determines sex remains unknown. A recent hypothesis proposes that conserved cellular mechanisms (calcium and redox; ‘CaRe’ status) sense temperature and identify genes and regulatory pathways likely to be involved in driving sexual development. We take advantage of the unique sex determining system of the model organism, Pogona vitticeps, to assess predictions of this hypothesis. P. vitticeps has ZZ male: ZW female sex chromosomes whose influence can be overridden in genetic males by high temperatures, causing male-to-female sex reversal. We compare a developmental transcriptome series of ZWf females and temperature sex reversed ZZf females. We demonstrate that early developmental cascades differ dramatically between genetically driven and thermally driven females, later converging to produce a common outcome (ovaries). We show that genes proposed as regulators of thermosensitive sex determination play a role in temperature sex reversal. Our study greatly advances the search for the mechanisms by which temperature determines sex.

Keywords:

Calcium signaling – Embryos – Gene expression – Gene regulation – Gonads – Heat shock response – Sex determination – Sexual differentiation

Introduction

Sex determination in vertebrates may be genetic or environmental. In genetic sex determination (GSD), offspring sex is determined by sex chromosomes inherited from each parent, which bear either a dominant gene on the heteromorphic sex chromosome (as with SRY in humans) [1,2], or a dosage sensitive gene on the homomorphic sex chromosome (as with DMRT1 in birds) [3]. However, some fish and many reptile species exhibit environmental sex determination (ESD), whereby a variety of external stimuli can determine sex, most commonly involving temperature (temperature dependent sex determination, TSD) [4,5]. While GSD and ESD are commonly viewed as a dichotomy, the reality is far more complex. Sex determination in vertebrates exists as a continuum of genetic and environmental influences [6] whereby genes and environment can interact to determine sex [79].

The genetic mechanisms that act in highly conserved pathways that ultimately yield testes or ovaries are quite well characterised [5,10,11]. Yet, despite decades of research on ESD systems, and TSD in particular, the upstream mechanisms by which an external signal is transduced to determine sex remains unknown [12]. Recent research led to the hypothesis that the cellular sensor initiating ESD is controlled by the balance of redox regulation and calcium (Ca2+) signalling (CaRe) [13]. The CaRe hypothesis proposes a link between CaRe sensitive cellular signalling and the highly conserved epigenetic processes that have been implicated in thermolabile sex (TSD and temperature sex reversal) [12,1417]. The CaRe hypothesis posits that in ESD systems a change in intracellular Ca2+ (probably mediated by thermosensitive transient receptor potential TRP channels) and increased reactive oxygen species (ROS) levels caused by high temperatures, alter the CaRe status of the cell, triggering cellular signalling cascades that drive differential sex-specific expression of genes to determine sex. The CaRe hypothesis makes several testable predictions for how an environmental signal is captured and transduced by the gonadal cells to deliver a male or a female phenotype.

Species in which genes and environment both influence sex determination provide unique opportunities to directly compare the regulatory and developmental processes involved in sex determination. By early gonad differentiation directed by genotype and temperature, it is possible to assess predictions of the CaRe hypothesis. In our model species, the central bearded dragon (Pogona vitticeps), we can compare female development via thermal and genetic cues because extreme temperatures (>32°C) override the male sex-determining signal from the ZZ sex micro-chromosomes to feminise embryos [8,18]. This makes it possible to distinguish between the previously confounded effects of thermal stress and phenotypic sex by comparing gene expression throughout embryonic development in sex reversed ZZf females with genetic ZWf females.

We can explore the predictions of the CaRe model, namely that under sex-reversing conditions, we will see differential regulation of: 1) genes involved in responding to Ca2+ influx and signalling; 2) genes involved in antioxidant and/or oxidative stress responses; 3) genes with known thermosensitivity, such as heat shock proteins; 4) candidate TSD genes, such as CIRBP and Jumonji family genes; 5) signal transduction pathways such as the JAK-STAT and NF-ĸB pathways.

We compared gene expression profiles in P. vitticeps embryonic gonads at three developmental stages (6, 12 and 15) [19,20] for ZWf and ZZf eggs incubated at 28°C and 36°C respectively (Fig 1). This allowed us to compare drivers of sex determination and differentiation under genetic or thermal influence. We found that very different regulatory processes are involved in temperature-driven regulation compared to gene-driven regulation, although both lead to a conserved outcome (ovaries, Fig 2). We discovered dramatic changes in cellular calcium homeostasis in the gonads of ZZf individuals incubated at high sex reversing temperatures, which fulfill predictions of the CaRe hypothesis that this is the key driver of temperature induced feminization. We argue that differential expression of calcium channels, and subsequent alterations of the intracellular environment combined with increased ROS production encode, then transduce, the thermal signal into altered gene expression, ultimately triggering male to female sex reversal in P. vitticeps.

Fig. 1. Schematic representation of experimental design used in this study to compare the differences between genetic sex determination and temperature dependent sex determination.
Schematic representation of experimental design used in this study to compare the differences between genetic sex determination and temperature dependent sex determination.
(A) Summary of experiment showing how the parental crosses were designed, and how eggs were allocated and incubated. Eggs from sex reversed females (ZZf) were initially incubated at 28°C for 10 days, then were switched to 36°C. Eggs were sampled at the same three developmental stages (6, 12, and 15) based on [19,20]. At stage 6 the gonad is bipotential, at stage 12 the gonad is in the early stages of differentiation, and it completely differentiated by stage 15. Eggs from concordant females (ZWf) were incubated at 28°C and sampled at the same three developmental stages as the ZZf eggs. (B) PCA plots showing the first and second principal components of read count per gene between ZZf (red) and ZWf (blue) at each stage of development.
Fig. 2.
Schematic overview of gene-driven (blue) and temperature-driven (red) female developmental pathways in Pogona vitticeps. The pathways are initially different (from stages 6 to 12), but they ultimately converge on highly similar expression profiles when ovarian differentiation has occurred by stage 15. Both pathways are characterised by repression of a male signal, however this signal is stronger in temperature-driven females and appears to require ongoing repression when compared with the gene-driven females.

Results

Gene-driven female determination in ZWf embryos

Comparisons between stages in ZWf embryos (Figs 3B and S1, S1 Data) showed that many genes were differentially expressed between stages 6 and 12 (210 genes downregulated and 627 genes upregulated at stage 12), but few genes were differentially expressed between stages 12 and 15 (2 genes upregulated at stage 15).

Fig. 3.
(A) Expression (transcripts per million, TPM) ± SE of three genes differentially expressed at all three developmental stages between ZZf and ZWf, with KDM6B and CIRPB (outlined in red) having consistently higher expression in ZZf embryos, and GCA having higher expression in ZWf. (B) Bar graphs representing the number of differentially expressed genes in all comparisons between ZZf and ZWf, and between developmental stages. MA plots of this data are available in S1 Fig. Differentially expressed genes were determined as having P values ≤ 0.01 and log2-fold changes of 1, -1.

SOX9 and GADD45G, genes strongly associated with male development in mammals, were downregulated from stage 6 to stage 12, whereas various female related genes were upregulated, such as PGR, ESR2, CYP19A1, and CYP17A1. BMP7, a regulator of germ cell proliferation was upregulated at stage 12 [21], as were components of the NOTCH signalling pathway (JAG2, DLL3, DLL4), which are required for the suppression of Leydig cell differentiation [22,23]. SRD5A2, whose product catalyses the 5-α reduction of steroid hormones such as testosterone and progesterone, was also upregulated [24,25].

Notably, there was little differential expression between stages 12 and 15, suggesting that genetically driven ovarian development is complete by stage 12 (S1 Data).

Temperature-driven female determination in ZZf embryos

Differential expression analysis of temperature-driven female development in ZZf embryos revealed many genes are differentially expressed between stages 6 and 12 (297 downregulated and 511 upregulated at stage 12) and no genes are differentially expressed between stage 12 and 15 (Figs 3 and S1, S1 Data), suggesting completion of the ovarian development by stage 12 also in ZZf females.

Upregulation of FZD1, a receptor for Wnt family proteins required for female development, suggests the activity of female pathways in ZZf embryos [26]. As was seen for ZWf females, canonical NOTCH ligands DLL3 and DLL4 were upregulated from stage 6 to stage 12 in ZZf females. However, this did not coincide with upregulation of JAG ligands or NOTCH genes, and the GO term “negative regulation of NOTCH signalling” was enriched within the group of genes upregulated from stage 6 to 12 in ZZf females (S2 Data). Further, PDGFB, which is required for Leydig cell differentiation, was upregulated [27]. Together, this suggests that the NOTCH signalling pathway may not be activated, and Leydig cell recruitment is not strongly repressed at stage 12 in ZZf. Alternatively, the absence of NOTCH signalling may indicate an important transition from progenitor cells to differentiated gonadal cell types in the early stages of the developing ovary [28]. These apparent differences in NOTCH signalling between ZZf and ZWf embryos suggests that ovarian development has progressed further in ZWf females.

Interestingly, genes typically associated with male development show diverse regulation in ZZf embryos, with some being downregulated and some being upregulated from stage 6 to 12. These included WNT5a and SFRP2, which are both involved in testicular development in mice [29,30], NCOA4, which enhances activity of various hormone receptors, and exhibits high expression in testes in mice during development [31], and HSD17B3, which catalyses androstenedione to testosterone [32]. Unlike what was observed in ZWf embryos, SOX9 and GADD45G were not differentially expressed between stages 6 and 12 in ZZf embryos. TGFBR3L, which is required for Leydig cell function in mouse testis [33], and NR5A1, SOX4, and AMHR2 [3437] were also differentially expressed between stages 6 and 12 (S1 Data).

A suite of genes typically associated with female development were upregulated from stage 6 to 12 [38], for example, FOXL2, CYP17A1, RSPO1, and ESRRG. As was also observed in stage 12 ZWfs, ESR2, BMP7, CYP19A1, and PGR, were more highly expressed at stage 12 in ZZfs. Notably, CYP19A1 was much more strongly upregulated at stage 12 in ZZfs compared with stage 12 in ZWfs (S1 Data). The increase in sex specific genes was also reflected in enriched GO terms at stage 12, which included “hormone binding”, “steroid hormone receptor activity”, and “female sex determination” (S2 Data).

Ovarian maintenance in sex reversed ZZf females

The maintenance of female gene expression and ovarian development at stage 12 in ZZf females may be centrally mediated by STAT4 (Fig 4). As a member of the JAK-STAT pathway, STAT4 is transduced by various signals, including reactive oxygen species, to undergo phosphorylation and translocate from the cytoplasm to nucleus [3941]. At stage 12, STAT4 is upregulated, alongside PDGFB compared to stage 6 in ZZf females. PDGFB is known to activate STAT4 [40]. Various STAT4 target genes, notably AMHR2, NR5A1, EGR1, and KDM6B [40] are also upregulated at stage 12 (S1 Data). Consistent with this link is the observation that a member of the same gene family, STAT3, is implicated in TSD in Trachemys scripta [17].

Fig. 4. Hypothesized pathway for the maintenance of the ovarian phenotype in stage 12 sex reversed ZZf Pogona vitticeps.
Hypothesized pathway for the maintenance of the ovarian phenotype in stage 12 sex reversed ZZf <i>Pogona vitticeps</i>.
Given the upregulation of these genes, it is likely that reactive oxygen species induce the phosphorylation, and subsequent activation and nuclear translocation of STAT4, likely mediated by PDGFB. Once in the nucleus, STAT4 is able to bind to promoter regions of known target genes, NR5A1, AMHR2 and EGR2 to regulate their expression and promote ovarian development.

Several targets of STAT4 are upregulated at stage 12, including AMHR2 and NR5A1. Though typically associated with male development, AMHR2 and NR5A1 may also have roles in ovarian development. Although it is the primary receptor for AMH, AMHR2 exhibits considerable evolutionary flexibility is sometimes associated with ovarian development (reviewed by [36]). NR5A1 is also often associated with male development, as it positively regulates expression of AMH and SOX9 in mammals [42]. However, NR5A1 can also interact with FOXL2 and bind to CYP19A1 promoter to promote female development [43,44]. The upregulation of FOXL2 and CYP19A1, but not AMH or SOX9, suggests that NR5A1 is involved in the establishment of the ovarian pathway in ZZf females.

EGR1 positively regulates DMRT1 expression through promoter binding in Sertoli cells, but knock-out of this gene can also cause female infertility in mice [4547]. EGR1 is also associated with female development in birds, likely controlling the production of steroid hormones [34]. As was observed for NR5A1, DMRT1 was also lowly expressed, suggesting it is not activated by EGR1 in ZZf females.

One explanation for these expression trends is that male-associated genes are not strongly repressed at this stage during the sex reversal process, and that more prolonged exposure to the sex reversing temperature is required to firmly establish the female phenotype. However, we argue that the results more strongly suggest ROS-induced activation of STAT4, and subsequent phosphorylation and translocation, probably mediated by PDGFB, allows for the transcriptional activation of NR5A1, AMHR2, and EGR1, which in the temperature driven process of sex reversal in the ZZf embryos serve to maintain the ovarian phenotype.

Differential regulation of female developmental pathways

To better understand differences in ovarian developmental pathways, we compared gene expression of ZZf with ZWf embryos at each developmental stage. There are large gene expression differences between normal ZWf females and ZZf sex reversed females early in development, before the bipotential gonad differentiates into an ovary. These differences are most pronounced early in development and diminish as development progresses. Stage 6 had the largest number of differentially expressed genes (DEGs) (281 genes higher expressed in ZWf embryos, 423 genes higher expressed in ZZf), with fewer DEGs at stage 12 (51 genes upregulated in ZWf, 63 genes upregulated in ZZf), and fewest at stage 15 (1 gene upregulated in ZWf, 2 genes upregulated in ZZf) (Figs 3B and S1, S3 Data). This suggests that the sex reversed embryos start out on a male developmental trajectory, which they pursue beyond the thermal cue (3 days when the eggs were switched to high incubation temperatures, Fig 1), but by stage 12 development has been taken over by female genes.

Gene ontology (GO) enrichment analysis showed important differences between ZZf and ZWf at stage 6, and provides independent support for the role of calcium and redox regulation in ZZf females as proposed by the CaRe model (Fig 5A–5D, S4 Data). GO processes enriched in the gene set higher expressed in ZZf at stage 6 included “oxidation-reduction processes”, “cytosolic calcium ion transport”, and “cellular homeostasis” (Fig 5A, S4 Data). GO function enrichment also included several terms related to oxidoreductase activities, as well as “active transmembrane transporter activity” (Fig 5C, S4 Data). No such GO terms were enriched in the gene set higher expressed in ZWf. Instead, enriched GO terms included “anatomical structure development”, and “positive regulation of developmental growth” (Fig 5B and 5D, S4 Data).

Fig. 5.
(A) A subset of GO processes and (C) GO functions enriched in stage 6 ZZf embryos compared with ZWf. (B) A subset of GO processes and (D) GO functions enriched in stage 6 ZWf embryos compared with ZZf. Complete results of GO analysis for all developmental stages in ZZf and ZWf for enriched GO processes and functions is provided in S2 Data. Differentially expressed chromatin modifier (E) and cellular stress (F) genes in Pogona vitticeps at stage 6 comparing ZWf and ZZf females.

Genes involved in female sex differentiation were higher expressed at stage 6 in normal ZWf embryos compared to sex reversed ZZf embryos (S3 Data). These included FOXL2, ESR2, PGR, and GATA6 [48,49]. Higher expression of LHX9, a gene with a role in bipotential gonad formation in mammals and birds, was more highly expressed in ZWf embryos [42,5052]. Two genes with well described roles in male development, SOX4 and ALDH1A2 [5355], were also higher expressed in ZWf embryos, suggesting they have an as yet unknown function in the early establishment of the ovarian trajectory in P. vitticeps.

Taken together, these results further suggest that ZWf females are committed to the female pathway earlier than ZZf females. This is not surprising, since ZWf females possess sex chromosomes from fertilisation, whereas ZZf individuals have had only 3 days of exposure to a sex reversal inducing incubation temperature (Fig 1A). This data is the first to demonstrate a difference in timing of genetic signals between gene and temperature driven development in the same species.

Three genes were constantly differentially expressed between ZWf and ZZf embryos at all three developmental stages (Fig 3A). GCA (grancalcin) was upregulated in ZWf embryos, and KDM6B and CIRBP were upregulated in ZZf embryos at all developmental stages. GCA is a calcium binding protein commonly found in neutrophils and is associated with the Nf-κB pathway [56,57]. It has no known roles in sex determination, but its consistent upregulation in ZWf embryos compared to ZZf embryos suggests GCA is associated with gene driven ovarian development, at least in P. vitticeps.

Further analysis of gene expression trends using K-means clustering analysis [58] was used to investigate genes associated with female development, and to determine to what extent these genes are shared between ZZf and ZWf embryos (Fig 6, S5 Data). Clusters with upward trends reflect genes likely to be associated with female development, so clusters 1 and 4 in ZWf (ZWC1 and ZWC4), and clusters 1 and 2 in ZZf (ZZC1 and ZZC2), were explored in greater detail (Fig 6, S5 Data).

Fig. 6.
K-means clustering analysis on normalized counts per million for ZZf (A) and ZWf (B) across all developmental stages. The colour depicts the correlation score of each gene in the cluster, where numbers approaching one (red) have the strongest correlation. All gene lists produced for each cluster are provided in S5 Data.

ZWC4 and ZZC2 shared 374 genes. Enriched GO terms included “germ cell development” and “reproductive processes” (S6 Data), consistent with a link with female development. Genes identified included FIGLA, a gene known to regulate oocyte-specific genes in the female mammalian sex determination pathway [59], and STRA8 which controls entry of oocytes into meiosis. Intriguingly, the GO term “spermatid development” was also enriched, encompassing many genes with known roles in testes function, including ADAD1 and UBE2J1 [60,61]. This suggests that genes involved in male sex determination in mammals may have been co-opted for use in the ovarian pathway in reptiles, so their roles require further investigation in other vertebrate groups, particularly given the complex nature of gene expression in sperm cell types.

ZWC1 and ZZC1 shared 998 genes. ZZC1 has about 700 unique genes and ZWC1 about 500. GO analysis on shared genes between these clusters (n = 998) revealed enrichment terms such as “kinase binding” and “intracellular signal transduction” (S5 Data). Genes unique to ZZC1 included members of heat shock protein families (HSPB11, HSPA4, HSP90AB1, HSPH1, HSPB1, HSPD1), heterogenous ribonucleoprotein particles (HNRNPUL1), mitogen activated proteins (including MAPK1, MAPK9, MAP3K8), and chromatin remodelling genes (KDM2B, KDM1A, KDM5B, KDM3B). GO enrichment for genes unique to ZZC1 included “mitochondrion organisation”, “cellular localisation”, and “ion binding”, while GO enrichment for genes unique to ZWC1 included “regulation of hormone levels” and numerous signalling related functions (S6 Data).

Taken together, our results show that although the same ovarian phenotype is produced in genetic and temperature induced females, this end is achieved via different gene expression networks. This is most pronounced at stage 6, after which the extent of the differences decreases through development. This reflects canalisation of the gonadal fate to a shared outcome (ovaries, Fig 2).

Signature of hormonal and cellular stress in ZWf females

Previous work on P. vitticeps has shown a more than 50-fold upregulation of a hormonal stress response gene, POMC, in sex reversed adult females, leading to the suggestion that induction of sex reversal is in response to temperature stress, or that it is an inherently stressful event, the effects of which persist into adult life [14]. We therefore investigated the expression of stress related genes in ZZf and ZWf embryos.

We found considerable evidence that ZZf embryos experience oxidative stress, likely resulting from increased ROS production (discussed in detail below). However, contrary to our expectations, we found that ZWf embryos showed higher expression than ZZf of hormonal stress genes and pathways that have been hypothesized to be involved in sex reversal (Fig 5E and 5F, S3 Data). Genes upregulated in ZWf embryos compared to ZZf embryos included STAT1, a component of the JAK-STAT pathway, with several roles in stress responses [62], and MAP3K1 and MAPK8, which are typically involved in mediating stress-related signal transduction cascades [6365]. TERF2IP is also upregulated; this gene is involved in telomere length maintenance and transcription regulation [66]. When cytoplasmic, TERF2IP associates with the l-kappa-B-kinase (IKK) complex and promotes IKK-mediated phosphorylation of RELA/p65, activating the NF-κB pathway and increasing expression of its target genes [67]. Notably two members of the IKK complex, IKBKG (also known as NEMO) and PRKCI, which are involved in NF-κB induction, were also upregulated in ZWf embryos compared to ZZf embryos (Fig 4B), implying activation of the NF-κB pathway [68]. This pathway is typically associated with transducing external environmental signals to a cellular response [69,70], but also has diverse roles in sex determination in mammals, fish, and invertebrate models (reviewed by Castelli et al., 2019).

CRH, another gene upregulated at stage 6 in ZWf females compared with ZZf females (Fig 5B, S3 Data), is best known for its role as a neuropeptide synthesised in the brain in response to stresses that trigger the hypothalamic-pituitary-adrenal (HPA) axis [71,72]. The role of CRH production in the gonads, particularly in ovaries is currently poorly understood [7376]. High CRH expression in ZWf gonads is the first observation of this in reptiles. The role of the hormonal stress response during embryonic development, and its apparent discordance with results observed in adults in P. vitticeps requires further investigation [14].

Cellular signalling cascades driving sex reversal

Results of this study provide considerable corroborative support for the CaRe model, which proposes a central role for calcium and redox in sensing and transducing environmental signals to determine sex. Many of the genes and pathways predicted by the CaRe model to be involved in sex reversal were shown to be upregulated in ZZf embryos at stage 6 compared to ZWf embryos. We use the CaRe model as a framework to understand the roles of each signalling participant in their cellular context during the initiation of sex reversal (Fig 7). This interpretation is also independently supported by GO analysis, showing enrichment of expected terms, such as “cytosolic calcium ion transport” (S4 Data), as well as k-means clustering analysis (S4 and S5 Data).

Fig. 7.
Hypothesised cellular environment (A) of a ZZf gonad at stage 6 in Pogona vitticeps based on differential expression analysis (B) using the CaRe model as a framework [13]. We used this approach to understand the cellular context responsible for driving sex reversal in ZZf samples. This reveals that calcium signalling likely plays a very important role in mediating the temperature signal to determine sex. Influx of intracellular calcium is likely mediated primarily by TRPV2, and may also be influenced by KCNN1 and CACNB3. This influx appears to trigger significant changes in the cell to maintain calcium homeostasis. MCU, ATP2B1, CALR and TRICB all play a role in this process by sequestering calcium and pumping it back out of the cell, in which KCNN1 and CACNB3 may have a role. Calcium signalling molecules C2CD2, C2CDL2, and S100Z are likely responsible for encoding and translating the calcium signal leading to changes in gene transcription. Changes in gene expression are likely mediated primarily by the two major Polycomb Repressive Complexes, PRC1 and PRC2. Members of these two complexes (PCGF1, PCGF6, KDM6B, and JARID2) transcriptionally regulate genes by controlling methylation dynamics of their targets, the latter two of which have been previously implicated in sex reversal [14,15]. ATF5 may also play a role in gene regulation, and alternative splicing, which has been implicated in sex reversal [14] may be mediated by CLK4. High temperatures necessarily increase cellular metabolism, which in turn increases the amount of reactive oxygen species (ROS) produced by the mitochondria. ROS can cause cellular damage at high levels, so trigger an antioxidant response, which is observed here in the upregulation of MGST1, PRDX3, TXNDC11 and FOXO3. Also of note is the upregulation of CIRBP, which has numerous functions in response to diverse cellular stresses, and has been implicated in TSD.

Cluster 6 in ZZf (Fig 6A) shows genes whose expression decreases after stage 6, so is likely to include genes responsible for the initial response to temperature and initiation of sex reversal, whose continuing action is not required once the ovarian trajectory has been established. Consistent with this assumption, as well as with predictions from the CaRe model, the 4050 genes in this cluster were enriched for GO terms that included “oxidation-reduction process” and various oxidoreductase activities (S5 and S6 Data).

Calcium transport, signalling, and homeostasis

Our data suggest that exposure to high temperatures may cause a rapid increase in cytosolic Ca2+ concentrations, as calcium influx is probably mediated by the thermosensitive calcium channel, TRPV2 [77,78]. TRPV2 was upregulated in stage 6 ZZf embryos compared to ZWf embryos (Fig 6). Transient receptor potential (TRP) ion channels, including TRPV2, have previously been implicated in TSD in Alligator sinensis and A. mississippiensis, as well as the turtle Trachemys scripta [7982].

TRPV2 mediated Ca2+ influx may trigger a cascade of changes within the gonadal cells of ZZf females, which restore calcium homeostasis, critical to avoid apoptosis [83,84]. We observed evidence of such a homeostatic response, with the upregulation of seven genes involved in Ca2+ transport and sequestration in ZZf females compared to ZWf females at stage 6 (Fig 6). Specifically, MCU, ATP2B1, ATP2B4, together regulate calcium homeostasis through active transport of calcium into the mitochondria and into the extracellular space [8587]. KCNN1 and CACNB3 encode proteins required for the formation of plasma membrane channels controlling the passage of Ca2+ [8890]. CACNB3 and KCNN1 have well characterised roles in the nervous system, and excitable cell types in muscle, but their association with TSD in embryonic gonads is novel [88,89]. Evidence is also building for a broader role for voltage-gated calcium channels, including CACNB3, in orchestrating Ca2+ signalling and gene regulation [91]. We suggest that KCNN1 and CACNB3 in gonads of TSD species play roles in mediating the homeostatic response to elevated cytosolic Ca2+ concentrations, and are involved in the subsequent modulation of Ca2+ signalling pathways.

TRPC4, another TRP family gene was, upregulated in stage 12 compared to stage 6 in both ZZf and ZWf embryos. TRPC4 is expressed in mouse sperm and inhibited by progesterone [92,93] but has no known association with sex determination. TRPC4 belongs to the TRPC superfamily, which all conduct calcium ions into the cell, typically through phospholipase C and calmodulin signalling pathways, G-protein-coupled receptors, and receptor tyrosine kinases [94,95]. Notably, PLCL2 a phospholipase gene, together with calmodulin genes CALM1 and CAMKK1, were upregulated alongside TRPC4 from stage 6 to stage 12 in ZZf embryos but not in ZWf embryos (S3 Data). Given TRPC4 is upregulated from stage 6 to 12 in both ZZf and ZWf females, it may play a more conserved role in ovarian development in P. vitticeps.

Several genes with functions in calcium metabolism were upregulated in stage 6 ZZf embryos compared to stage 6 ZWf embryos. CALR encodes a multifunctional protein that acts as a calcium binding storage protein in the lumen of the endoplasmic reticulum, so is also important for regulating Ca2+ homeostasis [84,96,97]. CALR is also present in the nucleus, where it may play a role in regulation of transcription factors, notably by interacting with DNA-binding domains of glucocorticoid and hormone receptors, inhibiting the action of androgens and retinoic acid [97101]. TMEM38B (commonly known as TRICB) is also found on the endoplasmic and sarcoplasmic reticula, where it is responsible for regulating the release of Ca2+ stores in response to changes in intracellular conditions [102].

MCU, ATP2B1, ATP2B4, KCNN1, CACNB3, CALR, and TMEM38B have no known roles in vertebrate sex determination, so their association with sex reversal in P. vitticeps is new. This upregulation during the early stage of sex reversal suggests that they are upstream modulators involved in the transduction of environmental cues that trigger sex determination cascades, which is consistent with predictions made by the CaRe hypothesis.

We hypothesize that intracellular Ca2+ increases in stage 6 ZZf gonads, and further observe that Ca2+ signalling related genes are also upregulated in ZZf females compared to stage 6 ZWf females (Fig 7). C2CD2 and C2CD2L are both thought to be involved in Ca2+ signalling, although there is no functional information about these genes. Of note is the significant upregulation of S100Z, which is a member of a large group of EF-hand Ca2+ binding proteins that play a role in mediating Ca2+ signalling [103]. The EF-hand domain is responsible for binding Ca2+, allowing proteins like that encoded by S100Z to ‘decode’ the Ca2+ biochemical signal and translate this to various targets involved in many cellular functions including Ca2+ buffering, transport, and enzyme activation [104,105]. PLCB1 also contains an EF-hand binding domain and behaves similarly, being activated by many extracellular stimuli and effecting numerous signalling cascades. It can translocate to the plasma membrane and nucleus, and release Ca2+ from intracellular stores [106]. Some Ca2+ related genes (GCA and CALM1) are also upregulated in ZWf embryos, but make only a small proportion of the overall response in differential gene expression (S3 Data).

Oxidative stress in response to high temperatures

The upregulation of antioxidant genes in ZZf compared to ZWf embryos suggests that the gonadal cells in the ZZf embryos are in a state of oxidative stress (Fig 7). As was proposed by the CaRe model, we see results consistent with the prediction that high incubation temperatures increase metabolism, which increases the production of reactive oxygen species (ROS) by the mitochondria, resulting in oxidative stress [13]. ROS are required for proper cellular function, but above an optimal threshold, they can cause cellular damage [107,108]. Crossing this threshold launches the antioxidant response, which causes the upregulation of antioxidant genes to produce protein products capable of neutralising ROS [109,110]. We observed upregulation of redox related genes, specifically of TXNDC11, PRDX3, MGST1 in ZZf embryos compared to ZWf embryos at stage 6. Also upregulated was FOXO3, which plays a role in oxidative stress responses, typically by mediating pro-apoptotic cascades [111,112]. Importantly, antioxidants play other cellular roles besides neutralisation of ROS. One of these is the alteration of cysteine resides through a process known as S-glutathionylation [113].

Various redox related genes were downregulated from stage 6 to 12 in ZZf embryos but not in ZW embryos, including GLRX and PRDX3 [114], as well as numerous genes involved in ROS induced DNA damage repair; LIG4, ENDOD1, and HERC2 [115]. This indicates a need for expression of these genes specifically in ZZf embryos in early stages that ceases in transition to stage 12. STAT4, a member of the ROS-induced JAK-STAT pathway [39], and DDIT4, which is involved stress responses to DNA damage [116], were both upregulated from stage 6 to stage 12 in ZZf embryos.

The vertebrate antioxidant response is typically initiated by NRF2, but we observed no differential expression of NRF2, only upregulation of some of its known targets in ZZf embryos [117]. This may mean that the action of NRF2 is depends more on its translocation from the cytoplasm to the nucleus to modulate transcription of target genes, a process that does not necessarily rely on increased expression of NRF2 [117]. Alternatively, NRF2 upregulation may have occurred prior to sampling.

Oxidative stress has previously been proposed to have a role in TSD, based on the upregulation of genes involved in oxidative stress response. One of these genes, UCP2, was upregulated at high male producing temperatures in A. mississippiensis [82]. UCP2, and others genes involved in oxidative stress responses, were also implicated in UV induced masculinisation in larvae of a thermosensitive fish species (Chirostoma estor) [118]. Notably, we found that UCP2 was upregulated between stages 6 to 12 in ZZf P. vitticeps embryos, suggesting a sustained response to thermal stress in the mitochondria (S3 Data).

Temperature response and cellular triage

We also observed upregulation of genes involved in response to more generalised environmental stress in ZZf compared to ZWf embryos, as expected since the embryos exposed to high temperature were experiencing a state of thermal stress (Fig 7). Notably, CIRBP a promising candidate for regulation of sex determination under thermal influence, is approximately 10-fold upregulated in ZZf compared to ZWf (S3 Data). CIRBP has a highly conserved role in generalised stress responses [119]. It has been suggested to be a putative sex determining gene in the TSD turtle Chelydra serpentina [120], and is differentially expressed at different incubation temperatures in Alligator sinensis [81]. We also observed the upregulation of CLK4 in ZZf compared to ZWf embryos, a gene that has been recently shown to be inherently thermosensitive, and to regulate splicing of temperature specific CIRBP isoforms [121].

We found that ATF5 is upregulated in ZZf embryos compared to ZWf embryos (Fig 7). ATF5 has diverse roles in stimulating gene expression or repression through binding of DNA regulatory elements. It is broadly involved in cell specific regulation of proliferation and differentiation, and may also be critical for activating the mitochondrial unfolded protein response [122]. This gene is induced in response to various external stressors, and is activated via phosphorylation by eukaryotic translation initiation factors, two of which (EIF1 and EIF4A2 [123] are also upregulated in ZZf embryos compared to ZWf embryos.

Though not well studied in the context of sex determination, heat shock factors and proteins have been implicated in female sex determination in mammals and fish, and may also play a conserved role in the ovarian pathway in P. vitticeps [79,124127]. Surprisingly, only one gene associated with canonical heat shock response (HSP40, also known as DNAJC28) was differentially expressed following exposure to high temperature in stage 6 ZZf females compared to ZWf embryos (S3 Data). This could mean either that a heat shock response occurs prior to sampling, or that P. vitticeps uses different mechanisms to cope with heat shock.

Chromatin remodelling

We observed upregulation of several components of two major chromatin remodelling complexes, polycomb repressive complexes PRC1 and PRC2, in both the genotype-directed ZWf and the temperature-directed ZZf female pathways in P. vitticeps (Fig 7). Chromatin modifier genes KDM6B and JARID2 are involved in regulation of gene expression during embryonic development and epigenetic modifications in response to environmental stimulus [128,129]. JARID2 and KDM6B were both upregulated in ZZf embryos compared to ZWf embryos in stages 6 and 12, and KDM6B was also upregulated at stage 15. These genes have recently been implicated in two TSD species (Alligator mississippiensis, and Trachemys scripta) and temperature sex reversed adult Pogona vitticeps [14,15,79,82].

We also found that two other members of the PRC1 complex, PCGF6 and PCGF1, were upregulated in ZZf embryos at stage 6 compared to ZWf embryos (Fig 7). PCGF6 is part of the non-canonical PRC1 complex (ncPRC1) that mediates histone H2A mono-ubiquitination at K119 (H2AK119ub) [130,131]. PCGF6 acts a master regulator for maintaining stem cell identity during embryonic development [132], and is known to bind to promoters of germ cell genes in developing mice [130]. PCFF1 exhibits similar functions by ensuring the proper differentiation of embryonic stem cells [133]. The ncPRC1 complex also promotes downstream recruitment of PRC2 and H3K27me3, so that complex synergistic interactions between PRC1 (both canonical and non-canonical) and PRC2 can occur [134,135].

We found that other components of both PRC1 and PRC2 complexes were also upregulated in ZWf embryos compared to ZZf embryos (Fig 5A, S3 Data). A member of the canonical PRC1 complex, PCGF2 (also known as MEL18), was upregulated in ZWf embryos compared to ZZf embryos [134]. This gene has previously been implicated in temperature induced male development in Dicentrachus labrax [136], and is required for coordinating the timing of sexual differential in female primordial germ cells in mammals [137]. KDM1A, a histone demethylase that is required for balancing cell differentiation and self-renewal [138], was upregulated in ZWf embryos compared to ZZf embryos. CHMP1A was upregulated in ZWf, and is likely to be involved in chromosome condensation, as well as targeting PcG proteins to regions with condensed chromatin [139].

Thus, we conclude that the initiation of sex reversal in ZZf P. vitticeps involves a complex cascade of cellular changes initiated by temperature. Our data are consistent with the predictions of the CaRe hypothesis that high temperatures are sensed by the cell via TRP channels, which causes an increase in intracellular increase of Ca2+. Coincident with this is an increase of ROS production in the mitochondria that causes a state of oxidative stress. Together, Ca2+ and ROS alter the CaRe status of the cell, trigger a suite of alternations in gene expression including chromatin remodelling, which drives sex reversal (Fig 7).

Discussion

We used the unique sex characteristics of our model reptile species, Pogona vitticeps, which determines sex genetically but sex reverses at high temperature, to assess predictions of the CaRe hypothesis [13]. By sequencing isolated embryonic gonads, we provide the first data to represent a suite of key developmental stages with comparable tissue types, and will be a valuable resource for this reptilian model system. There are few transcriptomes of GSD reptiles during embryonic development; the only dataset available prior to this study was a preliminary study of the spiny softshell turtle, Apalone spinifera [126], which was inadequate for the inter-stage comparisons required to explore genetic drivers of gonad differentiation.

Our analysis of expression data during embryogenesis of normal ZWf females and temperature sex reversed ZZf females revealed for the first time differences in gene-driven and temperature-driven female development in a single species. Early in development, prior to gonad differentiation, the initiation of the sex reversal trajectory differs from the genetic female pathway both in the timing and genes involved. As development proceeds, differences in expression patterns become less until the pathways converge on a conserved developmental outcome (ovaries). Our ability to compare two female types in P. vitticeps allowed us to avoid previously intractable confounding factors such as sex or species-specific differences, which provided unprecedented insight into parallel female pathways. We have provided new insight to the conserved evolutionary origins of the labile networks governing environmentally sensitive sex determination pathways. We have identified a suite of candidate genes, which now provide the necessary foundation for functional experiments in the future. This could include pharmacological manipulation of calcium signalling through alteration of intracellular calcium flux and concentration, such as interference with the TRPV4 channel, or via use of calcium chelators and ionophores, in an organ culture system [17,80]. Ongoing development of resources for P. vitticeps as an emerging model organism may also allow for RNA interference or gene editing experiments, whereby knock-down or knock-out of candidate genes like CIRBP, JARID2, and KDM6B, could be used to demonstrate their roles in sex reversal [15].

The maintenance of ovarian differentiation seems to require the operation of different pathways in gene and temperature driven female development. This may involve a pathway centrally mediated by STAT4 in sex reversed P. vitticeps, which has not been previously described, so requires additional confirmation with functional experiments. It will be interesting to determine if a role for these genes occurs in other species. Another STAT family gene, STAT3, has recently been demonstrated to play a critical role in the phosphorylation of KDM6B and subsequent demethylation of the DMRT1 promoter required for male development in T. scripta [17]. The involvement of different genes in the same family is intriguing in its implications; while different genes may be co-opted, natural selection may favour gene families with conserved functions even between evolutionarily disparate lineages.

Our data provided insight into the molecular landscape of the cell required to initiate temperature induced sex reversal. This is the first dataset to capture temperature-induced sex reversal in a reptile, and remarkably we have simultaneously implicated all functional candidates that have previously been identified to be involved in TSD across a range of other species (Table 1). Our results also identified novel genes involved with thermosensitive sex determination, and provide corroborative evidence for the CaRe hypothesis [13]. Importantly, our work highlights avenues for future studies to conduct functional experiments to definitively identify the genes and pathways implicated here in sex reversal. Observation and manipulation of intracellular calcium concentrations, as has been conducted in T. scripta [17], will also be crucial for fully understanding the role of calcium signalling in sex reversal.

Tab. 1. All genes, full gene names, functional categories and associations with either gene (ZWf) or temperature driven (ZZf) female development mentioned in the paper.
All genes, full gene names, functional categories and associations with either gene (ZWf) or temperature driven (ZZf) female development mentioned in the paper.
NA denotes a gene that was mentioned, but was not differentially expressed. Genes with an asterisk are those that have previously been implicated in thermosensitive sex determination cascades, either in Pogona vitticeps, or in another reptile species.

Our results highlight the complexity of initiating thermolabile systems. Indeed, it has been suggested that thermolabile sex determination involves system-wide displacement of gene regulation with multiple genes and gene products responding to temperature leading to the production of one sex or the other–a parliamentary system of sex determination [151]. We take an intermediate position, arguing for a central role for Calcium-Redox balance as the proximal cellular sensor for temperature, but interacting with other required thermosensitive genes or gene products (e.g. CLK4) to influence ubiquitous signalling pathways and downstream splicing regulation, epigenetic modification and sex gene expression. The level of interaction between each thermosensitive element remains to be explored. For example, if temperature can be sensed by both TRPV2 and CLK4, are both required to initiate sex reversal, or is the signal from only one sufficient? This raises the possibility that no single proximal sensor of the environmental exists, but that several thermosensitive elements early in development must come together to orchestrate alterations in gene expression.

It has been suggested that the products of TRP family genes act as mediators between the temperature signal and a cellular response through Ca2+ signalling and subsequent modulation of downstream gene targets [8082]. Notably, different TRP channels are implicated in two alligator species; TRVP4 in A. mississippiensis, but TRPV2, TRPC6, and TRPM6 in A. sinsensis. In T. scripta, TRPC3 and TRPV6 are upregulated at male producing temperatures (26°C), while TRPM4 and TRPV2 are upregulated at female producing temperatures (31°C), as is the case for TRPV2 in P. vitticeps [79]. The diversity of TRP channels recruited for roles in environmental sex determination hints at considerable evolutionary flexibility, perhaps the result of repeated and independent co-option of these channels in TSD species. As may be the case for STAT family genes, the evolution of environmentally sensitive sex determination pathways may involve the use of different genes within gene families that have conserved functions.

Our data also highlights the importance of chromatin remodelling genes in sex reversal in P. vitticeps. KDM6B and JARID2 have been previously implicated sex differentiation in adult P. vitticeps [14], embryonic T. scripta [15,17] and embryonic A. mississippiensis [14]. Sex-specific intron retention was observed in TSD alligators and turtle, and was exclusively associated with sex reversal in adult P. vitticeps [14]. Subsequently, knockdown of KDM6B in T. scripta caused male to female sex reversal by removing methylation marks on the promoter of DMRT1, a gene critical in the male sex determination pathway [15]. KDM6B and JARID2 have also been associated with TSD in another turtle species (Chrysemys picta) [126], female to male sex change in the bluehead wrasse, Thalassoma bifasciatum [140], and thermal responses in the European bass, Dicentrarchus labrax [141].

It is currently unknown if the unique splicing events in KDM6B and JARID2 in adult sex reversed P. vitticeps that cause intron retention and presumed gene inactivation, also occur in embryos. Given the high expression of these genes during embryonic development at sex reversing temperatures, it would be surprising if this pattern was observed. We also show a significant role for CIRBP as the only other gene, alongside KDM6B, to be consistently upregulated during sex reversal in all developmental stages assessed. CIRBP is a mRNA chaperone, which could be required to stabilise transcripts of crucial sex specific genes during oxidative, cellular and/or thermal stress. It has been proposed as a novel TSD candidate gene in the turtle, Chelydra serpentina [120]. This gene remains a promising candidate for mediating thermosensitive responses in TSD more broadly, and its role needs to be explored in more detail.

Conclusions

The alternative female pathways in P. vitticeps demonstrates that there is inherent flexibility in sex determination cascades even within the same species. This is consistent with the idea that, provided a functional gonad is produced, considerable variation in sex determining and differentiation processes at the early stages of development is tolerated under natural selection [151]. Perhaps this makes the astonishing variability in sex determination between diverse species less surprising. Our findings provide novel insights, and are a critical foundation for future studies of the mechanisms by which temperature determines sex.

Materials and methods

Ethics statement

All procedures were conducted in accordance with approved animal ethics protocols from the University of Canberra Animal Ethics Committee (AEC 17–17).

Animal breeding and egg incubations

Eggs were obtained during the 2017–18 breeding season from the research breeding colony at the University of Canberra. Breeding groups comprised three sex reversed females (ZZf) to one male (ZZ), and three concordant females (ZWf) to two males (Fig 1). Paternity was confirmed by SNP genotyping (S2 Fig). Females were allowed to lay naturally, and eggs were collected at lay or within two hours of lay. Eggs were inspected for viability as indicated by presence of vasculature in the egg, and viable eggs were incubated in temperature-controlled incubators (±1°C) on damp vermiculite (4 parts water to 5 parts vermiculate by weight). Clutches from sex reversed females (that is, ZZf x ZZm crosses) comprised eggs with only ZZ genotypes. These were initially incubated at 28°C (male producing temperature, MPT) to entrain and synchronise development. After 10 d of incubation, half of the eggs selected at random from each clutch was shifted to 36°C (female producing temperature, FPT). Clutches from ZWf x ZZm crosses were incubated at 28°C throughout the incubation period (Fig 1A). Sample sizes are given in Fig 1 and S7 Data.

Embryo sampling and genotyping

Eggs from both temperatures were sampled at times corresponding to three developmental stages (6, 12 and 15) [20], taking into account the differing developmental rates between 28°C and 36°C. These stages equate to the bipotential gonad, recently differentiated gonad, and differentiated gonad respectively [19]. Embryos were euthanized by intracranial injection of 0.1 ml sodium pentobarbitone (60mg/ml in isotonic saline). Individual gonads were dissected from the mesonephros under a dissection microscope and snap frozen in liquid nitrogen. Isolation of the gonad from the surrounding mesonephros was considered essential for studying transcriptional profiles within the gonad. Embryos from three different ZZf x ZZm clutches from each treatment class (temperature x stage) were selected for sequencing, and randomized across sequence runs to avoid batch effects. Embryos from concordant ZWf x ZZm crosses potentially yield both ZW and ZZ eggs, so these were genotyped using previously established protocols [8,20]. Briefly, this involved obtaining a blood sample from the vasculature on the inside of the eggshell on a FTA Elute micro card (Whatman). DNA was extracted from the card following the manufacturer protocols, and PCR was used to amplify a W specific region [8] so allowing the identification of ZW and ZZ samples.

RNA extraction and sequencing

RNA from isolated gonad samples was extracted in randomized batches using the Qiagen RNeasy Micro Kit (Cat. No. 74004) according to the manufacturer protocols. RNA was eluted in 14 μl of RNAase free water and frozen at -80°C prior to sequencing. Sequencing libraries were prepared in randomized batches using 50 ng RNA input and the Roche NimbleGen KAPA Stranded mRNA-Seq Kit (Cat. No. KK8420). Nine randomly selected samples were sequenced per lane using the Illumina HiSeq 2500 system, and 25 million read-pairs per sample were obtained on average. Read lengths of 2 x 150 bp were used. All samples were sequenced at the Kinghorn Centre for Clinical Genomics (Garvan Institute of Medical Research, Sydney). All sample RNA and library DNA was quantified using a Qubit Instrument (ThermoFisher Scientific, Scoresby, Australia), with fragment size and quality assessed using a Bioanalyzer (Agilent Technologies, Mulgrave, Australia).

Gene expression profiling

Paired-end RNA-seq libraries (.fastq format) were trimmed using trim_galore with default parameters (v0.4.1; https://www.bioinformatics.babraham.ac.uk/projects/trim_galore/, last access 21-Apr-2020). Trimmed reads were aligned to the Pogona vitticeps NCBI reference genome (pvi1.1, GenBank GCA_900067755.1; [142]) using STAR (v2.5.3; [143]), with splice-aware alignment guided by the accompanying NCBI gene (Pvi1.1) annotation (.gtf format). Likely PCR duplicates and non-unique alignments were removed using samtools (v1.5; [144]). Gene expression counts and normalised expression values (reported in TPM) were determined using RSEM (rsem-calculate-expression; v1.3.1; [145]).

Identification of non-sex reversed specimens

Normalised transcripts per million (TPM) for a panel of sex-specific genes (SOX9, AMH, DMRT1, FOXL2, CYP19A1, CYP17A1) were inspected across the three stages to identify if any samples showed aberrant expression patterns. This approach was also used to determine if any of the stage 12 and 15 samples from the 36°C treatment had not undergone sex reversal by comparing expression levels between ZWf and ZZf embryos; the rate of sex reversal is 96% at 36°C [8] (S3 Fig, S8 Data). The five samples from clutch 9 exhibited significantly higher expression values for SOX9, AMH, and DMRT1 and represented clear outliers. This was also supported by multidimensional scaling (MDS) plots, so the decision was made to regard the five samples from clutch 9 as aberrant and exclude them from subsequent analyses (S3 and S4 Figs, S9 Data). Any ZZf samples with male-like gene expression patterns (high expression for male-specific genes, and low expression for female-specific genes) were considered to have not been reversed (sex reversal is not 100% at 36°C) and were removed (two stage 15 samples).

Differential expression analysis

Differential expression analysis of ZZf and ZWf transcripts was conducted on raw counts using the EdgeR package (Bioconductor v 3.9 [146]) in R (v 1.2.1335, [147]), following standard procedures outlined in the EdgeR users guide [146,148]. Lowly expressed genes, which was applied to genes with fewer than ten counts across three samples, were removed from the raw counts (19,285 genes) so that the total number of genes retained was 17,075. Following conversion to a DGElist object in EdgeR, raw counts were normalised using the upper-quartile method (calcNormFactors function) [149]. Estimates for common negative binomial dispersion parameters were generated (estimateGLMCommonDisp function) [148], followed by generation of empirical Bayes dispersion estimates for each gene (estimateGLMTagwiseDisp function) [148,150]. A quasi-likelihood binomial generalised log-linear model was fitted (glmQLFit function) and the glmQFTest function was used to compare contrasts within the design matrix [151155]. A P-value cut-off of 0.01 and a log2-fold change threshold of 1 or -1 was applied to all contrasts (topTags function) [151]. Contrasts were used to assess differential expression between ZZf and ZWf samples across each developmental stage. Raw count (S10 Data) and expression files (S11 Data) from this analysis are supplied.

Gene ontology (GO) analysis was conducted for each set of differentially expressed genes using GOrilla [156,157]. The filtered count data file (17,075 genes) was used for the background gene set at a P-value threshold of 10−3.

K-means clustering analysis

K-means clustering analysis was performed on normalised counts per million extracted from the DGElist object produced by the initial process of the DGE analysis using edgeR (see above). Counts for each gene were averaged for each treatment group, and the number of clusters was selected using the sum of squared error approach, which was further validated by checking that each cluster centroid was poorly correlated with all other cluster centroids (maximum correlation 0.703 in ZWf clusters, and 0.65 in ZZf clusters). A total of 6 clusters was chosen, and clustering analysis was conducted using the kmeans function in R package stats v3.6.2. Resultant gene lists were sorted by unique and shared genes between clusters with similar trends between ZWf and ZZf (cluster 1 in ZWf, cluster 3 in ZZf, and cluster 3 in ZWf and cluster 5 in ZZf). Both unique and shared genes from each cluster and pairs of clusters (cluster 1 and 3, and clusters 3 and 5) were then analysed for gene ontology (GO) enrichment using GOrilla [156,157]. The filtered count data file (17,075 genes) was used for the background gene set at a P-value threshold of 10−3.

Supporting information

S1 Data [xlsx]
Differentially expressed genes between developmental stages 6 and 12, and 12 and 15 for ZZf and ZWf females generated from EdgeR’s “topTags” function.

S2 Data [xlsx]
Gene ontology (GO) enrichment for process and function for the stage 6 and 12 comparison for ZZf and ZWf ().

S3 Data [xlsx]
Differentially expressed genes between ZZf and ZWf for each developmental stage generated from EdgeR’s “topTags” function.

S4 Data [xlsx]
Gene ontology (GO) enrichment for process and function generated from GOrilla [,] at a significance threshold of ≤ 0.05 for differentially expressed genes at stages 6 and 12 for ZZf and ZWf samples ().

S5 Data [xlsx]
Gene list outputs for K-means clustering analysis (n = 6), and comparative information for matched clusters between ZZf and ZWf (ZZC1 and ZWC1, and ZZC2 and ZWC4) including genes that are unique and shared between each cluster.

S6 Data [xlsx]
Gene ontology (GO) process and function enrichment for genes in ZZf C1, and genes shared between ZZC5 and ZWC3 generated from GOrilla [,] at a significance threshold of ≤ 0.05.

S7 Data [xlsx]
Summary of all embryonic gonad samples sequenced for this study, including incubation temperature, genotype, parental cross, developmental stage, and clutch for each sample.

S8 Data [tpm]
Outputs from pairwise T-tests conducted between stage 15 (n = 7) and stage 15 ZZf samples suspected not to have undergone sex reversal (n = 2).

S9 Data [tpm]
Outputs from one way analysis of variance (ANOVAs) between all clutches (clutch 1, 2, 3, 6, and 9) across each developmental stage (6, 12 and 15) for a panel of sex specific genes (, and ) to determine whether clutch 9 exhibits aberrant expression levels.

S10 Data [csv]
Raw counts for all samples (n = 39) for all genes (n = 19,284) prior to any filtering or sample removal.

S11 Data [csv]
Raw expression values (TPM, transcripts per million) all samples (n = 39) for all genes (n = 19,284) prior to any filtering or sample removal.

S1 Fig [a]

S2 Fig [tif]
Network analysis of parental and offspring SNPs to confirm paternity of clutches used in this experiment.

S3 Fig [a]

S4 Fig [a]
Principal components analysis (PCA) plots performed on normalised counts per million for filtered genes following the EdgeR pipeline described in the materials and methods section.


Zdroje

1. Koopman P, Gubbay J, Vivian N, Goodfellow P, Lovell-Badge R. Male development of chromosomally female mice transgenic for Sry. Nature. 1991;351:117–21. doi: 10.1038/351117a0 2030730

2. Sinclair AH, Berta P, Palmer MS, Hawkins JR, Griffiths BL, Smith MJ, et al. A gene from the human sex-determining region encodes a protein with homology to a conserved DNA-binding motif. Nature. 1990;346:240–4. doi: 10.1038/346240a0 1695712

3. Smith CA, Roeszler KN, Ohnesorg T, Cummins DM, Farlie PG, Doran TJ, et al. The avian Z-linked gene DMRT1 is required for male sex determination in the chicken. Nature. 2009;461:267–71. doi: 10.1038/nature08298 19710650

4. Barske LA, Capel B. Blurring the edges in vertebrate sex determination. Curr Opin Genet Dev. 2008;18:499–505 doi: 10.1016/j.gde.2008.11.004 19152784

5. Capel B. Vertebrate sex determination: evolutionary plasticity of a fundamental switch. Nat Rev Genet. 2017;18:675–89. doi: 10.1038/nrg.2017.60 28804140

6. Sarre SD, Georges A, Quinn A. The ends of a continuum: Genetic and temperature-dependent sex determination in reptiles. BioEssays. 2004;26:639–45. doi: 10.1002/bies.20050 15170861

7. Holleley CE, Sarre SD, O’Meally D, Georges A. Sex reversal in reptiles: Reproductive oddity or powerful driver of evolutionary change? Sex Dev. 2016;10:279–87. doi: 10.1159/000450972 27794577

8. Holleley CE, O’Meally D, Sarre SD, Graves JAM, Ezaz T, Matsubara K, et al. Sex reversal triggers the rapid transition from genetic to temperature-dependent sex. Nature. 2015;523:79–82. doi: 10.1038/nature14574 26135451

9. Radder RS, Quinn AE, Georges A, Sarre SD, Shine R, Quinn AE, et al. Genetic evidence for co-occurrence of chromosomal and thermal sex-determining systems in a lizard. Biol Lett. 2008;4:176–8. doi: 10.1098/rsbl.2007.0583 18089519

10. Bachtrog D, Mank JE, Peichel CL, Kirkpatrick M, Otto SP, Ashman T-L, et al. Sex Determination: Why So Many Ways of Doing It? PLOS Biol. 2014;12:e1001899. doi: 10.1371/journal.pbio.1001899 24983465

11. Herpin A, Schartl M. Plasticity of gene-regulatory networks controlling sex determination: Of masters, slaves, usual suspects, newcomers, and usurpators. EMBO Rep. 2015;16:1260–74. doi: 10.15252/embr.201540667 26358957

12. Singh SK, Das D, Rhen T. Embryonic temperature programs phenotype in reptiles. Front Physiol. 2020; doi: 10.3389/fphys.2020.00035 32082193

13. Castelli MA, Whiteley SL, Georges A, Holleley CE. Cellular calcium and redox regulation: The mediator of vertebrate environmental sex determination? Biol Rev. 2020;95:680–95. doi: 10.1111/brv.12582 32027076

14. Deveson IW, Holleley CE, Blackburn J, Marshall Graves JAM, Mattick JS, Waters PD, et al. Differential intron retention in Jumonji chromatin modifier genes is implicated in reptile temperature-dependent sex determination. Sci Adv. 2017;3:e1700731. doi: 10.1126/sciadv.1700731 28630932

15. Ge C, Ye J, Weber C, Sun W, Zhang H, Zhou Y, et al. The histone demethylase KDM6B regulates temperature-dependent sex determination in a turtle species. Science. 2018;360:645–8. doi: 10.1126/science.aap8328 29748283

16. Georges A, Holleley CE. How does temperature determine sex? Science. 2018;360:601–2. doi: 10.1126/science.aat5993 29748270

17. Weber C, Zhou Y, Lee J, Looger L, Qian G, Ge C, et al. Temperature-dependent sex determination is mediated by pSTAT3 repression of Kdm6b. Science. 2020;3:303–6. doi: 10.1126/science.aaz4165 32299951

18. Quinn AE, Georges A, Sarre SD, Guarino F, Ezaz T, Graves JAM. Temperature sex reversal implies sex gene dosage in a reptile. Science. 2007;316:411. doi: 10.1126/science.1135925 17446395

19. Whiteley SL, Weisbecker V, Georges A, Gauthier ARG, Whitehead DL, Holleley CE. Developmental asynchrony and antagonism of sex determination pathways in a lizard with temperature-induced sex reversal. Sci Rep. 2018;8:1–9. doi: 10.1038/s41598-017-17765-5 29311619

20. Whiteley SL, Holleley CE, Ruscoe WA, Castelli MA, Whitehead DL, Lei J, et al. Sex determination mode does not affect body or genital development of the central bearded dragon (Pogona vitticeps). Evodevo. 2017; doi: 10.1186/s13227-017-0087-5 29225770

21. Ross A, Munger S, Capel B. Bmp7 regulates germ cell proliferation in mouse fetal gonads. Sex Dev. 2007;1:127–37. doi: 10.1159/000100034 18391523

22. Windley SP, Wilhelm D. Signaling pathways involved in mammalian sex determination and gonad development. Sex Dev. 2016;9:297–315.

23. Tang H, Brennan J, Karl J, Hamada Y, Raetzman L, Capel B. Notch signaling maintains Leydig progenitor cells in the mouse testis. Development. 2008;135:3745–53. doi: 10.1242/dev.024786 18927153

24. Krone N, Hanley NA, Arlt W. Age-specific changes in sex steroid biosynthesis and sex development. Best Pract Res Clin Endocrinol Metab. 2007;21:393–401. doi: 10.1016/j.beem.2007.06.001 17875487

25. Russell DW, Wilson JD. Steroid 5a-Rreducatse: Two genes/two enzymes. Annu Rev Biochem. 1994;63:25–61. doi: 10.1146/annurev.bi.63.070194.000325 7979239

26. Eid W, Opitz L, Biason-Lauber A. Genome-wide identification of CBX2 targets: Insights in the human sex development network. Mol Endocrinol. 2015;29:247–57. doi: 10.1210/me.2014-1339 25569159

27. Gnessi L, Basciani S, Mariani S, Arizzi M, Spera G, Wang C, et al. Leydig cell loss and spermatogenic arrest in platelet-derived growth factor (PDGF)-A-deficient mice. J Cell Biol. 2000;149:1019–25. doi: 10.1083/jcb.149.5.1019 10831606

28. Schmahl J, Rizzolo K, Soriano P. The PDGF signaling pathway controls multiple steroid-producing lineages. Genes Dev. 2008;22:3255–67. doi: 10.1101/gad.1723908 19056881

29. Chawengsaksophak K, Svingen T, Ng ET, Epp T, Spiller CM, Clark C, et al. Loss of Wnt5a disrupts primordial germ cell migration and male sexual development in mice. Biol Reprod. 2011;86:1–12.

30. Warr N, Siggers P, Bogani D, Brixey R, Pastorelli L, Yates L, et al. Sfrp1 and Sfrp2 are required for normal male sexual development in mice. Dev Biol. 2009;326:273–84. doi: 10.1016/j.ydbio.2008.11.023 19100252

31. Kollara A, Brown TJ. Variable expression of nuclear receptor coactivator 4 (NcoA4) during mouse embryonic development. J Histochem Cytochem. 2010;58:595–609. doi: 10.1369/jhc.2010.955294 20354146

32. Padua MB, Jiang T, Morse DA, Fox SC, Hatch HM, Tevosian SG. Combined loss of the GATA4 and GATA6 transcription factors in male mice disrupts testicular development and confers adrenal-like function in the testes. Endocrinology. 2015;156:1873–86. doi: 10.1210/en.2014-1907 25668066

33. Sarraj M, Chua HK, Umbers A, Loveland K, Findlay J, Stenvers KL. Differential expression of TGFBR3 (betaglycan) in mouse ovary and testis during gonadogenesis. Growth Factors. 2007;25:334–45. doi: 10.1080/08977190701833619 18236212

34. Carré GA, Couty I, Hennequet-Antier C, Govoroun MS. Gene expression profiling reveals new potential players of gonad differentiation in the chicken embryo. PLoS One. 2011;6:1–12. doi: 10.1371/journal.pone.0023959 21931629

35. Zhao L, Arsenault M, Ng ET, Longmuss E, Chau TCY, Hartwig S, et al. SOX4 regulates gonad morphogenesis and promotes male germ cell differentiation in mice. Dev Biol. 2017;423:46–56. doi: 10.1016/j.ydbio.2017.01.013 28118982

36. Adolfi MC, Nakajima RT, N RH, Schartl M. Intersex, hermaphroditism, and gonadal plasticity in vertebrates: Evolution of the Mullerian duct and Amh/Amhr2 signaling. Annu Rev Anim Biosci. 2019;7:149–72. doi: 10.1146/annurev-animal-020518-114955 30303691

37. Sekido R, Lovell-Badge R. Sex determination involves synergistic action of SRY and SF1 on a specific Sox9 enhancer. Nature. 2008;453:930–4. doi: 10.1038/nature06944 18454134

38. Rhen T, Schroeder A. Molecular mechanisms of sex determination in reptiles. Sex Dev. 2010;4:16–28. doi: 10.1159/000282495 20145384

39. Simon AR, Rai U, Fanburg BL, Cochran BH. Activation of the JAK-STAT pathway by reactive oxygen species. Am J Physiol Physiol. 1998;275:C1640–52. doi: 10.1152/ajpcell.1998.275.6.C1640 9843726

40. Good SR, Thieu VT, Mathur AN, Yu Q, Stritesky GL, Yeh N, et al. Temporal induction pattern of STAT4 target genes defines potential for Th1 lineage-specific programming. J Immunol. 2009;183:3839–47. doi: 10.4049/jimmunol.0901411 19710469

41. Rawlings JS, Rosler KM, Harrison DA. The JAK/STAT signaling pathway. J Cell Sci. 2004;117:1281–3. doi: 10.1242/jcs.00963 15020666

42. Eggers S, Ohnesorg T, Sinclair A. Genetic regulation of mammalian gonad development. Nat Rev Endocrinol. 2014;10:673–83. doi: 10.1038/nrendo.2014.163 25246082

43. Wang DS, Kobayashi T, Zhou LY, Paul-Prasanth B, Ijiri S, Sakai F, et al. Foxl2 up-regulates aromatase gene transcription in a female-specific manner by binding to the promoter as well as interacting with Ad4 binding protein/steroidogenic factor. Mol Endocrinol. 2007;21:712–25. doi: 10.1210/me.2006-0248 17192407

44. Park M, Shin E, Won M, Kim JH, Go H, Kim HL, et al. FOXL2 interacts with steroidogenic factor-1 (SF-1) and represses SF-1-induced CYP17 transcription in granulosa cells. Mol Endocrinol. 2010;24:1024–36. doi: 10.1210/me.2009-0375 20207836

45. Lei N, Heckert LL. Sp1 and Egr1 Regulate Transcription of the Dmrt1 Gene in Sertoli Cells. Biol Reprod. 2002;66:675–84. doi: 10.1095/biolreprod66.3.675 11870074

46. Topilko P, Schneider-Maunoury S, Levi G, Trembleau A, Gourdji D, Driancourt MA, et al. Multiple pituitary and ovarian defects in Krox-24 (NGFI-A, Egr-1)- targeted mice. Mol Endocrinol. 1998;12:107–22. doi: 10.1210/mend.12.1.0049 9440815

47. Lee SL, Sadovsky Y, Swirnoff AH, Polish JA. Luteinizing hormone deficiency and female infertility in mice lacking. Science. 1996;273:1219–21. doi: 10.1126/science.273.5279.1219 8703054

48. Cutting A, Chue J, Smith CA. Just how conserved is vertebrate sex determination? Dev Dyn. 2013;242:380–7. doi: 10.1002/dvdy.23944 23390004

49. Lydon JP, DeMayo FJ, Funk CR, Mani SK, Hughes AR, Montgomery CA, et al. Mice lacking progesterone receptor exhibit pleiotropic reproductive abnormalities. Genes Dev. 1995;9:2266–78. doi: 10.1101/gad.9.18.2266 7557380

50. Birk OS, Casiano DE, Wassif CA, Cogliati T, Zhao L, Zhao Y, et al. The LIM homeobox gene Lhx9 is essential for mouse gonad formation. Nature. 2000;403:909–13. doi: 10.1038/35002622 10706291

51. Oréal E, Mazaud S, Picard JY, Magre S, Carré-Eusébe D. Different patterns of anti-Müllerian hormone expression, as related to DMRT1, SF-1, WT1, GATA-4, Wnt-4, and Lhx9 expression, in the chick differentiating gonads. Dev Dyn. 2002;225:221–32. doi: 10.1002/dvdy.10153 12412004

52. Wilhelm D, Englert C. The Wilms tumor suppressor WT1 regulates early gonad development by activation of Sf1. Genes Dev. 2002;16:1839–51. doi: 10.1101/gad.220102 12130543

53. Duester G. Families of retinoid dehydrogenases regulating vitamin A function: Production of visual pigment and retinoic acid. Eur J Biochem. 2000;267:4315–24. doi: 10.1046/j.1432-1327.2000.01497.x 10880953

54. Bowles J, Knight D, Smith C, Wilhelm D, Richman J, Mamiya S, et al. Retinoid signaling determines germ cell fate in mice. Science. 2006;312:596–9. doi: 10.1126/science.1125691 16574820

55. Koubova J, Menke DB, Zhou Q, Cape B, Griswold MD, Page DC. Retinoic acid regulates sex-specific timing of meiotic initiation in mice. Proc Natl Acad Sci. 2006;103:2474–9. doi: 10.1073/pnas.0510813103 16461896

56. Roes J, Choi BK, Power D, Xu P, Segal AW. Granulocyte function in grancalcin-deficient mice. Mol Cell Biol. 2003;23:826–30. doi: 10.1128/mcb.23.3.826-830.2003 12529388

57. Maki M, Kitaura Y, Satoh H, Ohkouchi S, Shibata H. Structures, functions and molecular evolution of the penta-EF-hand Ca2+-binding proteins. Biochim Biophys Acta—Proteins Proteomics. 2002;1600(1–2):51–60. doi: 10.1016/s1570-9639(02)00444-2 12445459

58. Lloyd S. Least squares quantization in PCM. IEEE Trans Inf Theory. 1982;28:129–37.

59. Nef S, Schaad O, Stallings NR, Cederroth CR, Pitetti JL, Schaer G, et al. Gene expression during sex determination reveals a robust female genetic program at the onset of ovarian development. Dev Biol. 2005;287:361–77. doi: 10.1016/j.ydbio.2005.09.008 16214126

60. Greenlee AR, Shiao MS, Snyder E, Buaas FW, Gu T, Stearns TM, et al. Deregulated sex chromosome gene expression with male germ cell-specific loss of Dicer1. PLoS One. 2012;7:1–13. doi: 10.1371/journal.pone.0046359 23056286

61. Koenig PA, Nicholls PK, Schmidt FI, Hagiwara M, Maruyama T, Frydman GH, et al. The E2 ubiquitin-conjugating enzyme UBE2J1 is required for spermiogenesis in mice. J Biol Chem. 2014;289:34490–502. doi: 10.1074/jbc.M114.604132 25320092

62. La Fortezza M, Schenk M, Cosolo A, Kolybaba A, Grass I, Classen AK. JAK/STAT signalling mediates cell survival in response to tissue stress. Dev. 2016;143:2907–19. doi: 10.1242/dev.132340 27385008

63. Paul A, Wilson S, Belham CM, Robinson CJM, Scott PH, Gould GW, et al. Stress-activated protein kinases: Activation, regulation and function. Cell Signal. 1997;9:403–10. doi: 10.1016/s0898-6568(97)00042-9 9376221

64. Davis RJ. Signal transduction by the JNK group of MAP kinases. Cell. 2000;103:239–52. doi: 10.1016/s0092-8674(00)00116-1 11057897

65. Wang X, Martindale JL, Liu Y, Holbrook NJ. The cellular response to oxidative stress: influences of mitogen-activated protein kinase signalling pathways on cell survival. Biochem J. 1998;333:291–300. doi: 10.1042/bj3330291 9657968

66. Martinez P, Thanasoula M, Carlos AR, Gómez-López G, Tejera AM, Schoeftner S, et al. Mammalian Rap1 controls telomere function and gene expression through binding to telomeric and extratelomeric sites. Nat Cell Biol. 2010;12:768–80. doi: 10.1038/ncb2081 20622869

67. Teo H, Ghosh S, Luesch H, Ghosh A, Wong ET, Malik N, et al. Telomere-independent Rap1 is an IKK adaptor and regulates NF-κB-dependent gene expression. Nat Cell Biol. 2010;12:758–67. doi: 10.1038/ncb2080 20622870

68. Gilmore TD. Introduction to NF-κB: players, pathways, perspectives. Oncogene. 2006;25:6680. doi: 10.1038/sj.onc.1209954 17072321

69. Morgan MJ, Liu Z. Crosstalk of reactive oxygen species and NF-κB signaling. Cell Res. 2011;21:103–15. doi: 10.1038/cr.2010.178 21187859

70. Oeckinghaus A, Hayden MS, Ghosh S. Crosstalk in NF-κB signaling pathways. Nat Immunol. 2011;12:695–708. doi: 10.1038/ni.2065 21772278

71. McGuire NL, Bentley GE. Neuropeptides in the gonads: From evolution to pharmacology. Front Pharmacol. 2010;1:1–13. doi: 10.3389/fphar.2010.00001 21607058

72. Ulrich-Lai YM, Herman JP. Neural regulation of endocrine and autonomic stress responses. Nat Rev Neurosci. 2009;10:397–409. doi: 10.1038/nrn2647 19469025

73. Todd E V, Liu H, Muncaster S, Gemmell NJ. Bending genders: The biology of natural sex change in fish. Sex Dev. 2016 Mar;10:223–41. doi: 10.1159/000449297 27820936

74. Liu J, Liu X, Jin C, Du X, He Y, Zhang Q. Transcriptome profiling insights the feature of sex reversal induced by high temperature in tongue sole Cynoglossus semilaevis. Front Genet. 2019;10:1–15. doi: 10.3389/fgene.2019.00001 30804975

75. Wang Q, Liu K, Feng B, Zhang Z, Wang R, Tang L, et al. Gonad transcriptome analysis of high temperature induced sex reversal in Chinese Tongue Sole, Cynoglossus semilaevis. Front Genet. 2019;10:1–11. doi: 10.3389/fgene.2019.00001 30804975

76. Hattori R, Castaneda-Cortes D, Arias Padilla L, Strobl-Mazzulla P, Fernandino J. Activation of stress response axis as a key process in environment—induced sex plasticity in fish. Cell Mol Life Sci. 2020; doi: 10.1007/s00018-020-03532-9 32367192

77. Hilton JK, Rath P, Helsell CVM, Beckstein O, Van Horn WD. Understanding thermosensitive transient receptor potential channels as versatile polymodal cellular sensors. Biochemistry. 2015;54:2401–13. doi: 10.1021/acs.biochem.5b00071 25812016

78. Benham CD, Gunthorpe MJ, Davis JB. TRPV channels as temperature sensors. Cell Calcium. 2003;33:479–87. doi: 10.1016/s0143-4160(03)00063-0 12765693

79. Czerwinski M, Natarajan A, Barske L, Looger LL, Capel B. A timecourse analysis of systemic and gonadal effects of temperature on sexual development of the red-eared slider turtle Trachemys scripta elegans. Dev Biol. 2016;420:166–77. doi: 10.1016/j.ydbio.2016.09.018 27671871

80. Yatsu R, Miyagawa S, Kohno S, Saito S, Lowers RH, Ogino Y, et al. TRPV4 associates environmental temperature and sex determination in the American alligator. Sci Rep. 2015;5:1–10. doi: 10.1038/srep18581 26677944

81. Lin JQ, Zhou Q, Yang HQ, Fang LM, Tang KY, Sun L, et al. Molecular mechanism of temperature-dependent sex determination and differentiation in Chinese alligator revealed by developmental transcriptome profiling. Sci Bull. 2018;63:209–12.

82. Yatsu R, Miyagawa S, Kohno S, Parrott BB, Yamaguchi K, Ogino Y, et al. RNA-seq analysis of the gonadal transcriptome during Alligator mississippiensis temperature-dependent sex determination and differentiation. BMC Genomics. 2016;77:1–13. doi: 10.1186/s12864-016-2396-9 26810479

83. Maynard Case R, Eisner D, Gurney A, Jones O, Muallem S, Verkhratsky A. Evolution of calcium homeostasis: From birth of the first cell to an omnipresent signalling system. Cell Calcium. 2007;42:345–50. doi: 10.1016/j.ceca.2007.05.001 17574670

84. Carafoli E. The calcium-signalling saga: Tap water and protein crystals. Nat Rev Mol Cell Biol. 2003;4:326–32. doi: 10.1038/nrm1073 12671655

85. Penna E, Espino J, De Stefani D, Rizzuto R. The MCU complex in cell death. Cell Calcium. 2018;69:73–80. doi: 10.1016/j.ceca.2017.08.008 28867646

86. Hoffman NE, Zhang X, Gill DL, Shanmughapriya S, Rajan S, Jog NR, et al. Ca2+ signals regulate mitochondrial metabolism by stimulating CREB-mediated expression of the mitochondrial Ca2+ uniporter gene MCU. Sci Signal. 2015;8:73–80. doi: 10.1126/scisignal.2005673 25737585

87. Strehler E, Caride A, Filoteo A, Xiong Y, Penniston J, Enyedi A. Plasma membrane Ca2+-ATPases as dynamic regulators of cellular calcium handling. Ann N Y Acad Sci. 2013; doi: 10.1196 23631028

88. Stocker M. Ca2+-activated K+ channels: Molecular determinants and function of the SK family. Nat Rev Neurosci. 2004;5:758–70. doi: 10.1038/nrn1516 15378036

89. Faber ESL, Sah P. Functions of SK channels in central neurons. Clin Exp Pharmacol Physiol. 2007;34:1077–83. doi: 10.1111/j.1440-1681.2007.04725.x 17714097

90. Catterall W. Structure and regulation of voltage-gated Ca2+ channels. Annu Rev Immunol. 2004;22:485–501. doi: 10.1146/annurev.immunol.22.012703.104707 15032586

91. Campiglio M, Flucher BE. The role of auxiliary subunits for the functional diversity of voltage-gated calcium channels. J Cell Physiol. 2015;230:2019–31. doi: 10.1002/jcp.24998 25820299

92. Kumar A, Kumari S, Majhi RK, Swain N, Yadav M, Goswami C. Regulation of TRP channels by steroids: Implications in physiology and diseases. Gen Comp Endocrinol. 2015;220:23–32. doi: 10.1016/j.ygcen.2014.10.004 25449179

93. Miehe S, Crause P, Schmidt T, Löhn M, Kleemann HW, Licher T, et al. Inhibition of diacylglycerol-sensitive TRPC channels by synthetic and natural steroids. PLoS One. 2012;7:e35393. doi: 10.1371/journal.pone.0035393 22530015

94. Ramsey IS, Delling M, Clapham DE. An introduction to TRP channels. Annu Rev Physiol. 2006;68:619–47. doi: 10.1146/annurev.physiol.68.040204.100431 16460286

95. Schaefer M, Plant TD, Obukhov AG, Hofmann T, Gudermann T, Schultz G. Receptor-mediated regulation of the nonselective cation channels TRPC4 and TRPC5. J Biol Chem. 2000;275:17517–26. doi: 10.1074/jbc.275.23.17517 10837492

96. Brostrom MA, Brostrom CO. Calcium dynamics and endoplasmic reticular function in the regulation of protein synthesis: Implications for cell growth and adaptability. Cell Calcium. 2003;34:345–63. doi: 10.1016/s0143-4160(03)00127-1 12909081

97. Michalak M, Corbett EF, Mesaeli N, Nakamura K, Opas M. Calreticulin: One protein, one gene, many functions. Biochem J. 1999;344:281–92. 10567207

98. Dedhar S. Novel functions for calreticulin: interaction with integrins and modulation of gene expression? Trends Biochem Sci. 1994;19:269–71. doi: 10.1016/0968-0004(94)90001-9 8048166

99. Echevarria W, Leite MF, Guerra MT, Zipfel WR, Nathanson MH. Regulation of calcium signals in the nucleus by a nucleoplasmic reticulum. Nat Cell Biol. 2003;5:440–6. doi: 10.1038/ncb980 12717445

100. Burns K, Duggan B, Atkinson EA, Famulski KS, Nemer M, Bleackley RC, et al. Modulation of gene expression by calreticulin binding to the glucocorticoid receptor. Nature. 1994;367:476–80. doi: 10.1038/367476a0 8107808

101. Michalak M, Burns K, Andrin C, Mesaeli N, Jass GH, Busaan JL, et al. Endoplasmic reticulum form of calreticulin modulates glucocorticoid- sensitive gene expression. J Biol Chem. 1996;271:29436–45. doi: 10.1074/jbc.271.46.29436 8910610

102. Wang X, Su M, Gao F, Xie W, Zeng Y, Li D, et al. Structural basis for activity of TRIC counter-ion channels in calcium release. Proc Natl Acad Sci. 2019;116:4238–43. doi: 10.1073/pnas.1817271116 30770441

103. Santamaria-Kisiel L, Rintala-Dempsey A, Shaw G. Calcium-dependent and -independent interactions of the S100 protein family. Biochem J. 2006;396:201–14. doi: 10.1042/BJ20060195 16683912

104. Heizmann C. S-100 proteins. In: Offermanns S, Rosenthal W, editors. Encyclopedia of Molecular Pharmacology. Berlin: Springer; 2008. p. 123–45.

105. Heizmann CW. S100 proteins: structure, functions and pathology. Front Biosci. 2002;7:d1356. 11991838

106. Rebecchi MJ, Pentyala SN. Structure, function, and control of phosphoinositide-specific phospholipase C. Physiol Rev. 2000;80:1291–335. doi: 10.1152/physrev.2000.80.4.1291 11015615

107. Thannickal V, Fanburg B. Reactive oxygen species in cell signaling. Am J Physiol—Lung Cell Mol Physiol. 2000;279:1005–28. doi: 10.1152/ajplung.2000.279.6.L1005 11076791

108. Temple MD, Perrone GG, Dawes IW. Complex cellular responses to reactive oxygen species. Trends Cell Biol. 2005;15:319–26. doi: 10.1016/j.tcb.2005.04.003 15953550

109. Schenk H, Klein M, Erdbrügger W, Dröge W, Schulze-Osthoff K. Distinct effects of thioredoxin and antioxidants on the activation of transcription factors NF-kappa B and AP-1. Proc Natl Acad Sci. 1994;91:1672–6. doi: 10.1073/pnas.91.5.1672 8127864

110. Matsuzawa A. Thioredoxin and redox signaling: Roles of the thioredoxin system in control of cell fate. Arch Biochem Biophys. 2017;617:101–5. doi: 10.1016/j.abb.2016.09.011 27665998

111. Van Der Vos KE, Coffer PJ. FOXO-binding partners: It takes two to tango. Oncogene. 2008;27:2289–99. doi: 10.1038/onc.2008.22 18391971

112. Tao GZ, Lehwald N, Jang KY, Baek J, Xu B, Omary MB, et al. Wnt/β-catenin signaling protects mouse liver against oxidative stress-induced apoptosis through the inhibition of forkhead transcription factor FoxO3. J Biol Chem. 2013;288:17214–24. doi: 10.1074/jbc.M112.445965 23620592

113. Xiong Y, Uys JD, Tew KD, Townsend DM. S-Glutathionylation: From molecular mechanisms to health outcomes. Antioxid Redox Signal. 2011;15:233–70. doi: 10.1089/ars.2010.3540 21235352

114. Yankovskaya V, Horsefield R, Törnroth S, Luna-Chavez C, Miyoshi H, Léger C, et al. Architecture of succinate dehydrogenase and reactive oxygen species generation. Science. 2003;299:700–4. doi: 10.1126/science.1079605 12560550

115. Mikhed Y, Görlach A, Knaus UG, Daiber A. Redox regulation of genome stability by effects on gene expression, epigenetic pathways and DNA damage/repair. Redox Biol. 2015;5:275–89. doi: 10.1016/j.redox.2015.05.008 26079210

116. Cho SS, Kim KM, Yang JH, Kim JY, Park SJ, Kim SJ, et al. Induction of REDD1 via AP-1 prevents oxidative stress-mediated injury in hepatocytes. Free Radic Biol Med [Internet]. 2018;124:221–31. Available from: doi: 10.1016/j.freeradbiomed.2018.06.014 29909290

117. Tonelli C, Chio I, Tuveson D. Transcriptional Regulation by Nrf2. Antioxidants Redox Signal. 2018;29:1727–45. doi: 10.1089/ars.2017.7342 28899199

118. Corona-Herrera GA, Arranz SE, Martínez-Palacios CA, Navarrete-Ramírez P, Toledo-Cuevas EM, Valdez-Alarcón JJ, et al. Experimental evidence of masculinization by continuous illumination in a temperature sex determination teleost (Atherinopsidae) model: is oxidative stress involved? J Fish Biol. 2018;93:229–37. doi: 10.1111/jfb.13651 29931822

119. Zhong P, Huang H. Recent progress in the research of cold-inducible RNA-binding protein. Future Sci OA. 2017;3:FSO246. doi: 10.4155/fsoa-2017-0077 29134130

120. 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:557–71. doi: 10.1534/genetics.115.182840 26936926

121. Haltenhof T, Kotte A, de Bortoli F, Schiefer S, Meinke S, Emmerichs A-K, et al. A conserved kinase-based body temperature sensor globally controls alternative splicing and gene expression. Mol Cell. 2020;78:1–13. doi: 10.1016/j.molcel.2020.03.020 32243827

122. Wang YT, Lim Y, McCall MN, Huang K-T, Haynes CM, Nehrke K, et al. Cardioprotection by the mitochondrial unfolded protein response requires ATF5. Am J Physiol Circ Physiol. 2019;317:H472–8. doi: 10.1152/ajpheart.00244.2019 31274354

123. Zhou D, Palam LR, Jiang L, Narasimhan J, Staschke KA, Wek RC. Phosphorylation of eIF2 directs ATF5 translational control in response to diverse stress conditions. J Biol Chem. 2008;283:7064–73. doi: 10.1074/jbc.M708530200 18195013

124. Furukawa F, Hamasaki S, Hara S, Uchimura T, Shiraishi E, Osafune N, et al. Heat shock factor 1 protects germ cell proliferation during early ovarian differentiation in medaka. Sci Rep. 2019;6927:1–10. doi: 10.1038/s41598-019-43472-4 31061435

125. Metchat A, Akerfelt M, Bierkamp C, Delsinne V, Sistonen L, Alexandre H, et al. Mammalian heat shock factor 1 is essential for oocyte meiosis and directly regulates Hsp90α expression. J Biol Chem. 2009;284:9521–8. doi: 10.1074/jbc.M808819200 19158073

126. Radhakrishnan S, Literman R, Neuwald J, Severin A, Valenzuela N. Transcriptomic responses to environmental temperature by turtles with temperature-dependent and genotypic sex determination assessed by RNAseq inform the genetic architecture of embryonic gonadal development. PLoS One. 2017;12:e0172044. doi: 10.1371/journal.pone.0172044 28296881

127. Kohno S, Katsu Y, Urushitani H, Ohta Y, Iguchi T, Guillette JLJ. Potential contributions of heat shock proteins to temperature-dependent sex determination in the American Alligator. Sex Dev. 2010;4:73–87. doi: 10.1159/000260374 19940440

128. Aloia L, Di Stefano B, Di Croce L. Polycomb complexes in stem cells and embryonic development. Development. 2013;140:2525–34. doi: 10.1242/dev.091553 23715546

129. Marasca F, Bodega B, Orlando V. How polycomb-mediated cell memory deals with a changing environment: Variations in PcG complexes and proteins assortment convey plasticity to epigenetic regulation as a response to environment. Bioessays. 2018;40:1–13.

130. Endoh M, Endo TA, Shinga J, Hayashi K, Farcas A, Ma KW, et al. PCGF6-PRC1 suppresses premature differentiation of mouse embryonic stem cells by regulating germ cell-related genes. Elife. 2017;6:1–26.

131. Cohen I, Bar C, Ezhkova E. Activity of PRC1 and histone H2AK119 monoubiquitination: Revising popular misconceptions. Bioessays. 2020;1900192:1–8. doi: 10.1002/bies.201900192 32196702

132. Yang CS, Chang KY, Dang J, Rana TM. Polycomb group protein Pcgf6 acts as a master regulator to maintain embryonic stem cell identity. Sci Rep. 2016;6:1–12. doi: 10.1038/s41598-016-0001-8 28442746

133. Yan Y, Zhao W, Huang Y, Tong H, Xia Y, Jiang Q, et al. Loss of polycomb group protein Pcgf1 severely compromises proper differentiation of embryonic stem cells. Sci Rep. 2017;7:1–11. doi: 10.1038/s41598-016-0028-x 28127051

134. Fursova NA, Blackledge NP, Nakayama M, Ito S, Koseki Y, Farcas AM, et al. Synergy between variant PRC1 complexes defines polycomb-mediated gene repression. Mol Cell. 2019;74:1020–36 doi: 10.1016/j.molcel.2019.03.024 31029541

135. Blackledge NP, Farcas AM, Kondo T, King HW, McGouran JF, Hanssen LLP, et al. Variant PRC1 complex-dependent H2A ubiquitylation drives PRC2 recruitment and polycomb domain formation. Cell. 2014;157:1445–59. doi: 10.1016/j.cell.2014.05.004 24856970

136. Díaz N, Piferrer F. Lasting effects of early exposure to temperature on the gonadal transcriptome at the time of sex differentiation in the European sea bass, a fish with mixed genetic and environmental sex determination. BMC Genomics. 2015 Dec;16:2–16. doi: 10.1186/1471-2164-16-2 25555398

137. Yokobayashi S, Liang CY, Kohler H, Nestorov P, Liu Z, Vidal M, et al. PRC1 coordinates timing of sexual differentiation of female primordial germ cells. Nature. 2013;495:236–40. doi: 10.1038/nature11918 23486062

138. Shen H, Xu W, Lan F. Histone lysine demethylases in mammalian embryonic development. Exp Mol Med. 2017;49:e325–7. doi: 10.1038/emm.2017.57 28450736

139. Stauffer DR, Howard TL, Nyun T, Hollenberg SM. CHMP1 is a novel nuclear matrix protein affecting chromatin structure and cell-cycle progression. J Cell Sci. 2001;114:2383–93. 11559747

140. Todd E V, Ortega-Recalde O, Liu H, Lamm MS, Rutherford KM, Cross H, et al. Stress, novel sex genes and epigenetic reprogramming orchestrate socially-controlled sex change. Sci Adv. 2019;5:eaaw7006. doi: 10.1126/sciadv.aaw7006 31309157

141. Ribas L, Crespo B, Xavier D, Kuhl H, Rodríguez JM, Díaz N, et al. Characterization of the European Sea Bass (Dicentrarchus labrax) gonadal transcriptome during sexual development. Mar Biotechnol. 2019;21:359–73. doi: 10.1007/s10126-019-09886-x 30919121

142. Georges A, Li Q, Lian J, O’Meally D, Deakin J, Wang Z, et al. High-coverage sequencing and annotated assembly of the genome of the Australian dragon lizard Pogona vitticeps. Gigascience. 2015;45:1–10. doi: 10.1186/s13742-015-0085-2 26421146

143. Dobin A, Davis CA, Schlesinger F, Drenkow J, Zaleski C, Jha S, et al. STAR: Ultrafast universal RNA-seq aligner. Bioinformatics. 2013;29:15–21. doi: 10.1093/bioinformatics/bts635 23104886

144. Li H, Handsaker B, Wysoker A, Fennell T, Ruan J, Homer N, et al. The Sequence Alignment/Map format and SAMtools. Bioinformatics. 2009;25:2078–9. doi: 10.1093/bioinformatics/btp352 19505943

145. Li B, Dewey CN. RSEM: accurate transcript quantification from RNA-Seq data with or without a reference genome. BMC Bioinformatics. 2011;12:21–40. doi: 10.1186/1471-2105-12-21 21235786

146. Robinson MD, McCarthy DJ, Smyth GK. edgeR: A Bioconductor package for differential expression analysis of digital gene expression data. Bioinformatics. 2009;26:139–40. doi: 10.1093/bioinformatics/btp616 19910308

147. RStudio: Integrated development for R. Boston: RStudio Inc; 2015.

148. McCarthy DJ, Chen Y, Smyth GK. Differential expression analysis of multifactor RNA-Seq experiments with respect to biological variation. Nucleic Acids Res. 2012;40:4288–97. doi: 10.1093/nar/gks042 22287627

149. Anders S, Huber W. Differential expression analysis for sequence count data. Genome Biol. 2010;11:R106. doi: 10.1186/gb-2010-11-10-r106 20979621

150. Cox DR, Reid N. Parameter orthogonality and approximate conditional inference. J R Stat Soc B. 1987;49:1–39.

151. Chen Y, Lun ATL, Smyth GK. From reads to genes to pathways: Differential expression analysis of RNA-Seq experiments using Rsubread and the edgeR quasi-likelihood pipeline. F1000Research. 2016;5:1–49. doi: 10.12688/f1000research.8987.2 27508061

152. Lun A, Chen Y, Smyth G. It’s DE-licious: A recipe for differential expression analyses of RNA-seq experiments using quasi-likelihood methods in edgeR. In: Mathe E, Davis S, editors. Statistical Genomics. New York: Humana Press; 2016. p. 391–416.

153. Lund SP, Nettleton D, McCarthy DJ, Smyth GK. Detecting differential expression in RNA-sequence data using quasi-likelihood with shrunken dispersion estimates. Stat Appl Genet Mol Biol. 2012;11:1–42. doi: 10.1515/1544-6115.1826 23104842

154. Lun ATL, Smyth GK. No counts, no variance: allowing for loss of degrees of freedom when assessing biological variability from RNA-seq data. Stat Appl Genet Mol Biol. 2017;16:83–93. doi: 10.1515/sagmb-2017-0010 28599403

155. Phipson B, Lee S, Majewski IJ, Alexander WS, Smyth GK. Robust hyperparameter estimation protects against hypervariable genes and improves power to detect differential expression. Ann Appl Stat. 2016;10:946–63. doi: 10.1214/16-AOAS920 28367255

156. Eden E, Navon R, Steinfeld I, Lipson D, Yakhini Z. GOrilla: a tool for discovery and visualization of enriched GO terms in ranked gene lists. BMC Bioinformatics. 2009;10:48. doi: 10.1186/1471-2105-10-48 19192299

157. Eden E, Lipson D, Yogev S, Yakhini Z. Discovering motifs in ranked lists of DNA sequences. PLoS Comput Biol. 2007;3:0508–22. doi: 10.1371/journal.pcbi.0030039 17381235


Článek vyšel v časopise

PLOS Genetics


2021 Číslo 4
Nejčtenější tento týden
Nejčtenější v tomto čísle
Kurzy

Zvyšte si kvalifikaci online z pohodlí domova

Důležitost adherence při depresivním onemocnění
nový kurz
Autoři: MUDr. Eliška Bartečková, Ph.D.

Koncepce osteologické péče pro gynekology a praktické lékaře
Autoři: MUDr. František Šenk

Sekvenční léčba schizofrenie
Autoři: MUDr. Jana Hořínková, Ph.D.

Hypertenze a hypercholesterolémie – synergický efekt léčby
Autoři: prof. MUDr. Hana Rosolová, DrSc.

Multidisciplinární zkušenosti u pacientů s diabetem
Autoři: Prof. MUDr. Martin Haluzík, DrSc., prof. MUDr. Vojtěch Melenovský, CSc., prof. MUDr. Vladimír Tesař, DrSc.

Všechny kurzy
Přihlášení
Zapomenuté heslo

Zadejte e-mailovou adresu, se kterou jste vytvářel(a) účet, budou Vám na ni zaslány informace k nastavení nového hesla.

Přihlášení

Nemáte účet?  Registrujte se

#ADS_BOTTOM_SCRIPTS#