# Identification, expression, and endocrine-disruption of three ecdysone-responsive genes in the sentinel species Gammarus fossarum

D. Gouveia, F. Bonneton, C. Almunia, J. Armengaud, H. Queau, Davide

Degli Esposti, Olivier Geffard, Arnaud Chaumot

► To cite this version:

D. Gouveia, F. Bonneton, C. Almunia, J. Armengaud, H. Queau, et al.. Identification, expression, and endocrine-disruption of three ecdysone-responsive genes in the sentinel species Gammarus fossarum. Scientific Reports, 2018, 8 (3793), 10.1038/s41598-018-22235-7. hal-02073105

HAL Id: hal-02073105

https://hal.science/hal-02073105

Submitted on 19 Mar 2019

HAL is a multi-disciplinary open access
archive for the deposit and dissemination of sci-
entific research documents, whether they are pub-
lished or not. The documents may come from
teaching and research institutions in France or
abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est
destinée au dépôt et à la diffusion de documents
scientifiques de niveau recherche, publiés ou non,
émanant des établissements d'enseignement et de
recherche français ou étrangers, des laboratoires
publics ou privés.

<figure><img src="image_1.png" /></figure>

Distributed under a Creative Commons Attribution 4.0 International License
SCIENTIFIC REPORTS

<figure><img src="image_2.png" /></figure>

OPEN

# Identification, expression, and endocrine-disruption of three ecdysone-responsive genes in the sentinel species *Gammarus fossarum*

Received: 30 October 2017

Accepted: 14 February 2018

Published online: 28 February 2018

D. Gouveia<fnref n="1" />,², F. Bonneton<fnref n="1" />,³, C. Almunia², J. Armengaud², H. Quéau<fnref n="1" />, D. Degli-Esposti<fnref n="1" />, O. Geffard<fnref n="1" /> & A. Chaumot<fnref n="1" />

Taking advantage of a large transcriptomic dataset recently obtained in the sentinel crustacean amphipod *Gammarus fossarum*, we developed an approach based on sequence similarity and phylogenetic reconstruction to identify key players involved in the endocrine regulation of *G. fossarum*. Our work identified three genes of interest: the nuclear receptors RXR and E75, and the regulator broad-complex (BR). Their involvement in the regulation of molting and reproduction, along with their sensitivity to chemical contamination were experimentally assessed by studying gene expression during the female reproductive cycle, and after laboratory exposure to model endocrine disrupting compounds (EDCs): pyriproxyfen, tebufenozide and piperonyl butoxide. RXR expression suggested a role of this gene in ecdysis and post-molting processes. E75 presented two expression peaks that suggested a role in vitellogenesis, and molting. BR expression showed no variation during molting/ reproductive cycle. After exposure to the three EDCs, a strong inhibition of the inter-molt E75 peak was observed with tebufenozide, and an induction of RXR after exposure to pyriproxyfen and piperonyl butoxide. These results confirm the implication of RXR and E75 in hormonal regulation of female reproductive cycles in *G. fossarum* and their sensitivity towards EDCs opens the possibility of using them as specific endocrine disruption biomarkers.

Endocrine disrupting compounds (EDCs) are exogenous substances that interfere with hormone-regulated physiological processes and provoke adverse health effects in exposed organisms and/or their progeny<fnref n="1" />. EDCs typically interfere with hormone signaling, acting as agonist or antagonist, or with hormone synthesis, through anti-hormonal effects for example²,³. These compounds constitute a worldwide concern for potential human health implications⁴, especially due to the multiple developmental and reproductive disorders observed in wildlife⁵⁻⁹.

In aquatic environments, occurrence and biological effects of EDCs have been reported, such as masculinization and feminization events. This has led to the need of specific monitoring and the development of several robust, specific and reliable biomarkers of EDC exposure in vertebrates, such as the vitellogenin (Vtg) induction in male fishes¹⁰⁻¹². However, due to the molecular divergence acquired through animal evolution, vertebrate-derived biomarkers cannot be used in invertebrates, which represent the majority of animals and key species for ecosystem functioning. The need to develop relevant tools for identification and assessment of EDCs toxic effects and mode of action in invertebrates has been already highlighted by the scientific community. The discovery of imposex in marine gastropods due to Tributyltin (TBT)¹³, and intersex occurrence in natural populations of gammarids¹⁴,¹⁵ are two examples of the importance of developing more studies of EDCs effects in

<footnote n="1">¹Iirstea, UR RiverLy, Laboratoire d'écotoxicologie, centre de Lyon-Villeurbanne, 5 rue de la Doua CS 20244, F-69625, Villeurbanne, France. ²Laboratoire Innovations technologiques pour la Détection et le Diagnostic (Li2D), Service de Pharmacologie et Immunoanalyse (SPI), CEA, INRA, F-30207, Bagnols sur Cèze, France. ³IGFL, Université de Lyon, CNRS UMR5242, Ecole Normale Supérieure de Lyon, Université Claude Bernard Lyon 1, 46 allée d'Italie, F-69364, Lyon, France. Correspondence and requests for materials should be addressed to A.C. (email: arnaud.chaumot@irstea.fr)</footnote>
<figure><img src="image_3.png" /><figcaption>**Figure 1.** Putative ecdysone regulatory cascade in crustaceans. Ecdysone is secreted by the Y-organ, leading to its active form 20-hydroxyecdysone (20E). The heterodimer EcR/RXR receives 20E and initiates the transcription of the “early” response genes, which in turn regulates the transcription of “late” genes responsible for initiating molting and vitellogenesis. Modified from³⁰.</figcaption></figure>

reproductive function of invertebrate species. Nevertheless, due to a currently limited knowledge in invertebrate endocrinology, today there is still a scarce understanding of EDCs modes of action in the majority of invertebrate species.

Among crustacean species, the amphipod *Gammarus fossarum* is commonly used as a sentinel species in freshwater risk assessment¹⁶. Many toxicity markers are available in this species, mainly based on life-history traits such as reproductive features¹⁷, locomotor behavior/feeding rate¹⁸, neurotoxicity¹⁹ and genotoxicity²⁰. The recent advances in proteomics techniques led to the development of a specific methodology based on liquid chromatography coupled to tandem mass spectrometry for quantifying a Vtg protein in *G. fossarum*²¹. However, the proposal of this protein as an endocrine disruption (ED) biomarker for male feminization was discarded after laboratory and field experiments that yielded high variability and low Vtg inductions in organisms exposed to contaminants²². Similar conclusions were reached when analyzing expression levels of two Vtg genes in *Echinogammarus marinus*, highlighting the fact that invertebrate biomarkers must be developed and validated *de novo*, and not just derived from their eukaryotic counterparts²³. More recently, a proteogenomic study based on the alliance of transcriptomic and proteomic approaches²⁴ lead to the construction of a database (GFOSS database) containing 1,873 experimentally validated *Gammarus fossarum* specific-proteins²⁵. However, few proteins involved in hormonal regulation could be recognized from this database, and potential key candidates for ED biomarkers were not identified. Indeed, such proteins are expected to be present at very low levels and therefore are hence probably lost in the “background noise” of the proteome.

Taking this context into account, we hypothesized that a gene candidate approach considering a thorough bioinformatics mining of the transcriptomic database obtained in *G. fossarum* could be used to identify candidate biomarkers of endocrine disruptor effects in this species. As discussed in previous studies¹⁷,²⁶, the reproduction cycle of female amphipods such as gammarids (vitellogenesis/oogenesis, embryogenesis) is closely correlated with molting. Disruption of these processes by contamination results in reproductive impairment. In crustaceans, molting is regulated by a multi-hormonal system including specific neuropeptide hormones, juvenile hormones such as methyl farnesoate (MF), and biological active forms of ecdysteroids such as 20-hydroxyecdysone (20E) or ponasterone A (reviewed in²⁶ and³). Among arthropods, while the ecdysone pathway is well characterized in insects, much less is understood in crustacean species. In insects, 20E binds to and activates the ecdysone receptor/retinoid-X-receptor (EcR/RXR, mainly known as USP in insects) nuclear receptor complex, initiating a cascade of gene-regulatory events that will mediate molting and reproduction (Fig. 1). Early response genes from the insect ecdysone pathway include the nuclear receptors ecdysone induced protein 75B (E75), ecdysone induced protein 78 C (E78) and hormone receptor 3 (HR3), and the transcription factors broad-complex (BR) and ecdysone induced protein 74EF (E74)²⁷. In crustaceans, early response genes are thought to be similar to those of insects³,²⁸,²⁹.

The present study aimed at identifying gene candidate biomarkers of ED effects in the crustacean *G. fossarum*. Assuming a possible homology between crustacean and insect endocrine regulator genes, we mined the transcriptomic GFOSS database and performed sequence similarity searches and phylogenetic analysis in order to identify some key players potentially involved in endocrine regulation in *G. fossarum*. Since the ecdysone pathway is crucial for arthropod molting and reproduction, our work focused on ecdysone-responsive genes whose homologous sequences were identified in GFOSS. Through molecular phylogenetic analysis, we validated the annotations of the nuclear receptors E75 and RXR, and of the transcription factor BR in *G. fossarum* sequences. In order to obtain reliable nucleotide sequences of gene candidates and primers for subsequent gene expression studies, these genes were cloned and sequenced. Their involvement in the endocrine control of molting and reproductive process and sensitivity to EDCs were evaluated by studying their profile during the female reproductive cycle and following laboratory exposure to three chemicals acting on the hormonal regulation of target arthropods: pyriproxyfen, tebufenozide and piperonyl butoxide.
<table><caption>**Table 1.** Sequence similarity search results after blastp of RXR, EcR and BR query sequences against the GFOSS database.</caption><thead><tr><th>Query sequence (Accesion n° - name - [species])</th><th>Top 3 blast results</th><th>Score</th><th>Identitites</th><th>Percentage</th><th>Expect</th><th>GFOSS annotation</th></tr></thead><tbody><tr><td rowspan="3">XP_018017941.1 PREDICTED: ecdysone receptor-like isoform X1 [Hyalella azteca]</td><td>GFOSS_83504_fr5</td><td>94</td><td>39/86</td><td>45</td><td>1,0E-18</td><td>ultraspiracle</td></tr><tr><td>GFOSS_2900_fr4</td><td>70</td><td>95/389</td><td>24</td><td>1,0E-11</td><td>Nuclear hormone receptor E75</td></tr><tr><td>seq. 406923_fr4</td><td>70</td><td>30/47</td><td>63</td><td>2,0E-11</td><td>—</td></tr><tr><td rowspan="3">XP_018008501.1 PREDICTED: broad-complex core protein isoforms 1/2/3/4/5-like isoform X4 [Hyalella azteca]</td><td>GFOSS_1452_fr3</td><td>281</td><td>162/230</td><td>70</td><td>7,00E-75</td><td>broad-complex, isoform Z2</td></tr><tr><td>seq. 986730_fr4</td><td>276</td><td>132/139</td><td>94</td><td>2,00E-73</td><td>—</td></tr><tr><td>seq. 1061041_fr1</td><td>274</td><td>132/139</td><td>94</td><td>7,00E-73</td><td>—</td></tr><tr><td rowspan="3">AHF65151.1 retinoid X receptor [Eriocheir sinensis]</td><td>seq. 307974_fr2</td><td>228</td><td>113/198</td><td>57</td><td>3,0E-59</td><td>—</td></tr><tr><td>GFOSS_83504_fr5</td><td>159</td><td>71/90</td><td>78</td><td>2,0E-38</td><td>ultraspiracle</td></tr><tr><td>GFOSS_2900_fr4</td><td>114</td><td>81/307</td><td>26</td><td>1,0E-24</td><td>Nuclear hormone receptor E75</td></tr></tbody></table>

## Results and Discussion

**Phylogenetic sequence validation.** In the present study, based on the GFOSS RNAseq-derived transcriptome database, we aimed to identify *G. fossarum* specific sequences homologue to key genes known to be involved in endocrine regulation of arthropods. Seeing the importance of ecdysone pathway in arthropod molting and reproduction, we focused on three ecdysone-responsive genes: RXR, EcR and BR. After blast sequence similarity searches against GFOSS, the top three results were kept for phylogenetic validation of their annotations. As nuclear receptors are paralog genes with a highly evolutionary conserved domain called DNA-binding domain (DBD) and other conserved amino-acyl sequences in functional areas in the ligand-binding domain (LBD), the BLAST searches for RXR and EcR sequences yielded common results (Table 1). Sequences GFOSS_83504 (annotated as ultraspiracle - same as RXR - in GFOSS) and GFOSS_2900 (annotated as nuclear hormone receptor E75 in GFOSS database) were found among the top three blast results of both searches using the query sequences of RXR and EcR from the decapod *Eriocheir sinensis* and the amphipod *Hyalella azteca*, respectively. Concerning BR, blast results yielded one sequence GFOSS_1452 (annotated as broad-complex isoform Z2) and two unassembled reads seq. 986730_fr4 and seq. 1061041_fr1.

The sequence annotations available in GFOSS (reported in Table 1) were previously generated by an automatic bioinformatic pipeline based on reciprocal blast procedure against public sequence databases²⁵. Such approach can be flawed due to affection of functional annotations from a paralog sequence to the gene identified in the species of interest. Then, in order to select the *G. fossarum* ortholog genes of RXR, EcR, and BR, phylogenetic analyses were conducted to validate the annotation of the candidate sequences. Each *G. fossarum* sequence from Table 1 was thus aligned with ortholog sequences from different crustacean/insect species, and paralog sequences from different members of the corresponding multigenic family, i.e. nuclear receptor family for RXR and EcR, and BTB-ZF family for BR (sequence accession numbers are listed in Supplementary Table 1). The phylogenetic trees for candidate sequences are represented in Fig. 2 (EcR), Fig. 3 (BR) and Fig. 4 (RXR).

As presented in Fig. 2, sequence GFOSS_2900 did not cluster with EcR sequences, but rather with sequences from the nuclear receptor E75, confirming the bioinformatic functional annotation on GFOSS. The unassembled read GFOSS_seq. 406923 (47 amino acids long) also clustered with E75 sequences (data not shown) and therefore it also does not correspond to a *G. fossarum* EcR sequence. As expected, both these *G. fossarum* sequences showed the closest relationship with E75 from *Hyalella azteca*, the only genome-sequenced crustacean amphipod. These results show that no EcR transcripts were sequenced in the GFOSS RNAseq-derived transcriptome database. Because seq. 406923 is a short non-assembled read prone to have sequence errors, it was not used for primer design and sequence amplification. However, the E75 annotated sequence GFOSS_2900 was kept for subsequent studies. E75 is in fact one of the primary target genes of the EcR/RXR heterodimer, and regulates further downstream transcriptional cascades for molting and reproductive processes in both insects and decapod crustaceans³⁰⁻³².

Figure 3 shows the phylogenetic tree obtained from the sequence dataset constructed for BR. CHINMO, TTK and lola genes from crustaceans and insects were considered to build this tree for the BTB-Zf family. As expected, the *G. fossarum* sequence GFOSS_1452 clustered together with the BR sequence from *H. azteca*, along with the other BR branches for crustaceans and insects. Similarly to the EcR results, the unassembled reads seq. 986730 and seq. 1061041 were not considered for the subsequent studies. Of note, this tree evidences that available annotation of crustacean sequences from this BTB-Zf family could be questioned due to potential duplication events unfound in the insects, and which apparently may have led to misannotations. For example, lola genes from *Caligus clemensi* and *Lepeophteirus salmonis* clustered with the crustacean BR branch. This means that this specific family of genes is still relatively unknown and poorly described in crustacean species, and some of the annotations available by successive sequence similarity derivation when annotating new genomes could be erroneous. However, given the good blast score obtained and the posterior phylogenetic analysis for sequence GFOSS_1452, we validated the BR annotation for this sequence.

In Fig. 4 are represented two phylogenetic trees, one for each RXR *G. fossarum* sequence identified by blast. Two sequences were retained for this gene in order to obtain a higher alignment length. In fact, sequence GFOSS_83504 covers only the DBD domain of the RXR gene, while the unassembled read GFOSS_seq. 307974 covers another part of the gene, the ligand-binding domain (LBD). Analyzing the trees, both sequences clustered with the RXR branch, and were thus validated as being RXR sequences.
<figure><img src="image_4.png" /><figcaption>**Figure 2.** Phylogenetic tree of the GFOSS_2900 sequence. Tree construction was based on a homologous sequence dataset constructed using ortholog and paralog sequences from the same multigenic family as EcR (nuclear receptors E78, E75 and HR3).</figcaption></figure>

RXR, E75 and BR were successfully cloned and sequenced using the primers from Table 2. For E75, a 419 nucleotide sequence was obtained, covering a portion of both the DBD and LBD of this nuclear receptor. Smaller sequences were obtained for BR and RXR (278 and 220 nucleotides, respectively). The 220 nucleotide RXR sequence corresponds to primers designed in the DBD domain, using only the GFOSS_83504 sequence. Other primer pairs were tested using forward sequences on the GFOSS_83504 and reverse sequences on seq. 307974, in order to obtain a sequence covering both DBD and LBD. However, we failed to amplify the RXR sequence with the pairs of primers tested.

**Expression during female reproductive cycle.** In order to assess the possible involvement of the three ecdysone-inducible genes in reproductive processes, we quantified their expression during the female reproductive/molt cycle. The bar graphs in Fig. 5 represent the expression levels of the three genes through the different stages. E75 and RXR transcripts showed important modulations during the female reproductive cycle, while BR remained relatively stable.

**E75.** The expression pattern of E75 presents a strong peak at the C2/D1 stage, with a fold change of 12 relative to stage B. The C2 developmental stage is characterized by the quick increase in yolk formation and oocyte size, namely the start of secondary vitellogenesis<sup>17,33</sup>. This continuous increase from the C1/C2 until the C2/D1 stage suggests that the expression of E75 could be related to the onset of secondary vitellogenesis and oogenesis. This expression peak was validated through the analysis of individual replicates, where statistically significant difference was observed between C2D1 and D1D2 stages (p value < 0.0001, Supplementary Figure 2). A second slight increase in the levels of this transcript in the final D2 stage (Fig. 5) could indicate the existence of a second short peak of expression, as observed by the individual-specific increases obtained in biological replicates (Supplementary Figure 2). However, no statistically significant results corroborate the existence of this peak. It has been demonstrated that E75 plays a role in molting, oogenesis and vitellogenesis of several insects<sup>34–37</sup>. It has been also suggested that E75 plays an important role in the interplay between the ecdysone and JH pathways<sup>38,39</sup>. In *Litopenaeus vannamei*, a crustacean decapod that is grouped with amphipods in Malacostraca crustaceans, E75 is positively regulated by ecdysone, regulates subsequent chitinase activities, and has molt stage-specific expression patterns<sup>32</sup>. The levels of EcR and E75 elevated significantly during vitellogenic stages I and II in the hepatopancreas of the crab *Oziothelphusa senex senex*<sup>31</sup>. Priya et al.<sup>40</sup> showed that E75 is modulated through the molt process of *Fenneropenaeus chinensis*, and silencing the gene leads to a complete arrest of molting. However, other studies<sup>41</sup> showed that E75 is marginally induced by ecdysteroids and presents constant levels of expression
<figure><img src="image_5.png" /><figcaption>**Figure 3.** Phylogenetic tree of the GFOSS_1452 sequence. Tree construction was based on a homologous sequence dataset constructed using ortholog and paralog sequences from the same multigenic family as BR (transcription factors tramtrack TTK and longitudinals-like protein lola). CHINMO was used as an outgroup.</figcaption></figure>

during the molt cycle and embryo development of the branchiopod *Daphnia magna*. The results obtained herein suggest that, similarly to other malacostraca crustaceans, E75 may play a role in female reproductive-related processes of *G. fossarum*, notably in vitellogenesis.

**BR.** BR presented an up-and-down expression pattern, tending to a slight decrease in the later stages of the molt cycle. However, no deep modulations can be observed. BR is known to play a role in ecdysone and JH responsive processes such as metamorphosis42,43, vitellogenesis and oogenesis42,44 in several insects. Little is known about this gene in crustacean species, however it has been suggested that BR acts as a regulator of the action of methyl farnesoate (MF)28. Another study highlighted the importance of BR during ovarian development of *P. monodon*45. In that study, BR expression was more abundantly expressed in the ovaries than other tissues and injections of serotonin, progesterone and ecdysone all increased BR expression levels. This little information on BR function in crustacean species is complicated by the unknown evolutionary origin of this gene. Because there was no evidence of a BR gene ortholog in the non-insect and non-crustacean arthropod genomes, it has been suggested that BR was gained in the Pancrustacea46. Our results also highlight the limitations of using molecular information obtained from insects (or other phylogenetically distant species) for functional annotation of crustacean sequences (cross-species homology). The fact that the GFOSS putative BR sequence was not modulated during the reproductive processes opens several questions. The fact that this gene is a complex locus with many isoforms can be a limiting factor for the analysis. In order to investigate the role of different possible isoforms in *G. fossarum*, it would be necessary to test other BR-annotated sequences obtained by sequence homology and to check their expression patterns. The fact that the whole body of gammarids was used could have also diluted the response of this gene, and tissue specific modulations were therefore not detected.

**RXR.** An intense peak at the post-molt stage A was observed for the RXR transcript. This peak was validated through the analysis of individual replicates, where a significant difference was observed between stage A and C2D1 (*p* value < 0.0001, Supplementary Figure 2). This stage lasts for only one to two days and is the shortest of the *G. fossarum* molting cycle. After this expression peak, RXR expression gradually decreased until the C1/C2 stage, stabilizing until the end of the cycle. This expression pattern suggests that RXR may be involved in ecdysis, since the expression peak takes place during this process. There are several studies on both RXR tissue-specific expression, and during the reproductive cycle of crustaceans. However, different trends are observed between
<figure><img src="image_6.png" /><figcaption>**Figure 4.** Phylogenetic tree of the (A) GFOSS_83504 sequence and (B) GFOSS_seq. 307974 read. Tree construction was based on a homologous sequence dataset constructed using ortholog and paralog sequences from the same multigenic family as RXR (nuclear receptor HNF4). EcR was used as an outgroup.</figcaption></figure>

studies, due to inter-species variability and to the existence of several RXR isoforms. For example, in the lobster *Homarus americanus*, RXR was more expressed in the mandibular organ (MO) than in other tissues⁴⁷. In MO and in ovaries, an expression peak occurred in the intermediate stage of the reproductive cycle (secondary vitellogenesis), while for hepatopancreas the peak occurred only during the early stages⁴⁷. In the decapod shrimp *Fenneropenaeus chinensis*, two isoforms of RXR were more expressed in both the gills and nerve cord, and their expression during the molt cycle using whole-body organisms showed a peak at late premolt and postmolt stages, suggesting a role in ecdysis⁴⁸. Another decapod, *Portunus trituberculatus*, presented higher expressions in both
<table><caption>**Table 2.** Nucleotide sequences of the primers used in the PCR and qPCR experiments.</caption><thead><tr><th rowspan="2"></th><th rowspan="2">Primer ID</th><th rowspan="2">5' - &gt; 3' nucleotide sequence</th><th rowspan="2">Fragment (bp)</th><th colspan="2">Ct standard deviation</th></tr><tr><th>Experiment 1</th><th>Experiment 2</th></tr></thead><tbody><tr><td rowspan="6">Target genes</td><td>RXR FW1</td><td>AAACCTTTGTGCCATCTGTG</td><td rowspan="3">220</td><td rowspan="3">—</td><td rowspan="3">—</td></tr><tr><td>RXR RV3</td><td>GAACAGCCTCCCTCTTCATG</td></tr><tr><td>BRC FW1</td><td>GGGAGCAGACCAACAGTTCT</td></tr><tr><td>BRC RV2</td><td>TGACCTCGCCATGGTACATA</td><td rowspan="3">278</td><td rowspan="3">—</td><td rowspan="3">—</td></tr><tr><td>E75 FW1</td><td>TCGTCAATGCAGCAAGAATC</td></tr><tr><td>E75 RV2</td><td>GGATGTTCTTGGCGAAGGT</td></tr><tr><td rowspan="6">Reference genes</td><td>EF FW</td><td>TTCAAGTATGCCTGGGTGCT</td><td rowspan="3">82</td><td rowspan="3">0,42</td><td rowspan="3">0,34</td></tr><tr><td>EF RV</td><td>CGAACTTCCAGAGAGCAATGTC</td></tr><tr><td>GADPH FW</td><td>CGCTGGCCAGAACATCATTC</td></tr><tr><td>GADPH RV</td><td>CGGCCTTGATGTCGTCGTAA</td><td rowspan="3">—</td><td rowspan="3">0,84</td><td rowspan="3">0,25</td></tr><tr><td>18 S FW</td><td>TGGGGGAGGTAGTGACGAAATC</td></tr><tr><td>18 S RV</td><td>CCTGCGCTCGATACAGACATTC</td></tr></tbody></table>

the Y-organ (YO) and the ovaries, as well as a further expression peak in the YO during the post-molt stages⁴⁹. The authors highlight that RXR has an effect on molting but its expression does not correspond to the ecdysone titers observed in their previous study⁵⁰. Some studies already proposed that MF and RXR may have a specific signaling pathway independent of the ecdysone pathway in decapods. This pathway is thought to participate in Vtg gene expression and ovarian development⁴⁷,⁵¹,⁵². Nevertheless, in line with most crustacean literature published, the results presented herein suggest that RXR might be involved in the ecdysone-mediated ecdysis process of *G. fossarum*.

**Expression following exposure to endocrine disrupting chemicals.** In order to understand the susceptibility of E75, RXR and BR to EDCs, gammarids were exposed to three chemicals known to act on endocrine control of insect development. The chemicals were used at sub-lethal concentrations (PYR 5 µg/L, TEB 0.5 µg/L and PBO 150 µg/L) chosen based on previous experiments made in our laboratory, and therefore known to provoke reproductive impairments in female *G. fossarum*, particularly in vitellogenesis and oogenesis (Supplementary Figure 3).

Pyriproxyfen (PYR) is a potent juvenile hormone analog with insecticidal activity and little mammalian toxicity⁵³, being very effective in inhibiting metamorphosis and embryogenesis in insects. Tebufenozide (TEB) is an analog for the molting hormone 20E. In lepidopteran species, the action of TEB results in the initiation of a precocious and incomplete molt (the new cuticle is formed before the shed of the old cuticle). This pesticide is not cleared by the organism, consequently inducing a continuous perturbation that maintains the expression of “early” response genes and the repression of genes normally expressed after the ecdysone peak for completing the molt⁵⁴. Piperonyl butoxide (PBO) has no direct pesticidal activity of its own, but acts as a synergist of other insecticides by inhibiting some detoxification enzymes, namely the cytochrome P450 monooxygenases and some esterases⁵⁵.

The durations of exposure were chosen based on the results presented above (Fig. 5). Our aim was to target the peaks of maximum expression for RXR in the beginning of the molting cycle (4 days, stage B) and for E75 in the onset of vitellogenesis (14 days, stage C2/D1). For RXR, since the stage A is very short, the stage B was selected for the first termination date. The nine day exposure (C1 molt stage) was intended to be an intermediate stage between the latter two in order to verify the occurrence of possible RXR/E75 inductions, as well as to confirm the trends between stages observed in the previous experiment (comparing the controls).

The values for the normalized expression in each condition, for the three target genes and at every exposure time, are shown in Fig. 6. Analyzing solely the control samples we can observe the same trends as before between the B, C1 and C2/D1 stages for the three genes analyzed, reinforcing the repeatability of the described expression pattern during female cycle. The more evident example is E75 expression, very similar between day-4 and day-9 experiments, with the peak appearing at day-14 when females are at their C2/D1 reproductive stage (fold change of 13 between day-4 and day-14). For both RXR and BR control samples, the expression pattern throughout the reproductive stage is also similar to the previous experiment, as shown in Table 3. Comparing control and exposed samples, some effects of contamination were observed, suggesting that the female endocrine system was disturbed upon exposure to the pesticides.

To this respect, the E75 transcript was the more responsive to contamination. At day-4 of exposure there was an increase in the order of 1.5 to 2 fold of E75 expression in all three contaminated conditions, a trend that repeated itself at day-9 (fold changes between 1.6 and 2.2). At day-14, transcript levels returned to normal in PYR and PBO exposures, while the TEB exposure led to a 2.8-fold decrease of E75 expression. The first increase in E75 expression means that the contaminants may induce, in a first phase, an over-expression of the “early” ecdysone-responsive genes in *G. fossarum*. These effects are concomitant with the role of this insecticide as an analog of the molt hormone 20E. However, the inhibition observed for the TEB exposed organisms at day-14 suggests that the ecdysteroid agonist effect of TEB is no longer taking place. Although in other species TEB is known to induce precocious molts through the continuous expression of the “early” ecdysone-responsive genes, in the
<figure><img src="image_7.png" /><figcaption>**Figure 5.** Expression profiles of E75, BR and RXR during the female reproductive cycle. Expression values of target genes were normalized to the expression of the reference gene EF. Each bar is the mean obtained from three technical replicates. Each condition (A, B, C1, C1/C2, C2, C2/D1, D1, D1/D2, D2 and D2 F) contains mRNAs from five individuals pooled together. The represented standard deviations derive from the individual analysis of the five biological replicates (Supplementary Figure 2) from the selected stages.</figcaption></figure>

<figure><img src="image_8.png" /><figcaption>**Figure 6.** mRNA levels of E75, RXR and BR in control and exposed organisms for the three exposure times studied. Expression values of target genes were normalized to the expression of the reference gene EF. Each bar is the mean ± SD obtained from three experimental replicates. Each condition (control, PYR, TEB, PBO) contains mRNAs from five individuals pooled together. Statistical significance was accepted with fold changes > 1.5 and p-value < 0.05.</figcaption></figure>

<table><caption>**Table 3.** Fold change expressions observed between different reproductive stages (B, C1 and C2/D1) in the two experiments performed.</caption><thead><tr><th rowspan="2"></th><th rowspan="2"></th><th>Experiment 1<br>Reproductive cycle</th><th>Experiment 2<br>Contaminant exposures</th></tr><tr><th>fold change</th><th>fold change</th></tr></thead><tbody><tr><td rowspan="3">E75</td><td>B/C1</td><td>1,27</td><td>1,07</td></tr><tr><td>B/C2D1</td><td>11,68</td><td>13,30</td></tr><tr><td>C1/C2D1</td><td>9,22</td><td>12,38</td></tr><tr><td rowspan="3">BR</td><td>B/C1</td><td>-1,19</td><td>-1,35</td></tr><tr><td>B/C2D1</td><td>-1,32</td><td>-1,76</td></tr><tr><td>C1/C2D1</td><td>-1,11</td><td>-1,31</td></tr><tr><td rowspan="3">RXR</td><td>B/C1</td><td>-1,43</td><td>-1,61</td></tr><tr><td>B/C2D1</td><td>-3,18</td><td>-3,81</td></tr><tr><td>C1/C2D1</td><td>-2,23</td><td>-2,36</td></tr></tbody></table>

studied species we have rather observed an inhibition of the early response-gene E75 after a 14-day exposure, which could reflect a general disruption of the normal ongoing of the female reproductive molting cycle.

E75 is a key molecule in the reproductive processes of arthropods, and since is an early ecdysone-responsive gene, is prone to be affected by interactions of chemical compounds with the ecdysteroid receptor which is at the upstream of the ecdysone regulatory cascade. In a recent study⁵⁶, an Adverse Outcome Pathway (AOP) was developed for explaining the molecular events that lead to lethal molting due to ecdysone receptor agonism in arthropods. The mechanism proposed was based on the literature available for both insects and crustaceans. The expression of the E75 gene was considered the first key event of the molecular chain leading to the adverse outcome, reinforcing its potential to be used as a biomarker of endocrine disturbances.
Some effects were also seen in the expression of the **RXR** transcript at day-14 of exposure. PYR and PBO led to an over-expression of RXR at the C2/D1 stage (fold changes of 2.8 and 2.1 respectively). As discussed before, RXR and MF may have an ecdysone-independent regulatory pathway, which can explain the over-expression of RXR following exposure to the MF analog PYR. However, the possible mechanism that led to an induction of RXR provoked by PBO is not clear. It is known that this molecule does not have any insecticidal activity, but it inhibits the cytochrome P450 oxydases⁵⁷, which participate in the ecdysone synthesis and inactivation pathways. If this occurs, we would rather expect an inhibition of RXR expression due to an inhibition of the 20E expression.

Concerning **BR** expression, some statistically significant modulations were observed, but with fold changes <1.5. This slight over-expression of this transcript at day-14 of exposure with PYR and PBO (fold changes of 1.43 and 1.47 respectively and p value < 0.001) follows the modulations observed for RXR, but establishing an eventual connection is difficult since the other early gene E75 was not modulated by the same exposure.

## Conclusion

Based on a RNAseq-derived transcriptomic database and phylogenetic analyses, we identified potential key players of the hormonal system of *G. fossarum*: the nuclear receptors E75, RXR and the regulator BR. Our analyses, in line with most of the information available in the literature, provided some relevant insights into the possible roles of E75 and RXR in vitellogenesis, oogenesis and ecdysis. Exposure to environmental concentrations of three EDCs provoked significant modulations in the levels of these transcripts. The role of BR in reproductive *G. fossarum* female processes was not elucidated in this study. Allied to the fact that this gene revealed less sensitivity to contamination than the others, BR is not proposed as a biomarker of EDC exposure in female *G. fossarum*.

These results provide the first elements for understanding the molecular mechanisms of endocrine regulation in gammarids. Despite calling for more mechanistic studies in order to unravel the roles of these genes in specific physiological processes, their susceptibility to being modulated by an exposure to EDCs was verified in this study. The rupture of the expression cycle of these genes is therefore proposed as a biomarker of EDCs effects in female *G. fossarum*. To our knowledge, our study is the first to analyze the impact of endocrine disruptors in ecdysone-responsive genes of *G. fossarum*.

## Methods

**Candidate selection and phylogenetic validation.** *Selection of Gammarus fossarum specific sequences.*
Given the importance of ecdysone in arthropod molting and reproduction, we focused on three ecdysone-responsive genes: RXR, EcR and BR. For each candidate, a BLASTP homology search was performed against our RNAseq-derived transcriptome database GFOSS²⁵, using query sequences from phylogenetically close species obtained from the NCBI protein database. The EcR and BR query sequences used for blast were from the recently genome-sequenced amphipod *Hyallela azteca*, while the RXR sequence belonged to the crab *Eriocheir sinensis*. For this gene, decapods were the closest order to have sequenced RXR sequences in the NCBIr protein database. The BLAST search was conducted against both assembled contigs and unassembled reads databases. The top 3 blast scores with E-values inferior to $E^{-10}$ and alignment length superior to 25% were chosen for the following analysis (Table 1).

*Sequence alignment and phylogenetic analysis.* For each candidate, a homologous sequence dataset was constructed using ortholog and paralog sequences from the same multigenic family as the candidate gene. The deduced *G. fossarum* sequences were aligned with other known sequences from a diversity of crustacean and insect species, obtained from GenBank (National Center for Biotechnology Information - NCBI) (listed in Supplementary Table 1).

Multiple sequence alignments and phylogenetic trees were performed in the SeaView software version 4.6.1⁵⁸. Alignments were made using the Multiple Sequence Comparison by Log-Expectation (MUSCLE)⁵⁹ program using default parameters. Phylogenetic tree model was determined using the IQ-TREE web server model selection⁶⁰. Trees were built also in the SeaView, using the PhyML program based on the maximum-likelihood principle⁶¹, with the LG model with 4 substitution rate categories to estimate the gamma parameter shape, and 100 bootstrap replicates for branch support. Default settings were used for all other parameters.

**Gene expression studies.** *Sampling and maintenance of organisms.* Gammarids were collected from the River Le Pollon in France (45°57′25.8″N 5°15′43.6″E) and acclimatized to laboratory conditions as previously described³³. Organisms were collected by kick sampling using a net, and quickly transported to the laboratory in plastic buckets containing freshwater from the station. In the laboratory, organisms were kept in 30 L tanks continuously supplied with drilled groundwater (500µS/cm) and under constant aeration for at least 10 days. A 16/8 h light/dark photoperiod was maintained and the temperature was kept at 12 ± 1 °C. Organisms were fed *ad libitum* with alder leaves (*Alnus glutinosa*), previously conditioned for 6 ± 1 days in water. Sexually mature male and female organisms (selection based on the developmental stage of the embryos; well-developed juveniles in the female marsupium) were sampled for RNA extraction and gene sequencing protocols.

*Selection of females at different phases of reproductive cycle.* For gene expression studies, five females from each main and intermediate stages of the molt cycle (A, B, C1, C1/C2, C2, C2/D1, D1, D1/D2, D2, and final D2) were sampled. The selection of females at specific stages was determined through observation of characteristic integumental changes on the dactylopodite and protopodite during the molt cycle (described in¹⁷). Briefly, in *Gammarus fossarum* females, oogenesis and embryogenesis occurs simultaneously with the molting cycle. Six molt stages were previously defined¹⁷, A and B for postmolt, C1 and C2 for intermolt, and D1 and D2 for premolt. Organisms were incubated overnight at 4 °C in 300 µL of RNAlater, frozen in liquid nitrogen and immediately stored at -80 °C until total RNA extraction.
**Laboratory exposure to endocrine disrupting chemicals.** Coupled male and female organisms were sampled and isolated in a separated aquarium. Immediately after spawning, females were isolated and placed in 500 mL glass beakers (10 females per beaker, 3 beakers per condition). There were three conditions corresponding to three different exposure times: 4, 9 and 14 days. Solvent controls (uncontaminated drilled groundwater spiked with 0.005% acetone) were used to take into account any potential solvent effect in the interpretation of the subsequent gene expression analysis. Organisms were thus exposed to three EDCs: pyriproxyfen (PYR 5 µg/L), tebufenozide (TEB 0.5 µg/L) and piperonyl butoxide (PBO 150 µg/L). Stock solutions were prepared in acetone at concentrations of 100, 10 and 3000 mg/L for PYR, TEB and PBO, respectively, adding small aliquots of these stock solutions to the dilution water. Renewal of the media was performed manually every two days, and temperature was maintained at 16 °C throughout the experiment. At the end of the exposure (days 4, 9, and 14), five females from each condition were individually sampled, rapidly weighed, and incubated overnight at 4 °C with 300 µL of RNA later. After incubation, samples were frozen in liquid nitrogen and stored at -80 °C until the RNA extraction protocol.

**Total RNA extraction and cDNA synthesis.** Whole-body organisms were disrupted with a TissueRuptor (Qiagen). For total RNA extraction, the RNeasy® Fibrous Tissue Mini Kit was used according to the manufacturer's instructions (Qiagen). DNase I was used during the protocol to remove possible genomic contaminations. Total RNA was subsequently quantified using a NanoDrop 2000 spectrophotometer (Thermo Scientific), followed by a qualitative analysis using an Agilent Bioanalyzer 2100 RNA Nanochip (Agilent technologies, Santa Clara, CA). A quantity of 1 µg of RNA from each sample was converted into cDNA using the SuperScript® III First-Strand Synthesis System (Thermo Fisher Scientific) following the manufacturer's protocol, and conserved at -20 °C until further analysis.

**Gradient PCR amplification, cloning and sequencing.** PCR was performed with a 50 µL mix of Millipore water containing 1 µL of cDNA, 2 µL of 10 mM dNTP, 2 µL of each primer at 10 µM, 5 µL of RT buffer, 3 µL of 25 mM MgCl₂, and 0.5 µL of native Taq DNA polymerase (Thermo Fisher Scientific) (200U/µL). Primers (listed in Table 2) were manually designed, based on sequence alignments of *G. fossarum* deducted sequences with other arthropod sequences of the same gene.

The following PCR program was used (Biometra TGradient): initial denaturation step at 94 °C for 3 min, followed by 39 amplification cycles (denaturation 30 s at 94 °C; annealing 30 s at 56, 58.3, 59.7, and 62 °C; elongation 2 min at 72 °C). A final extension step was made with 10 min incubation at 72 °C. The PCR products were run on a 1% agarose gel (Supplementary Figure 1).

PCR amplicons were purified using the QIAquick® PCR Purification Kit (Qiagen). Purified fragments were inserted into pGEM®-T Vectors (pGEM®-T Vector System I, Promega), and transformed into DH5α competent cells. Clones containing target fragments were isolated and submitted for direct colony PCR verification using the following program: initial denaturation step at 94 °C for 10 min, followed by 40 amplification cycles (denaturation 30 s at 94 °C; annealing 30 s at 55 °C; elongation 1 min at 72 °C). A final extension step was made with 5 min incubation at 72 °C. The PCR products were run on a 1% agarose gel. To obtain high-purity plasmids, positive clones were purified with the QIAprep® Spin Miniprep Kit (Qiagen) following the manufacturer's protocol. Samples were sent to Beckman Coulter Genomics (GENEWIZ®, UK) for sequencing.

**Quantitative real-time polymerase chain reaction (qPCR).** For target genes, the same primers were used as for gene cloning (Table 2). Primers for reference genes, also listed in Table 2, were selected from previous studies that performed qPCR studies in *Gammarus* species<sup>6,63</sup>. Three reference genes were tested: Glyceraldehyde 3-phosphate dehydrogenase (GAPDH), 18S RNA, and elongation factor 1 α (EF). In order to select a suitable reference gene, the stability of these three genes was assessed in cDNA samples from the physiological and exposure experiments. EF was considered as the most stable, presenting the lowest standard deviation for all Ct values obtained in the different samples, as shown in Table 2. The EF gene was therefore chosen to normalize expression data in both experiments.

For the quantitative PCR reactions, the iTaq™ Universal SYBR® Green Supermix (Biorad) kit was used. For the reactions, 25 ng (and 8.5 ng in the contamination experiment) of template cDNA and 300 mM of each primer were used, and the procedure was according the manufacturer's instructions. PCR was performed in a CFX96 Touch™ Real-Time PCR Detection System (Biorad), using the following program: 30 sec at 95 °C, 40 cycles of 95 °C for 5 s, 60 °C for 30 s then a temperature increment was programmed for the melting curve (65–95 °C with 0.5 °C increments at 5 s/step). Samples for each condition were a pool of five biological replicates (RNAs extracted from the individual females were pooled together for the reverse transcription and subsequent qPCR). For the physiological experiment, we also analyzed the individual biological replicates in specific stages of the female reproductive cycle. The stages were selected based on the results from the pools, with the purpose of validating the observed peaks of expression. Therefore, for E75 we analyzed individuals from the C2D1, D1D2, and D2F conditions, while for RXR we analyzed individuals from the A, C2D1, and D2F conditions. All samples were analyzed in triplicate, and controls without DNA were run in every plate. Gene expression levels were analyzed using the relative quantification method (ΔCt)<sup>64</sup>.

**Statistical analysis.** Statistical analyses for testing differences between conditions were performed with the GraphPad Prism Version 7.02 software, using unpaired t-tests corrected for multiple comparisons, according to the the Holm-Sidak method. Fold change was calculated using the 2<sup>-ΔΔCt</sup> method normalized to the control condition. Significant differences were accepted when meeting both *p* < 0.05 and fold change > 1.5.
References

1. Birnbaum, L. S. State of the Science of Endocrine Disruptors. *Environmental health perspectives* **121**, a107-a107, https://doi.org/10.1289/ehp.1306695 (2013).

2. Rodriguez, E. M., Medesani, D. A. & Fingerman, M. Endocrine disruption in crustaceans due to pollutants: a review. *Comp Biochem Physiol A Mol Integr Physiol* **146**, 661-671 (2007).

3. LeBlanc, G. A. Crustacean endocrine toxicology: A review. *Ecotoxicology* **16**, 61-81, https://doi.org/10.1007/s10646-006-0115-z (2007).

4. Cravedi, J.-P., Zalko, D., Savouret, J.-F., Menuet, A. & Jégou, B. Le concept de perturbation endocrinienne et la santé humaine. *Med Sci (Paris)* **23**, 198-204 (2007).

5. Kloas, W. et al. In *Annals of the New York Academy of Sciences* Vol. 1163, 187-200 (2009).

6. Jobling, S. & Tyler, C. R. In *Pure and Applied Chemistry* Vol. 75, 2219 (2003).

7. Oetken, M., Bachmann, J., Schulte-Oehlmann, U. & Oehlmann, J. Evidence for endocrine disruption in invertebrates. *International review of cytology* **236**, 1-44, https://doi.org/10.1016/s0074-7696(04)36001-8 (2004).

8. Höss, S. & Weltje, L. Endocrine disruption in nematodes: Effects and mechanisms. *Ecotoxicology* **16**, 15-28, https://doi.org/10.1007/s10646-006-0108-y (2007).

9. Oehlmann, J. et al. Endocrine disruption in prosobranch molluscs: Evidence and ecological relevance. *Ecotoxicology* **16**, 29-43, https://doi.org/10.1007/s10646-006-0109-x (2007).

10. Sumpter, J. P. & Jobling, S. Vitellogenesis as a biomarker for estrogenic contamination of the aquatic environment. *Environ Health Perspect* **7**, 173-178 (1995).

11. Kirby, M. F. et al. Surveys of plasma vitellogenin and intersex in male flounder (Platichthys flesus) as measures of endocrine disruption by estrogenic contamination in United Kingdom estuaries: temporal trends, 1996 to 2001. *Environ Toxicol Chem* **23**, 748-758 (2004).

12. Scott, A. P., Sanders, M., Stentiford, G. D., Reese, R. A. & Katsiadaki, I. Evidence for estrogenic endocrine disruption in an offshore flatfish, the dab (Limanda limanda L.). *Mar Environ Res* **64**, 128-148 (2007).

13. Smith, B. S. Sexuality in the american mud snail, nassarius obsoletus say. *Journal of Molluscan Studies* **39**, 377-378, https://doi.org/10.1093/oxfordjournals.mollus.a065117 (1971).

14. Jungmann, D., Ladewig, V., Ludwichowski, K.-U., Petzsch, P. & Nagel, R. Intersexuality in Gammarus fossarum KOCH - a common inducible phenomenon? *Archiv für Hydrobiologie* **159**, 511-529, https://doi.org/10.1127/0003-9136/2004/0159-0511 (2004).

15. Ford, A. T. et al. Endocrine disruption in a marine amphipod? Field observations of intersexuality and de-masculinisation. *Marine environmental research* **58**, 169-173, https://doi.org/10.1016/j.marenvres.2004.03.013 (2004).

16. Kunz, P. Y., Kienle, C. & Gerhardt, A. Gammarus spp. in aquatic ecotoxicology and water quality assessment: toward integrated multilevel tests. *Reviews of environmental contamination and toxicology* **205**, 1-76, https://doi.org/10.1007/978-1-4419-5623-1_1 (2010).

17. Geffard, O. et al. Ovarian cycle and embryonic development in Gammarus fossarum: application for reproductive toxicity assessment. *Environmental toxicology and chemistry/SETAC* **29**, 2249-2259, https://doi.org/10.1002/etc.268 (2010).

18. Coulaud, R. et al. *In situ* feeding assay with Gammarus fossarum (Crustacea): Modelling the influence of confounding factors to improve water quality biomonitoring. *Water Research* **45**, 6417-6429, https://doi.org/10.1016/j.watres.2011.09.035 (2011).

19. Xuereb, Bt, Chaumot, A., Mons, R., Garric, J. & Geffard, O. Acetylcholinesterase activity in Gammarus fossarum (Crustacea Amphipoda). Intrinsic variability, reference levels, and a reliable tool for field surveys. *Aquatic Toxicology* **93**, 225-233, https://doi.org/10.1016/j.aquatox.2009.05.006 (2009).

20. Lacaze, E. et al. DNA damage in caged Gammarus fossarum amphipods: a tool for freshwater genotoxicity assessment. *Environ Pollut* **159**, 1682-1691 (2011).

21. Simon, R. et al. Mass spectrometry assay as an alternative to the enzyme-linked immunosorbent assay test for biomarker quantitation in ecotoxicology: Application to vitellogenin in Crustacea (Gammarus fossarum). *Journal of Chromatography A* **1217**, 5109-5115, https://doi.org/10.1016/j.chroma.2010.06.015 (2010).

22. Jubeaux, G. et al. Vitellogenin-like protein measurement in caged Gammarus fossarum males as a biomarker of endocrine disruptor exposure: Inconclusive experience. *Aquatic Toxicology* **122-123**, 9-18, https://doi.org/10.1016/j.aquatox.2012.05.007 (2012).

23. Short, S., Yang, G., Kille, P. & Ford, A. T. Vitellogenin is not an appropriate biomarker of feminisation in a crustacean. *Aquatic toxicology (Amsterdam, Netherlands)* **153**, 89-97, https://doi.org/10.1016/j.aquatox.2013.11.014 (2014).

24. Armengaud, J. et al. Non-model organisms, a species endangered by proteogenomics. *J Proteomics* **105**, 5-18 (2014).

25. Trapp, J. et al. Proteogenomics of gammarus fossarum to document the reproductive system of amphipods. *Molecular and Cellular Proteomics* **13**, 3612-3625, https://doi.org/10.1074/mcp.M114.038851 (2014).

26. Hyne, R. V. Review of the reproductive biology of amphipods and their endocrine regulation: identification of mechanistic pathways for reproductive toxicants. *Environmental toxicology and chemistry / SETAC* **30**, 2647-2657, https://doi.org/10.1002/etc.673 (2011).

27. King-Jones, K. & Thummel, C. S. Nuclear receptors-a perspective from Drosophila. *Nature reviews. Genetics* **6**, 311-323 (2005).

28. Sin, Y. W. et al. Identification of putative ecdysteroid and juvenile hormone pathway genes in the shrimp Neocaridina denticulata. *General and comparative endocrinology*, https://doi.org/10.1016/j.ygcen.2014.07.018 (2014).

29. Regier, J. C., Shultz, J. W. & Kambic, R. E. Pancrustacean phylogeny: hexapods are terrestrial crustaceans and maxillopods are not monophyletic. *Proc R Soc B* **272**, https://doi.org/10.1098/rspb.2004.2917 (2005).

30. Bonneton, F., Chaumot, A. & Laudet, V. Annotation of Tribolium nuclear receptors reveals an increase in evolutionary rate of a network controlling the ecdysone cascade. *Insect biochemistry and molecular biology* **38**, 416-429, https://doi.org/10.1016/j.ibmb.2007.10.006 (2008).

31. Girish, B. P., Swetha, C. H. & Reddy, P. S. Expression of RXR, EcR, E75 and VtG mRNA levels in the hepatopancreas and ovary of the freshwater edible crab, Oziothelphusa senex senex (Fabricius, 1798) during different vitelogenic stages. *Science of Nature* **102**, https://doi.org/10.1007/s00114-015-1272-9 (2015).

32. Qian, Z. et al. Identification of ecdysteroid signaling late-response genes from different tissues of the Pacific white shrimp, Litopenaeus vannamei. *Comparative Biochemistry and Physiology - A Molecular and Integrative Physiology* **172**, 10-30, https://doi.org/10.1016/j.cbpa.2014.02.011 (2014).

33. Jubeaux, G. et al. Vitellogenin-like proteins in the freshwater amphipod Gammarus fossarum (Koch, 1835): functional characterization throughout reproductive process, potential for use as an indicator of oocyte quality and endocrine disruption biomarker in males. *Aquatic toxicology (Amsterdam, Netherlands)* **112-113**, 72-82, https://doi.org/10.1016/j.aquatox.2012.01.011 (2012).

34. Bialecki, M., Shilton, A., Fichtenberg, C., Segraves, W. A. & Thummel, C. S. Loss of the ecdysteroid-inducible E75A orphan nuclear receptor uncouples molting from metamorphosis in Drosophila. *Developmental Cell* **3**, 209-220, https://doi.org/10.1016/s1534-5807(02)00204-6 (2002).

35. Swevers, L., Eystathioy, T. & Iatrou, K. The orphan nuclear receptors BmE75A and BmE75C of the silkmoth Bombyx mori: Hornmonal control and ovarian expression. *Insect biochemistry and molecular biology* **32**, 1643-1652, https://doi.org/10.1016/s0965-1748(02)00104-2 (2002).

36. Pierceall, W. E. et al. E75 expression in mosquito ovary and fat body suggests reiterative use of ecdysone-regulated hierarchies in development and reproduction. *Molecular and Cellular Endocrinology* **150**, 73-89, https://doi.org/10.1016/s0303-7207(99)00022-2 (1999).
37. Buszczak, M. et al. Ecdysone response genes govern egg chamber development during mid-oogenesis in Drosophila. *Development* **126**, 4581–4589 (1999).

38. Dubrovskaya, V. A., Berger, E. M. & Dubrovsky, E. B. Juvenile hormone regulation of the E75 nuclear receptor is conserved in Diptera and Lepidoptera. *Gene* **340**, 171–177, https://doi.org/10.1016/j.gene.2004.07.022 (2004).

39. Guo, W. C. et al. Nuclear receptor ecdysone-induced protein 75 is required for larval-pupal metamorphosis in the Colorado potato beetle Leptinotarsa decemlineata (Say). *Insect Molecular Biology* **25**, 44–57, https://doi.org/10.1111/imb.12197 (2016).

40. Priya, T. A. J., Li, F., Zhang, J., Yang, C. & Xiang, J. Molecular characterization of an ecdysone inducible gene E75 of Chinese shrimp Fenneropenaeus chinensis and elucidation of its role in molting by RNA interference. *Comparative Biochemistry and Physiology - B Biochemistry and Molecular Biology* **156**, 149–157, https://doi.org/10.1016/j.cbpb.2010.02.004 (2010).

41. Hannas, B. R. & LeBlanc, G. A. Expression and ecdyroid responsiveness of the nuclear receptors HR3 and E75 in the crustacean Daphnia magna. *Molecular and Cellular Endocrinology* **315**, 208–218, https://doi.org/10.1016/j.mce.2009.07.013 (2010).

42. Reza, A. M. S. et al. Hormonal control of a metamorphosis-specific transcriptional factor Broad-Complex in silkworm. *Comparative Biochemistry and Physiology Part B: Biochemistry and Molecular Biology* **139**, 753–761, https://doi.org/10.1016/j.cbpc.2004.09.009 (2004).

43. Zhou, B. & Riddiford, L. M. Hormonal Regulation and Patterning of the Broad-Complex in the Epidermis and Wing Discs of the Tobacco Hornworm, Manduca sexta. *Developmental Biology* **231**, 125–137, https://doi.org/10.1006/dbio.2000.0143 (2001).

44. Yang, C. et al. The Broad Complex isoform 2 (BrC-Z2) transcriptional factor plays a critical role in vitellogenin transcription in the silkworm Bombyx mori. *Biochimica et Biophysica Acta (BBA) - General Subjects* **1840**, 2674–2684, https://doi.org/10.1016/j.bagen.2014.05.013 (2014).

45. Buaklin, A., Sittikankaew, K., Khamamtong, B., Menasveta, P. & Klinbunga, S. Characterization and expression analysis of the Broad-complex (Br-c) gene of the giant tiger shrimp Penaeus monodon. *Comparative Biochemistry and Physiology - B Biochemistry and Molecular Biology* **164**, 280–289, https://doi.org/10.1016/j.cbpb.2013.02.004 (2013).

46. Qu, Z. et al. How Did Arthropod Sesquiterpenoids and Ecdysteroids Arise? Comparison of Hormonal Pathway Genes in Noninsect Arthropod Genomes. *Genome Biology and Evolution* **7**, 1951–1959, https://doi.org/10.1093/gbe/evv120 (2015).

47. Tiu, S. H. K., Hult, E. F., Yagi, K. J. & Tobe, S. S. Farnesoic acid and methyl farnesoate production during lobster reproduction: Possible functional correlation with retinoid X receptor expression. *General and comparative endocrinology* **175**, 259–269, https://doi.org/10.1016/j.ygcen.2011.11.016 (2012).

48. Priya, T. A. J. et al. Molecular characterization and effect of RNA interference of retinoid Xreceptor (RXR) on E75 and chitinase gene expression in Chinese shrimp Fenneropenaeus chinensis. *Comparative Biochemistry and Physiology - B Biochemistry and Molecular Biology* **153**, 121–129, https://doi.org/10.1016/j.cbpb.2009.02.009 (2009).

49. Tang, J., Zhu, D. F., Cui, X. Y., Xie, X. & Qiu, X. E. Molecular cloning, characterization and expression analysis of the retinoid X receptor in the swimming crab, portunus trituberculatus (Miers, 1876) (Decapoda, Portunidae). *Crustaceana* **87**, 312–327, https://doi.org/10.1163/15685403-00003286 (2014).

50. WANG Chun-Jian, Z. D.-F., Yi-Zhou, Q. I., Ze-Hui, H. U. & Xi, X. I. E. Shen jie. molt-inhibiting hormone levels and ecdysteroid titer during a molt cycle of portunus trituberculatus. *Acta Hydrobiologica Sinica* **37**, 22–28, https://doi.org/10.7541/2013.22 (2013).

51. Gong, J. et al. The retinoid X receptor from mud crab: New insights into its roles in ovarian development and related signaling pathway. *Scientific reports* **6**, https://doi.org/10.1038/srep23654 (2016).

52. Nagaraju, G. P. C., Rajitha, B. & Borst, D. W. Molecular cloning and sequence of retinoid X receptor in the green crab Carcinus maenas: A possible role in female reproduction. *Journal of Endocrinology* **210**, 379–390, https://doi.org/10.1530/joe-11-0154 (2011).

53. Miyamoto, J., Hirano, M., Takimoto, Y. & Hatakoshi, M. In *Pest Control with Enhanced Environmental Safety* Vol. 524 ACS Symposium Series Ch. 11, 144–168 (American Chemical Society, 1993).

54. Retnakaran, A. et al. Mode of action of the ecdysone agonist tebufenozide (RH-5992), and an exclusion mechanism to explain resistance to it. *Pest Management Science* **57**, 951–957, https://doi.org/10.1002/ps.377 (2001).

55. Zimmer, C. T. et al. Use of the synergist piperonyl butoxide can slow the development of alpha-cypermethrin resistance in the whitefly Bemisia tabaci. *Insect Molecular Biology* **26**, 152–163, https://doi.org/10.1111/imb.12276 (2017).

56. Song, Y., Villeneuve, D. L., Toyota, K., Iguchi, T. & Tollefsen, K. E. Ecdysone Receptor Agonism Leading to Lethal Molting Disruption in Arthropods: Review and Adverse Outcome Pathway Development. *Environmental Science and Technology* **51**, 4142–4157, https://doi.org/10.1021/acs.est.7b00480 (2017).

57. Hardstone, M. C. et al. Development of multifunctional metabolic synergists to suppress the evolution of resistance against pyrethroids in insects that blood feed on humans. *Pest Management Science* **71**, 842–849, https://doi.org/10.1002/ps.3856 (2015).

58. Gouy, M., Guindon, S. & Gascuel, O. *SeaView version 4: A multiplatform graphical user interface for sequence alignment and phylogenetic tree building*. (Mol Biol Evol. 2010 Feb;27(2):221-4. doi: 10.1093/molbev/msp259. Epub 2009 Oct 23., 2010).

59. Edgar, R. C. MUSCLE: multiple sequence alignment with high accuracy and high throughput. *Nucleic Acids Research* **32**, 1792–1797, https://doi.org/10.1093/nar/gkh340 (2004).

60. Trifinopoulos, J., Nguyen, L.-T., von Haeseler, A. & Minh, B. Q. W-IQ-TREE: a fast online phylogenetic tool for maximum likelihood analysis. *Nucleic Acids Research* **44**, W232–W235, https://doi.org/10.1093/nar/gkw256 (2016).

61. Guindon, S. et al. New algorithms and methods to estimate maximum-likelihood phylogenies: assessing the performance of PhyML 3.0. *Syst Biol* **59**, 307–321 (2010).

62. Issartel, J., Boulo, V., Wallon, S., Geffard, O. & Charmantier, G. Cellular and molecular osmoregulatory responses to cadmium exposure in Gammarus fossarum (Crustacea, Amphipoda). *Chemosphere* **81**, 701–710 (2010).

63. Cottin, D., Foucreau, N., Hervant, F. & Piscart, C. Differential regulation of hsp70 genes in the freshwater key species Gammarus pulex (Crustacea, Amphipoda) exposed to thermal stress: effects of latitude and ontogeny. *Journal of Comparative Physiology B* **185**, 303–313, https://doi.org/10.1007/s00360-014-0885-1 (2015).

64. Livak, K. J. & Schmittgen, T. D. Analysis of relative gene expression data using real-time quantitative PCR and the 2(-Delta Delta C(T)) Method. *Methods* **25**, 402–408 (2001).

## Acknowledgements

We thank the “Institut National de Recherche en Sciences et Technologies pour l’Environnement et l’Agriculture” (Irstea, France), the “Institut de Génomique Fonctionnelle” (IGFL, France), the “Commissariat pour l’Energie Atomique et aux Energies Alternatives” (CEA, France), and the ANR program “ProteoGam” (ANR-14-CE21-0006-02) for financial support. We also want to thank Abderrahman Khila, Fabienne Gas, Virginie Nouvel, Séverine Viala and Florent Decugis for providing useful insights and technical support in the experiments performed.

## Author Contributions

D.G., A.C., F.B. and O.G. conceptualized the project and designed the experiments. D.G., H.Q., and C.A. performed the experiments and collected data. D.G., A.C., F.B., D.D.E., J.A., C.A. and O.G. contributed to data interpretation. D.G. wrote the main manuscript text which all authors have reviewed.
**Additional Information**

**Supplementary information** accompanies this paper at https://doi.org/10.1038/s41598-018-22235-7.

**Competing Interests:** The authors declare no competing interests.

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

<figure><img src="image_9.png" /></figure>

**Open Access** This article is licensed under a Creative Commons Attribution 4.0 International License, which permits use, sharing, adaptation, distribution and reproduction in any medium or format, as long as you give appropriate credit to the original author(s) and the source, provide a link to the Creative Commons license, and indicate if changes were made. The images or other third party material in this article are included in the article's Creative Commons license, unless indicated otherwise in a credit line to the material. If material is not included in the article's Creative Commons license and your intended use is not permitted by statutory regulation or exceeds the permitted use, you will need to obtain permission directly from the copyright holder. To view a copy of this license, visit http://creativecommons.org/licenses/by/4.0/.

© The Author(s) 2018