Extraordinary variability in gene activation and repression programs during gonadal sex differentiation across vertebrates.

Genes involved in gonadal sex differentiation have been traditionally thought to be fairly conserved across vertebrates, but this has been lately questioned. Here, we performed the first comparative analysis of gonadal transcriptomes across vertebrates, from fish to mammals. Our results unambiguously show an extraordinary overall variability in gene activation and repression programs without a phylogenetic pattern. During sex differentiation, genes such as dmrt1, sox9, amh, cyp19a and foxl2 were consistently either male- or female-enriched across species while many genes with the greatest expression change within each sex were not. We also found that downregulation in the opposite sex, which had only been quantified in the mouse model, was also prominent in the rest of vertebrates. Finally, we report 16 novel conserved markers (e.g., fshr and dazl) and 11 signaling pathways. We propose viewing vertebrate gonadal sex differentiation as a hierarchical network, with conserved hub genes such as sox9 and amh alongside less connected and less conserved nodes. This proposed framework implies that evolutionary pressures may impact genes based on their level of connectivity.

Exploring the Effects of Rearing Densities on Epigenetic Modifications in the Zebrafish Gonads.

Rearing density directly impacts fish welfare, which, in turn, affects productivity in aquaculture. Previous studies have indicated that high-density rearing during sexual development in fish can induce stress, resulting in a tendency towards male-biased sex ratios in the populations. In recent years, research has defined the relevance of the interactions between the environment and epigenetics playing a key role in the final phenotype. However, the underlying epigenetic mechanisms of individuals exposed to confinement remain elucidated. By using zebrafish (Danio rerio), the DNA methylation promotor region and the gene expression patterns of six genes, namely dnmt1, cyp19a1a, dmrt1, cyp11c1, hsd17b1, and hsd11b2, involved in the DNA maintenance methylation, reproduction, and stress were assessed. Zebrafish larvae were subjected to two high-density conditions (9 and 66 fish/L) during two periods of overlapping sex differentiation of this species (7 to 18 and 18 to 45 days post-fertilization, dpf). Results showed a significant masculinization in the populations of fish subjected to high densities from 18 to 45 dpf. In adulthood, the dnmt1 gene was differentially hypomethylated in ovaries and its expression was significantly downregulated in the testes of fish exposed to high-density. Further, the cyp19a1a gene showed downregulation of gene expression in the ovaries of fish subjected to elevated density, as previously observed in other studies. We proposed dnmt1 as a potential testicular epimarker and the expression of ovarian cyp19a1a as a potential biomarker for predicting stress originated from high densities during the early stages of development. These findings highlight the importance of rearing densities by long-lasting effects in adulthood conveying cautions for stocking protocols in fish hatcheries.

Single-cell-resolution transcriptome map revealed novel genes involved in testicular germ cell progression and somatic cells specification in Chinese tongue sole with sex reversal.

Female-to-male sex reversals (pseudomales) are common in lower vertebrates and have been found in natural populations, which is a concern under rapid changes in environmental conditions. Pseudomales can exhibit altered spermatogenesis. However, the regulatory mechanisms underlying pseudomale spermatogenesis remain unclear. Here, we characterized spermatogenesis in Chinese tongue sole (Cynoglossus semilaevis), a species with genetic and environmental sex determination, based on a high-resolution single-cell RNA-seq atlas of cells derived from the testes of genotypic males and pseudomales. We identified five germ cell types and six somatic cell types and obtained a single-cell atlas of dynamic changes in gene expression during spermatogenesis in Chinese tongue sole, including alterations in pseudomales. We detected decreased levels of Ca2+ signaling pathway-related genes in spermatogonia, insufficient meiotic initiation in spermatocytes, and a malfunction of somatic niche cells in pseudomales. However, a cluster of CaSR genes and MAPK signaling factors were upregulated in undifferentiated spermatogonia of pseudomales. Additionally, we revealed that Z chromosome-specific genes, such as piwil2, dhx37, and ehmt1, were important for spermatogenesis. These results improve our understanding of reproduction after female-to-male sex-reversal and provide new insights into the adaptability of reproductive strategies in lower vertebrates.

Development of epigenetic biomarkers for the identification of sex and thermal stress in fish using DNA methylation analysis and machine learning procedures.

The sex ratio is a key ecological demographic parameter crucial for population viability. However, the epigenetic mechanisms operating during gonadal development regulating gene expression and the sex ratio remain poorly understood. Moreover, there is interest in the development of epigenetic markers associated with a particular phenotype or as sentinels of environmental effects. Here, we profiled DNA methylation and gene expression of 10 key genes related to sex development and stress, including steroidogenic enzymes, and growth and transcription factors. We provide novel information on the sex-related differences and on the influence of elevated temperature on these genes in zebrafish, a species with mixed genetic and environmental influences on sex ratios. We identified both positive (e.g., amh, cyp11c and hsd11b2) and negative (e.g., cyp11a1 and dmrt1) correlations in unexposed males, and negative correlation (amh) in exposed females between DNA methylation and gene expression levels. Further, we combined DNA methylation analysis with machine learning procedures and found a series of informative CpGs capable not only of correctly identifying sex (based on cyp19a1a DNA methylation levels) but also of identifying whether males and females had been exposed to abnormally elevated temperature when young (based on amh and foxl2a DNA methylation levels, respectively). This was achieved in the absence of conspicuous morphological alterations of the gonads. These DNA methylation-based epigenetic biomarkers represent molecular resources that can correctly recapitulate past thermal history and pave the way for similar findings in other species to assess potential ecological effects of environmental disturbances in the context of climate change.

Deciphering sex-specific miRNAs as heat-recorders in zebrafish.

In the last decade, a plethora of microRNAs (miRNAs) has been reported in a wide variety of physiological processes, including reproduction, in many aquatic organisms. However, miRNAome alterations occurred by environmental cues due to water temperature increment have not yet been elucidated. With the aim to identify epigenetic regulations mediated by miRNAs in the gonads in a climate change scenario, the animal model zebrafish (Danio rerio) were subjected to high temperatures during sex differentiation, a treatment that results in male-skewed sex ratios in the adulthood. Once the fish reached adulthood, gonads were sequenced by high-throughput technologies and a total of 23 and 1 differentially expressed miRNAs in ovaries and testes, respectively, were identified two months after the heat treatment. Most of these heat-recorder miRNAs were involved in human sex-related cancer and about 400 predicted-target genes were obtained, some with reproduction-related functions. Their synteny in the zebrafish genome was, for more than half of the predicted target genes, in the chromosomes 7, 2, 4, 3 and 11 in the ovaries, chromosome 4 being the place where the sex-associated-region (sar) is localized in wild zebrafish. Further, spatial localization in the gonads of two selected heat-recorder miRNAs (miR-122-5p and miR-146-5p) showed exclusive expression in the ovarian germ cells. The present study expands the catalog of sex-specific miRNAs and deciphers, for the first time, thermosensitive miRNAs in the zebrafish gonads that might be used as potential epimarkers to predict environmental past events.

Improved biomarker discovery through a plot twist in transcriptomic data analysis.

BACKGROUND: Transcriptomic analysis is crucial for understanding the functional elements of the genome, with the classic method consisting of screening transcriptomics datasets for differentially expressed genes (DEGs). Additionally, since 2005, weighted gene co-expression network analysis (WGCNA) has emerged as a powerful method to explore relationships between genes. However, an approach combining both methods, i.e., filtering the transcriptome dataset by DEGs or other criteria, followed by WGCNA (DEGs + WGCNA), has become common. This is of concern because such approach can affect the resulting underlying architecture of the network under analysis and lead to wrong conclusions. Here, we explore a plot twist to transcriptome data analysis: applying WGCNA to exploit entire datasets without affecting the topology of the network, followed with the strength and relative simplicity of DEG analysis (WGCNA + DEGs). We tested WGCNA + DEGs against DEGs + WGCNA to publicly available transcriptomics data in one of the most transcriptomically complex tissues and delicate processes: vertebrate gonads undergoing sex differentiation. We further validate the general applicability of our approach through analysis of datasets from three distinct model systems: European sea bass, mouse, and human. RESULTS: In all cases, WGCNA + DEGs clearly outperformed DEGs + WGCNA. First, the network model fit and node connectivity measures and other network statistics improved. The gene lists filtered by each method were different, the number of modules associated with the trait of interest and key genes retained increased, and GO terms of biological processes provided a more nuanced representation of the biological question under consideration. Lastly, WGCNA + DEGs facilitated biomarker discovery. CONCLUSIONS: We propose that building a co-expression network from an entire dataset, and only thereafter filtering by DEGs, should be the method to use in transcriptomic studies, regardless of biological system, species, or question being considered.

Fish as Model Systems to Study Epigenetic Drivers in Human Self-Domestication and Neurodevelopmental Cognitive Disorders.

Modern humans exhibit phenotypic traits and molecular events shared with other domesticates that are thought to be by-products of selection for reduced aggression. This is the human self-domestication hypothesis. As one of the first types of responses to a novel environment, epigenetic changes may have also facilitated early self-domestication in humans. Here, we argue that fish species, which have been recently domesticated, can provide model systems to study epigenetic drivers in human self-domestication. To test this, we used in silico approaches to compare genes with epigenetic changes in early domesticates of European sea bass with genes exhibiting methylation changes in anatomically modern humans (comparison 1), and neurodevelopmental cognitive disorders considered to exhibit abnormal self-domestication traits, i.e., schizophrenia, Williams syndrome, and autism spectrum disorders (comparison 2). Overlapping genes in comparison 1 were involved in processes like limb morphogenesis and phenotypes like abnormal jaw morphology and hypopigmentation. Overlapping genes in comparison 2 affected paralogue genes involved in processes such as neural crest differentiation and ectoderm differentiation. These findings pave the way for future studies using fish species as models to investigate epigenetic changes as drivers of human self-domestication and as triggers of cognitive disorders.

Environmentally-induced sex reversal in fish with chromosomal vs. polygenic sex determination.

Sex ratio depends on sex determination mechanisms and is a key demographic parameter determining population viability and resilience to natural and anthropogenic stressors. There is increasing evidence that the environment can alter sex ratio even in genetically sex-determined species (GSD), as elevated temperature can cause female-to-male sex reversal (neomales). Alarmingly, neomales are being discovered in natural populations of several fish, amphibian and reptile species worldwide. Understanding the basis of neomale development is important for conservation biology. Among GSD species, it is unknown whether those with chromosomal sex determination (CSD), the most common system, will better resist the influence of high temperature than those with polygenic sex determination (PSD). Here, we compared the effects of elevated temperature in two wild zebrafish strains, Nadia (NA) and Ekkwill (EKW), which have CSD with a ZZ/ZW system, against the AB laboratory strain, which has PSD. First, we uncovered novel sex genotypes and the results showed that, at control temperature, the masculinization rate roughly doubled with the addition of each Z chromosome, while some ZW and WW fish of the wild strains became neomales. Surprisingly, we found that at elevated temperatures WW fish were just as likely as ZW fish to become neomales and that all strains were equally susceptible to masculinization. These results demonstrate that the Z chromosome is not essential for male development and that the dose of W buffers masculinization at the control temperature but not at elevated temperature. Furthermore, at the elevated temperature the testes of neomales, but not of normal males, contained more spermatozoa than at the control temperature. Our results show in an unprecedented way that, in a global warming scenario, CSD species may not necessarily be better protected against the masculinizing effect of elevated temperature than PSD species, and reveal genotype-by-temperature interactions in male sex determination and spermatogenesis.

Switches, stability and reversals in the evolutionary history of sexual systems in fish.

Sexual systems are highly diverse and have profound consequences for population dynamics and resilience. Yet, little is known about how they evolved. Using phylogenetic Bayesian modelling and a sample of 4614 species, we show that gonochorism is the likely ancestral condition in teleost fish. While all hermaphroditic forms revert quickly to gonochorism, protogyny and simultaneous hermaphroditism are evolutionarily more stable than protandry. In line with theoretical expectations, simultaneous hermaphroditism does not evolve directly from gonochorism but can evolve slowly from sequential hermaphroditism, particularly protandry. We find support for the predictions from life history theory that protogynous, but not protandrous, species live longer than gonochoristic species and invest the least in male gonad mass. The distribution of teleosts' sexual systems on the tree of life does not seem to reflect just adaptive predictions, suggesting that adaptations alone may not fully explain why some sexual forms evolve in some taxa but not others (Williams' paradox). We propose that future studies should incorporate mating systems, spawning behaviours, and the diversity of sex determining mechanisms. Some of the latter might constrain the evolution of hermaphroditism, while the non-duality of the embryological origin of teleost gonads might explain why protogyny predominates over protandry in teleosts.

Single-Cell Atlas of the Chinese Tongue Sole (Cynoglossus semilaevis) Ovary Reveals Transcriptional Programs of Oogenesis in Fish.

Oogenesis is a highly orchestrated process that depends on regulation by autocrine/paracrine hormones and growth factors. However, many details of the molecular mechanisms that regulate fish oogenesis remain elusive. Here, we performed a single-cell RNA sequencing (scRNA-seq) analysis of the molecular signatures of distinct ovarian cell categories in adult Chinese tongue sole (Cynoglossus semilaevis). We characterized the successive stepwise development of three germ cell subtypes. Notably, we identified the cellular composition of fish follicle walls, including four granulosa cell types and one theca cell type, and we proposed important transcription factors (TFs) showing high activity in the regulation of cell identity. Moreover, we found that the extensive niche-germline bidirectional communications regulate fish oogenesis, whereas ovulation in fish is accompanied by the coordination of simultaneous and tightly sequential processes across different granulosa cells. Additionally, a systems biology analysis of the homologous genes shared by Chinese tongue sole and macaques revealed remarkably conserved biological processes in germ cells and granulosa cells across vertebrates. Our results provide key insights into the cell-type-specific mechanisms underlying fish oogenesis at a single-cell resolution, which offers important clues for exploring fish breeding mechanisms and the evolution of vertebrate reproductive systems.

Transposon-induced epigenetic silencing in the X chromosome as a novel form of dmrt1 expression regulation during sex determination in the fighting fish.

BACKGROUND: Fishes are the one of the most diverse groups of animals with respect to their modes of sex determination, providing unique models for uncovering the evolutionary and molecular mechanisms underlying sex determination and reversal. Here, we have investigated how sex is determined in a species of both commercial and ecological importance, the Siamese fighting fish Betta splendens. RESULTS: We conducted association mapping on four commercial and two wild populations of B. splendens. In three of the four commercial populations, the master sex determining (MSD) locus was found to be located in a region of ~ 80 kb on LG2 which harbours five protein coding genes, including dmrt1, a gene involved in male sex determination in different animal taxa. In these fish, dmrt1 shows a male-biased gonadal expression from undifferentiated stages to adult organs and the knockout of this gene resulted in ovarian development in XY genotypes. Genome sequencing of XX and YY genotypes identified a transposon, drbx1, inserted into the fourth intron of the X-linked dmrt1 allele. Methylation assays revealed that epigenetic changes induced by drbx1 spread out to the promoter region of dmrt1. In addition, drbx1 being inserted between two closely linked cis-regulatory elements reduced their enhancer activities. Thus, epigenetic changes, induced by drbx1, contribute to the reduced expression of the X-linked dmrt1 allele, leading to female development. This represents a previously undescribed solution in animals relying on dmrt1 function for sex determination. Differentiation between the X and Y chromosomes is limited to a small region of ~ 200 kb surrounding the MSD gene. Recombination suppression spread slightly out of the SD locus. However, this mechanism was not found in the fourth commercial stock we studied, or in the two wild populations analysed, suggesting that it originated recently during domestication. CONCLUSIONS: Taken together, our data provide novel insights into the role of epigenetic regulation of dmrt1 in sex determination and turnover of SD systems and suggest that fighting fish are a suitable model to study the initial stages of sex chromosome evolution.

Unraveling the genotype by environment interaction in a thermosensitive fish with a polygenic sex determination system.

In most animals, sex determination occurs at conception, when sex chromosomes are segregated following Mendelian laws. However, in multiple reptiles and fishes, this genetic sex can be overridden by external factors after fertilization or birth. In some species, the genetic sex may also be governed by multiple genes, further limiting our understanding of sex determination in such species. We used the European sea bass (Dicentrarchus labrax) as a model and combined genomic (using a single nucleotide polymorphism chip) and transcriptomic (RNA-Sequencing) approaches to thoroughly depict this polygenic sex determination system and its interaction with temperature. We estimated genetic sex tendency (eGST), defined as the estimated genetic liability to become a given sex under a liability threshold model for sex determination, which accurately predicts the future phenotypic sex. We found evidence that energetic pathways, concerning the regulation of lipids and glucose, are involved in sex determination and could explain why females tend to exhibit higher energy levels and improved growth compared to males. Besides, early exposure to high-temperature up-regulated sox3, followed by sox9a in individuals with intermediate eGST, but not in individuals showing highly female-biased eGST, providing the most parsimonious explanation for temperature-induced masculinization. This gonadal state was maintained likely by DNA methylation and the up-regulation of several genes involved in histone modifications, including jmjd1c Overall, we describe a sex determination system resulting from continuous genetic and environmental influences in an animal. Our results provide significant progress in our understanding of the mechanisms underlying temperature-induced masculinization in fish.

Epigenetic mapping of the somatotropic axis in Nile tilapia reveals differential DNA hydroxymethylation marks associated with growth.

In vertebrates, the somatotropic axis comprising the pituitary gland, liver and muscle plays a major role in myogenesis. Its output in terms of muscle growth is highly affected by nutritional and environmental cues, and thus likely epigenetically regulated. Hydroxymethylation is emerging as a DNA modification that modulates gene expression but a holistic characterization of the hydroxymethylome of the somatotropic axis has not been investigated to date. Using reduced representation 5-hydroxymethylcytosine profiling we demonstrate tissue-specific localization of 5-hydroxymethylcytosines at single nucleotide resolution. Their abundance within gene bodies and promoters of several growth-related genes supports their pertinent role in gene regulation. We propose that cytosine hydroxymethylation may contribute to the phenotypic plasticity of growth through epigenetic regulation of the somatotropic axis.

Epigenetic mechanisms in sex determination and in the evolutionary transitions between sexual systems.

The hypothesis that epigenetic mechanisms of gene expression regulation have two main roles in vertebrate sex is presented. First, and within a given generation, by contributing to the acquisition and maintenance of (i) the male or female function once during the lifetime in individuals of gonochoristic species; and (ii) the male and female function in the same individual, either at the same time in simultaneous hermaphrodites, or first as one sex and then as the other in sequential hermaphrodites. Second, if environmental conditions change, epigenetic mechanisms may have also a role across generations, by providing the necessary phenotypic plasticity to facilitate the transition: (i) from one sexual system to another, or (ii) from one sex-determining mechanism to another. Furthermore, if the environmental change lasts enough time, epimutations could facilitate assimilation into genetic changes that stabilize the new sexual system or sex-determining mechanism. Examples supporting these assertions are presented, caveats or difficulties and knowledge gaps identified, and possible ways to test this hypothesis suggested. This article is part of the theme issue 'Challenging the paradigm in sex chromosome evolution: empirical and theoretical insights with a focus on vertebrates (Part I)'.

The first mitochondrial 5-methylcytosine map in a non-model teleost (Oreochromis niloticus) reveals extensive strand-specific and non-CpG methylation.

DNA methylation is one of the main epigenetic mechanisms that regulate gene expression in a manner that depends on the genomic context and varies considerably across taxa. This DNA modification was first found in nuclear genomes of eukaryote several decades ago and it has also been described in mitochondrial DNA. It has recently been shown that mitochondrial DNA is extensively methylated in mammals and other vertebrates. Our current knowledge of mitochondrial DNA methylation in fish is very limited, especially in non-model teleosts. In this study, using whole-genome bisulfite sequencing, we determined methylation patterns within non-CpG (CH) and CpG (CG) contexts in the mitochondrial genome of Nile tilapia, a non-model teleost of high economic importance. Our results demonstrate the presence of mitochondrial DNA methylation in this species predominantly within a non-CpG context, similarly to mammals. We found a strand-specific distribution of methylation, in which highly methylated cytosines were located on the minus strand. The D-loop region had the highest mean methylation level among all mitochondrial loci. Our data provide new insights into the potential role of epigenetic mechanisms in regulating metabolic flexibility of mitochondria in fish, with implications in various biological processes, such as growth and development.

Genomic Basis of Striking Fin Shapes and Colors in the Fighting Fish.

Resolving the genomic basis underlying phenotypic variations is a question of great importance in evolutionary biology. However, understanding how genotypes determine the phenotypes is still challenging. Centuries of artificial selective breeding for beauty and aggression resulted in a plethora of colors, long-fin varieties, and hyper-aggressive behavior in the air-breathing Siamese fighting fish (Betta splendens), supplying an excellent system for studying the genomic basis of phenotypic variations. Combining whole-genome sequencing, quantitative trait loci mapping, genome-wide association studies, and genome editing, we investigated the genomic basis of huge morphological variation in fins and striking differences in coloration in the fighting fish. Results revealed that the double tail, elephant ear, albino, and fin spot mutants each were determined by single major-effect loci. The elephant ear phenotype was likely related to differential expression of a potassium ion channel gene, kcnh8. The albinotic phenotype was likely linked to a cis-regulatory element acting on the mitfa gene and the double-tail mutant was suggested to be caused by a deletion in a zic1/zic4 coenhancer. Our data highlight that major loci and cis-regulatory elements play important roles in bringing about phenotypic innovations and establish Bettas as new powerful model to study the genomic basis of evolved changes.

A genome-wide association study, supported by a new chromosome-level genome assembly, suggests sox2 as a main driver of the undifferentiatiated ZZ/ZW sex determination of turbot (Scophthalmus maximus).

BACKGROUND: Understanding sex determination (SD) across taxa is a major challenge for evolutionary biology. The new genomic tools are paving the way to identify genomic features underlying SD in fish, a group frequently showing limited sex chromosome differentiation and high SD evolutionary turnover. Turbot (Scophthalmus maximus) is a commercially important flatfish with an undifferentiated ZW/ZZ SD system and remarkable sexual dimorphism. Here we describe a new long-read turbot genome assembly used to disentangle the genetic architecture of turbot SD by combining genomics and classical genetics approaches. RESULTS: The new turbot genome assembly consists of 145 contigs (N50 = 22.9 Mb), 27 of them representing >95% of its estimated genome size. A genome wide association study (GWAS) identified a ~ 6.8 Mb region on chromosome 12 associated with sex in 69.4% of the 36 families analyzed. The highest associated markers flanked sox2, the only gene in the region showing differential expression between sexes before gonad differentiation. A single SNP showed consistent differences between Z and W chromosomes. The analysis of a broad sample of families suggested the presence of additional genetic and/or environmental factors on turbot SD. CONCLUSIONS: The new chromosome-level turbot genome assembly, one of the most contiguous fish assemblies to date, facilitated the identification of sox2 as a consistent candidate gene putatively driving SD in this species. This chromosome SD system barely showed any signs of differentiation, and other factors beyond the main QTL seem to control SD in a certain proportion of families.

Do the Offspring of Sex Reversals Have Higher Sensitivity to Environmental Perturbations?

Sex determination systems in vertebrates vary along a continuum from genetic (GSD) to environmental sex determination (ESD). Individuals that show a sexual phenotype opposite to their genotypic sex are called sex reversals. Aside from genetic elements, temperature, sex steroids, and exogenous chemicals are common factors triggering sex reversal, a phenomenon that may occur even in strict GSD species. In this paper, we review the literature on instances of sex reversal in fish, amphibians, reptiles, birds, and mammals. We focus on the offspring of sex-reversed parents in the instances that they can be produced, and show that in all cases studied the offspring of these sex-reversed parents exhibit a higher sensitivity to environmental perturbations than the offspring of non-sex-reversed parents. We suggest that the inheritance of this sensitivity, aside from possible genetic factors, is likely to be mediated by epigenetic mechanisms such as DNA methylation, since these mechanisms are responsive to environmental cues, and epigenetic modifications can be transmitted to the subsequent generations. Species with a chromosomal GSD system with environmental sensitivity and availability of genetic sex markers should be employed to further test whether offspring of sex-reversed parents have greater sensitivity to environmental perturbations. Future studies could also benefit from detailed whole-genome data in order to elucidate the underlying molecular mechanisms. Finally, we discuss the consequences of such higher sensitivity in the context of global climate change.

Footprints of global change in marine life: Inferring past environment based on DNA methylation and gene expression marks.

Ocean global warming affects the distribution, life history and physiology of marine life. Extreme events, like marine heatwaves, are increasing in frequency and intensity. During sensitive stages of early fish development, the consequences may be long-lasting and mediated by epigenetic mechanisms. Here, we used European sea bass as a model to study the possible long-lasting effects of a marine heatwave during early development. We measured DNA methylation and gene expression in four tissues (brain, muscle, liver and testis) and detected differentially methylated regions (DMRs). Six genes were differentially expressed and contained DMRs three years after exposure to increased temperature, indicating direct phenotypic consequences and representing persistent changes. Interestingly, nine genes contained DMRs around the same genomic regions across tissues, therefore consisting of common footprints of developmental temperature in environmentally responsive loci. These loci are, to our knowledge, the first metastable epialleles (MEs) described in fish. MEs may serve as biomarkers to infer past life history events linked with persistent consequences. These results highlight the importance of subtle phenotypic changes mediated by epigenetics to extreme weather events during sensitive life stages. Also, to our knowledge, it is the first time the molecular effects of a marine heatwave during the lifetime of individuals are assessed. MEs could be used in surveillance programs aimed at determining the footprints of climate change on marine life. Our study paves the way for the identification of conserved MEs that respond equally to environmental perturbations across species. Conserved MEs would constitute a tool of assessment of global change effects in marine life at a large scale.