Impact of GC content on de novo gene birth

Domain: biology
Published: 2026-01-06
Editorial outcome: major_revision

Abstract
Noncoding regions in eukaryotes are extensively expressed and represent a significant source of novel microproteins, some of which become fixed as de novo genes. However, the structural properties of these unevolved products and the features driving their fixation remain poorly understood. Particularly, the influence of nucleotide composition (GC content) on their structural properties and evolutionary trajectories is still unclear. Here, we predict the foldability and sequence properties of millions of microproteins potentially encoded in the noncoding open reading frames (ORFs) of 3,379 eukaryotic genomes with GC contents ranging from 18% to 79%. Depending on GC content, these microproteins exhibit distinct structural properties, suggesting different cellular impacts if non-genic regions are pervasively expressed. Using phylostratigraphy, de novo gene search, and ancestral sequence reconstruction, we trace the evolution of several hundred de novo proteins across 22 organisms with varying GC contents. We show that de novo genes preferentially emerge from GC-rich ORFs with folding potential, revealing that the interplay between GC content and foldability-rooted in the structure of the genetic code-shapes the emergence of novel genes.

Introduction
Genomes can be broadly divided into coding and noncoding regions, conceptually considered as distinct entities. Noncoding regions, which consist of untranslated regions that may play regulatory roles, are primarily constrained at the nucleotide level, such as GC content (i.e., percentage of GC bases), and evolve independently of protein-related constraints. In contrast, coding regions, being translated into functional proteins, must simultaneously satisfy the constraints imposed by both DNA and proteins. Recent studies showed that (i) noncoding regions harbor vast reservoirs of small Open Reading Frames (ORFs) that can be expressed at the protein level, and (ii) the transition from noncoding to coding sequences (CDSs) is possible1–9. Notably, advances in OMIC technologies revealed that regions annotated as noncoding can be transcribed and translated into small proteins, with thousands of microproteins reported across multiple eukaryotic species3,5,10–13. Furthermore, several studies have demonstrated that transitions from noncoding to CDSs are far more common than previously thought6,14,15. Notably, hundreds of newly expressed products have been shown to be functional based on CRISPR-based genetic screens, or shown to undergo selection, thereby establishing them as part of the coding repertoire3,15–19.

De novo gene birth, which refers to the emergence of novel genes from ancestrally noncoding DNA, has been demonstrated to significantly contribute to genome evolution in various eukaryotic species14,15,20–22. Initially considered unlikely due to the stark differences between coding and noncoding regions, de novo genes have garnered increasing attention over the past decade, prompting researchers to elucidate their properties, functions, and mechanisms of emergence. These genes display peculiar characteristics: they have unique amino-acid compositions, are expressed in specific tissues or conditions, and are longer than random ORFs, while shorter than canonical genes6,23–27. Although some features of de novo genes (e.g., short length, weak expression) might also reflect nonfunctional, noisily expressed ORFs, the fact that several display tissue-specific expression28–30, undergo selection15,16,18,19, or have been experimentally shown to be functional19,31–33, supports the view that they can reasonably be classified as coding.

Taken together, these observations establish de novo gene birth as a significant driver of genetic novelty and raise a fundamental question: how does a noncoding sequence become coding? The mechanism of pseudogenization34,35, whereby protein-coding genes revert to noncoding DNA, is easier to grasp since breaking appears simpler than creating. However, de novo gene birth from scratch is harder to envision, especially given the unevolved nature (from a protein perspective) of noncoding sequences. This intriguing question can be reframed in terms of constraints faced by coding and noncoding sequences: how does a noncoding sequence governed only by DNA-level constraints transition to a CDS subject to both DNA and protein-level constraints?

Noncoding regions display distinct nucleotide distributions from coding ones, resulting in different amino-acid patterns if translated23,36. The expression of noncoding regions is, therefore, expected to generate aberrant or harmful proteins rather than functional templates for novel genes. This apparent paradox—how compositionally biased sequences give rise to functional proteins—has stimulated numerous studies, leading to various models of de novo gene birth. The “preadaptation” model posits that only sequences “preadapted” to avoid deleterious effects, such as aggregation, can transition into CDSs37. These de novo proteins, characterized by higher disorder, are therefore well-tolerated by the cell and evolve towards greater structural order over time. However, Tretyachenko et al. showed that random sequences can display secondary structure content and physicochemical properties similar to natural proteins, making them well-tolerated by the cell38. Heames et al. reported hundreds of putative de novo proteins in humans and flies that resemble random sequences in structural potential, though with a higher proportion of soluble and protease-sensitive sequences. Interestingly, their findings suggest that structure-forming potential can emerge from random sequences without prior selection39. In contrast to the preadaptation model, the proto-gene model proposes a gradual transition, where de novo genes evolve from transitory proto-genes arising from pervasive expression of non-genic sequences6. In yeast, de novo proteins were predicted to be less disordered than old ones and expected to acquire greater solubility during evolution. Recently, Tassios et al. reported a widespread potential for transmembrane (TM) domains encoded in intergenic regions of several hundred budding yeast genomes40. Interestingly, Vakirlis et al. found that ORFs with a high TM propensity could be beneficial to Saccharomyces cerevisiae41. The authors proposed a TM-first model, whereby newly translated TM peptides benefit from the membrane environment that prevents them from aggregation, and may later evolve towards greater solubility and broader cellular functions. On the other hand, using computational methods, we showed that yeast intergenic ORFs encode a wide diversity of structural elements, including aggregation-prone peptides, and a majority of ORFs with folding potential23. The latter are not necessarily expected to form well-defined 3D structures on their own but contain elementary structural blocks that may promote folding in specific contexts (e.g., upon binding to a partner or membrane, or through chaperone assistance39) and could support the evolution of more stable structures. Using ancestral sequence reconstruction, we found that the ancestral sequences of yeast de novo genes were enriched in sequences with folding potential. Similar studies in Drosophila melanogaster revealed that a majority of de novo proteins are partially or well-folded according to Alphafold2 predictions25. Interestingly, ancestral sequence analyses indicated that most of the well-folded proteins were likely born well-folded as well. While computational studies on the structural properties of de novo proteins have accumulated23,26,42,43, experimental analyses remain scarce, with only a few de novo proteins structurally characterized44–48. Notably, the yeast Bsc4 and the fly Goddard de novo proteins were shown to fold into partially ordered 3D structures containing secondary structure elements but lacking the stable tertiary structure typical of globular proteins44,47. Collectively, these findings underscore the importance of understanding how protein foldability contributes to the establishment of novel protein-coding genes, while raising questions about their generalizability across species. In particular, protein foldability depends on amino-acid composition, which in turn is influenced by DNA GC content.

In this work, we investigate how GC content impacts (i) the ORF reservoir of an organism and (ii) the transition between noncoding and CDSs. Across 3379 eukaryotic species with various GC contents, we analyze the sequence and structural properties of presumed noncoding ORFs and CDSs. Then, focusing on 22 species, we employ phylostratigraphy, de novo gene search, and ancestral sequence reconstruction to investigate the transitory stages leading to the establishment of canonical proteins. This work reveals (i) GC-content ranges that seem more favorable to de novo gene birth, and (ii) features promoting the emergence of novel genes. Our results support a unified model of de novo gene birth, integrating previous observations through the lens of GC content and specific protein features coordinated by the structure of the genetic code.

Results
Impact of GC content on CDS and iORF codons
We established a dataset of 3379 eukaryotic genomes with GC contents ranging from 17.5 to 78.5% (see Supplementary Data File 1 for the genome list). We extracted their CDSs and all stop-to-stop ORFs of at least 60 nucleotides located in inter-CDS regions (iORFs), since large fractions of microproteins do not start with methionine10,23,49. Consistent with previous studies50, GC content is strongly correlated between CDSs and noncoding regions (Pearson’s r=0.85, p < 1 × 10−16), yet remains systematically higher in CDSs (Fig. S1A). To assess GC effects on the codon frequencies of CDSs, we analyzed the distribution of various codon groups in CDSs and iORFs across the 3379 eukaryotic genomes (Fig. 1A). While some codons display comparable frequencies in CDSs and iORFs, lying on the diagonal (e.g., GTN, CAN, AAN), others appear enriched or depleted in CDSs, lying either above or below the diagonal, respectively (e.g., GCN, GAN, TTN). To discern whether codon frequencies are simply driven by the higher GC content of CDSs or reflect constraints operating at the protein level, we generated scrambled CDSs preserving nucleotide composition but disrupting codons and, consequently, protein-related constraints. Then, to quantify whether differences between real and scrambled CDSs stem solely from protein-level constraints or are also influenced by adaptive codon preferences (i.e., translational selection), we calculated the correlation between GC3 and the GC content of noncoding regions across the 3379 genomes. In line with Bénitière et al.51, we observed a strong correlation (Pearson R2=0.64, r=0.79, p < 1 × 10−16), supporting that codon usage is primarily driven by GC content (Fig. S1B). Figure 1B shows the log fold change (LFC) of each codon group representing the deviation of its frequency in real CDSs compared to scrambled ones. Figure 1C illustrates the variation in its frequency across genomes for iORFs (x-axis) and CDSs (y-axis). Variation is quantified by the standard deviation (SD) of a codon group’s frequency distribution, with high SDs indicating high variability according to GC content, and low SDs indicating more stable frequencies. Codons near the diagonal vary similarly in CDSs and iORFs, whereas those in the middle-right region of the plot (e.g., NTN, NGN) vary widely in iORFs but remain constrained in CDSs (Fig. 1C).

Several codon groups exhibit small LFCs, reflecting similar frequencies in real and scrambled CDSs (Fig. 1B). Notably, mixed AT/GC codons (e.g., TCN, CAN, or GTN) typically show little variation in both iORFs and CDSs (bottom-left corner of Fig. 1C), reflecting low sensitivity to GC content. Overall, most codons have absolute LFCs below 0.5, but four groups show stronger deviations (|LFCs|>0.75): TAN, TGN, CGN, and GAN. Their median frequencies differ significantly from those of scrambled CDSs (Mann–Whitney U test (two-sided), all p < 1 × 10−16) (Fig. S1C). TAN and TGN comprise the three stop codons and the rare amino acids tryptophan and cysteine, explaining their strong depletion in CDSs relative to expectations (Figs. 1B and S1C). CGN codons, which encode arginine, are also markedly depleted (Figs. 1B and S1C). CpG sites, consisting of cytosine followed by guanine, are prone to methylation and are hypermutable because methylated cytosines frequently undergo deamination to thymine52. Such mutations can generate TGN codons, which may introduce deleterious changes and thereby contribute to the observed depletion of CGN codons. GAN codons, encoding aspartate and glutamate, display the highest LFC, and a low frequency variation in CDSs (0.116 ± 0.005) across genomes (Fig. 1B, C). Previous studies showed that high proportions of negative charges on protein surfaces prevent promiscuous interactions with the negatively charged, abundant ribosomes53, and stretches of such residues optimize translation by facilitating the nascent peptide’s passage through the negatively charged ribosomal exit tunnel54. In yeast, we also found a correlation between negative charge frequency and protein abundance, suggesting that a higher frequency of negative charges in abundant proteins may favor their translation efficiency and reduce nonspecific interactions with ribosomes23. Finally, the NTN group, encoding most hydrophobic residues, stands out in Fig. 1C, showing the widest variation in iORFs but constrained frequencies in CDSs (0.27 ± 0.01). This suggests that selection acts to maintain optimal hydrophobicity for protein folding, preserving enough hydrophobic residues in GC-rich species for proper folding while reducing their frequency in AT-rich species to prevent aggregation. This aligns with Dill’s classic theory, which proposed that protein stability depends on an optimal fraction of hydrophobic residues55. Overall, our findings reveal codon groups whose frequencies are tightly regulated to meet protein requirements at any GC content. Interestingly, these properties, ensuring proper foldability, expression, and interactability, together constitute a hallmark of CDSs.

Intermediate GC-content genomes encode microproteins akin to CDSs
Typically, hydrophobic residues are encoded by AT-rich codons, whereas GC-rich codons predominantly encode charged and polar residues. In noncoding regions, where no selection acts on codons, AT-rich and GC-rich sequences are expected to encode preferentially hydrophobic and hydrophilic sequences, respectively. We assessed the foldability of the proteins and microproteins encoded in the CDSs and iORFs of the 3379 genomes using the Hydrophobic Cluster Analysis (HCA) method56–58 (see “Methods”). HCA is fast and does not require multiple sequence alignments, making it ideal for characterizing millions of fast-evolving iORFs. Based on the sole information of the amino acid sequence, HCA identifies (i) clusters of hydrophobic residues, which have been shown to constitute signatures of folded domains, and (ii) hydrophilic linkers connecting clusters and corresponding to loops or disordered regions. From these features, HCA computes a foldability score. Low scores indicate disordered proteins enriched in hydrophilic residues, while high scores correspond to sequences with long hydrophobic clusters, prone to aggregation or forming TM segments. Intermediate HCA scores reflect sequences with balanced proportions of hydrophobic clusters and linkers expected to be associated with folding potential23,58. The latter may either fold into well-defined 3D structures or harbor latent folding potential, i.e., the capacity to adopt structured conformations under specific conditions, and/or evolve more stable structures. Figure 2A represents the HCA scores of iORFs and CDSs according to GC content. We show that GC content has a substantial impact on iORF foldability. AT-rich and GC-rich species predominantly harbor iORFs encoding aggregation-prone or disordered microproteins, respectively, while iORFs of intermediate GC-content genomes tend to display intermediate HCA scores, reflecting a greater folding potential. As a control, we predicted the disorder and aggregation propensities of iORFs using IUPred59 and TANGO60 (Fig. S2A, B). Results align with HCA predictions, showing that disorder propensity positively correlates with GC content, while aggregation propensity anticorrelates with it. This supports the hypothesis that AT-rich genomes predominantly encode aggregation-prone iORFs, whereas GC-rich ones encode disordered ones. Consequently, depending on GC content, iORFs-encoded microproteins may sample distinct structural spaces, potentially generating ensembles with various structural properties and cellular impact if translated. Conversely, although the GC content slightly influences CDS foldability, all genomes generally encode foldable CDSs and display tighter HCA score distributions than iORFs, suggesting foldability as a key constraint in CDSs. This observation highlights the dual pressures imposed by both DNA composition and protein structure and suggests that foldability has played a pivotal role in the establishment of the genetic code, allowing CDSs to maintain foldability regardless of GC content. Finally, we show that in intermediate GC-content species, iORFs and CDSs exhibit similar folding potentials with significant overlaps in their HCA score distributions, whereas in extreme GC-content species, their foldability distributions differ markedly. Previously, we showed that yeast iORFs had HCA cluster sizes similar to CDSs23. To assess the impact of GC content on the size of HCA clusters, we analyzed their size distributions across species. Figure S2C shows that intermediate GC-content iORFs display HCA cluster sizes comparable to CDSs, while extreme AT-rich and GC-rich iORFs deviate in opposite directions (larger and smaller clusters, respectively). This suggests that intermediate GC-content species harbor a reservoir of iORFs potentially encoding the structural building blocks of proteins.

Similar trends are observed for other sequence properties, with CDSs and iORFs in extreme GC-content species exhibiting greater differences in some codon frequencies. In intermediate GC-content species, iORFs display frequencies of codons coding for positively charged residues similar to CDSs, whereas in GC-rich species, these frequencies significantly exceed those of CDSs (Fig. 2B). Although GAN codon distributions in iORFs and CDSs never overlap, the difference between iORF and CDS distributions is much greater in extreme GC-content species (Fig. 2C). Nonetheless, the consistent difference in GAN frequencies between iORFs and CDSs across all genomes suggests that while CDSs display high GAN frequencies, matching these frequencies is not crucial to enter the coding repertoire, and this criterion might be adjusted later in evolution. Finally, as expected, iORF length increases with GC content, due to the AT-richness of stop codons (Fig. S2D).

Overall, these results suggest that, aside from ORF length, iORFs in intermediate GC-content species exhibit sequence and structural characteristics more akin to CDSs than those in extreme GC-content species.

Transition from noncoding to coding sequences
We next investigated the features that could determine de novo gene birth and characterized the evolutionary intermediates leading to canonical proteins. We used DENSE, a fully automated workflow that identifies genes with strong evidence of de novo emergence61. DENSE screens the nr database for taxonomically restricted genes (TRGs) and tests whether they were ancestrally noncoding by searching for homology traces within syntenic orthologous noncoding regions of the outgroup species lacking the gene (see “Methods”). We applied DENSE to 22 well-annotated species (focal species) with well-annotated neighboring species, spanning GC contents from 19% (P. falciparum) to 64.5% (C. reinhardtii) (see Supplementary Data File 2 for the genome list). TRGs were defined at the genus level, except for H. sapiens and M. musculus, where the level was set to “Hominoidea” and “Murinae,” respectively, to ensure enough neighboring species. Table 1 summarizes the number of detected de novo genes and associated CDS isoforms across the 22 genomes. Interestingly, no de novo genes were found in H. opuntiae, and no more than two were identified in extremely AT-rich species (P. falciparum, E. histolytica, and Rotaria sp.).

Ancestral sequence reconstruction is a powerful approach to investigate the properties of the ancestral ORFs that gave rise to de novo genes23,47,62. We reconstructed the ancestral sequences of the identified de novo genes using a protocol inspired by our previous work23, as well as an alternative protocol described in Vakirlis et al.62 to ensure robustness. For each species, Fig. 3A shows the GC content distribution of its iORFs and the ancestral ORFs of its de novo genes. Only species with at least 20 reconstructed ancestral ORFs are shown (see Fig. S3A for the others species). Ancestral ORFs consistently exhibit higher GC content than the other iORFs of the same species (one-sided Mann–Whitney U test, all p < 1 × 10−3, Table S3). This trend persists even in GC-rich species, where ancestral ORFs exceed the already high GC levels of their noncoding regions. Moreover, this pattern also holds for the ancestral ORFs predicted by the PREQUEL-based protocol of Vakirlis et al.62, and for the de novo genes’ orthologous noncoding regions in the first outgroup species where the de novo gene is absent (Figs. S4 and S5). These orthologous noncoding loci also show higher GC contents, further supporting the trend predicted for ancestral sequences. Using GenEra63, a program predicting gene age through phylostratigraphy, we monitored the evolution of GC content over time. For all species, GC content appears to reach a plateau from the first phylostrata and remains strongly constrained throughout all evolutionary intermediates (Figs. 4 and S6). These findings suggest that GC-rich iORFs are more likely to give rise to novel genes, supporting the hypothesis that increased GC content is a key feature of the transition from noncoding to CDSs. However, the observed increase in GC content may also impact the properties of the resulting products, potentially influencing the de novo gene emergence process.

We thus investigated how the GC content of ancestral ORFs affects their structural properties, GAN content, and length. Figure 3B shows the HCA score distributions of the iORFs and ancestral ORFs for each species with more than 20 reconstructed ancestral ORFs (see Fig. S3B for the others species). In many species, ancestral ORFs exhibit distinct foldability from other iORFs, a trend also observed for the ancestral sequences predicted by the PREQUEL-based protocol of Vakirlis et al.62, and the de novo genes’ orthologous noncoding loci in the neighboring species (Figs. S4 and S5). The GC-content increase in ancestral ORFs seems generally associated with a decrease in the HCA scores, except for the most GC-rich species, where the iORF HCA scores are already low (Fig. 3B). This pattern reflects that GC-rich codons predominantly encode hydrophilic or disorder-promoting residues, such as proline and glycine. Nonetheless, ancestral ORFs typically remain within the foldable range of HCA scores, except in H. sapiens and M. musculus, where distributions span disordered and foldable regions due to their substantial GC-content increases. Remarkably, AT-rich species feature a high discrepancy in foldability between their ancestral ORFs and current iORFs, with ancestral ORFs displaying significantly lower HCA scores (Fig. 3B). In highly GC-rich species, ancestral iORFs remain predicted to harbor folding potential despite their GC-content increase, normally expected to promote local disorder. To determine whether ancestral ORF foldability stems solely from their GC content or is a selected trait, we scrambled their nucleotides, preserving GC content while disrupting codons and, therefore, protein-level constraints. The foldability of scrambled sequences, thus, represent the expectation based only on nucleotide composition. In several species with high levels of AT or GC (e.g., A. thaliana, S. cerevisiae, and L. infantum), ancestral ORFs show significantly different HCA scores from their scrambled counterparts. Specifically, L. infantum’s ancestral ORFs have higher HCA scores than their scrambled counterparts, while most AT-rich species shift toward higher local disorder, though remaining within the foldable bin (Figs. 3B and S3B). This suggests that selection retained iORFs with a certain level of foldability, specifically under extreme GC contexts, where foldable iORFs arise less frequently. Depending on their nucleotide composition, ancestral ORFs exhibit specific codon distributions that either further amplify the GC-content effect in AT-rich species or buffer it in GC-rich ones (Fig. S7). Specifically, ancestral ORFs from GC-rich species are usually enriched in cysteines and consistently depleted in disorder-promoting residues (glycine and proline) compared to their scrambled counterparts, potentially enhancing their folding potential. Conversely, ancestral ORFs from AT-rich species are generally depleted in hydrophobic residues and enriched in glycines and prolines, potentially promoting local disorder and solubility. This is clearly illustrated by the principal component analysis of their amino-acid propensities, which clusters species by GC content, except Plasmodium vivax, which stands apart and, interestingly, is associated with 42% GC but was ancestrally AT-rich (Fig. S7). Finally, HCA scores across phylostrata reveal that although young phylostrata (TRGs), including de novo genes, are generally predicted with folding potential, their foldability varies significantly depending on the species’ GC content (Figs. 4 and S8B, C). In contrast, older phylostrata (Eukaryotes and Cellular organisms) display lower variations in their HCA scores irrespective of species’ GC contents (F-test, p = 6.2 × 10−6) (Fig. S8C).

No clear pattern emerges in GAN codon frequency distributions of ancestral ORFs. In most species, their GAN frequency is similar to that of their scrambled sequences and that of iORFs, supporting our hypothesis that GAN frequency is not critical in the early stages of de novo gene birth (Fig. S3C). Although GAN codons may influence translation efficiency and protein interactability, de novo microproteins are typically expressed at low levels21, minimizing their immediate cellular impact. Nonetheless, our evolutionary kinetics reveal a rapid adjustment of GAN frequency, as shown with the red arrows indicating transitions to distinct evolutionary regimes (see “Methods”) (Figs. 4 and S9). The latest transition generally coincides with the onset of a plateau in GAN frequencies, usually reached at the ancestral ORFs or de novo gene levels (8/19) or at the TRG level (7/19), suggesting that GAN frequency stabilizes early in evolution. Finally, consistent with prior studies20,23,64,65, ancestral ORFs are significantly longer than other iORFs, regardless of GC content, and longer than their scrambled counterparts (Figs. 4, S10 and S11), suggesting that their increased size reflects a selected criterion rather than resulting solely from higher GC content. Moreover, in line with previous studies6,22,66,67, we confirm that CDS size evolution is a continuous process, as evidenced by the increase in ORF size across phylostrata (Figs. 4 and S10).

GC content in relation to de novo gene counts
AT-rich and GC-rich species appear to harbor fewer de novo genes than intermediate-GC species. We first examined whether technical factors could bias de novo gene detection. We searched for the divergence time for each focal species and its nearest neighbor from TimeTree68 or the literature69 (Table 1). For most species, including P. falciparum, P. vivax, and E. histolytica, which have very few detected de novo genes, divergence times are short (<6–7 Mya). Interestingly, H. opuntiae harbors only 27 orphans, despite a long divergence time with its neighbor (41 Mya), showing that interpreting patterns of molecular innovation from divergence time is complex and may be inaccurate when comparing species with different generation times and evolutionary rates. To further assess the genomic distance between focal species and their neighbors, we analyzed the homology profiles of noncoding segments across each focal species’ neighbors to see whether species with few de novo genes were associated with poor homology signals61 (Fig. S12). Species with few de novo genes were not consistently linked to poor homology profiles: E. histolytica and Rotaria sp. have only one de novo gene despite homology profiles similar to A. thaliana (99 de novo genes), whereas C. javanica harbors 63 de novo genes despite a poor homology detectability. Additionally, low-complexity sequences, frequent in extreme GC-content genomes, may also weaken homology signals between de novo candidates and their orthologous noncoding regions in the neighbors. However, Fig. S13 shows that amino-acid sequences of extreme GC-content species display complexity comparable to those of other species. These findings suggest that although technical factors may influence de novo gene detection, they alone cannot explain the low de novo gene counts observed for AT- and GC-rich species.

We, thus, investigated biological features that might underlie the discrepancies in the number of de novo genes across the 22 species. Most AT- and GC-rich species are parasitic or opportunistic and, as typical for parasites, have small genomes (<10,000 genes). The three exceptions are Rotaria sp., T. versicolor, and C. reinhardtii, which are free-living organisms with large gene repertoires despite few de novo genes. AT-rich species seem to contain even fewer de novo genes than GC-rich ones, potentially due to their more extreme GC content (20–30%) compared to the most GC-rich species, C. reinhardtii (64.5%). Furthermore, aggregation, being more detrimental than disorder, may also influence innovation strategies in extreme AT-rich species. P. falciparum displays high recombination rates in subtelomeric antigen families, generating antigenic diversity70,71. E. histolytica shows evidence of bacterial gene acquisition through horizontal gene transfer (HGT)72, and harbors abundant retrotransposons implicated in genome rearrangements73. Finally, Rotaria sp. exhibits frequent HGTs from non-metazoan sources74.

We next examined whether a phylogenetic signal exists across the 3379 species based on their GC content (Table S6). GC-rich species mostly consist of Fungi (82%), with the remaining species belonging to SAR (8%), Chlorophytes (4%), and Discoba (3.5%). In contrast, AT-rich species are more widely distributed across the tree of life, with 67% Opisthokontes (mainly Arthropods, Nematodes, and Fungi), 22% Streptophytes, and 7% SAR (almost all Apicomplexa). Intermediate GC-content species display a distinct pattern, dominated by Opisthokontes (88%)—split between Metazoa (57%) and Fungi (43%), and small proportions of Streptophytes (9%) and SAR (3%). Finally, we investigated the relationship between genome size and GC content, revealing two genome-size populations across the 3379 genomes (Fig. S14). Large genomes (>1000 Mbp) are mostly found in intermediate GC-content species, whereas GC-rich species typically have small genomes (<100 Mbp). Unlike GC-rich species, AT-rich species feature a broader range of genome sizes, with Arthropods and Nematodes associated with larger genomes, while Fungi and Apicomplexa display smaller genomes. Overall, although we observe a relationship between species phylogenetic distribution and GC content, and extreme GC-content species tend to have fewer de novo candidates, these patterns do not directly match genome size and lifestyle. Some parasites (e.g., C. javanica, L. infantum) possess many de novo genes, whereas some large-genome species have few (e.g., Rotaria sp., T. versicolor, C. reinhardtii). Conversely, small genomes, such as that of S. cerevisiae, still contain many de novo genes, highlighting a complex relationship between these factors. Altogether, these results provide a basis for future studies on how lifestyle, phylogeny, and GC content influence de novo gene emergence and whether GC content limits the ability of some species to evolve through de novo gene birth.

Methods
Collection of CDSs and iORFs
From the NCBI genomes datasets (https://www.ncbi.nlm.nih.gov/datasets/genome/), we selected all the annotated genomes with an assembly level of at least “scaffold” (contig < scaffold < chromosome < complete) as available on the 6th of March 2023. We kept only one genome per organism name (see Supplementary Data File 1 for the list of genomes). To that aim, we prioritize genomes sharing the same organism’s name as follows (decreasing priority): (1) reference/representative genome, (2) highest assembly level, (3) most recent publication date, (4) highest N50. The CDSs were directly extracted from the CDS FASTA files available on the NCBI. The iORFs were defined from stop-to-stop using the ORFtrack program from the ORFmine package88 (default parameters: 60 nucleotides, iORFs do not overlap any annotated feature (e.g., tRNA, rRNA, CDS…) by more than 70% of their length). To limit computational storage and calculations, Fig. 2 and S2 were generated using randomly selected subsets of 1000 CDSs and 1000 iORFs per genome. Similarly, the iORF distributions of the different features investigated in Figs. 3, 4, S3–S6, S8–S10 and S15 were calculated for 1000 randomly selected iORFs.

Detection of de novo genes
The prediction of de novo genes was realized using DENSE61 version 1.0.8. DENSE was run with default parameters, except for the number of outgroup species, which was set to “2.” DENSE first identifies TRGs using GenEra63. The TRG rank was set to “genus” for all species except for H. sapiens and M. musculus, where it was set to “Hominoidea” and “Murinae,” respectively, to ensure a minimum number of neighbors. For each TRG, it then tests whether the gene was ancestrally noncoding. Therefore, using TBLASTN, DENSE searches for homology traces in noncoding regions of outgroup species among the neighbors of the focal species. These species are those (i) lacking the candidate gene and (ii) branching after the last common ancestor node of the gene family. To validate TBLASTN hits, DENSE uses microsynteny to confirm that the detected noncoding hits in outgroup species are located at the same locus as the focal de novo gene candidate. Candidates with one syntenic noncoding hit in two or more outgroup species are then retained for the next filtering step. DIAMONDv281 may miss distant homologs that have diverged beyond recognition. To mitigate this risk, we screened the retained candidates against the Pfam database using HMMER89 (e value ≤ 1 × 10−3). Candidates with significant matches outside the genus of the focal species were removed from the de novo candidate set. For genes with multiple isoforms, all isoforms were screened. If any isoform yielded a significant match, the gene and its isoforms were conservatively excluded. This step eliminated 15 candidates. We also predicted the 3D structure of each de novo protein candidate using ESMFold90 and searched for structural homologs among the 15 million AlphafoldDB cluster representatives91 with Foldseek92 (e value ≤ 1 × 10−2). Of 10 matches, five overlapped the HMMER-filtered set. The remaining five had Foldseek TM-scores < 0.5, indicating weak structural alignment. TM-scores were reassessed with TM-align93, which is more accurate. One of them yielded a significant TM-score (0.64) and was conservatively removed as well. Links to the used genomes and their annotation are provided in Supplementary Data File 2.

Phylogenetic trees
Genome trees for the reconstruction of ancestral sequences and the calculation of ORF conservation profiles were constructed with OrthoFinder94 version 2.5.5. First, for every focal genome and its neighbors, we extracted the longest isoform of each gene to produce an amino acid FASTA file per genome. Orthofinder identifies orthogroups and constructs one unrooted gene tree for each orthogroup using MAFFT95 for alignments and FastTree96 for tree inference. It then produces a single unrooted species tree with STAG, subsequently rooted with STRIDE through gene tree reconciliation. The command used was “docker run davidemms/orthofinder:2.5.5.2 orthofinder -f <directory_with_input_fasta_files> -t 90 -a 4 -M msa -T fasttree -A mafft.”

Reconstructing the ancestral sequences that gave rise to de novo genes
For each de novo CDS, the DENSE workflow retrieves its best matches across the neighbor genomes that can be either coding or noncoding if no homologous CDS is found. De novo CDSs were then aligned with these homologous sequences. All aligned sequences were previously extended by 99 nucleotides upstream and downstream to include the genomic context, particularly in the case of start or stop displacement. These alignments were generated using MACSE97 version 2.07. Finally, the ancestral sequences were reconstructed using the PRANK98 program version 170427 (e.g., «prank -d=aligned.fna -t=phyml_tree.nwk -o=ouput -once -showanc -prunetree -keep”). As a control, we also applied the ancestral sequence reconstruction method PREQUEL from the PHAST99 package following the same protocol as published in Vakirlis et al.62. For both protocols, we then extracted the ancestral sequence corresponding to the most recent common ancestor between the focal sequence and the first noncoding outgroup. The presence of the second outgroup, therefore, strengthens the confidence in the inference of the ancestral state. From the three frames of the reconstructed ancestral sequences, we then extracted all ORFs defined stop-to-stop of at least 60 nucleotides. Finally, we established an optimal Smith-Waterman100 search of each ancestral ORF with the focal de novo CDS to determine those that actually gave rise to the gene of interest. This search was performed using the ssearch36 program101 (e value ≤ 1 × 10−3).

Scrambled and random sequences
Scrambled CDSs were generated by shuffling their nucleotides. To prevent any internal stop codons, when a stop codon was generated by blending, the codon and the four surrounding codons (five upstream and five downstream) were blended together to minimize bias towards codons similar to stop codons (i.e., sharing stop codons’ nucleotides). To increase statistical power, 1000 scrambled sequences were generated for each species. To compare the size of ancestral ORFs with random sequences of the same GC content (Fig. S11), we generated for each GC percentage bin (0.5%) a random DNA sequence of one million nucleotides long, spanning a range of 0% to 100% GC content. We then extracted the ORFs (≥ 20 codons) from each of these sequences. For each set of ancestral ORFs, we used, as a control, the corresponding set of ORFs with the closest GC content.

GC and HCA distributions of the noncoding regions of the first outgroup species orthologous to the predicted de novo genes
One thousand iORFs were drawn from the species that constitute the first outgroup of each detected de novo gene. These species are represented in proportion to the frequency with which they serve as the first outgroup (e.g., 200 iORFs from species X if species X is the first outgroup in 20% of cases). Orthologous intergenic regions were provided by DENSE as the intergenic regions in the first outgroup species associated with the best TBLASTN102 alignment with a given de novo CDS. We then identified the ORFs within these orthologous intergenic regions that align with the focal CDS (e value ≤ 1 × 10−3) and calculated their HCA scores. The calculation of the GC content of the orthologous noncoding region was performed on the region that aligned through TBLASTN with the nucleotide sequence of the input de novo genes.

Phylostratigraphy and subsequent analyses
The age of each CDS of every focal genome was extracted from the GenEra63 results produced during the execution of DENSE. GenEra63 screens each focal species’ proteome with DIAMONDv281 against the NCBI non-redundant protein database (downloaded on March 28, 2022). Based on the presence/absence of each CDS in the NCBI phylogeny, the CDS is then assigned an evolutionary age. The latter corresponds to the most distant node in the phylogeny where DIAMONDv2 has detected the CDS. If a gene has multiple isoforms, the isoforms are assigned the age of the oldest isoform.

To detect breakpoints (i.e., transitions) in GAN phylostratigraphies, we computed for each species the median GAN codon frequency of each phylostratum. We then fitted three Bayesian piecewise-linear models (null: no change, one breakpoint, and two breakpoints) using the mcp package103. Model performance was evaluated by leave-one-out cross-validation104, extracting the expected log-predictive density (ELPD) and its standard error. Pairwise ELPD differences were calculated, and a more complex model (e.g., two breakpoints over one breakpoint, or one breakpoint over the null model) was accepted only if its improvement (ΔELPD) exceeded the standard error. Posterior mean breakpoints were then extracted from the model estimates and mapped back onto the corresponding phylostrata.

Inter-CDS and intergenic regions
Inter-CDS and intergenic regions were respectively obtained by subtracting all the “CDS” and “genes” intervals from the whole genome GFF, using BEDTOOLS version 2.30.0105 (bedtools subtract).

Conservation profiles
Conservation profiles are made of three columns that reflect the conservation across the neighbors of each focal species: (1) the CDS of de novo genes, (2) a random subset of 1000 intergenic segments, and (3) 1000 randomly selected CDSs that were not predicted as de novo. CDSs of de novo genes: Percentage of de novo CDSs with a coding hit in each neighbor species (BLASTP, e value ≤ 1 × 10−3, query coverage by subject ≥ 50%). Noncoding segments: 1000 intergenic segments of 100 codons were randomly selected from each focal species. The percentage of intergenic segments with a hit in each neighboring species was then estimated using TBLASTN (e value ≤ 1 × 10−3, query coverage by subject ≥ 50%). CDSs: 1000 CDSs were randomly selected from each focal species. The percentage of CDSs with a coding hit in each neighboring species was then estimated using BLASTP (e value ≤ 1 × 10−3, query coverage by subject ≥ 50%).

Calculation of sequence descriptors
Sequence size and GC composition were determined using the faSize and faCount programs from the UCSC utilities (https://hgdownload.cse.ucsc.edu/admin/exe/linux.x86_64/). Codon composition was determined using the R package coRdon (https://www.bioconductor.org/packages/release/bioc/html/coRdon.html). The HCA score was calculated using the ORFold program from the ORFmine suite57,58,88. Aggregation propensity was assessed with TANGO60. Following Linding et al.60, residues within segments of more than five consecutive residues predicted to adopt a β-aggregated conformation (probability > 5%) were considered aggregation-prone. Sequence aggregation propensity was defined as the fraction of residues in such segments. Disorder propensity was evaluated with IUPred59 using the short prediction option. Residues in segments of more than 5 consecutive residues with disorder probability higher than 0.5 were considered disordered. Disorder propensity was calculated as the fraction of residues in these segments.

Identification of young iORFs under selection in S. cerevisiae
We screened all S. cerevisiae iORFs against the full set of iORFs of its five neighbors using DIAMONDv2 (e value ≤ 1 × 10−3 and coverage ≥ 70%) and retained those with significant matches in at least four neighbors. These sequences were aligned with MAFFT and tested for signature of selection using multiple CODEML branch models (from PAML106 package), following Zhang et al.107, but with a distinct subset of trees for the two-ratio models. This subset included all three versions where a single branch or clade was labeled with a distinct ratio. Chi-squared test p values were Bonferroni-corrected, with significance set at 0.05. Candidates showing signatures of selection were then screened with BLASTP and TBLASTN against the 332 yeast species from Shen et al.108. We removed those with evidence of likely ancient origin, defined as either significant BLASTP matches (e value ≤ 1 × 10−3 and coverage ≥ 70%) or TBLASTN hits with ORF intactness ≥ 70%.

Statistical analyses
Statistical analysis and graphics (ggplot2 (Wickham 2016)) were produced with R version 4.3.2 (2023-10-31)109. Figure S2 involves median comparisons between bins comprising several sets of 1000 ORFs. To address p value issues with large samples, we empirically estimated p values using an iterative sampling approach. For each comparison, we computed the observed median difference (Dobs). We then drew 10,000 random samples of 1000 individuals from each bin being compared and combined them. In each iteration, the pooled sample was split into two equal subsamples to calculate the associated median difference (Dsamp). The p value was defined as the proportion of Dsamp values exceeding Dobs. Boxplots in Figs. 3, 4, S1, S3–S6, S8–S10, S13 and S15 show the median (center line), first and third quartiles (bounds of box), and whiskers extending to the most extreme data points within 1.5 × the interquartile range (IQR): points beyond are outliers. Notches in Figs. 2 and S2 indicate an approximate 95% confidence interval for the median (1.58 × IQR/√n).

Peer review rounds
[{"round": 1, "reviewer_comments": {"reviewer_1": "In this manuscript Roginski, Papadopoulos and Lopes present a broad computational analysis of intergenic ORFs across thousands of genomes and then de novo genes and their evolution across 22 eukaryotic lineages. Multiple findings are noteworthy, but what stands out is how GC% seems to influence and constrain the potential for de novo gene emergence and the structural properties of the resulting de novo proteins. The manuscript encompasses an impressive amount of methodologically sound and careful work and is overall very well written. This work will have substantial impact on the field of novel gene evolution, and will no doubt leave a significant mark in the field of genome biology more generally. The authors beautifully build on their previous work in budding yeast and use cutting edge tools that they have developed to present a comprehensive comparison across many species. The scope and depth of the comparison and the fact that it focuses, among others, on the link between genetic code and protein structure makes this work particularly novel. I've listed my comments and some concerns below which I am sure the authors can address to make an already strong manuscript even stronger.\n1) It's important to back up as much as possible the claim that very GC-rich and very GC-poor outliers have less de novo genes. Several points to do that:\n- I understand that in Table 1, the CDS column contains all isoforms, since we see a value of 143,864 for human. Does the de novo CDS column also contain all isoforms? This should in any case be clarified in the text. It is important because it can bias downstream analyses.\n- It should be clear in the text that the TRGs are defined at the genus-level (if this is the case, and unless it's there already and I missed it) and that de novo genes are always species-specific. It also seems that the TRGs do not contain all isoforms, which seems a bit strange but maybe there's a reason, it would be good to have some clarification.\n- To further support the link between GC and de novo gene emergence, the authors could also include the numbers of all species-specific genes, not just the ones with de novo evidence. One might expect the trend to hold, especially for those species with close phylogenetic neighbors, where we don't expect many novel genes due to extreme divergence. This could also further show that detectability of orthologous noncoding regions is not an issue.\n- The authors could control for other factors that could explain/confound the difference in number of de novo genes. For example, free-living versus parasitic lifestyle (and potentially other organismal traits), evolutionary distance to closest neighbor, genome size and total number of genes. The list of the 22 species contains very different organisms with a range of genome properties. For example, do insects with a similar size genome as Dmel but extreme GC% have lower de novo genes, or even just species-specific genes? (maybe useful https://pmc.ncbi.nlm.nih.gov/articles/PMC10978632/) . I leave it up to the authors how much they would like to do on this point, as this could also be future work. In that case, it would be useful to include such potentially confounding factors in the discussion.\n2) The finding that ancestral iORFs of de novo genes have higher GC (and that they are longer) than average is a key finding. However there is an important parameter involved, and that is the process of gene finding and annotation, by which these CDSs have been identified as CDSs. Could it be that the way that genes are identified in these genomes biases them towards resembling known, conserved coding sequence GC%, i.e. ORFs that have too low GC% compared to known, conserved genes have a very low chance of being included in gene catalogues in the first place. This can artificially result in the following: out of, say, 100 de novo species-specific CDSs that truly exist in a genome, only the top 50 in terms of GC% will be called by the annotation pipeline, the rest will not be identified as genes. This in turn will mean that when the ancestral ORFs of these CDSs are examined, they will have higher GC% on average than all the other iORFs because they were detected based on that. The same bias could explain the Vakirlis et al. 2018 results as well. I would like an argument by the authors against this, or some analysis that shows that this is not the case. Alternatively this possibility could be presented in the discussion.\nOn the same point, the authors write on l. 404 \"iORFs, which gave rise to novel genes, are distinct from typical iORFs. They form a unique population characterized by higher GC content...\"\nI think the adjectives \"distinct\" and \"unique\" do not reflect the data here. Sure the distributions show clear differences, but there is still overlap. I would suggest using something milder.\n3) l. 66 \"These observations clearly demonstrate that despite their singularities, these newly emerged genes can unequivocally be classified as coding\"\nIt is important to also note that some of the features of de novo genes (i.e. short length, weak expression) can also be interpreted as them being spurious, \"noisily\" expressed ORFs, and so this needs to be considered.\n4) Figure 1 and associated analysis and comparison of codon frequencies and GCs. I have no doubt about the conclusions here, but some statistics are missing to back up the comments the authors make about the various comparisons, since some of the differences are not that clearly visible by just looking at the figure.\n5) l. 365: \"Accordingly, older phylostrata converge toward more homogeneous HCA values despite significant variations in GC contents in the corresponding CDSs\"\nI'm not sure I see that on the plot. Can this evolution of variability be statistically tested so that it's more convincing?\n6) l. 374 \"Nonetheless, our evolutionary kinetics revealed a rapid adjustment of GAN codon frequency from the youngest phylostrata\"\nAgain, statistical comparison between younger and older phylostrata is needed to support this, even if it seems to be true by looking at Figure S9.\n7) Some methodological questions, to be clarified in the text if needed.\n- For de novo gene detection, are the outgroups always 2, or at least 2?\n- Some additional details in the section \"Phylogenetic trees\" would be helpful.\n- Phylostratigraphy section: please provide a brief description about how GenEra works for the non-specialist reader.\n- \"Inter-CDS and intergenic regions\": Are other genomic elements, apart from CDSs, eg. transposable elements, ncRNA genes etc. included in the sequences that remain after the subtraction? If so, could this be affecting the iORFs in some species more than others?\n8) In Figure 5A, I'm not sure I understand the arrows. There is a pink arrow named \"codon selection\" that points to ancestral iORFs, which I find confusing. Some details in the legend about how to interpret this would be very helpful.", "reviewer_2": "The manuscript by Roginski et al. examines the impact of GC content on de novo gene birth, an important topic in gene origination, given that species in different clades may have widely varying GC contents, yet all likely possess de novo genes in their genomes. The authors identified stop-to-stop intergenic ORFs from more than 3,000 species, analyzed their hydrophobicity, and used ancestral sequence reconstruction to study the evolution of ORF hydrophobicity. The topic is both interesting and timely, though it is unclear why the authors focused solely on HCA scores, given the wide range of ordered and disordered genes in any given species, and that there is no clear agreement and empirical evidence on whether young de novo genes are ordered. Below are my comments related to the methods and results.\nMajor Comments\n1) While it is understandable that a straightforward approach to identifying ORFs is through STOP-to-STOP analysis, these are not truly iORFs. I wonder how the results would differ if the authors used a start-to-stop criterion for ORF identification. Ideally, having start-to-stop data for all species would be the best approach, but this may not be feasible within the revision timeframe. I suggest that the authors compare STOP-to-STOP ORFs with start-to-stop ORFs in a subset of species across different clades to provide insights into the impact of the methodological approach.\n2) In addition to concerns regarding STOP-to-STOP ORFs, I am curious about how GC content influences the origin and occurrence of ATG—the primary start codon in the vast majority of species.\n3) The authors should provide more details on the detection of de novo genes using DENSE. One major challenge in identifying de novo genes is distinguishing them from rapidly evolving genes. It may be beneficial to incorporate methods such as HMMER or foldseek to further exclude potential artefacts.\n4) The presentation of Figure 1 is somewhat unclear, though I believe I have understood its intended meaning. Comparisons between CDSs, iORFs, and scrambled sequences may all be significantly different—if so, what is the significance of these differences? Lines 148–181 may contain results related to Figure 1. However, without clear panel separations within the figure, it is difficult to interpret these paragraphs in relation to the data presented fully. Providing clearer figure annotations or structured panel divisions would enhance readability.\n5) Figure 2A is very interesting. It shows that there is a strong constraint on HCA scores in the long-term evolution of genes —though not the outliers. I wonder, for the de novo CDS identified in Table 1, what are the HCA scores for de novo CDS with very high and very low GC contents? Comparing the extremes to the annotated CDSs may provide insights into whether there are HCA constraints for de novo genes during birth or whether the HCAs change gradually after gene birth.\n6) In Figure 3A, the authors found that ancestral iORFs consistently exhibit higher GC content than other iORFs in the analyzed species, which is consistent with a previously published study. However, this is an unusual result that is difficult to explain biologically and may instead be a computational artefact. I strongly recommend that the authors use at least one additional method for ancestral sequence reconstruction to assess the robustness of these findings.\n7) In the paper the authors used HCA scores as a major parameter important for evolution as well as gene origination, as indicated in several figures as well as in figure 5A. However, there is no clear evidence that a functional novel gene needs to be folded, especially at the early stage.\nMinor Comments\nLine 70: In the section titled “How do we pass from one world to another”, the intended meaning is unclear. The use of the term “world” is somewhat confusing, as there is no distinct separation between different worlds in the context of gene birth, particularly concerning the reservoir of genes.\nLines 86-94: Many proto-genes tend to be AT-rich, especially given that transmembrane domains often originate in T-rich or AT-rich regions, as these regions are more likely to generate amino acids with high hydrophobicity. How does this reconcile with the authors' observations?\nLine 214-217, if the raised question not answered, it should belong to discussion.\nLines 561-564: The authors' approach of blending surrounding codons to prevent internal stop codons in scrambled sequences is commendable. However, does blending with four codons introduce any bias in codon composition within these regions?\nThe readability of the manuscript can be improved. Specifically, figure and figure panels can be cited more rigorously. Some result sections contain guesses and open questions which belong to the discussion.", "reviewer_3": "In their manuscript, “Impact of GC content on de novo gene birth,” Roginski et al. explore the influence of GC content in non- coding open reading frames on the structure of potentially encoded microproteins. Notably, as a major merit, their analysis spans an impressive dataset of 3,379 eukaryotic genomes—far exceeding the scope of previous studies, which have been limited either to specific evolutionary lineages or to smaller genome sets (e.g., 187 genomes in Basile et al., 2017). For a subset of 22 organisms, they employ a combination of ancestral sequence reconstruction (ASR)—used here more extensively than in prior work—and a streamlined, custom-built pipeline for de novo gene identification based on their prior developed tool DENSE. They compare ancestral ORFs with scrambled ORFs. This is the second major merit of the study, however I would suggest to compare this also with degenerated codons. Additionally, the authors consciously avoid reliance on 3D structural predictions, which have recently come under scrutiny for their (lack of) reliability in the context of de novo or random sequences. However, the use of hydrophobic cluster analysis (HCA) seems somehow questionable for the short ORFs (60 nucleotides, i.e., 20 amino acids) analysed because it implicietely assumes the presence of secondary structural elements, or at least the principal capacity of the sequence to form such.\nThe manuscript itself lacks engagement with foundational literature on the structural characteristics of de novo proteins and biophysics. In a few cases, the citations do not support what is written, either because they are reviews or simply not matching. In some instances, the language is overly colorful for a scientific article, and certain metaphors or illustrative phrases do not add clarity but instead may detract from the manuscript’s overall tone.\nIn the following, I am suggesting changes line-by-line:\nL1: Specify non-coding regions in eukaryotes.\nL24: Change “investigated” to “predicted” to emphasize that the study is limited to computational work.\nL29: What is actually meant by “the nature of the associated evolutionary stages”?\nL33: Please highlight in 1–2 sentences the actual results of the study.\nL43: Use “conceptually classified” instead of “long viewed” to highlight the conceptual aspect.\nL42–55: The use of words such as “materials” and “worlds” is not scientific and is uncommon.\nIt impedes readability and does not seem as conceptually intriguing as it was possibly intended to be.\nL49: Refer to Grandchamp et al. 2023 & 2024, and Iyengar & Bornberg-Bauer 2023 regarding the transition from non-coding to coding.\nL53: “Domains” can be misunderstood. Is it in the sense of Pfam?\nL55: None of these citations found positive selection. Citation 14 is a review.\nBy following the citations within it, the authors can find studies supporting their statement.\nAnother, more recent, review that deals extensively with selection can also be cited:\nZhao, L., Svetec, N., Begun, D.J. De novo genes. Annual Review of Genetics 58 (2024).\nThis paper is explicit on selection and methodological limits to detect selection in the context.\nOther primary studies that investigate selection include:\n• Zhao, L., Saelao, P., Jones, C.D., Begun, D.J. Origin and spread of de novo genes in Drosophila melanogaster\npopulations. Science 343(6172), 769–772 (2014)\n• Schmitz, J.F., Chain, F.J., Bornberg-Bauer, E. Evolution of novel genes in three-spined stickleback populations. Heredity\n125(1–2), 50–59 (2020)\n• Durand et al. Turnover of ribosome-associated transcripts from de novo ORFs produces gene-like characteristics available for\nde novo gene emergence in wild yeast populations. Genome Research (2019)\n• Vara, C., Montañés, J.C., Albà, M.M. High polymorphism levels of de novo ORFs in a Yoruba human population. Genome\nBiology and Evolution 16(7), 126 (2024)\nL64: Not all of these citations address the issues described in this sentence, particularly regarding length. A references in\nwhich this issue is\ncentrally adresses should be added:\n• Lebherz, M.K., Iyengar, B.R., Bornberg-Bauer, E. Modeling Length Changes in De Novo Open Reading Frames during\nNeutral Evolution, Genome Biology and Evolution, 16(7), July 2024, evae129\nL66: Double-check your references regarding purifying selection. “Singularities” is not descriptive, and “unequivocally” is too\nstrong.\nL69: Which results? What does “permeability” refer to? The phrase “passing from one world to another” is not scientific; see\nalso L75, which is likewise incorrect.\nL80: Replace “Their” with a specific subject. Be precise.\nL82: Same as above—what does “This” refer to?\nL83: The preadaptation model is generally flawed, as it does not explain how adaptation can occur without selection.\nNeutrality would be a more appropriate null model, especially since random sequences have been experimentally\nshown to be foldable and soluble (see Tretyachenko et al., 2017, Scientific Reports). Random sequences have\nalso been experimentally compared for disorder, solubility, and foldability\n(Heames et al. 2023, Nature Ecology & Evolution)\nL88: The authors have computationally predicted in yeast.\nL90: The transition from “less-disordered and soluble” to “transmembrane” is not smooth. If this is meant as a summary of\nstructural studies, the structure should follow the sequence:\ndisorder → structure → higher-order (TM, multimer).\nExperimental studies should be discussed here:\n• Bungard et al. 2017\n• Baalsrud et al. 2018\n• Lange et al. 2021\n• Matsuo et al. 2021\nRecent computational studies to cite:\n• Middendorf, L., Iyengar, B.R., & Eicholt, L.A. (2024). Sequence, Structure, and Functional Space of Drosophila De Novo\nProteins. GBE\n• Chen, J., Li, Q., Xia, S., Arsala, D., Sosa, D., Wang, D., Long, M. The Rapid Evolution of De Novo Proteins in Structure and\nComplex. GBE 16(6), June 2024, evae107\nL99: What is meant by “folding potential,” and how is it assessed?\nL101: Lange et al. 2021 and Vakirlis et al. 2024 also used ASR and should be cited here.\nL121: What is the length distribution? Is length the only selection criterion? This would likely bias for certain age groups.\nL123: What is being correlated with what, and why “although”? This sentence is unclear.\nL125–6: Uncommon phrasing—should be more precise. Who are “they”?\nL128: Selection by what? Evolution?\nL131: “Inherited” is confusing here.\nL132: Possibly meant: “due to evolutionary pressure to prefer certain amino acids”?\nTesting for codon bias would be useful here—note that codon preferences are generally\nstronger in lower eukaryotes (see recent papers by Laurent Duret and colleagues).\nL133: While the scrambled CDS are a great idea, comparisons should also consider degenerate codons.\nL163–4: See again:\nIyengar, B.R., Grandchamp, A., Bornberg-Bauer, E. How antisense transcripts can evolve to encode novel proteins. Nature\nCommunications 15(1), 6187 (2024).\nAlso, add statistics to support “strong.” Why is it called “modular”? What is modular?\nL167: Background in biophysics should be cited here. Foldability also depends on protein length; see work by Ken Dill and\nothers, e.g.\nSen et al. Chemical Physics Letters 405, 224–228 (2005)\nL173: Ribosomes are not freely floating in the cell.\nL178: What causes the enhancement of diffusion in the cytosol?\nL181: Avoid “world” metaphors in general. Use “a hallmark” instead of “the hallmark.”\nL187: Change to “encode preferentially…”\nL189: While the use of HCA over flawed structure predictions is appreciated, is HCA applicable to proteins of only ~20\namino acids?\nThese lack a hydrophobic core and likely do not fold in any case (see above).\nAnything shorter than 60 amino acids may give random results, unless they are extensions of existing proteins with a folding core, see work by Isabelle Callebaut and collaborators.\nL208: What is the cause and what is the consequence in the relationship between DNA and protein material?\nL215: It would be interesting to know what “these” organisms are. Is there truly no phylogenetic signal?\nEven if not, this should be demonstrated.\nL382: “Modularity” is again confusing here, since no extensions, domains, or similar features are analyzed.\nL411: See again : Iyengar, B.R., et al Nature Communications (2024).\nL414: Add substantial references here.\nL417: What exactly did they show, and where? This is unclear.\nL424: Why not try running some aggregation prediction tools?\nL432: Please explain which model of de novo gene emergence is being referred to.\nL443: Add in brackets what this species actually is.\nL460 & L478: Definitions of “effort zone,” “minimal effort,” and “GC space” seem unnecessary and may not aid in explaining\nthe results.\nL472: Strange reference formatting—check reference No. 58.\nL496: What was at odds with what in these five publications?", "reviewer_author_response": "Response Le)er\nImpact of GC content on de novo gene birth\nRoginski et al.\nFirst of all, we would like to thank the three Reviewers for their posi6ve feedback on our work and for their\nvery insigh;ul comments, which we believe have helped us strengthen our conclusions and improve the\nclarity and content of the manuscript.\nWe have made every effort to address all points raised and hope the Reviewers will find the revised version\nsubstan6ally improved as a result of their valuable comments/sugges6ons.\nBy conven6on, and to avoid any confusion, in this Response leCer, all figures from the manuscript are\nreferred to with Arabic numerals, while the figures generated specifically for this Response leCer are labeled\nwith alphabe6cal leCers and highlighted in grey. Reviewer comments are shown in grey font, our responses\nin black, and sentences quoted from the revised manuscript in blue.\nReviewer #1 (Remarks to the Author):\nIn this manuscript Roginski, Papadopoulos and Lopes present a broad computational analysis of intergenic\nORFs across thousands of genomes and then de novo genes and their evolution across 22 eukaryotic\nlineages. Multiple findings are noteworthy, but what stands out is how GC% seems to influence and\nconstrain the potential for de novo gene emergence and the structural properties of the resulting de novo\nproteins. The manuscript encompasses an impressive amount of methodologically sound and careful work\nand is overall very well written. This work will have substantial impact on the field of novel gene evolution,\nand will no doubt leave a significant mark in the field of genome biology more generally. The authors\nbeautifully build on their previous work in budding yeast and use cutting edge tools that they have\ndeveloped to present a comprehensive comparison across many species. The scope and depth of the\ncomparison and the fact that it focuses, among others, on the link between genetic code and protein\nstructure makes this work particularly novel. I've listed my comments and some concerns below which I am\nsure the authors can address to make an already strong manuscript even stronger.\nWe thank the Reviewer for their very positive and encouraging comments regarding the novelty and impact\nof our work. We have carefully considered all suggestions and concerns and revised the manuscript\naccordingly. We believe these revisions have strengthened the study, and we hope the updated version\nmeets the Reviewer’s expectations. Below, we provide a point-by-point response to each comment.\n1\n1) It's important to back up as much as possible the claim that very GC-rich and very GC-poor outliers have\nless de novo genes. Several points to do that:\n- I understand that in Table 1, the CDS column contains all isoforms, since we see a value of 143,864 for\nhuman. Does the de novo CDS column also contain all isoforms? This should in any case be clarified in the\ntext. It is important because it can bias downstream analyses.\nWe apologize for the lack of clarity. Yes, the de novo CDSs include all isoforms. To avoid any confusion, we\nhave now also added the number of de novo genes and clarified this point in the main text.\n“Table 1 summarizes the number of detected de novo genes and associated CDS isoforms across the 22\ngenomes.”\n- It should be clear in the text that the TRGs are defined at the genus-level (if this is the case, and unless it's\nthere already and I missed it) and that de novo genes are always species-specific. It also seems that the\nTRGs do not contain all isoforms, which seems a bit strange but maybe there's a reason, it would be good\nto have some clarification.\nTRGs are indeed defined at the genus level, except for H. sapiens and M. musculus, for which the level was\nset to ‘Hominoidea’ and 'Murinae', respectively, to ensure a minimum number of neighbors. This is now\nexplicitly stated in the Results and Methods sections as follows:\n“TRGs were defined at the genus level, except for H. sapiens and M. musculus, where the level was set to\n‘Hominoidea’ and ‘Murinae’, respecGvely, to ensure enough neighboring species. “\n“The TRG rank was set to ‘genus’ for all species except for H. sapiens and M. musculus, where it was set to\n'Hominoidea' and 'Murinae', respecGvely, to ensure a minimum number of neighbors.”\nTRGs, when referred to as TRGs (CDSs), do include all isoforms. We have now added to Table 1 the number\nof genes that are taxonomically restricted, referred to as TRGs. For clarity, we now report both the number\nof genes and the number of CDSs (including all isoforms) for each gene category: genes/CDSs, TRGS/TRGs\n(CDSs), Orphan genes and their associated CDSs, and de novo genes with their associated CDSs.\nDe novo genes are not always species-specific. We also consider de novo genes that are genus-specific. As\na control, we recalculated the GC content and HCA scores using only the ancestral ORFs of species-specific\nde novo genes (Figure A). The trends remain consistent when focusing solely on this subset. Although most\nde novo genes are species-specific, in D. melanogaster, the short divergence time with its neighbors would\nresult in the loss of a substantial fraction of de novo genes and associated CDSs (66 de novo CDSs from\nDrosophila de novo genes vs. 26 de novo CDSs from Drosophila melanogaster de novo genes) if only orphan\nde novo genes were considered (see Table A for orphan counts). Since the manuscript is already dense and\nthe trends are similar, we chose to consider all genus-specific de novo genes and not include this additional\nFigure in the manuscript.\nGC% 19% 29% 31% 35.5% 35.5% 38% 40.5% 41.5% 42% 42% 43.5% 51% 53% 55.5% 56.5% 57.5% 58.5% 59.5% 60.5% 64%\n#ancORFs 4 7 1 103 24 208 158 191 29 21 2173 58 118 9 3 17 14 55 4 17\nFigure A: GC content and HCA scores of the ancestral ORFs of species-specific de novo genes. A For each species with at\nleast one detected orphan de novo gene, GC content distribu9ons are shown as boxplots for 1,000 randomly selected iORFs\n(purple) and for the ancestral ORFs of the CDSs associated with species-specific de novo genes (light purple)(see Table 1 for\nthe full scien9fic names). B For the same species as in panel A, HCA score distribu9ons are shown as boxplots for 1,000\nrandomly selected iORFs (purple), the ancestral ORFs of the CDSs associated with species-specific de novo genes (light purple),\nand their corresponding scrambled version (grey). The number of reconstructed ancestral iORFs is shown at the top of each\nplot. For species with at least 30 ancestral ORFs, we compared the medians of their associated distribu9ons with those of\niORFs or scrambled sequences using a one-sided Mann-Whitney U test applying Bonferroni correc9on for mul9ple\ncomparisons. Asterisks denote significance levels: *p < 5 × 10-2, **p < 1 × 10-2, ***p < 1 × 10-3. “ns” indicates p-values > 0.05.\n“NA” indicates species with less than 30 ancestral ORFs, for which the Mann-Whitney test was thus not applicable.\n3\n%CG\nFRO\nA\n75\n50\n25\nGC% 19% 29% 31% 35.5% 35.5% 38% 40.5% 41.5% 42% 42% 43.5% 51% 53% 55.5% 56.5% 57.5% 58.5% 59.5% 60.5% 64%\n#ancORFs 4 7 1 103 24 208 158 191 29 21 2173 58 118 9 3 17 14 55 4 17\nerocs\nACH\nNA NA NA *** NA *** *** *** NA NA *** ns ns NA NA NA NA ns NA NA\n5\n0\n−5\n−10\naggreg.\nprop.\nfolding\npotential\ndisord.\nprop.\nB\nName de novo CDS orphan de novo CDS\nAtha 104 81\nCcav 13 10\nCele 28 19\nCjav 63 63\nCrei 11 11\nDmel 66 26\nEhis 1 0\nHopu 0 0\nHsap 97 91\nLinf 52 44\nMmus 77 72\nMres 4 4\nOsat 1892 1620\nPfal 2 2\nPmac 5 5\nPost 33 27\nPviv 9 9\nRccf 4 4\nRsia 1 1\nScer 150 137\nTcar 2 1\nTver 14 11\nTable A: Number of CDSs of de novo genes and species-specific de novo genes. The first column lists the abbreviated\nscientific names of each species (see Table 1 for the correspondence with full scientific names). The two remaining columns\nindicate, respectively, the total number of de novo CDSs and the number of CDSs associated with species-specific de novo\ngenes (including all isoforms in each case).\n- To further support the link between GC and de novo gene emergence, the authors could also include the\nnumbers of all species-specific genes, not just the ones with de novo evidence. One might expect the trend\nto hold, especially for those species with close phylogenetic neighbors, where we don't expect many novel\ngenes due to extreme divergence. This could also further show that detectability of orthologous noncoding\nregions is not an issue.\nThe relationship between orphan gene counts and potentially undetected de novo genes is not\nstraightforward, as it is influenced by multiple factors. For instance, the number of orphan genes is expected\nto depend on the phylogenetic distance between each focal species and its closest relative. As suggested\nby the Reviewer, we searched for the divergence times for each focal species and its nearest neighbor from\nTimeTree1 or the literature2, and added this information when available to Table 1. The observed patterns\nvary, reflecting the difficulty in interpreting the number of orphan genes.\nFor instance, Drosophila melanogaster exhibits few orphan genes (28) but many taxonomically restricted\ngenes (966), likely due to its closest relative, which diverged only ~4.6 Mya. This short evolutionary distance\nis also supported by a remarkable detectable homology across their respective noncoding regions, as shown\nby the noncoding segment homology profile calculated for D. melanogaster and its neighboring species\n(Figure S12).\nOn the other hand, some AT- or GC-rich species with few de novo candidates (such as P. falciparum, E.\nhistolytica, M. restricta, T. caries, and Rhodotorula sp.) also have relatively few orphan genes, supporting\nthe hypothesis that de novo genes are indeed rare in these species. In contrast, other AT- or GC-rich species\nlike P. macacae, Rotaria sp., T. versicolor, and C. reinhardtii harbor many orphans (Table 1). However, the\n4\nhigh orphan counts do not necessarily indicate undetected de novo genes, but may instead point to\nalternative mechanisms such as horizontal gene transfer (HGT), retrotransposition, or structural genome\nrearrangements, that can generate orphan gene gains independently of de novo emergence.\nFor example:\nPneumocystis species have expanded several gene families (e.g., msg-A1, A3, D), with strong host-specific\ndiversification, reflecting rapid adaptation to mammalian hosts through gene duplication and diversification\nrather than de novo emergence3.\nRotaria species (bdelloid rotifers) show extensive HGT from bacteria, fungi, and protists, with estimates\nsuggesting that 10–15% of their genes are of non-metazoan origin, making HGT a major source of their\norphan gene content4,5.\nC. reinhardtii exhibits a highly dynamic genome characterized by frequent transposon activity, structural\nvariation, and numerous retrogenes, all of which offer alternative mechanisms for generating novel genes\nbeyond de novo birth6–8.\nWhile the respective contribution of these aforementioned mechanisms remains challenging to quantify,\nwe caution against interpreting orphan counts as direct proxies for undetected de novo gene count\nestimates, especially in taxa where alternative gene innovation mechanisms exist. This important question\nis now examined and presented along with our response to the next Reviewer’s comment in a new Results\nsection of the manuscript (see next comment).\n- The authors could control for other factors that could explain/confound the difference in number of de\nnovo genes. For example, free-living versus parasitic lifestyle (and potentially other organismal traits),\nevolutionary distance to closest neighbor, genome size and total number of genes. The list of the 22 species\ncontains very different organisms with a range of genome properties. For example, do insects with a similar\nsize genome as Dmel but extreme GC% have lower de novo genes, or even just species-specific genes?\n(maybe useful https://pmc.ncbi.nlm.nih.gov/articles/PMC10978632/) . I leave it up to the authors how\nmuch they would like to do on this point, as this could also be future work. In that case, it would be useful\nto include such potentially confounding factors in the discussion.\nWe thank the Reviewer for these insightful suggestions, all of which we have fully implemented and\nincorporated into the revised manuscript, except for the analysis of other insects – see below. Although our\ndataset is substantial for the analyses conducted on de novo genes and their ancestral sequences in relation\nto GC content, it remains limited to draw definitive conclusions about how these factors influence species’\nmolecular innovation strategies, including de novo gene birth. Nevertheless, the requested analyses are\npresented and discussed in the Results section, as they raise interesting hypotheses and important\nquestions that, altogether, provide the basis for future research\nWe first examined whether technical factors could bias de novo gene detection. Specifically, we searched\nfor the divergence time of each focal species and its nearest neighbor from TimeTree1 or the literature2\nwhen available. Except for H. opuntiae and C. elegans, divergence times with the focal species’ closest\nneighbors are typically short (< 6-7 Mya). This also applies to P. falciparum, P. vivax, and E. histolytica, for\nwhich we were able to detect only very few de novo genes. The divergence time for H. opuntiae is much\nlonger (41 Mya). Although this could have explained our difficulty in detecting de novo genes in H. opuntiae,\nthe species has only 27 orphan genes, showing that interpreting patterns of molecular innovation from\ndivergence time is complex and may be inaccurate when comparing species with different generation times\nand evolutionary rates.\nhomology profiles of noncoding segments across the neighbors of each focal species9 to determine whether\nspecies with few de novo genes were associated with poor homology signals, which could explain our\ninability to retrieve the homologous noncoding regions of de novo gene candidates (Figure S12). Yet the\nsignal is not straightforward: species with few de novo genes are not consistently linked to poor homology\nprofiles. For example, E. histolytica and Rotaria sp. have only one de novo gene each, despite homology\nprofiles comparable to those of A. thaliana (99 de novo genes), revealing significant differences in the\nnumber of detected candidates despite similar detection conditions. Conversely, C. javanica harbors 63 de\nnovo genes despite a poor homology profile. Additionally, low-complexity sequences, frequent in extreme\nGC-content genomes, may also weaken homology signals between de novo candidates and their\northologous noncoding regions in the neighbors. However, Figure S13 shows that amino-acid sequences of\nextreme GC-content species display complexity comparable to those of other species. Overall, these\nanalyses show that the relationship between the number of detected de novo genes and the proximity of\nthe focal’s neighbors is complex, and although technical factors may influence de novo gene detection,\nthey alone cannot explain the low de novo gene counts observed for AT- and GC-rich species.\nAs suggested by the Reviewer, we searched for other insect species with high AT or GC content to assess\nwhether they harbor fewer de novo genes than D. melanogaster. However, at the time of revision (June\n2025), we were unable to identify additional AT- or GC-rich insects in NCBI with comparable genome size,\nlifestyle, and a sufficient number of sequenced and annotated neighbors. The only potentially “suitable”\nspecies was Apis mellifera (32% GC content). Although A. mellifera has fewer de novo genes than D.\nmelanogaster (28 vs. 54, respectively), this comparison is not conclusive because A. mellifera has a smaller\ngenome (9,935 vs. 13,900 genes, respectively) and a drastically different effective population size, as only\nthe queen transmits genes while most individuals are sterile workers. Robust conclusions would require\nclades containing closely related species with contrasting GC contents and enough examples for statistical\npower. Therefore, we chose not to elaborate on this case in the manuscript, as this very interesting question\ncalls for more comprehensive comparative analyses for proper interpretation.\nWe next investigated biological features that might underlie the discrepancies in the number of de novo\ngenes across the 22 focal species. Most AT- and GC-rich species are parasitic or opportunistic and, as typical\nfor parasites, have small genomes with fewer than 10,000 genes. The three exceptions are Rotaria sp., T.\nversicolor, and C. reinhardtii, which are free-living organisms with remarkably large gene repertoires despite\nlow counts of de novo genes. Although our dataset is too small for definitive conclusions, AT-rich species\nseem to contain even fewer de novo genes than GC-rich ones. One should note that AT-rich species under\nstudy are more extreme (20-30% GC content) than GC-rich ones, where the most GC-rich species in the\ndataset, C. reinhardtii, reaches “only” 64.5% GC content. Furthermore, aggregation is known to be more\ndetrimental than disorder, which may influence innovation strategies in extreme AT-rich species. Notably,\nP. falciparum displays high recombination rates in subtelomeric antigen families, generating antigenic\ndiversity10–12. E. histolytica shows evidence of bacterial gene acquisitions through horizontal gene\ntransfer13\nand harbors abundant retrotransposons implicated in genome rearrangements14. Finally, Rotaria sp.\nexhibits frequent HGTs from non-metazoan sources4. While it remains difficult to draw general conclusions,\nthese studies raise the hypothesis that these species rely on alternative mechanisms (e.g., recombination,\nhorizontal gene transfer…) to drive innovation.\nWe next examined whether a phylogenetic signal exists across the 3,379 species based on their GC content.\nNew Table S1 shows that GC-rich species mostly consist of Fungi (82%), with the remaining species\nbelonging to SAR (8%), Chlorophytes (4%), and Discoba (3.5%). In contrast, AT-rich species are more widely\ndistributed across the tree of life, with 67% Opisthokontes (mainly Arthropods, Nematodes, and Fungi),\n22% Streptophytes, and 7% SAR (almost all Apicomplexa). Intermediate GC-content species display a\ndistinct pattern, dominated by Opisthokontes (88%) - split between Metazoa (57%) and Fungi (43%), and\nsmall proportions of Streptophytes (9%) and SAR (3%). Finally, we investigated the relationship between\ngenome size and GC content, revealing two genome-size populations across the 3,379 genomes (new Figure\nS14). Large genomes (> 1,000 Mbp) are mostly found in intermediate GC-content species, whereas GC-rich\nspecies typically have small genomes (< 100 Mbp). Unlike GC-rich species, AT-rich species feature a broader\n6\nrange of genome sizes, with Arthropods and Nematodes associated with larger genomes, while Fungi and\nApicomplexa display smaller genomes.\n10000\n1000\n100\n10\n1\n0 25 50 75 100\nGenome GC content (%)\nFigure S14: Genome size with respect to genome GC content. Scatter plot showing the size (Mbp) of each of the 3,379\ngenomes plotted against their GC content. Genomes are colored by main taxonomical groups (see Table S1 for a detailed list\nof the represented groups). The two vertical lines mark 35% and 50% GC content.\nOverall, although we observe a relationship between species phylogenetic distribution and GC content, and\nextreme GC-content species tend to have fewer de novo candidates, these patterns do not directly match\ngenome size and lifestyle. Some parasites (e.g., C. javanica, L. infantum) possess many de novo genes,\nwhereas some large-genome species have few (e.g., Rotaria sp., T. versicolor, C. reinhardtii). Conversely,\nsmall genomes, such as that of S. cerevisiae, still contain many de novo genes, highlighting a complex\nrelationship between these factors.\nThese observations, nevertheless, raise exciting questions and provide a basis for future studies to further\nexplore how species’ lifestyle, phylogeny, and GC content influence de novo gene emergence and whether\nGC content limits the ability of some species to evolve through de novo gene birth.\nFinally, for each focal species, we have added genome size and divergent time from its closest neighbor\n(when available) to Table 1. We have also included the distribution of the 3,379 species across the tree of\nlife according to their GC-content bin (Table S1), their associated genome sizes in relation to GC content\n(Figure S14), and added a new subsection in the Results section describing these findings.\n“AT-rich and GC-rich species appear to harbor fewer de novo genes than intermediate-GC species. We first\nexamined whether technical factors could bias de novo gene detection. We searched for the divergence time\nfor each focal species and its nearest neighbor from TimeTree68 or the literature69 (Table 1). For most species,\nincluding P. falciparum, P. vivax, and E. histolytica, which have very few detected de novo genes, divergence\ntimes are short (<6–7 Mya). Interestingly, H. opuntiae harbors only 27 orphans, despite a long divergence\ntime with its neighbor (41 Mya), showing that interpreting patterns of molecular innovation from divergence\n7\nelacs\ngol\n−\n)pbM(\nezis\nemoneG\nAmoebozoa Discoba Metazoa Sar\nClass\nChlorophyta Fungi Other Streptophyta\ntime is complex and may be inaccurate when comparing species with different generation times and\nevolutionary rates. To further assess the genomic distance between focal species and their neighbors, we\nanalyzed the homology profiles of noncoding segments across each focal species' neighbors to see whether\nspecies with few de novo genes were associated with poor homology signals60 (Figure S12). Species with few\nde novo genes were not consistently linked to poor homology profiles: E. histolytica and Rotaria sp. have\nonly one de novo gene despite homology profiles similar to A. thaliana (99 de novo genes), whereas C.\njavanica harbors 63 de novo genes despite a poor homology detectability. Additionally, low-complexity\nsequences, frequent in extreme GC-content genomes, may also weaken homology signals between de novo\ncandidates and their orthologous noncoding regions in the neighbors. However, Figure S13 shows that\namino-acid sequences of extreme GC-content species display complexity comparable to those of other\nspecies. These findings suggest that although technical factors may influence de novo gene detection, they\nalone cannot explain the low de novo gene counts observed for AT- and GC-rich species.\nWe, thus, investigated biological features that might underlie the discrepancies in the number of de novo\ngenes across the 22 species. Most AT- and GC-rich species are parasitic or opportunistic and, as typical\nfor parasites, have small genomes (<10,000 genes). The three exceptions are Rotaria sp., T. versicolor, and C.\nreinhardtii, which are free-living organisms with large gene repertoires despite few de novo genes. AT-rich\nspecies seem to contain even fewer de novo genes than GC-rich ones, potentially due to their more extreme\nGC-content (20-30%) compared to the most GC-rich species, C. reinhardtii (64.5%). Furthermore, aggregation,\nbeing more detrimental than disorder, may also influence innovation strategies in extreme AT-\nrich species. P. falciparum displays high recombination rates in subtelomeric antigen families, generating\nantigenic diversity70,71. E. histolytica shows evidence of bacterial gene acquisitions through horizontal gene\ntransfer (HGT)72, and harbors abundant retrotransposons implicated in genome rearrangements73. Finally,\nRotaria sp. exhibits frequent HGTs from non-metazoan sources74.\nWe next examined whether a phylogenetic signal exists across the 3,379 species based on their GC content\n(Table S1). GC-rich species mostly consist of Fungi (82%), with the remaining species belonging to SAR (8%),\nChlorophytes (4%), and Discoba (3.5%). In contrast, AT-rich species are more widely distributed across the\ntree of life, with 67% Opisthokontes (mainly Arthropods, Nematodes, and Fungi), 22% Streptophytes, and\n7% SAR (almost all Apicomplexa). Intermediate GC-content species display a distinct pattern, dominated by\nOpisthokontes (88%) - split between Metazoa (57%) and Fungi (43%), and small proportions of\nStreptophytes (9%) and SAR (3%). Finally, we investigated the relationship between genome size and GC\ncontent, revealing two genome-size populations across the 3,379 genomes (Figure S14). Large genomes (>\n1,000 Mbp) are mostly found in intermediate GC-content species, whereas GC-rich species typically have\nsmall genomes (< 100 Mbp). Unlike GC-rich species, AT-rich species feature a broader range of genome sizes,\nwith Arthropods and Nematodes associated with larger genomes, while Fungi and Apicomplexa display\nsmaller genomes. Overall, although we observe a relationship between species phylogenetic distribution\nand GC content, and extreme GC-content species tend to have fewer de novo candidates, these patterns do\nnot directly match genome size and lifestyle. Some parasites (e.g., C. javanica, L. infantum) possess many\nde novo genes, whereas some large-genome species have few (e.g., Rotaria sp., T. versicolor, C. reinhardtii).\nConversely, small genomes, such as that of S. cerevisiae, still contain many de novo genes, highlighting a\ncomplex relationship between these factors. Altogether, these results provide a basis for future studies on\nhow lifestyle, phylogeny, and GC content influence de novo gene emergence and whether GC content limits\nthe ability of some species to evolve through de novo gene birth. “\n2) The finding that ancestral iORFs of de novo genes have higher GC (and that they are longer) than average\nis a key finding. However there is an important parameter involved, and that is the process of gene finding\nand annotation, by which these CDSs have been identified as CDSs. Could it be that the way that genes are\nidentified in these genomes biases them towards resembling known, conserved coding sequence GC%, i.e.\nORFs that have too low GC% compared to known, conserved genes have a very low chance of being included\nin gene catalogues in the first place. This can artificially result in the following: out of, say, 100 de novo\nspecies-specific CDSs that truly exist in a genome, only the top 50 in terms of GC% will be called by the\nannotation pipeline, the rest will not be identified as genes. This in turn will mean that when the ancestral\nORFs of these CDSs are examined, they will have higher GC% on average than all the other iORFs because\n8\nthey were detected based on that. The same bias could explain the Vakirlis et al. 2018 results as well. I\nwould like an argument by the authors against this, or some analysis that shows that this is not the case.\nAlternatively this possibility could be presented in the discussion.\nWe thank the Reviewer for this very insightful comment, which raises an important point regarding\npotential biases introduced by genome annotation protocols, especially in the context of de novo gene\ndiscovery. We fully agree that annotation remains a non-trivial task, particularly for distinguishing between\nyoung coding sequences from noncoding ORFs.\nTo explore this issue, we focused on Saccharomyces cerevisiae, whose genome is arguably the best\nannotated among eukaryotes. To minimize annotation biases arising from heterogeneous annotation\npipelines - an issue we previously discussed in Roginski et al.9 - we prioritized genomes annotated using\nstandardized approaches, such as those from NCBI or specific genome consortia when possible (e.g.,\nCaenorhabditis, Oryza sativa). In the case of S. cerevisiae, we deliberately used the NCBI annotation rather\nthan the manually curated SGD version to maintain as much consistency as possible with its related species.\nFor neighboring genomes, we relied on the NCBI S. paradoxus genome, which was available at the start of\nthe project. For the other Saccharomyces species, we used the published annotations at that time: S.\narboricola from Yue et al. 201715, and S. mikatae, S. kudriavzevii, and S. uvarum from Scannell et al. 201116.\nOne should note that the manual annotation of S. cerevisiae, available in the Saccharomyces Genome\nDatabase (SGD, version R64-5-1), includes 723 more annotated ORFs than the NCBI version that we used.\nHowever, most of them (683) are labeled \"Dubious\", meaning, according to SGD, that they are “unlikely to\nencode an expressed protein” and, therefore, may not correspond to true protein-coding genes.\nInterestingly, in the work of Vakirlis et al.17, 87% of the ORFs classified as Dubious were predicted as\n“Emerging” based on the annotation version used in their study (Supplemental Data 1 file from Vakirlis et\nal.17). The authors also note that “there is no evidence that emerging ORFs correspond to canonical protein-\ncoding genes”. In any case, Figure B shows that Dubious ORFs are longer and GC-richer than iORFs in general\n(Mann-Whitney U test (one-sided), all P < 1x10-16) and exhibit GC content similar to de novo and canonical\nCDSs. While we chose not to consider them for the aforementioned reasons (unlikely to represent actual\nprotein-coding genes and to maintain consistency with neighboring species’ annotations), their properties\nsuggest that including them is unlikely to affect our conclusions based on NCBI annotation. Nevertheless,\nwe cannot exclude the possibility that shorter ORFs with atypical nucleotide composition were also missed\nin the SGD annotation.\nFigure B: GC content and length of Dubious ORFs. A GC-content distributions for iORFs (purple), Dubious ORFs from SGD\nannotation (grey), de novo CDSs predicted from NCBI annotation (orange), and CDSs from NCBI annotation (dark orange). B\nSame as panel A for ORF length.\nNCBI annotation. Therefore, we screened all S. cerevisiae iORFs against the full set of iORFs of its five\nneighbors using DIAMONDv2 (e-value threshold £ 1x10-3 and coverage threshold ³ 70%) and retained those\nwith significant matches in at least four neighbors. We next tested these ORFs for signature of selection and\nidentified 29 ORFs (new Table S2). We screened them with BLASTP and TBLASTN against the 332 yeast\nspecies from Shen et al.18. We removed six ORFs that showed evidence of likely ancient origin, defined as\neither significant BLASTP matches (e-value £ 1x10-3 and coverage ³ 70%) or TBLASTN hits with ORF\nintactness ³ 70% (new Table S2). We then reconstructed the ancestral sequences of the remaining 23\ncandidates and evaluated their GC content and size. Although their GC enrichment was less pronounced\nthan that of the 227 ancestral ORFs reconstructed based on the NCBI annotation, these sequences still\ndiffer from iORFs in both GC content and length, remaining longer and GC-richer (Mann-Whitney U test (one-\nsided), P = 2.7x10-7 and 2.1x10-2, respectively) (new Figure S15). While we chose to rely on automated\nannotation pipelines for consistency with the analyses performed on the other 21 focal species (since\nmining all their iORFs under selection, as done here for yeast, would be computationally intractable), these\npotentially overlooked candidates are unlikely to affect our general conclusion. They are few in number and\nremain overall longer and GC-richer than iORFs. However, it is worth noting that longer sequences are\ninherently more likely to yield significant matches when searching for homologs and are, therefore, more\nlikely to be retained when searching for candidates under selection\n1.00\n0.75\n*\n0.50\n0.25\n0.00\nFigure S15: GC content and ORF length of iORFs and the ancestral ORFs of iORFs under selection. A GC-content distributions\nfor iORFs (purple) and ancestral ORFs of the 23 iORFs predicted under selection and of recent origin (i.e., no homolog where\nthe ORF is intact beyond Saccharomyces – see Table S2 for details) (grey). B Length (number of codons) distributions for iORFs\n(purple) and the ancestral ORFs of the 23 iORFs predicted under selection (grey). Medians of the distribu9ons were compared\nusing a one-sided Mann-Whitney U test. Asterisks denote significance levels: *p < 5 × 10-2, **p < 1 × 10-2, ***p < 1 × 10-3.\nAs a further control, we considered the recent study by Wacholder et al.19, which compiled 412 ribosome\nprofiling datasets from S. cerevisiae across diverse conditions, including stress, as well as various wild-type\nand mutant strains. Their analysis framework and the unprecedented compiled dataset greatly improved\nthe detection sensitivity of translated ORFs and provided the most comprehensive translational landscape\nof yeast to date. By doing so, the authors identified about 19,000 translated noncanonical ORFs. A key\nstrength of this study is that it does not rely on ORF length as a filtering criterion, instead considering all\npotential ORFs genome-wide, regardless of their size. Although longer ORFs are statistically more likely to\nbe detected due to greater read coverage, the identified translated noncanonical ORFs spanned a broad\n10\n%CG\nA\n150\n***\n100\n50\n0\n)snodoc(\nhtgneL\nB\niORFs ancestral ORFs of iORFs under selection\nrange of lengths (see Figure 2C in Wacholder et al.19). They then established a comparative genomics\nframework across different evolutionary scales (e.g., within and between species) to identify noncanonical\nORFs under purifying selection. Interestingly, only 14 showed such signatures, potentially representing\noverlooked de novo gene candidates. These 14 ORFs are, on average, longer and GC-richer than typical\niORFs (Figure C). Among them, nine have significant BLASTP matches outside the Saccharomyces genus,\nsuggesting an ancient origin (see Table S5 in Wacholder et al.19). We, therefore, focused on the five\nremaining candidates of recent origin and successfully reconstructed the ancestral sequences for two of\nthem. Ancestral sequence reconstruction requires at least three sequences to provide an outgroup, and the\nother three candidates lacked sufficient homologs. The reconstructed ancestral sequences were longer\nthan the median iORF size (61 and 130 residues) and were associated with GC contents (47% and 42%)\ncomparable to the ancestral ORFs reconstructed in our study.\n100\n75\n**\n50\n25\nfor iORFs (purple) and the 14 ORFs predicted under selection in Wacholder et al. (dark purple). B Length (number of codons)\ndistributions for iORFs (purple) and the 14 ORFs predicted under selection in Wacholder et al. (dark purple). Medians of the\ndistribu9ons were compared using a one-sided Mann-Whitney U test. Asterisks denote significance levels: *p < 5 × 10-2, **p\n< 1 × 10-2, ***p < 1 × 10-3.\nThese two additional analyses suggest that, regardless of their exact status (which remains challenging to\ninfer), whether emerging, coding, or relics of pseudogenes whose ancestral origin is no longer detectable,\nthese potentially overlooked ORFs are few in number and, overall, would not alter the main trends if\nincluded in our analyses. Yet, it should be noted that neither annotation pipelines nor experimental\nmethods can fully capture the complexity of the ORF repertoire. Notably, Vakirlis et al. identified 28\nunannotated ORFs with adaptive potential in yeast through large-scale overexpression screens17. These\nthymine-rich ORFs often encode transmembrane domains and show no observable fitness impact upon\ndisruption, raising questions about their evolutionary fate. The authors interestingly discussed the\npossibility that beneficial effects observed under laboratory conditions may not have the same effect in\nnatural environments, particularly since their screen involved artificially high expression levels. They further\npoint out that adaptive potential (i) does not guarantee fixation, since some ORFs may be beneficial under\ncertain conditions but deleterious in others, and (ii) may involve tradeoffs that hinder long-term retention.\nTheir work underscores the complexity of annotating such ORFs and highlights the possibility that some\nrelevant and/or adaptive ORFs may escape detection/annotation. It also raises important questions about\nthe evolutionary fate and genomic status of these potentially overlooked ORFs, especially in light of two\nkey observations: CDSs across nearly all eukaryotic genomes are significantly more GC-rich than noncoding\nregions (Figure S1A), and transmembrane domains remain underrepresented in proteomes20.\n11\n)%(\ntnetnoc\nCG\n100\n**\n75\n50\n25\n0\n)snodoc\nfo\nrebmun(\nhtgneL\niORFs conserved nORFs from Wacholder et al. 2023\nIn summary, while annotation strategies can introduce biases, we used S. cerevisiae, the best-annotated\neukaryotic genome, to assess their potential influence on our conclusions. Our analyses suggest that the\neffect of these possible annotation biases is likely limited. We now address this point explicitly in the\nDiscussion, as we believe it raises important and timely questions for future work.\n“Investigating newly emerged genes can be affected by annotation, since annotation protocols often rely\non prior knowledge, thereby potentially overlooking ORFs with unusual patterns. To assess this, we focused\non S. cerevisiae, one of the best-annotated genomes, and searched for potential additional coding ORFs that\nmay have been missed during annotation. We screened yeast iORFs against those of its five neighboring\nspecies, retained those with DIAMONDv281 matches in at least four neighbors, and tested them for\nsignatures of purifying selection. We identified 23 ORFs under selection, specific to Saccharomyces (Table\nS2), and reconstructed their ancestral sequences. Although their GC enrichment was less pronounced than\nthat of the 227 ancestral ORFs reconstructed based on the S. cerevisiae annotation, these sequences still\ndiffer from iORFs in both GC content and length, remaining longer and GC-richer (Mann-Whitney U test\n(one-sided), P = 2.7x10-7 and 2.1x10-2, respectively)(Figure S15). As an additional control, we considered\nthe study by Wacholder et al.4, who compiled 412 ribosome profiling datasets from S. cerevisiae, producing\nthe most comprehensive yeast translational landscape to date. Their annotation-independent approach\nidentified approximately 19,000 translated unannotated ORFs. However, only 14 showed signatures of\nselection, nine of which had homologs outside Saccharomyces and were likely ancient4. Among the five\nremaining candidates, three were lacking homologs for ancestral reconstruction. Interestingly, the two\nreconstructed ancestral sequences are longer (61 and 130 residues) and GC-richer (47% and 42%) than\ntypical iORFs.\nThese results suggest that potentially overlooked de novo genes are not numerous and, overall, would not\ncontradict our conclusions if included in our analyses. Yet, it should be noted that neither annotation\npipelines nor experimental methods can fully capture ORF diversity. Notably, Vakirlis et al. identified 28\nunannotated ORFs with adaptive potential in yeast through large-scale overexpression screens41. These\nthymine-rich ORFs often encode transmembrane domains and show no fitness impact upon disruption. The\nauthors interestingly discussed the possibility that beneficial effects observed under laboratory conditions\nmay not apply in natural environments, and that adaptive potential does not guarantee fixation, as trade-\noffs may hinder long-term retention. Their work underscores the complexity of ORF annotation and the\npossibility that some adaptive ORFs escape detection, and raises important questions about the evolutionary\nfate of these ORFs regarding two key observations: CDSs in nearly all eukaryotes are GC-richer than\nnoncoding regions (Figure S1A), and transmembrane domains remain underrepresented in proteomes82.”\nOn the same point, the authors write on l. 404 \"iORFs, which gave rise to novel genes, are distinct from\ntypical iORFs. They form a unique population characterized by higher GC content...\"\nI think the adjectives \"distinct\" and \"unique\" do not reflect the data here. Sure the distributions show clear\ndifferences, but there is still overlap. I would suggest using something milder.\nWe agree with the reviewer that the adjectives \"distinct\" and \"unique\" were too strong, as the distributions,\nwhile clearly different, still show some overlap. During our revision of the Discussion to incorporate\nadditional points raised by the reviewers, we needed to reduce the overall length of the section. As a result,\nwe removed this sentence, thereby also addressing the concern about wording.\n3) l. 66 \"These observations clearly demonstrate that despite their singularities, these newly emerged genes\ncan unequivocally be classified as coding\"\nIt is important to also note that some of the features of de novo genes (i.e. short length, weak expression)\ncan also be interpreted as them being spurious, \"noisily\" expressed ORFs, and so this needs to be\nconsidered.\naccordingly, which now reads:\n“Although some features of de novo genes (e.g., short length, weak expression) might also reflect\nnonfunctional, noisily expressed ORFs, the evidence that several display tissue-specific expression28–30,\nundergo selection15,16,18,19, or have been experimentally shown to be functional19,31–33, supports the view\nthat they can reasonably be classified as coding.”\n4) Figure 1 and associated analysis and comparison of codon frequencies and GCs. I have no doubt about\nthe conclusions here, but some statistics are missing to back up the comments the authors make about the\nvarious comparisons, since some of the differences are not that clearly visible by just looking at the figure.\nWe apologize for the lack of clarity and statistical analyses supporting the interpretations, and thank the\nReviewer for this helpful feedback. In response, we now provide statistical analyses to support our claims,\nwhich also allowed us to better describe the trends observed in the initial plot. Specifically, we calculated\nthe log fold change (LFC) of each codon group to quantify its frequency deviation in real CDSs compared to\nscrambled ones. These results are now presented in a new panel (Figure 1B – see below), enabling a clearer\nclassification of codon groups that was more difficult to infer from the qualitative analysis of the initial plot.\nIn addition, we quantified the variation in codon group frequencies across genomes for iORFs (x-axis) and\nCDSs (y-axis). This analysis, now presented by a new panel (Figure 1C – see below), explicitly measures the\nextent to which the frequency of a codon group varies across genomes in iORFs and CDSs. Variation was\nassessed by the standard deviations (SDs) of each codon group’s frequency distribution, with high SDs\nindicating high variability according to GC content, and low SDs indicating more stable frequencies. Codons\nnear the diagonal vary similarly in CDSs and iORFs, whereas those in the middle-right region of the plot (e.g.,\nNTN, NGN) vary widely in iORFs but remain constrained in CDSs (Figure 1C).\nFinally, we now report the distribution of codon group frequencies in real and scrambled CDSs across the\n3,379 genomes as boxplots (new panel C of Figure S1) and include statistical comparisons in the text where\nappropriate. The revised text now reads:\n“Figure 1B shows the log fold change (LFC) of each codon group representing the deviation of its frequency\nin real CDSs compared to scrambled ones. Figure 1C illustrates the variation in its frequency across genomes\nfor iORFs (x-axis) and CDSs (y-axis). Variation is quantified by the standard deviation (SD) of a codon group’s\nfrequency distribution, with high SDs indicating high variability according to GC content, and low SDs\nindicating more stable frequencies. Codons near the diagonal vary similarly in CDSs and iORFs, whereas\nthose in the middle-right region of the plot (e.g., NTN, NGN) vary widely in iORFs but remain constrained in\nCDSs (Figure 1C).\nSeveral codon groups exhibit small LFCs, reflecting similar frequencies in real and scrambled CDSs (Figure\n1B). Notably, mixed AT/GC codons (e.g., TCN, CAN, or GTN) typically show little variation in both iORFs and\nCDSs (bottom-left corner of Figure 1C), reflecting low sensitivity to GC content. Overall, most codons have\nabsolute LFCs below 0.5, but four groups show stronger deviations (LFCs > 0.75): TAN, TGN, CGN, and GAN.\nTheir median frequencies differ significantly from those of scrambled CDSs (Mann-Whitney U test, (two-\nsided), all P < 1 × 10-16)(Figure S1C). TAN and TGN comprise the three stop codons and the rare amino acids\ntryptophan and cysteine, explaining their strong depletion in CDSs relative to expectations (Figure 1B, Figure\nS1C). CGN codons, which encode arginine, are also markedly depleted (Figure 1B, Figure S1C). CpG sites,\nconsisting of cytosine followed by guanine, are prone to methylation and are hypermutable because\nmethylated cytosines frequently undergo deamination to thymine51. Such mutations can generate TGN\ncodons, which may introduce deleterious changes and thereby contribute to the observed depletion of CGN\ncodons. GAN codons, encoding aspartate and glutamate, display the highest LFC, and a low frequency\nvariation in CDSs (0.116 ± 0.005) across genomes (Figure 1B-C). Previous studies showed that high\nproportions of negative charges on protein surfaces prevent promiscuous interactions with the negatively\ncharged, abundant ribosomes52, and stretches of such residues optimize translation by facilitating the\nnascent peptide’s passage through the negatively charged ribosomal exit tunnel53. In yeast, we also found\na correlation between negative charge frequency and protein abundance, suggesting that a higher\nfrequency of negative charges in abundant proteins may favor their translation efficiency and reduce\nnonspecific interactions with ribosomes23. Finally, the NTN group, encoding most hydrophobic residues,\nstands out in Figure 1C, showing the widest variation in iORFs but constrained frequencies in CDSs (0.27 ±\n0.01). This suggests that selection acts to maintain optimal hydrophobicity for protein folding, preserving\nenough hydrophobic residues in GC-rich species for proper folding while reducing their frequency in AT-rich\nspecies to prevent aggregation. This aligns with Dill’s classic theory, which proposed that protein stability\ndepends on an optimal fraction of hydrophobic residues54.”\nacross the 3,379 genomes. Posi9ve values indicate enrichment in real CDSs, nega9ve values indicate deple9on in real CDSS.\nError bars represent standard devia9on. Detailed distribu9ons for each codon group are shown in Supplementary Figure S1C.\nC The SDs on the x-axis and y-axis reflect the spread of codon group frequencies across the 3,379 genomes in iORFs and CDSs,\nrespec9vely. Codon groups are colored according to their absolute LFC value, with light dots corresponding to codon groups\nwith small frequency devia9ons compared to scrambled CDSs, and dark dots indica9ng codon groups whose CDS frequency\ndeviates from that of scrambled CDSs.\ndistribu9ons were compared using a one-sided Mann-Whitney U test. Asterisks denote significance levels: *p < 5 × 10-2, **p\n< 1 × 10-2, ***p < 1 × 10-3.\n5) l. 365: \"Accordingly, older phylostrata converge toward more homogeneous HCA values despite\nsignificant variations in GC contents in the corresponding CDSs\"\nI'm not sure I see that on the plot. Can this evolution of variability be statistically tested so that it's more\nconvincing?\nWe agree that the trend was not clear in Figure S8 because it required comparing all plots simultaneously\nand matching corresponding phylostrata across species. In response to the Reviewer’s comment, we added\na new panel (C) to Figure S8 showing the distribution of median HCA scores for young phylostrata (TRGs)\nacross all species, along with that of the median HCA scores for the oldest phylostrata shared by all studied\nspecies (Eukaryotes and Cellular organisms). This representation highlights that HCA score medians vary\nmore across young phylostrata than across the two oldest ones, with TRGs showing greater variance (F-\ntest, P = 5.6×10⁻⁶).\nWe have introduced the variance comparison (F-test) in the main text and incorporated the corresponding\ndistributions as new panel C of Figure S8.\n“Finally, HCA scores across phylostrata reveal that although young phylostrata (TRGs), including de novo\ngenes, are generally predicted with folding potential, their foldability varies significantly depending on the\nspecies’ GC content (Figure 4, Figure S8B-C). In contrast, older phylostrata (Eukaryotes and Cellular\norganisms) display lower variations in their HCA scores irrespective of species’ GC contents (F-test, P =\n5.6x10-6)(Figure S8C).”\n10\n5\n0\n−5\n−10\nPanel C of Figure S8: Median HCA score distributions for young and old phylostrata. Distributions of the median HCA scores\nfor young phylostrata from all species and for the two oldest phylostrata. Young phylostrata include levels from de novo\ngenes to genus (except for H. sapiens and M. musculus, where they extend to ‘Hominoidea’ and ‘Murinae’, respec9vely). The\nvariances of the two distributions were compared using an F-test (P = 5.6x10-6).\n6) l. 374 \"Nonetheless, our evolutionary kinetics revealed a rapid adjustment of GAN codon frequency from\nthe youngest phylostrata\"\nAgain, statistical comparison between younger and older phylostrata is needed to support this, even if it\nseems to be true by looking at Figure S9.\nWe thank the Reviewer for raising this point. We have now introduced a new test to identify breakpoints\nwithin a sequence of values for each species – specifically, the medians of the GAN distribution of each\nphylostratum. To do so, we fitted three Bayesian piecewise-linear models (null: no breakpoint, one\nbreakpoint, and two breakpoints) using the mcp package21. Model performance was evaluated by leave-\none-out cross-validation22, extracting the expected log-predictive density (ELPD) and its standard error.\nPairwise ELPD differences were calculated, and a more complex model (e.g., two breakpoints over one\nbreakpoint, or one breakpoint over null model) was accepted only when its improvement (DELPD) exceeded\nthe combined standard error. Posterior mean breakpoints were then extracted from the model estimates\nand mapped back onto the corresponding phylostrata for the 19 species in which at least one breakpoint\nwas detected.\nThis analysis allowed us to pinpoint breakpoints across GAN phylostratigraphies, i.e., the phylostrata\nmarking the transition to a distinct evolutionary regime. The latest breakpoint generally coincides with the\nonset of a plateau in GAN frequencies. In most cases, this occurs at the ancestral ORFs or de novo levels\n(8/19) or at the TRG level (7/19), suggesting that GAN frequency plateaus early in evolution. The exceptions\nconsist of H. sapiens and the two plant species (A. thaliana and O. sativa), where the latest breakpoint\noccurs slightly later, at Boreoeutheria, Brassicaceae, and Commelinids, respectively.\nFigure 4 and S9 have been updated to display the localization of these breakpoints along the GAN kinetics.\nThe text now reads as:\n“Although GAN codons may influence translaGon efficiency and protein interactability, de novo\nmicroproteins are typically expressed at low levels21, minimizing their immediate cellular impact.\nNonetheless, our evoluGonary kineGcs reveal a rapid adjustment of GAN frequency, as shown with the red\narrows indicaGng transiGons to disGnct evoluGonary regimes (see Methods)(Figure 4 and S9). The latest\ntransiGon generally coincides with the onset of a plateau in GAN frequencies, usually reached at the\nancestral ORFs or de novo gene levels (8/19) or at the TRG level (7/19), suggesGng that GAN frequency\nstabilizes early in evoluGon.”\n16\nerocs\nACH\nP. falciparum | 19 GC%\nNA\n10\n5\n0\n−5\n−10\nerocs\nACH\nE. histolytica | 24 GC%\nNA\n10\n5\n0\n−5\n−10\nerocs\nACH\nPneumocystis sp. | 29 GC%\nNA 10\n5\n0\n−5\n−10\nerocs\nACH\nRotaria sp. | 31 GC%\n10 ***\n5\n0\n−5\n−10\nerocs\nACH\nA. thaliana | 35.5 GC%\n10 ***\n5\n0\n−5\n−10\nerocs\nACH\nC. elegans | 35.5 GC%\n10 ***\n5\n0\n−5\n−10\nerocs\nACH\nS. cerevisiae | 38 GC%\n10 ***\n5\n0\n−5\n−10\nerocs\nACH\nH. sapiens | 40.5 GC%\n10 ***\n5\n0\n−5\n−10\nerocs\nACH\nM. musculus | 41.5 GC%\n10 ***\n5\n0\n−5\n−10\nerocs\nACH\nD. melanogaster | 42 GC%\nNA 10\n5\n0\n−5\n−10\nerocs\nACH\nP. vivax | 42 GC%\n10 ***\n5\n0\n−5\n−10\nerocs\nACH\nO. sativa | 43.5 GC%\nns\n10\n5\n0\n−5\n−10\nerocs\nACH\nP. ostreatus | 51 GC%\nns\n10\n5\n0\n−5\n−10\nerocs\nACH\nC. javanica | 53 GC%\nNA\n10\n5\n0\n−5\n−10\nerocs\nACH\nM. restricta | 55.5 GC%\nNA\n10\n5\n0\n−5\n−10\nerocs\nACH\nT. caries | 56.5 GC%\nNA\n10\n5\n0\n−5\n−10\nerocs\nACH\nT. versicolor | 57.5 GC%\nNA\n10\n5\n0\n−5\n−10\nerocs\nACH\nC. cavernicola | 58.5 GC%\n10 *\n5\n0\n−5\n−10\nerocs\nACH\nL. infantum | 59.5 GC%\nNA 10\n5\n0\n−5\n−10\nerocs\nACH\nRhodotorula sp. | 60.5 GC%\nNA 10\n5\n0\n−5\n−10\nerocs\nACH\nA\nC. reinhardtii | 64 GC%\n5\n0\n−5\n−10 o\nv o d\ne\nn ol\nd\nerocs\nACH\nB 10\n*** 5\n0\n−5\n−10 g n d\no\nu ol\ny\nerocs\nACH\nC\n0.3\nNA\n0.2\n0.1\n0.0\nFigure S9: Phylostratigraphies of GAN codon frequencies. For each species, GAN codon frequency distributions are shown\nfor iORFs (1,000 randomly selected iORFs) (purple), ancestral ORFs (light purple), and CDSs binned by their corresponding\nphylostratum (shades of orange). Only species with at least one detected de novo gene are represented (H. opuntiae was\nthus excluded). The first phylostratum corresponds to de novo CDSs (light yellow), while the last one corresponds to CDSs\nshared by all cellular organisms (dark orange). Phylostrata are defined according to NCBI taxonomy and vary in number\ndepending on the taxonomic resolution around each focal species. Red arrows mark transitions in GAN frequencies (see\nMethods). Medians of the distributions of GAN codon frequencies were compared between iORFs and ancestral ORFs using\na one-sided Mann-Whitney U test with Bonferroni correction. “NA” indicates species with less than 30 ancestral iORFs.\n17\nycneuqerf\nsnodoc\nNAG\nP. falciparum | 19 GC%\n0.3\nNA\n0.2\n0.1\n0.0\nycneuqerf\nsnodoc\nNAG\nE. histolytica | 24 GC%\n0.3\nNA\n0.2\n0.1\n0.0\nycneuqerf\nsnodoc\nNAG\nPneumocystis sp. | 29 GC%\n0.3\nNA 0.2\n0.1\n0.0\nycneuqerf\nsnodoc\nNAG\nRotaria sp. | 31 GC%\n0.3\n0.2 ***\n0.1\n0.0\nycneuqerf\nsnodoc\nNAG\nA. thaliana | 35.5 GC%\n0.3\n0.2 *\n0.1\n0.0\nycneuqerf\nsnodoc\nNAG\nC. elegans | 35.5 GC%\n0.3\nns\n0.2\n0.1\n0.0\nycneuqerf\nsnodoc\nNAG\nS. cerevisiae | 38 GC%\n0.3\nns\n0.2\n0.1\n0.0\nycneuqerf\nsnodoc\nNAG\nH. sapiens | 40.5 GC%\n0.3\nns\n0.2\n0.1\n0.0\nycneuqerf\nsnodoc\nNAG\nM. musculus | 41.5 GC%\n0.3\n0.2 **\n0.1\n0.0\nycneuqerf\nsnodoc\nNAG\nD. melanogaster | 42 GC%\n0.3\nNA 0.2\n0.1\n0.0\nycneuqerf\nsnodoc\nNAG\nP. vivax | 42 GC%\n0.3\n0.2 ***\n0.1\n0.0\nycneuqerf\nsnodoc\nNAG\nO. sativa | 43.5 GC%\n0.3\nns\n0.2\n0.1\n0.0\nycneuqerf\nsnodoc\nNAG\nP. ostreatus | 51 GC%\n0.3\nns\n0.2\n0.1\n0.0\nycneuqerf\nsnodoc\nNAG\nC. javanica | 53 GC%\n0.3\nNA\n0.2\n0.1\n0.0\nycneuqerf\nsnodoc\nNAG\nM. restricta | 55.5 GC%\n0.3\nNA\n0.2\n0.1\n0.0\nycneuqerf\nsnodoc\nNAG\nT. caries | 56.5 GC%\n0.3\nNA\n0.2\n0.1\n0.0\nycneuqerf\nsnodoc\nNAG\nT. versicolor | 57.5 GC%\n0.3\nNA\n0.2\n0.1\n0.0\nycneuqerf\nsnodoc\nNAG\nC. cavernicola | 58.5 GC%\n0.3\nns\n0.2\n0.1\n0.0\nycneuqerf\nsnodoc\nNAG\nL. infantum | 59.5 GC%\n0.3\nNA\n0.2\n0.1\n0.0\nycneuqerf\nsnodoc\nNAG\nRhodotorula sp. | 60.5 GC%\n0.3\nNA\n0.2\n0.1\n0.0\nycneuqerf\nsnodoc\nNAG\nC. reinhardtii | 64 GC%\n7) Some methodological questions, to be clarified in the text if needed.\n- For de novo gene detection, are the outgroups always 2, or at least 2?\nThe number of outgroup species with syntenic noncoding matches is set to at least two. This criterion was\nstated in the Methods section, which we have now re-edited to make this clearer:\n“Candidates with one syntenic noncoding hit in two or more outgroup species are then retained for the next\nfiltering step.”\n- Some additional details in the section \"Phylogenetic trees\" would be helpful.\nWe apologize for this section being too brief and have now revised it as follows:\n“Genome trees for the reconstruction of ancestral sequences and the calculation of ORF conservation\nprofiles were constructed with OrthoFinder93 version 2.5.5. First, for every focal genome and its neighbors,\nwe extracted the longest isoform of each gene to produce an amino acid FASTA file per genome. Orthofinder\nidentifies orthogroups and constructs one unrooted gene tree for each orthogroup using MAFFT94 for\nalignments and FastTree95 for tree inference. It then produces a single unrooted species tree with STAG,\nsubsequently rooted with STRIDE through gene tree reconciliation. The command used was `docker run\ndavidemms/orthofinder:2.5.5.2 orthofinder -f <directory_with_input_fasta_files> -t 90 -a 4 -M msa -T\nfasttree -A mafft`.”\n- Phylostratigraphy section: please provide a brief description about how GenEra works for the non-\nspecialist reader.\nWe agree with the Reviewer that this section was too brief and have revised it to better introduce GenEra.\n“The age of each CDS of every focal genome was extracted from the GenEra62 results produced during the\nexecution of DENSE. GenEra62 screens each focal species’ proteome with DIAMONDv281 against the NCBI\nnon-redundant protein database (downloaded on March 28, 2022). Based on the presence/absence of each\nCDS in the NCBI phylogeny, the CDS is then assigned an evolutionary age. The latter corresponds to the most\ndistant node in the phylogeny where DIAMONDv2 has detected the CDS. If a gene has multiple isoforms, the\nisoforms are assigned the age of the oldest isoform. ”\n- \"Inter-CDS and intergenic regions\": Are other genomic elements, apart from CDSs, eg. transposable\nelements, ncRNA genes etc. included in the sequences that remain after the subtraction? If so, could this\nbe affecting the iORFs in some species more than others?\nWe defined iORFs with no overlap with any annotated genomic feature, in order to avoid confounding\neffects from overlapping elements (e.g., noncoding ORFs in alternative reading frames potentially\ninfluenced by CDSs). This is stated in the Methods section:\n“The iORFs were defined from stop-to-stop using the ORFtrack program from the ORFmine package87\n(default parameters: 60 nucleotides, iORFs do not overlap any annotated feature (e.g., tRNA, rRNA, CDS…)\nby more than 70% of their length).”\n8) In Figure 5A, I'm not sure I understand the arrows. There is a pink arrow named \"codon selection\" that\npoints to ancestral iORFs, which I find confusing. Some details in the legend about how to interpret this\nwould be very helpful.\nas follows:\n“Figure 5: A Foldability in relation to GC content. GC content is divided into two main zones: extreme (pink)\nand intermediate (green). iORF foldability (purple) is strongly impacted by GC content. The increase in GC\ncontent observed for ancestral ORFs is consistently associated with a decrease in HCA scores (grey). In\nextreme GC-content species, ancestral ORFs retained by selection exhibit specific codon distributions,\ncompared to their scrambled versions, which either amplify the effect of increased GC content in AT-rich\nspecies or counterbalance it in GC-rich species, thereby maintaining a certain level of foldability (pink). CDSs\noptimize their amino acid composition to maintain foldability irrespective of GC content. This schematic\nillustrates how the genetic code integrates both DNA- and protein-level constraints. B Two reversals\noperating at different timescales. The immediate reversal between DNA and protein levels is mediated by\nthe genetic code, while the evolutionary reversal from noncoding to coding sequences is driven by de novo\ngene birth. Red dashed arrows indicate DNA- and protein-level constraints acting on noncoding and coding\nregions.“\nReviewer #2 (Remarks to the Author):\nThe manuscript by Roginski et al. examines the impact of GC content on de novo gene birth, an important\ntopic in gene origination, given that species in different clades may have widely varying GC contents, yet all\nlikely possess de novo genes in their genomes. The authors identified stop-to-stop intergenic ORFs from\nmore than 3,000 species, analyzed their hydrophobicity, and used ancestral sequence reconstruction to\nstudy the evolution of ORF hydrophobicity. The topic is both interesting and timely, though it is unclear why\nthe authors focused solely on HCA scores, given the wide range of ordered and disordered genes in any\ngiven species, and that there is no clear agreement and empirical evidence on whether young de novo\ngenes are ordered. Below are my comments related to the methods and results.\nWe would first like to thank the Reviewer for their positive comments and constructive feedback, which we\nhave fully addressed in our response letter and, where appropriate, in the revised manuscript. Before\naddressing the specific comments raised by the Reviewer, we would like to clarify our rationale for using\nHCA.\nWe chose HCA because it is fast and alignment-free, which was crucial for studying millions of fast-evolving\nORFs. Moreover, HCA provides a single score that spans the structural spectrum (from disorder to\naggregation propensities), allowing us to describe the structural properties of a sequence with one unified\nmetric, rather than relying on multiple tools with distinct outputs (e.g., IUPred23, TANGO24, etc).\naggregation propensities of iORFs with IUPred and TANGO. These results, now presented in Figure S2, are\nconsistent with the trends predicted by HCA. We observe a positive correlation between disorder\npropensity and GC content, and an anticorrelation with aggregation propensity, supporting the hypothesis\nthat AT-rich genomes contain hydrophobic, aggregation-prone iORFs, while GC-rich ones mostly encode\nhydrophilic, disordered iORFs. While these independent tools help corroborate the structural trends\nidentified with HCA, we remain aware that all of these approaches are computational. However, this\ncomputational framework was essential for exploring such questions across 3,379 eukaryotic genomes,\nencompassing broad ranges of GC content, phylogenetic diversity, and lifestyles - something that would be\nintractable experimentally. These findings now provide a basis for targeted experimental validation in\nselected species, akin to a recent study performed on Drosophila melanogaster and Homo sapiens25.\nWe included this new analysis in the manuscript:\n“As a control, we predicted the disorder and aggregation propensities of iORFs using IUPred58 and TANGO59\n(Figure S2A-B). Results align with HCA predictions, showing that disorder propensity positively correlates\nwith GC content, while aggregation propensity anticorrelates with it. This supports the hypothesis that AT-\nrich genomes predominantly encode aggregation-prone iORFs, whereas GC-rich ones encode disordered\nones.”\n1.00 *** *** *** *** ns ns *** *** *** *** *** *** ***\n0.75\n0.50\n0.25\n0.00\n15−20 20−25 25−30 30−35 35−40 40−45 45−50 50−55 55−60 60−65 65−70 70−75 75−80\nGenome GC content (%)\nPanels A and B of Figure S2: A Boxplots showing the distribution of disorder propensity with respect to GC content, as\npredicted by IUPred, for CDSs (orange) and iORFs (purple) across the 3,379 eukaryotic genomes. B Boxplots showing the\ndistribution of aggregation propensity with respect to GC content, as predicted by TANGO, for CDSs (orange) and iORFs\n(purple) across the 3,379 eukaryotic genomes. For all panels, each species is represented by 1,000 CDSs and 1,000 iORFs that\nwere randomly selected. Due to the large sample sizes, medians of the distributions were compared using an empirical\npermutation test (see Methods). Asterisks denote significance levels: *p < 5 × 10-2, **p < 1 × 10-2, ***p < 1 × 10-3.\nMajor Comments\n20\nni\nseudiser\nfo\nnoitroporP\n)derPUI(\nsnoiger\nderedrosid\nA\n1.00 *** *** *** *** *** *** *** *** *** *** *** *** ***\n0.75\n0.50\n0.25\n0.00\n15−20 20−25 25−30 30−35 35−40 40−45 45−50 50−55 55−60 60−65 65−70 70−75 75−80\nGenome GC content (%)\nni\nseudiser\nfo\nnoitroporP\n)OGNAT(\nsnoiger\nenorp−noitagergga\nB\n25 *** *** *** * * ns ns *** *** *** *** *** ***\n20\n15\n10\n5\n0\n15−20 20−25 25−30 30−35 35−40 40−45 45−50 50−55 55−60 60−65 65−70 70−75 75−80\nGenome GC content (%)\nhtgnel\nretsulc\nACH\n)sAA\nfo\nrebmun(\nC\n*** *** *** *** *** *** *** *** *** *** *** *** ***\n1000\n500\n0\n15−20 20−25 25−30 30−35 35−40 40−45 45−50 50−55 55−60 60−65 65−70 70−75 75−80\nGenome GC content (%)\nhtgneL\n)snodoc\nfo\nrebmun(\nD\nCDSs iORFs\n1) While it is understandable that a straightforward approach to identifying ORFs is through STOP-to-STOP\nanalysis, these are not truly iORFs. I wonder how the results would differ if the authors used a start-to-stop\ncriterion for ORF identification. Ideally, having start-to-stop data for all species would be the best approach,\nbut this may not be feasible within the revision timeframe. I suggest that the authors compare STOP-to-\nSTOP ORFs with start-to-stop ORFs in a subset of species across different clades to provide insights into the\nimpact of the methodological approach.\nWe define iORFs from stop-to-stop, not because it is the most straightforward computationally, but\nbecause growing evidence shows that translation frequently initiates at non-AUG codons, making ATG-\nbased definitions potentially too restrictive. Multiple studies have demonstrated that translation does not\nrely exclusively on canonical ATG start codons (e.g., Slavoff et al.26, Cuevas et al.27, Ingolia et al.28, Na et\nal.29, Papadopoulos et al.30, Zheng & Zhao31).\nFor instance, Slavoff et al.26 reported that in H. sapiens, only 43% of small ORF-encoded peptides start with\nATG, 21% with near-cognate codons, and the remaining 36% with other codons (see Figure 2D in Slavoff et\nal.). Similarly, in our previous analyses of Ribo-seq data from S. cerevisiae30, we found that only ~50% of\ndetected translation events initiate at ATGs. Recently, Zheng & Zhao identified 993 peptides from\nunannotated ORFs in D. melanogaster via mass spectrometry, of which only 331 began with canonical start\ncodons31. Notably, those with and without ATG start shared comparable properties in amino acid\ncomposition and genomic context (e.g., intergenic, UTRs, sense, antisense).\nThus, defining iORFs from stop-to-stop ensures we do not systematically exclude many ORFs that could\ninitiate translation with non-canonical starts. While this definition may add a short N-terminal extension\nthat in fact would not be translated (if the ORF were to be translated), unlike canonical CDSs - whose\nupstream regions differ markedly from their evolved sequence, iORFs are unevolved and their putative\nupstream regions (if translated) are not expected to differ significantly in amino acid composition and\npatterning from the rest of the sequence.\nNevertheless, to address the Reviewer’s concern, we re-analyzed the dataset used for the original Figure 2\n(based on a representative subset of iORFs from the 3,379 genomes), this time defining iORFs from start-\nto-stop (with start = ATG). As shown in Figure D, the overall trends remain consistent: iORFs in AT-rich\nspecies tend to display higher HCA scores, while those from intermediate and GC-rich species show\nintermediate and lower scores, respectively. Similarly, the frequencies of positively and negatively charged\nresidues across genomes are in line with our initial findings using the stop-to-stop definition.\nWhile this complementary analysis confirms the robustness of our findings, we chose not to include it in\nthe main manuscript to maintain focus, and since the manuscript is already quite lengthy. However, we\nremain fully open to providing it as supplementary material should the Editors or Reviewers find it helpful.\nor iORFs defined from start-to-stop (purple) across the 3,379 eukaryo9c genomes as a func9on of genome GC content.\nHorizontal lines delineate the boundaries of the low, intermediate, and high HCA score categories, reflec9ng the three\ncategories of fold poten9al (i.e., disorder-prone, folding poten9al, or aggrega9on-prone). B Distribu9ons of the frequencies\nof the codons encoding posi9vely charged residues of the CDSs (orange) or iORFs defined from start-to-stop (purple) of the\n3,379 eukaryo9c genomes according to genome GC content. C Distribu9ons of the GAN codon frequencies of the CDSs\n(orange) or iORFs defined from start-to-stop purple) of the 3,379 eukaryo9c genomes according to genome GC content. For\nall panels, each species is represented by a representa9ve subset of 1,000 CDSs and 1,000 iORFs selected randomly.\n2) In addition to concerns regarding STOP-to-STOP ORFs, I am curious about how GC content influences the\norigin and occurrence of ATG—the primary start codon in the vast majority of species.\nFigure E shows the frequency of ATGs in iORFs with respect to GC content. We observe that GC content\nimpacts the frequency of ATGs in iORFs, with ATGs being more frequent in iORFs from AT-rich genomes. In\nfact, this is simply a result of ATGs being AT-rich. However, for the aforementioned reasons, we prefer not\nto include this analysis as well, for the manuscript’s conciseness, and since our study focuses on stop-to-\nstop iORFs. That said, if the Reviewers and the Editor find this result informative and relevant to the\nmessage of the manuscript, we will add it as a supplemental figure.\n22\naggreg.\nfoldable\ndisorder\n10\n5\n0\n−5\n−10\n15−20 20−25 25−30 30−35 35−40 40−45 45−50 50−55 55−60 60−65 65−70 70−75 75−80\nGenome GC content (%)\nerocs\nACH\nA\n0.3\n0.2\n0.1\n0.0\n15−20 20−25 25−30 30−35 35−40 40−45 45−50 50−55 55−60 60−65 65−70 70−75 75−80\nGenome GC content (%)\nycneuqerf\nsnodoc\nAA\nevitisoP\nB\n0.15\n0.10\n0.05\n0.00\n15−20 20−25 25−30 30−35 35−40 40−45 45−50 50−55 55−60 60−65 65−70 70−75 75−80\nGenome GC content (%)\nycneuqerf\nsnodoc\nNAG\nC\nCDSs iORFs\ny = −1.27e−04x + 0.02348, R² = 0.166\nr = −0.408, p = 2.27e−135\n0.03\n0.02\n0.01\n0.00\n20 40 60 80\nGenome GC content (%)\nFigure E: Frequency of ATGs as a function of GC content. Mean frequency of ATGs in iORFs calculated for each genome (n =\n3,379) with respect to GC content.\n3) The authors should provide more details on the detection of de novo genes using DENSE. One major\nchallenge in identifying de novo genes is distinguishing them from rapidly evolving genes. It may be\nbeneficial to incorporate methods such as HMMER or foldseek to further exclude potential artefacts.\nWe thank the Reviewer for this insightful comment regarding fast-evolving genes and acknowledge that\nour description of the DENSE pipeline lacked sufficient detail. We have expanded the method section to\nbetter introduce the methodology.\nDENSE9 is a two-step procedure. It first identifies gene candidates that are taxonomically restricted, by\ndefault at the genus level. To do so, it relies on GenEra32, which screens with DIAMONDv233 the focal\nspecies’ proteome against the nr database. Each CDS is then assigned an evolutionary age, corresponding\nto the most distant node in the NCBI phylogeny where DIAMONDv2 has detected the CDS. In cases where\na gene has multiple isoforms, the age of the oldest one is assigned to all isoforms. For each TRG, it then\ntests whether the gene was ancestrally noncoding. Therefore, using TBLASTN, DENSE searches for\nhomology traces in noncoding regions of outgroup species among the neighbors of the focal species. These\nspecies are those (i) lacking the candidate gene and (ii) branching after the last common ancestor node of\nthe gene family. To validate TBLASTN hits, DENSE uses microsynteny to confirm that the detected\nnoncoding hits in outgroup species are located at the same locus as the focal de novo gene candidate.\nCandidates with at least one syntenic noncoding hit in two or more outgroup species are retained as\npotential candidates and then filtered based on remote or structural homology according to the Reviewer’s\nsuggestion. Indeed, we fully agree that DIAMONDv2 may fail to detect distant homologs that have diverged\nbeyond recognition. To mitigate this risk, we implemented a filtering strategy to exclude candidates\npotentially corresponding to fast-evolving genes:\nOn the one hand, we screened all de novo CDS candidates against the Pfam database using HMMER34 (e-\nvalue £ 1x10-3), without applying a coverage threshold, in order to remain conservative - that is, we chose\nto eliminate any candidate showing a detectable domain match, regardless of the extent of the match, to\nreduce the risk of retaining false positives (i.e., fast-evolving genes that may have retained only partial\ndomain similarity). Any candidate with a significant match outside the genus of the focal species was\nremoved from the de novo CDS set. For genes with multiple isoforms, all isoforms were screened. If any\n23\nsFROi\nni\nGTA\nfo\nycneuqerf\nnaem\nisoform yielded a significant match, the gene and its isoforms were conservatively excluded. This step\neliminated 15 candidates.\nOn the other hand, we predicted the 3D structure of each de novo protein candidate using ESMFold35 and\nsearched for structural homologs among the 15 million AlphafoldDB cluster representatives36 with\nFoldseek37 (e-value £ 1x10-2). One could argue that structural homologs can reflect convergence instead of\ndivergent remote homologs, and distinguishing between these two scenarios is complex. However, we\nprefer to be conservative and remove all potential remote homologs that would have been missed by\nDIAMONDv2. Foldseek returned 10 significant matches. Among them, five overlapped the HMMER-filtered\nset, confirming their evolutionary relationship (HMMER relies on sequence information). The remaining five\ncandidates had Foldseek TM-scores < 0.5, indicating weak structural alignment. However, Foldseek TM-\nscores have been debated by users (e.g., GitHub issues). TM-scores were thus reassessed with TM-align38,\nwhich is more accurate but also much slower, and, therefore, could not be employed to screen the 15\nmillion AFDB representatives. One of the five Foldseek hits yielded a significant TM-score (0.64) with TM-\nalign and was conservatively removed as well.\nIn total, 16 candidates were filtered out of the original 2,551 de novo genes.\nAll figures, tables, results, and downstream analyses were updated accordingly to account for the curated\ndataset.\nThe Methods section now provides a more detailed description of DENSE along with an explanation of the\nfiltering procedure.\n“DENSE was run with default parameters, except for the number of outgroup species, which was set to ‘2’.\nDENSE first identifies TRGs using GenEra62. The TRG rank was set to ‘genus’ for all species except for H.\nsapiens and M. musculus, where it was set to 'Hominoidea' and 'Murinae', respectively, to ensure a minimum\nnumber of neighbors. For each TRG, it then tests whether the gene was ancestrally noncoding. Therefore,\nusing TBLASTN, DENSE searches for homology traces in noncoding regions of outgroup species among the\nneighbors of the focal species. These species are those (i) lacking the candidate gene and (ii) branching after\nthe last common ancestor node of the gene family. To validate TBLASTN hits, DENSE uses microsynteny to\nconfirm that the detected noncoding hits in outgroup species are located at the same locus as the focal de\nnovo gene candidate. Candidates with one syntenic noncoding hit in two or more outgroup species are then\nretained for the next filtering step. DIAMONDv281 may miss distant homologs that have diverged beyond\nrecognition. To mitigate this risk, we screened the retained candidates against the Pfam database using\nHMMER88 (e-value £ 1x10-3). Candidates with significant matches outside the genus of the focal species\nwere removed from the de novo candidate set. For genes with multiple isoforms, all isoforms were screened.\nIf any isoform yielded a significant match, the gene and its isoforms were conservatively excluded. This step\neliminated 15 candidates. We also predicted the 3D structure of each de novo protein candidate using\nESMFold89 and searched for structural homologs among the 15 million AlphafoldDB cluster representatives90\nwith Foldseek91 (e-value £ 1x10-2). Of 10 matches, five overlapped the HMMER-filtered set. The remaining\nfive had Foldseek TM-scores < 0.5, indicating weak structural alignment. TM-scores were reassessed with\nTM-align92, which is more accurate. One of them yielded a significant TM-score (0.64) and was\nconservatively removed as well.”\n4) The presentation of Figure 1 is somewhat unclear, though I believe I have understood its intended\nmeaning. Comparisons between CDSs, iORFs, and scrambled sequences may all be significantly different—\nif so, what is the significance of these differences? Lines 148–181 may contain results related to Figure 1.\nHowever, without clear panel separations within the figure, it is difficult to interpret these paragraphs in\nrelation to the data presented fully. Providing clearer figure annotations or structured panel divisions would\nenhance readability.\nscrambled ones differ significantly for all codon group frequencies, as now shown with the new panel C of\nFigure S1 (Mann-Whitney U test (two-sided), all P < 1 × 10-3). To better quantify these differences, we\ncalculated the log fold change (LFC) of each codon group, representing its deviation in real CDSs compared\nto scrambled ones. These results are now presented in the new panel B of Figure 1, allowing a clearer\nclassification of codon groups that was more difficult to infer from the qualitative analysis of the initial plot.\nBy doing so, we highlight some codon groups displaying higher deviations, which we relate in the text to\nspecific properties (stop codons, foldability, negatively charged residues, and CpGs).\nIn addition, we quantified the variation in codon group frequencies across genomes for iORFs (x-axis) and\nCDSs (y-axis). This analysis, now presented by the new panel C of Figure 1, explicitly measures the extent to\nwhich the frequency of a codon group varies across genomes in iORFs and CDSs. Variation was assessed by\nthe standard deviations (SDs) of each codon group’s frequency distribution, with high SDs indicating high\nvariability according to GC content, and low SDs indicating more stable frequencies. Codons near the\ndiagonal vary similarly in CDSs and iORFs, whereas those in the middle-right region of the plot (e.g.,\nNTN, NGN) vary widely in iORFs but remain constrained in CDSs (Figure 1C).\nWe hope the statistical analyses, together with the new Figures (Figure 1B-C and Figure S1C) and revised\ntext, have improved the clarity and readability of this section. The text now reads:\n“Figure 1B shows the log fold change (LFC) of each codon group representing the deviation of its frequency\nin real CDSs compared to scrambled ones. Figure 1C illustrates the variation in its frequency across genomes\nfor iORFs (x-axis) and CDSs (y-axis). Variation is quantified by the standard deviation (SD) of a codon group’s\nfrequency distribution, with high SDs indicating high variability according to GC content, and low SDs\nindicating more stable frequencies. Codons near the diagonal vary similarly in CDSs and iORFs, whereas\nthose in the middle-right region of the plot (e.g., NTN, NGN) vary widely in iORFs but remain constrained in\nCDSs (Figure 1C).\nSeveral codon groups exhibit small LFCs, reflecting similar frequencies in real and scrambled CDSs (Figure\n1B). Notably, mixed AT/GC codons (e.g., TCN, CAN, or GTN) typically show little variation in both iORFs and\nCDSs (bottom-left corner of Figure 1C), reflecting low sensitivity to GC content. Overall, most codons have\nabsolute LFCs below 0.5, but four groups show stronger deviations (LFCs > 0.75): TAN, TGN, CGN, and GAN.\nTheir median frequencies differ significantly from those of scrambled CDSs (Mann-Whitney U test, (two-\nsided), all P < 1 × 10-16)(Figure S1C). TAN and TGN comprise the three stop codons and the rare amino acids\ntryptophan and cysteine, explaining their strong depletion in CDSs relative to expectations (Figure 1B, Figure\nS1C). CGN codons, which encode arginine, are also markedly depleted (Figure 1B, Figure S1C). CpG sites,\nconsisting of cytosine followed by guanine, are prone to methylation and are hypermutable because\nmethylated cytosines frequently undergo deamination to thymine51. Such mutations can generate TGN\ncodons, which may introduce deleterious changes and thereby contribute to the observed depletion of CGN\ncodons. GAN codons, encoding aspartate and glutamate, display the highest LFC, and a low frequency\nvariation in CDSs (0.116 ± 0.005) across genomes (Figure 1B-C). Previous studies showed that high\nproportions of negative charges on protein surfaces prevent promiscuous interactions with the negatively\ncharged, abundant ribosomes52, and stretches of such residues optimize translation by facilitating the\nnascent peptide’s passage through the negatively charged ribosomal exit tunnel53. In yeast, we also found\na correlation between negative charge frequency and protein abundance, suggesting that a higher\nfrequency of negative charges in abundant proteins may favor their translation efficiency and reduce\nnonspecific interactions with ribosomes23. Finally, the NTN group, encoding most hydrophobic residues,\nstands out in Figure 1C, showing the widest variation in iORFs but constrained frequencies in CDSs (0.27 ±\n0.01). This suggests that selection acts to maintain optimal hydrophobicity for protein folding, preserving\nenough hydrophobic residues in GC-rich species for proper folding while reducing their frequency in AT-rich\nspecies to prevent aggregation. This aligns with Dill’s classic theory, which proposed that protein stability\ndepends on an optimal fraction of hydrophobic residues54.”\nacross the 3,379 genomes. Posi9ve values indicate enrichment in real CDSs, nega9ve values indicate deple9on in real CDSs.\nError bars represent standard devia9on. Detailed distribu9ons for each codon group are shown in Supplementary Figure S1C.\nC The SDs on the x-axis and y-axis reflect the spread of codon group frequencies across the 3,379 genomes in iORFs and CDSs,\nrespec9vely. Codon groups are colored according to their absolute LFC value, with light dots corresponding to codon groups\nwith small frequency devia9ons compared to scrambled CDSs, and dark dots indica9ng codon groups whose CDS frequency\ndeviates from that of scrambled CDSs.\n`\nPanel C of Figure S1. Distributions of codon group frequencies for scrambled (grey) and real (orange) CDSs. Medians of the\ndistributions were compared using a one-sided Mann-Whitney U test. Asterisks denote significance levels: *p < 5 × 10-2, **p\n< 1 × 10-2, ***p < 1 × 10-3.\n5) Figure 2A is very interesting. It shows that there is a strong constraint on HCA scores in the long-term\nevolution of genes —though not the outliers. I wonder, for the de novo CDS identified in Table 1, what are\n26\nthe HCA scores for de novo CDS with very high and very low GC contents? Comparing the extremes to the\nannotated CDSs may provide insights into whether there are HCA constraints for de novo genes during birth\nor whether the HCAs change gradually after gene birth.\nWe apologize if this was not clearly conveyed in the text. De novo gene HCA scores are presented in Figure\nS8, as indicated by the legend “The first phylostratum corresponds to de novo CDSs (light yellow), while\nthe last one corresponds to CDSs shared by all cellular organisms (dark orange).”. However, we\nacknowledge that this requires examining all individual plots. To better illustrate the overall trend and\naddress the Reviewer’s concern, we added a new panel (B) showing the distribution of the median HCA\nscores for each species’ de novo CDSs (21 medians, as H. opuntae has no de novo genes) together with that\nfor the two oldest phylostrata, i.e., Eukaryotes and Cellular Organisms (phylostrata shared by all species in\nthe study). We represented the medians rather than pooling all de novo genes to account for the large\ndifferences in de novo gene counts across species. As shown in Figure S8B, de novo CDSs overall fall in the\nfoldable bin but display higher variability and lower HCA scores, potentially reflecting greater local disorder.\nThis trend holds for young phylostrata (TRGs), now also represented in the new panel C of Figure S8. We\nadded two sentences that summarize this:\n“Finally, HCA scores across phylostrata reveal that although young phylostrata (TRGs), including de novo\ngenes, are generally predicted with folding potential, their foldability varies significantly depending on the\nspecies’ GC content (Figure 4, Figure S8B-C). In contrast, older phylostrata (Eukaryotes and Cellular\norganisms) display lower variations in their HCA scores irrespective of species’ GC contents (F-test, P =\n5.6x10-6)(Figure S8C).”\nPanels B and C of Figure S8: Median HCA score distributions for young and old phylostrata. B Distributions of median HCA\nscores for the de novo gene phylostratum across all species and for the two oldest phylostrata (Eukaryotes and Cellular\norganisms). The Mann-Whitney test for the median comparison was not applicable due to the small size of the de novo bin\n(one median per species, i.e., 21 values). C Distributions of the median HCA scores for young phylostrata from all species and\nfor the two oldest phylostrata. Young phylostrata include levels from de novo genes to genus (except for H. sapiens and M.\nmusculus, where they extend to ‘Hominoidea’ and ‘Murinae’, respectively). The variances of the two distributions were\ncompared using an F-test.\n6) In Figure 3A, the authors found that ancestral iORFs consistently exhibit higher GC content than other\niORFs in the analyzed species, which is consistent with a previously published study. However, this is an\nunusual result that is difficult to explain biologically and may instead be a computational artefact. I strongly\nrecommend that the authors use at least one additional method for ancestral sequence reconstruction to\nassess the robustness of these findings.\nPRANK39. As a control, we also applied the ancestral sequence reconstruction (ASR) method PREQUEL from\nthe PHAST package40, following the protocol described in Vakirlis et al.41, who compared different ASR\nprotocols. The patterns observed with PREQUEL are consistent with those obtained with PRANK, which\naligns with the findings of Vakirlis et al., who reported good overall agreement between ASR protocols41.\nSpecifically, in our analysis, PREQUEL ancestral ORFs remain predicted to have higher GC contents than\niORFs and generally display intermediate HCA scores (Figure S5), reflecting folding potential, except for the\ntwo mammals, which were also predicted by PRANK to be enriched in disordered sequences.\nTwo differences are nevertheless noted (i) in Rotaria sp. (only one de novo CDS), PREQUEL did not predict\nany ancestral ORF, and (ii) in Rhodotorula sp. (only four de novo CDSs), ancestral ORFs reconstructed by\nPREQUEL were split between folding potential and disorder, whereas those reconstructed by PRANK were\npredicted to have folding potential. However, in this latter case, the small number of reconstructed\nsequences (n = 4) limits the accuracy of the comparison between the two protocols.\nWe have now included a new Figure S5 showing the GC content and HCA scores of the ancestral sequences\nreconstructed by PREQUEL and have commented on this analysis in the Results section.\n“Moreover, this pauern also holds for the ancestral ORFs predicted by the PREQUEL-based protocol of\nVakirlis et al.61”\n“In many species, ancestral ORFs exhibit disGnct foldability from other iORFs, a trend also observed for the\nancestral sequences predicted by the PREQUEL-based protocol of Vakirlis et al.61”\nGC% 19% 24% 29% 35.5% 35.5% 38% 40.5% 41.5% 42% 42% 43.5% 51% 53% 55.5% 56.5% 57.5% 58.5% 59.5% 60.5% 64%\n#ancORFs 3 1 8 135 30 249 149 129 78 24 2319 54 99 5 5 18 17 59 4 11\nFigure S5: ProperSes of ancestral ORFs reconstructed using PREQUEL. A For each species with at least one detected de novo\ngene (H. opun7ae and Rotaria sp. were thus excluded), GC content distribu9ons are shown as boxplots for 1,000 randomly\nselected iORFs (purple) and for the ancestral ORFs of all de novo CDSs (light purple)(see Table 1 for the full scien9fic names).\nB For the same species as in panel A, HCA score distribu9ons are shown as boxplots for 1,000 randomly selected iORFs\n(purple), the ancestral ORFs of all de novo CDSs (light purple), and their corresponding scrambled version (grey). The number\nof reconstructed ancestral iORFs is shown at the top of each plot. For species with at least 30 ancestral ORFs, we compared\nthe medians of their associated distribu9ons with those of iORFs or scrambled sequences using a one-sided Mann-Whitney\nU test applying Bonferroni correc9on for mul9ple comparisons. Asterisks denote significance levels: *p < 5 × 10-2, **p < 1 ×\n10-2, ***p < 1 × 10-3. “ns” indicates p-values > 0.05. “NA” indicates species with less than 30 ancestral ORFs, for which the\nMann-Whitney test was thus not applicable.\n7) In the paper the authors used HCA scores as a major parameter important for evolution as well as gene\norigination, as indicated in several figures as well as in figure 5A. However, there is no clear evidence that\na functional novel gene needs to be folded, especially at the early stage.\nWe approached de novo gene emergence from a structural perspective for the following reasons:\n(i) Protein structural properties (including disorder) are closely linked to function. Investigating\nthe structural features of novel CDSs and their ancestral sequences is therefore key to\ndetermining whether these ORFs exhibit specific structural properties that may influence their\n29\n%CG\nA\n75\n50\n25\n0\niORFs PREQUEL ancORFs\nPfal Ehis Pmac Atha Cele Scer HsapMmusDmel Pviv Osat Post Cjav Mres Tcar Tver Ccav Linf Rccf Crei\nGC% 19% 24% 29% 35.5% 35.5% 38% 40.5% 41.5% 42% 42% 43.5% 51% 53% 55.5% 56.5% 57.5% 58.5% 59.5% 60.5% 64%\n#ancORFs 3 1 8 135 30 249 149 129 78 24 2319 54 99 5 5 18 17 59 4 11\nerocs\nACH\nNA NA NA *** ** *** *** *** *** NA *** ns ns NA NA NA NA ns NA NA\n5\n0\n−5\n−10\niORFs PREQUEL ancORFs scrambled PREQUEL ancORFs\naggreg.\nprop.\nfolding\npotential\ndisord.\nprop.\nB\nfixation (whether disordered, with folding potential, or aggregation prone), and to better\nunderstanding the relationship between structure and function in emerging proteins.\n(ii) Several computational and experimental studies have investigated this question. Some\nreported AT-rich protogenes enriched in transmembrane domains17. Other experimental\nstudies demonstrated that de novo proteins could fold into partially ordered 3D structures or\nbe associated with some folding potential25,42–45. These experimental studies are consistent\nwith other computational analyses, including our own work on yeast, which, although limited\nto a small number of species (e.g., S. cerevisiae, D. melanogaster, O. sativa…), showed that de\nnovo proteins may harbor folding potential, with some even predicted as well-folded46–49.\nMotivated by these analyses, we extended the analysis to a larger scale. Notably, because structural\nproperties are directly influenced by amino acid composition, which itself is shaped by GC content (at least\nin unevolved sequences), we analyzed a comprehensive dataset of eukaryotic genomes (all annotated\ngenomes in NCBI with an assembly level of at least “scaffold” as available on March 6, 2023), encompassing\na broad range of GC contents.\nAs explained earlier, we used HCA as our primary metric, since it is fast and alignment-free (which was\ncritical for large-scale analyses) and provides a single score spanning the structural spectrum from disorder\nto aggregation. However, to strengthen our predictions, we have now included two complementary and\nindependent analyses (TANGO and IUPred), which have corroborated the predictions obtained by HCA.\nMinor Comments\nLine 70: In the section titled “How do we pass from one world to another”, the intended meaning is unclear.\nThe use of the term “world” is somewhat confusing, as there is no distinct separation between different\nworlds in the context of gene birth, particularly concerning the reservoir of genes.\nWe acknowledge that the term “word” was not appropriate and have removed all its occurrences from the\ntext.\nLines 86-94: Many proto-genes tend to be AT-rich, especially given that transmembrane domains often\noriginate in T-rich or AT-rich regions, as these regions are more likely to generate amino acids with high\nhydrophobicity. How does this reconcile with the authors' observations?\nWe have now included a new paragraph in the Discussion section regarding the TM peptides with adaptive\npotential reported in Vakirlis et al.17.\n“Yet, it should be noted that neither annotation pipelines nor experimental methods can fully capture ORF\ndiversity. Notably, Vakirlis et al. identified 28 unannotated ORFs with adaptive potential in yeast through\nlarge-scale overexpression screens41. These thymine-rich ORFs often encode transmembrane domains and\nshow no fitness impact upon disruption. The authors interestingly discussed the possibility that beneficial\neffects observed under laboratory conditions may not apply in natural environments, and that adaptive\npotential does not guarantee fixation, as trade-offs may hinder long-term retention. Their work underscores\nthe complexity of ORF annotation and the possibility that some adaptive ORFs escape detection, and raises\nimportant questions about the evolutionary fate of these ORFs regarding two key observations: CDSs in\nnearly all eukaryotes are GC-richer than noncoding regions (Figure S1A), and transmembrane domains\nremain underrepresented in proteomes82.”\nLine 214-217, if the raised question not answered, it should belong to discussion.\nWe thank the Reviewer for pointing this out. This sentence has now been removed from the text.\nscrambled sequences is commendable. However, does blending with four codons introduce any bias in\ncodon composition within these regions?\nWe thank the Reviewer for catching this. Indeed, blending more codons reduces bias but also increases\ncomputational cost. We compared our approach with scrambling using five upstream and downstream codons, which\nprovides a good balance between smoother blending and computational efficiency, and led to the same\nresults. As this approach offers an optimal compromise, we retained it in the present version of the\nmanuscript.\nThe readability of the manuscript can be improved. Specifically, figure and figure panels can be cited more\nrigorously. Some result sections contain guesses and open questions which belong to the discussion.\nWe are grateful to the Reviewer for raising these important points. We have thoroughly revised the\nmanuscript and updated Figure 1 accordingly, which we believe has significantly improved its clarity and\nreadability.\nReviewer #3 (Remarks to the Author):\nIn their manuscript, “Impact of GC content on de novo gene birth,” Roginski et al. explore the influence of\nGC content in non-coding open reading frames on the structure of potentially encoded microproteins.\nNotably, as a major merit, their analysis spans an impressive dataset of 3,379 eukaryotic genomes—far\nexceeding the scope of previous studies, which have been limited either to specific evolutionary lineages\nor to smaller genome sets (e.g., 187 genomes in Basile et al., 2017). For a subset of 22 organisms, they\nemploy a combination of ancestral sequence reconstruction (ASR)—used here more extensively than in\nprior work—and a streamlined, custom-built pipeline for de novo gene identification based on their prior\ndeveloped tool DENSE. They compare ancestral ORFs with scrambled ORFs. This is the second major merit\nof the study, however I would suggest to compare this also with degenerated codons. Additionally, the\nauthors consciously avoid reliance on 3D structural predictions, which have recently come under scrutiny\nfor their (lack of) reliability in the context of de novo or random sequences. However, the use of\nhydrophobic cluster analysis (HCA) seems somehow questionable for the short ORFs (60 nucleotides, i.e.,\n20 amino acids) analysed because it implicietely assumes the presence of secondary structural elements,\nor at least the principal capacity of the sequence to form such.\nThe manuscript itself lacks engagement with foundational literature on the structural characteristics of de\nnovo proteins and biophysics. In a few cases, the citations do not support what is written, either because\nthey are reviews or simply not matching. In some instances, the language is overly colorful for a scientific\narticle, and certain metaphors or illustrative phrases do not add clarity but instead may detract from the\nmanuscript’s overall tone.\nWe are grateful to the Reviewer for the positive feedback on our work. We acknowledge that the initial\nversion lacked references to key literature on the structural analyses of de novo proteins and biophysics,\nand that some metaphors were not appropriate for a scientific article. We have thoroughly revised the\nmanuscript to address these issues. In the following point-by-point responses, we also address the\nReviewer’s comments on the use of HCA or degenerate codons.\nL1: Specify non-coding regions in eukaryotes.\nDone.\n“Noncoding regions in eukaryotes are extensively expressed and represent a significant source of novel\nmicroproteins, some of which become fixed as de novo genes.”\nL24: Change “investigated” to “predicted” to emphasize that the study is limited to computational work.\nDone.\n“Here, we predicted the foldability and sequence properties of millions of microproteins potentially encoded\nin the noncoding open reading frames (ORFs) of 3,379 eukaryotic genomes with GC contents ranging from\n18% to 79%.”\nL29: What is actually meant by “the nature of the associated evolutionary stages”?\nWe thank the Reviewer for this question. By “stages”, we actually meant \"intermediates,\" which is more\nappropriate as it refers to the products (i.e., ancestral ORFs), rather than the process. All occurrences of\n“stages” have therefore been replaced with “intermediates” in the text. Regarding this sentence (L29), it\nhas been removed from the abstract to reduce word count.\nL33: Please highlight in 1–2 sentences the actual results of the study.\nWe thank the Reviewer for this relevant comment and have revised the abstract accordingly, clarifying the\nfinal sentence to better present our main results.\n“We show that de novo genes preferentially emerge from GC-rich ORFs with folding potential, revealing that\nGC content and foldability, coordinated by the genetic code structure, together shape the emergence of\nnovel genes.”\nL43: Use “conceptually classified” instead of “long viewed” to highlight the conceptual aspect.\nDone.\n“Genomes can be broadly subdivided into coding and noncoding regions, conceptually classified as distinct\nand impermeable entities”\nL42–55: The use of words such as “materials” and “worlds” is not scientific and is uncommon.\nIt impedes readability and does not seem as conceptually intriguing as it was possibly intended to be.\nWe agree with the Reviewer and thank them for pointing this out. All occurrences of these words have been\nremoved from the entire text. A few examples below:\n32\n“Recent studies showed that (i) noncoding regions harbor vast reservoirs of small Open Reading Frames\n(ORFs) that can be expressed at the protein level, and (ii) the transiFon from noncoding to coding\nsequences (CDSs) is possible1-9”\n“Hundreds of newly expressed products have been shown to be functional based on CRISPR-based genetic\nscreens, or evidence of selection, thereby establishing them as part of the coding repertoire3,15–22.”\nL49: Refer to Grandchamp et al. 2023 & 2024, and Iyengar & Bornberg-Bauer 2023 regarding the transition\nfrom non-coding to coding.\nDone.\nL53: “Domains” can be misunderstood. Is it in the sense of Pfam?\nWe agree and replaced “transitions between these presumed impermeable domains” with “transitions from\nnoncoding to coding sequences”. The sentence now reads:\n“Interestingly, several studies have demonstrated that transitions from noncoding to coding sequences are\nfar more common than previously thought6,14,15”\nL55: None of these citations found positive selection. Citation 14 is a review.\nWe thank the Reviewer for catching this. We revised the sentence to explicitly mention CRISPR-based\ngenetic screens, thereby justifying the decision to retain the citation to Chen et al., 202052. We also\nincorporated several references suggested by the reviewer (see point-by-point responses below) and\nremoved those that were not appropriate. Below is the revised sentence:\n“Notably, hundreds of newly expressed products have been shown to be functional based on CRISPR-based\ngenetic screens, or shown to undergo selection, thereby establishing them as part of the coding\nrepertoire3,15–19.”\nBy following the citations within it, the authors can find studies supporting their statement.\nAnother, more recent, review that deals extensively with selection can also be cited:\nZhao, L., Svetec, N., Begun, D.J. De novo genes. Annual Review of Genetics 58 (2024).\nThis paper is explicit on selection and methodological limits to detect selection in the context.\nDone.\nOther primary studies that investigate selection include:\n• Zhao, L., Saelao, P., Jones, C.D., Begun, D.J. Origin and spread of de novo genes in Drosophila melanogaster\npopulations. Science 343(6172), 769–772 (2014)\nDone.\n• Schmitz, J.F., Chain, F.J., Bornberg-Bauer, E. Evolution of novel genes in three-spined stickleback\npopulations. Heredity 125(1–2), 50–59 (2020)\nDone.\n33\n• Durand et al. Turnover of ribosome-associated transcripts from de novo ORFs produces gene-like\ncharacteristics available for de novo gene emergence in wild yeast populations. Genome Research (2019)\nThis reference is not relevant to our sentence, as it reports that the ORFs detected as translated do not\nshow evidence of significant purifying selection. -> \"Altogether, tORFs do not display significant purifying\nselection, but it appears that as a neutral pool, they provide raw material with gene-like characteristics for\nselection to act.”\n• Vara, C., Montañés, J.C., Albà, M.M. High polymorphism levels of de novo ORFs in a Yoruba human\npopulation. Genome Biology and Evolution 16(7), 126 (2024)\nThis study does not explicitly test de novo ORFs for selection and, therefore, was not included as a\nreference.\nWe instead added these references:\n*** Zhang, H. et al. Genome-wide maps of ribosomal occupancy provide insights into adap6ve evolu6on\nand regulatory roles of uORFs during Drosophila development. PLOS Biol. 16, e2003903 (2018).\nhttps://journals.plos.org/plosbiology/article?id=10.1371/journal.pbio.2003903\nIn this ar6cle, the authors showed that “the majority of the newly fixed uORFs in D. melanogaster were\ndriven by posiGve Darwinian selecGon.”\n*** Gubala, A. M. et al. The Goddard and Saturn Genes Are Essen6al for Drosophila Male Fer6lity and\nMay Have Arisen De Novo. Mol. Biol. Evol. 34, 1066–1082 (2017).\nhttps://pmc.ncbi.nlm.nih.gov/articles/PMC5400382/pdf/msx057.pdf\nIn this article, the authors showed that “Goddard is consistently found in single-copy and evolves under\npurifying selection. In contrast, saturn has diversified through gene duplication and positive selection.”\nL64: Not all of these citations address the issues described in this sentence, particularly regarding length.\nA references in which this issue is\ncentrally adresses should be added:\n• Lebherz, M.K., Iyengar, B.R., Bornberg-Bauer, E. Modeling Length Changes in De Novo Open Reading\nFrames during Neutral Evolution, Genome Biology and Evolution, 16(7), July 2024, evae129\nWe thank the Reviewer for pointing this out. We have carefully re-examined the associated references and\nremoved the citation to Patraquim et al., 2022, as it was not appropriate for this sentence. The remaining\nreferences were retained since each addresses specific aspects of the sentence. The retained references\nand their associated points are listed below. We have also incorporated Lebherz et al., GBE 2024, with\nrespect to ORF length.\n*** Amino acid compositions:\n– Papadopoulos et al., Genome Res., 2021\n– Montañés et al., Mol. Biol. Evol., 2023\n– Chen et al., Genome Biol. Evol., 2024\n34\n*** ORF length:\n– Carvunis et al., Nature, 2012\n– Papadopoulos et al., Genome Res., 2021\n– Montañés et al., Mol. Biol. Evol., 2023\n– Lebherz et al., Genome Biol. Evol., 2024\n*** Expression in specific tissues or conditions:\n– Peng & Zhao, Nat. Commun., 2024\nL66: Double-check your references regarding purifying selection. “Singularities” is not descriptive, and\n“unequivocally” is too strong.\nThe Reviewer is correct regarding the use of “singularities” and “unequivocally,” which have now been\nreplaced with “distinct features” and “confidently,” respectively.\nMoreover, we carefully re-examined all related references. The original sentence has been removed and\nreplaced with the following (to also address a comment from Reviewer 1), along with the relevant\nreferences:\n“Although some features of de novo genes (e.g., short length, weak expression) might also reflect\nnonfunctional, noisily expressed ORFs, the evidence that several display tissue-specific expression28–30,\nundergo selection15,16,18,19, or have been experimentally shown to be functional19,31–33, supports the view\nthat they can reasonably be classified as coding.”\nL69: Which results? What does “permeability” refer to? The phrase “passing from one world to another” is\nnot scientific; see also L75, which is likewise incorrect.\nWe thank the Reviewer again for pointing this out and have removed all occurrences of “word” and revised\nthese sentences for better clarity, as follows:\n“These results, which establish de novo gene birth as a significant driver of geneGc novelty, pose a\nfundamental quesGon: how does a noncoding sequence become coding? “\n“This intriguing quesGon can be reframed in terms of constraints faced by coding and noncoding sequences:\nhow does a noncoding sequence governed only by DNA-level constraints transiGon to a coding sequence\nsubject to both DNA and protein-level constraints?”\nL80: Replace “Their” with a specific subject. Be precise.\nWe apologize for this ambiguity and have edited the sentence for improved clarity.\n“The expression of noncoding regions is, therefore, expected to generate aberrant or harmful proteins\nrather than funcGonal templates for novel genes.”\nL82: Same as above—what does “This” refer to?\nWe acknowledge that the sentence was unclear and have revised it to improve clarity.\n35\n“This apparent paradox - how composiGonally biased noncoding regions give rise to funcGonal proteins -\nhas sGmulated numerous studies, leading to various models of de novo gene birth.”\nL83: The preadaptation model is generally flawed, as it does not explain how adaptation can occur without\nselection. Neutrality would be a more appropriate null model, especially since random sequences have\nbeen experimentally\nshown to be foldable and soluble (see Tretyachenko et al., 2017, Scientific Reports). Random sequences have\nalso been experimentally compared for disorder, solubility, and foldability\n(Heames et al. 2023, Nature Ecology & Evolution)\nWe fully agree with the Reviewer and apologize for omitting these important works. We have now\nintroduced them in the Introduction.\n“However, Tretyachenko et al. showed that random sequences can display secondary structure content and\nphysicochemical properties similar to natural proteins, making them well tolerated by the cell38. Heames et\nal. reported hundreds of putative de novo proteins in humans and flies that resemble random sequences in\nstructural potential, though with a higher proportion of soluble and protease-sensitive sequences.\nInterestingly, their findings suggest that structure-forming potential can emerge from random sequences\nwithout prior selection39.”\nL88: The authors have computationally predicted in yeast.\nThe sentence has now been edited and reads:\n“In yeast, de novo proteins were predicted to be less disordered than old ones and expected to acquire\ngreater solubility during evoluGon.”\nL90: The transition from “less-disordered and soluble” to “transmembrane” is not smooth. If this is meant\nas a summary of structural studies, the structure should follow the sequence:\ndisorder → structure → higher-order (TM, multimer).\nWe may have misunderstood the Reviewer’s point, but we do not intend to directly link the proto-gene\nmodel of Carvunis et al.53, which suggests a transition from less disordered sequences to higher solubility,\nto the TM-first model, which posits that transmembrane de novo peptides can provide beneficial effects\nand later evolve broader cellular functions.\nIn this paragraph, our aim was to present different models, including the non-exclusive proto-gene and TM-\nfirst models, to highlight different findings that attribute an important role to foldability in the transition\nfrom noncoding to coding sequences. We have also added the works of Tretyachenko et al.57 and Heames\net al.25, which were missing and are highly relevant to these questions.\nExperimental studies should be discussed here:\n• Bungard et al. 2017, Done.\n• Baalsrud et al. 2018, Unless we have missed it, this study does not report any experimental\ncharacterization.\n• Lange et al. 2021, Done.\n• Matsuo et al. 2021 Done.\nRecent computational studies to cite:\n• Middendorf, L., Iyengar, B.R., & Eicholt, L.A. (2024). Sequence, Structure, and Functional Space of\nDrosophila De Novo Proteins. GBE Done.\n36\n• Chen, J., Li, Q., Xia, S., Arsala, D., Sosa, D., Wang, D., Long, M. The Rapid Evolution of De Novo Proteins in\nStructure and Complex. GBE 16(6), June 2024, evae107 Done.\nThe Reviewer is correct that these references were missing. We have now included a paragraph to introduce\nthem.\n“While computational studies on the structural properties of de novo proteins have accumulated23,26,42,43,\nexperimental analyses remain scarce, with only a few de novo proteins structurally characterized44–47.”\n“Notably, the yeast Bsc4 and the fly Goddard de novo proteins were shown to fold into partially ordered 3D\nstructures containing secondary structure elements but lacking the stable tertiary structure typical of\nglobular proteins44,47.”\nL99: What is meant by “folding potential,” and how is it assessed?\nSince this question is closely related to the one referring to line 189 (see below), we have provided a\ncommon answer addressing both. Please refer to our response at L189 (below - highlighted in yellow).\nNonetheless, we have added a clarification in the Introduction and have also specified that this was\nassessed computationally.\n“On the other hand, using computational methods, we showed that yeast intergenic ORFs encode a wide\ndiversity of structural potentials, including aggregation-prone peptides, and a majority of ORFs with\nfolding potential23. The latter are not necessarily expected to form well-defined 3D structures but contain\nelementary structural blocks that may promote folding in specific contexts (e.g., binding with a partner,\nchaperone assistance39) or the evolution of more stable structures.”\nL101: Lange et al. 2021 and Vakirlis et al. 2024 also used ASR and should be cited here.\nOur inten6on with this sentence was not to detail ASR methods or literature, but to report on the structural\nfeatures of ancestral sequences. However, we have now cited these two references at the relevant point in\nthe updated manuscript:\n“Ancestral sequence reconstrucGon is a powerful approach to invesGgate the properGes of the ancestral\niORFs that gave rise to de novo genes23,47,61”\nL121: What is the length distribution? Is length the only selection criterion? This would likely bias for certain\nage groups.\nWe defined iORFs based on two criteria: (i) a minimum length of 60 nucleotides and (ii) the absence of\noverlap with any annotated genomic feature, to avoid confounding effects from overlapping elements (e.g.,\nnoncoding ORFs in CDS alternative reading frames potentially influenced by the CDS codons).\nWe focused on ORFs longer than 60 nucleotides for two main reasons:\n(i) To ensure the robustness of downstream analyses such as HCA, TANGO, and IUPred, codon and amino\nacid frequencies, GC content assessment, homologous searches (DIAMONDv2), and selection tests (dN/dS)\nfor yeast iORFs. Very short sequences tend to yield unstable or uninformative predictions in these contexts.\n(ii) To maintain consistency with our previous work30,48 and related ongoing studies, facilitating\nreproducibility and interoperability across studies.\nvery short ORFs (< 20 codons). Calculating the length distribution for all real iORFs - or even subsets of iORFs\n- without a minimal length threshold across all genomes would be computationally prohibitive, as the\nnumber of ORFs increases dramatically when the threshold is lowered.\nWe acknowledge that this length threshold may exclude some potentially relevant short ORFs. However,\nORFs shorter than 20 codons pose limitations for the structural, evolutionary, and compositional analyses\ncentral to our study, making this threshold a necessary trade-off between inclusiveness and analytical\nreliability.\n200\n100\n0\n0 25 50 75 100\nMedian GC content\nFigure F: Length of random sequences as a func5on of GC content. ORF length as a func0on of GC content for\nORFs from random DNA sequences. Therefore, we randomly generated DNA sequences (one million nucleo0des\neach) ranging from 0% to 100% GC content in 0.5% increments (201 bins). ORFs (stop-to-stop) of at least one\ncodon were extracted from each random DNA sequence. The plot shows the median GC content and median\nlength (number of codons). The x-axis represents median GC content, and the y-axis represents the median length\n(number of codons).\nL123: What is being correlated with what, and why “although”? This sentence is unclear.\nWe apologize for the unclear wording in the original sentence. It has been rephrased as follows:\n“Consistent with previous studies49, GC content is strongly correlated between CDSs and noncoding regions\n(Pearson’s r = 0.85, P < 1x10-16), yet remains systemaGcally higher in CDSs (Figure S1A).”\nL125–6: Uncommon phrasing—should be more precise. Who are “they”?\nIn the revised version, we removed this sentence to reduce the word count.\n38\n)snodoc\nfo\nrebmun(\nhtgnel\nnaideM\nORFs (>= 1 codon) from random DNA sequences\nL128: Selection by what? Evolution?\nWe apologize for the lack of clarity in our previous formulation, which has now been revised as follows:\n“To assess GC effects on the codon frequencies of CDSs, we analyzed the distribuGon of various codon groups\nin CDSs and iORFs across the 3,379 eukaryoGc genomes (Figure 1A).”\nL131: “Inherited” is confusing here.\nThe reviewer is correct, and the sentence has been edited accordingly:\n“To discern whether codon frequencies are simply driven by the higher GC content of CDSs or reflect\nconstraints operaGng at the protein level, we generated scrambled CDSs preserving nucleoGde composiGon\nbut disrupGng codons and, consequently, protein-related constraints.”\nL132: Possibly meant: “due to evolutionary pressure to prefer certain amino acids”?\nTesting for codon bias would be useful here—note that codon preferences are generally\nstronger in lower eukaryotes (see recent papers by Laurent Duret and colleagues).\nWe combined the answers related to L132 and L133, since they are closely connected. Our response is\nprovided below (in the next comment).\nL133: While the scrambled CDS are a great idea, comparisons should also consider degenerate codons.\nWe are not entirely sure that we have accurately understood the reviewer's intention. Based on the\nReviewer’s comments and the recent article of Bénitière et al.58, we understood the question to be whether\nthe observed differences between real CDSs and scrambled CDSs can be attributed solely to protein-level\nconstraints or whether they might also be influenced by adaptive codon preferences (i.e., translational\nselection). Bénitière et al. reported, on a set of 223 metazoans, a strong correlation between average GC3\nand the GC content of introns (Figure 1B of Bénitière et al.58). This suggests, as mentioned by the authors,\nthat “nonadaptive processes are the primary factor driving the observed variation in codon usage across\nspecies”, at least in metazoans. To assess whether the effect holds for the 3,379 eukaryotes analyzed in the\nstudy, we computed the correlation between GC3 and the GC content of noncoding regions (Figure S1B).\nWe observe a similar trend, with a strong positive correlation (Pearson R2 = 0.64, r = 0.79, P < 1x10-16),\nsupporting the view that codon usage is predominantly determined by GC content.\nFigure G shows, for each amino acid coded by at least two codons, the correlation between the observed\nfrequencies of its associated codons in a genome and their expected frequencies based on the nucleotide\ncomposition of the genome’s noncoding regions. While most codons exhibit significant positive correlations\nacross the 3,379 genomes, some are associated with weaker correlations. As expected, these include\ncodons belonging to sixfold degenerate families (Leucine, Serine, and Arginine). Moreover, in fourfold\ndegenerate families, codons ending with G (NNG) generally show much weaker correlations between their\nobserved and expected frequencies than the other three synonymous codons (see Valine, for example).\nNotably, they often display low R2 values, suggesting that their variance is largely driven by factors other\nthan GC content. To our knowledge, this pattern has not been explicitly reported; its cause remains unclear\nand warrants further investigation beyond the scope of this study.\nenrichment/depletion signals independent of third-base usage. Thus, the unusual behavior of NNG codons\ndoes not affect our conclusions regarding codon group enrichment.\nGiven that the manuscript is already dense and because this observation, though interesting, opens\nquestions beyond the focus of our study, we chose not to include it in the manuscript. Instead, we\nintegrated the analysis on the relationship between GC3 and the GC content of noncoding regions (Figure\nS1B). We added the following sentence to the results section:\n“Then, to quantify whether differences between real and scrambled CDSs stem solely from protein-level\nconstraints or are also influenced by adaptive codon preferences (i.e., translational selection), we calculated\nthe correlation between GC3 and the GC content of noncoding regions across the 3,379 genomes. In line\nwith Bénitière et al.50, we observed a strong correlation (Pearson R2 = 0.64, r = 0.79, P < 1x10-16), supporting\nthat codon usage is primarily driven by GC content (Figure S1B).”\nPanel B of Figure S1: GC3 as a funcSon of interCDS GC content. For each genome, the mean GC3 content of its CDSs is ploked\nagainst the mean GC content of its interCDS regions. Each genome is colored according to its global GC content, ranging from\ngreen (AT-rich) to blue (GC-rich) (n = 3,379).\n40\nobserved frequencies of its codons in each genome (n = 3,379) rela9ve to their expected frequencies, calculated from the\nnucleo9de composi9on of the genome’s noncoding regions. The Pearson correla9on coefficients are given for each codon.\nL163–4: See again:\nIyengar, B.R., Grandchamp, A., Bornberg-Bauer, E. How antisense transcripts can evolve to encode novel\nproteins. Nature Communications 15(1), 6187 (2024).\nThis reference concerns antisense transcripts and cannot be directly related to iORFs, as these are\nconstrained by the evolution of the overlapping CDS.\nDue to the restructuring of this section, we have now removed the part of this sentence referring to the\nevolution of longer sequences. We now cite this reference in the Discussion, where we discuss ORF length\nand note that, since stop codons are AT-rich, they are depleted in GC-rich sequences.\n“Moreover, as stop codons are AT-rich, higher GC content likely produces longer ORFs79”\nAlso, add statistics to support “strong.” Why is it called “modular”? What is modular?\nWe apologize for not providing statistical support for the observations described in Figure 1. We now\ninclude the distribution of codon group frequencies across the 3,379 genomes for real and scrambled CDSs,\npresented as boxplots (Figure S1C), along with the associated mean log fold-changes (LFCs) (Figure 1B).\nThe TAN and TGN codon groups exhibit significantly lower frequencies in real CDSs compared to scrambled\nsequences (Mann-Whitney U test (one-sided), all P < 1×10-16), and display among the strongest negative\nLFCs, consistent with counter-selection in coding regions.\nThis is now explicitly mentioned in the main text:\n“Overall, most codons have absolute LFCs below 0.5, but four groups show stronger deviaGons (LFCs > 0.75):\nTAN, TGN, CGN, and GAN. Their median frequencies differ significantly from those of scrambled CDSs (Mann-\nWhitney U test, (two-sided), all P < 1×10-16)(Figure S1C). TAN and TGN comprise the three stop codons and\nthe rare amino acids tryptophan and cysteine, explaining their strong depleGon in CDSs relaGve to\nexpectaGons (Figure 1B, Figure S1C).”\nPanel B of Figure 1. Mean log2 fold-change (LFC) of codon group frequencies in real CDSs relative to scrambled CDSs across\nthe 3,379 genomes. Positive values indicate enrichment in real CDSs, negative values indicate depletion in real CDSs. Error\nbars represent standard deviation. Detailed distributions for each codon group are shown in Supplementary Figure S1C.\ndistributions were compared using a one-sided Mann-Whitney U test. Asterisks denote significance levels: *p < 5 × 10-2, **p\n< 1 × 10-2, ***p < 1 × 10-3.\nL167: Background in biophysics should be cited here. Foldability also depends on protein length; see work\nby Ken Dill and others, e.g.\nSen et al. Chemical Physics Letters 405, 224–228 (2005)\nWe thank the Reviewer for bringing to our attention the beautiful paper by K. Dill (1985)59, which we now\nexplicitly cite when discussing the optimal balance between hydrophobic and hydrophilic residues required\nfor proper protein folding.\n“This aligns with Dill’s classic theory, which proposed that protein stability depends on an optimal fraction\nof hydrophobic residues54.”\nWe also cite Dill’s article together with that of Shen et al. (2005)60 in the Discussion when addressing the\nrelationship between protein length and foldability.\n“Moreover, as stop codons are AT-rich, higher GC content likely produces longer ORFs79, which may enhance\nfoldability by facilitating burial of hydrophobic residues and exposure of hydrophilic ones, thereby resulting\nin higher stability54,80.”\nL173: Ribosomes are not freely floating in the cell.\nWe combined the responses to the two comments from the Reviewer because they are closely related. See\nour reply below.\nL178: What causes the enhancement of diffusion in the cytosol?\ninfluences diffusion in living cells. Using surface-modified fluorescent proteins (GFP variants) with varying\nnet charges, they measured cytoplasmic diffusion coefficients via FRAP in E. coli, Lactococcus lactis, and\nHaloferax volcanii. They demonstrated that in E. coli, positively charged proteins diffuse over 100× slower\nthan their negatively charged counterparts. They showed that this effect was weaker in L. lactis and H.\nvolcanii due to electrostatic screening.\nThe observed slowing of diffusion for positively charged proteins is attributed to nonspecific binding to\nnegatively charged ribosomes. Supporting this, sucrose gradient experiments showed that positively\ncharged +25 GFP co-localizes with ribosomes, whereas negatively charged −7 GFP does not.\nThese findings are highly relevant to our study because they provide a mechanistic explanation for how a\nnet negative charge on proteins can limit promiscuous ribosome interactions and thereby favor more\nefficient cytosolic mobility. While this work was conducted in prokaryotes, its biophysical principles are\nlikely applicable to eukaryotes and support our broader interpretation.\nRegarding our original wording, we agree that the phrase “enhance cytosolic diffusion” was inappropriate,\nas the effect is relative rather than absolute. We have revised the sentence to:\n“In yeast, we also found a correlation between negative charge frequency and protein abundance,\nsuggesting that a higher frequency of negative charges in abundant proteins may favor their translation\nefficiency and reduce nonspecific interactions with ribosomes23.”\nL181: Avoid “world” metaphors in general. Use “a hallmark” instead of “the hallmark.”\nWe thank the Reviewer again for bringing this to our attention. We removed the word “world” and revised\nthe sentence as follows:\n“Interestingly, these properties, ensuring proper foldability, expression, and interactability, together\nconstitute a hallmark of CDSs.”\nL187: Change to “encode preferentially…”\nDone.\n“In noncoding regions, where no selection acts on codons, AT-rich and GC-rich sequences are expected to\nencode preferentially hydrophobic and hydrophilic sequences, respectively.”\nL189: While the use of HCA over flawed structure predictions is appreciated, is HCA applicable to proteins\nof only ~20 amino acids?\nThese lack a hydrophobic core and likely do not fold in any case (see above). Anything shorter than 60\namino acids may give random results, unless they are extensions of existing proteins with a folding core, see work\nby Isabelle Callebaut and collaborators.\nWe thank the reviewer for raising this important point. We fully agree that protein length is a key\ndeterminant of folding, and we have already acknowledged this in the Discussion section (“Moreover, as\nstop codons are AT-rich, higher GC content likely produces longer ORFs79, which may enhance foldability by\nfacilitating burial of hydrophobic residues and exposure of hydrophilic ones, thereby resulting in higher\nstability54,80.”). In our study, we use the term folding potential to refer to intrinsic sequence features that\nmay favor folding or the evolution of a foldable protein, rather than to claim that these iORFs currently\n44\nadopt stable, well-packed, and well-defined 3D structures. [That said, we fully agree that this needed to be\nstated explicitly, and it is now clarified in the revised version - see below.]\nNotably, regarding the folding potential, we show that iORFs (depending on GC content) can encode\nstructural building blocks that resemble those found in folded proteins (see next paragraph and new Figure\nS2C). These structural building blocks are identified as HCA clusters. The latter have been shown to\ncorrespond to the elementary structural units of folded proteins62,63. The number, size, and composition of\nthese clusters, along with the intervening linker regions, are integrated into the HCA score, which serves as\na proxy for predicting an ORF’s folding potential62–64. These structural building blocks may not be sufficient\nto drive autonomous folding today, but they could facilitate folding in specific contexts (e.g., upon\ninteraction with partners, oligomerization, or chaperone assistance25), or evolve toward more stable\nstructures over time through length extension, refinement of residue composition and patterning, or\ndomain fusion. In this sense, these iORFs host latent structural potential that may be revealed through\nmolecular interactions, ORF extension and/or amino acid sequence optimization.\nNotably, we have now included a new analysis of HCA cluster sizes in iORFs and CDSs across the 3,379\neukaryotic genomes (new Figure S2C). We find that iORFs in intermediate GC-content species encode HCA\nclusters with size distributions comparable to those of CDSs, while extreme AT-rich and GC-rich species\ndeviate in opposite directions (larger and smaller clusters, respectively)(new Figure S2C). This suggests that\nintermediate GC-content species harbor a reservoir of iORFs whose sequence composition and organization\nconfer a higher propensity to adopt or evolve into what might be described as protofolds.\nOverall, our results suggest that such proto-structural features are more common in intermediate GC-\ncontent species, which encode iORFs with a better balance of hydrophobic and hydrophilic residues (Figures\n2A and S2C). In contrast, iORFs in GC-rich or AT-rich species often skew toward disorder-prone (high\nhydrophilic residue content) or aggregation-prone sequences (high hydrophobic residue content),\nrespectively.\nFurthermore, while we agree that very short proteins are unlikely to fold independently, in intermediate\nGC-content species (35-50% GC), the median iORF length is approximately 30-40 amino acids (depending\non GC content), and ancestral iORFs are even longer, with lengths often exceeding 60 residues (Figure S11),\nthereby likely further enhancing their folding potential. Overall, our results suggest that iORFs from\nintermediate GC-content species share more characteristics with CDSs (e.g., hydrophobic residue content,\nelementary structural units of proteins) than iORFs from species with extreme GC content. These findings\nsupport the hypothesis that this reservoir may include microproteins with the capacity to evolve more\nstable structures. Indeed, if selected, their folding potential could further increase through ORF extension,\nand/or optimization of amino acid composition and patterning via positive selection, as discussed in Aubel\net al.45.\nAltogether, our findings suggest that intermediate GC-content species offer a favorable background for the\nemergence of structurally competent novel microproteins, with more balanced amino acid compositions\nand appropriate structural motifs.\nWe have revised the Results and Discussion sections to incorporate these new analyses.\nL212 “Intermediate HCA scores reflect sequences with balanced proportions of hydrophobic clusters and\nlinkers expected to be associated with folding potential23,57. The latter may either fold into well-defined 3D\nstructures or harbor latent folding potential, i.e., the capacity to adopt structured conformations under\nspecific conditions, and/or evolve more stable structures. Figure 2A represents the HCA scores of iORFs and\nCDSs according to GC content. We show that GC content has a substantial impact on iORF foldability. AT-\nrich and GC-rich species predominantly harbor iORFs encoding aggregation-prone or disordered\nmicroproteins, respectively, while iORFs of intermediate GC-content genomes tend to display intermediate\nHCA scores, reflecting a greater folding potential.”\n45\n1.00 *** *** *** *** ns ns *** *** *** *** *** *** ***\n0.75\nL235 “Previously, we showed that yeast iORFs had HCA cluster sizes similar to CDSs23. To assess the impact\n0.50\nof GC content on the size of HCA clusters, we analyzed their size distributions across species. Figure S2C\nshows that in0.t2e5rmediate GC-content iORFs display HCA cluster sizes comparable to CDSs, while extreme AT-\nrich and GC-rich iORFs deviate in opposite directions (larger and smaller clusters, respectively). This suggests\n0.00\nthat intermediate G 15−C20 -con 20t−2e 5 nt s 2 p5−e30 cies 3 h0−3a 5 rbo 3 r5− a40 res 4 e0−r45 voir 4 o5−f50 iORF 50−s5 5 pote 55n−60 tially 60 −e65 ncod 65−i7n 0 g th 7 e0−7 5 struc 75t−8u 0 ral building\nblocks of proteins.”\nGenome GC content (%)\nL461 “Our results suggest that, depending on GC content, iORFs may encode structural building blocks (HCA\nclusters) resembling those found in folded proteins (Figure S2C). While these structural building blocks may\nnot be sufficient to drive autonomous folding today, they carry latent folding potential that could facilitate\nfolding in specific contexts or promote the evolution of more stable structures through length extension,\nand/or optimization of amino acid composition and patterning via positive selection, as discussed in Aubel\net al42.”\nPanel C of Figure S2: HCA cluster length of CDSs and iORFs relaSve to GC content. Boxplots showing the distribu9on of HCA\ncluster length (number of amino acids) with respect to GC content for CDSs (orange) and iORFs (purple) across the 3,379\neukaryo9c genomes. Each species is represented by 1,000 CDSs and 1,000 iORFs that were randomly selected. Due to the\nlarge sample sizes, medians of the distribu9ons were compared using an empirical permuta9on test (see Methods). Asterisks denote significance levels: *p < 5 × 10-2, **p < 1 × 10-2, ***p < 1 × 10-3.\nL208: What is the cause and what is the consequence in the relationship between DNA and protein\nmaterial?\nWe are not entirely sure we have correctly understood the reviewer’s question. If the question concerns\nthe direction of causality between DNA and protein constraints - that is, who shaped whom - our study does\nnot aim to resolve this. Instead, our results indicate that protein foldability imposes a strong constraint on\nCDSs in eukaryotic species, regardless of their underlying nucleotide composition.\nThe phrase “reconciliation between DNA and protein constraints” was not intended to imply a\nunidirectional causal relationship, but to emphasize the dual requirements that coding sequences must\nsatisfy: encoding foldable proteins while being composed of nucleotides that vary in GC content.\nNotably, our findings suggest that the need to maintain protein foldability may have influenced the\nestablishment and structure of the genetic code itself. This suggests a co-evolutionary scenario in which\nfoldability acted as a key constraint, shaping the correspondence between nucleotide triplets and amino\nacids, thereby enabling CDSs to meet both DNA and protein-level constraints.\nWe have revised the corresponding paragraph to clarify this point:\n46\nni\nseudiser\nfo\nnoitroporP\n)derPUI(\nsnoiger\nderedrosid\nA\n1.00 *** *** *** *** *** *** *** *** *** *** *** *** ***\n0.75\n0.50\n0.25\n0.00\n15−20 20−25 25−30 30−35 35−40 40−45 45−50 50−55 55−60 60−65 65−70 70−75 75−80\nGenome GC content (%)\nni\nseudiser\nfo\nnoitroporP\n)OGNAT(\nsnoiger\nenorp−noitagergga\nB\n25 *** *** *** * * ns ns *** *** *** *** *** ***\n20\n15\n10\n5\n0\n15−20 20−25 25−30 30−35 35−40 40−45 45−50 50−55 55−60 60−65 65−70 70−75 75−80\nGenome GC content (%)\nhtgnel\nretsulc\nACH\n)sAA\nfo\nrebmun(\nC\n*** *** *** *** *** *** *** *** *** *** *** *** ***\n1000\n500\n0\n15−20 20−25 25−30 30−35 35−40 40−45 45−50 50−55 55−60 60−65 65−70 70−75 75−80\nGenome GC content (%)\nhtgneL\n)snodoc\nfo\nrebmun(\nD\nCDSs iORFs\n“Conversely, although the GC content slightly influences CDS foldability, all genomes generally encode\nfoldable CDSs and display tighter HCA score distributions than iORFs, suggesting foldability as a key\nconstraint in CDSs. This observation highlights the dual pressures imposed by both DNA composition and\nprotein structure and suggests that foldability has played a pivotal role in the establishment of the genetic\ncode, allowing CDSs to maintain foldability regardless of GC content.”\nL215: It would be interesting to know what “these” organisms are. Is there truly no phylogenetic signal?\nEven if not, this should be demonstrated.\nWe thank the Reviewer for this insightful question. We investigated the taxonomic distribution of the 3,379\nspecies with respect to GC content. We indeed observe distinct distributions across the tree of life\ndepending on GC content. These results are now compiled within a new table (Table S1) and presented in\nthe Results section.\n“We next examined whether a phylogenetic signal exists across the 3,379 species based on their GC content\n(Table S1). GC-rich species mostly consist of Fungi (82%), with the remaining species belonging to SAR (8%),\nChlorophytes (4%), and Discoba (3.5%). In contrast, AT-rich species are more widely distributed across the\ntree of life, with 67% Opisthokontes (mainly Arthropods, Nematodes, and Fungi), 22% Streptophytes, and\n7% SAR (almost all Apicomplexa). Intermediate GC-content species display a distinct pattern, dominated by\nOpisthokontes (88%) - split between Metazoa (57%) and Fungi (43%), and small proportions of\nStreptophytes (9%) and SAR (3%).”\nL382: “Modularity” is again confusing here, since no extensions, domains, or similar features are analyzed.\nWhile proteins are known to evolve through domain shuffling, ORF extension, etc, we agree with the\nReviewer that this aspect has not been explicitly tackled in the present manuscript. We have therefore\nremoved this statement, and the sentence now sticks to the current data and reads as follows:\n“Moreover, in line with previous studies6,22,66,67, we confirm that CDS size evolution is a continuous process,\nas evidenced by the increase in ORF size across the phylostrata (Figure 4, S10).”\nL411: See again : Iyengar, B.R., et al Nature Communications (2024).\nDone. This reference is now cited as follows:\n“Moreover, as stop codons are AT-rich, higher GC content likely produces longer ORFs79”\nL414: Add substantial references here.\nWe thank the Reviewer again for highlighting the missing references here and have added the following\ncitations:\nDill, Biochemistry, 1985\nShen et al., Chem. Phys. Lett., 2005\nAubel et al., Genome Biol. Evol., 2024\n“Moreover, as stop codons are AT-rich, higher GC content likely produces longer ORFs79, which may enhance\nfoldability by facilitating burial of hydrophobic residues and exposure of hydrophilic ones, thereby resulting\nin higher stability .”\n42,54,80\nSince we needed to reduce the length of the revised manuscript, we removed this sentence.\nL424: Why not try running some aggregation prediction tools?\nWe thank the Reviewer for this insightful comment. We therefore used TANGO and IUPred to assess the\naggregation and disorder propensities of the ORFs (1,000 iORFs and CDSs per genome) across the 3,379\ngenomes. This analysis, now presented in Figure S2 (panels A and B) and discussed in the Results section,\nfurther supports the trends observed with HCA and strengthens our conclusions.\nL220 “As a control, we predicted the disorder and aggregaGon propensiGes of iORFs using IUPred58 and\nTANGO59 (Figure S2A-B). Results align with HCA predicGons, showing that disorder propensity posiGvely\ncorrelates with GC content, while aggregaGon propensity anGcorrelates with it.“\nL506 “In contrast, in extreme GC-content species, iORFs rarely meet CDS-like criteria, with AT-rich and GC-\nrich genomes skewed toward aggregaGon-prone or disorder-prone sequences, respecGvely (Figure S2A-B).”\n1.00 *** *** *** *** ns ns *** *** *** *** *** *** ***\n0.75\n0.50\n0.25\n0.00\n15−20 20−25 25−30 30−35 35−40 40−45 45−50 50−55 55−60 60−65 65−70 70−75 75−80\nGenome GC content (%)\nPanel A of Figure S2: AggregaSon propensity of CDSs and iORFs with respect to GC content. A Boxplots showing the\ndistribu9on of disorder propensity with respect to GC content, as predicted by IUPred, for CDSs (orange) and iORFs (purple)\nacross the 3,379 eukaryo9c genomes. Each species is represented by 1,000 CDSs and 1,000 iORFs that were randomly selected.\nDue to the large sample sizes, medians of the distribu9ons were compared using an empirical permuta9on test (see Methods).\nAsterisks denote significance levels: *p < 5 × 10-2, **p < 1 × 10-2, ***p < 1 × 10-3.\nL432: Please explain which model of de novo gene emergence is being referred to.\nWe apologize, but the Reviewer’s comment is unclear to us, as the sentence in question (see L432 of the\ninitial version) does not refer to any model of de novo gene emergence. The original sentence was:\n“Nevertheless, we cannot rule out that since genomes with extreme GC contents often carry many\nnucleotide regions with low complexity, this may affect the detection of de novo genes when searching for\ntheir homologous noncoding regions in the neighboring species.”\nIn any case, the sentence has been revised in the current version as follows:\n“Additionally, low-complexity sequences, frequent in extreme GC-content genomes, may also weaken\nhomology signals between de novo candidates and their orthologous noncoding regions in the neighbors.”\n48\nni\nseudiser\nfo\nnoitroporP\n)derPUI(\nsnoiger\nderedrosid\nA\n1.00 *** *** *** *** *** *** *** *** *** *** *** *** ***\n0.75\n0.50\n0.25\n0.00\n15−20 20−25 25−30 30−35 35−40 40−45 45−50 50−55 55−60 60−65 65−70 70−75 75−80\nGenome GC content (%)\nni\nseudiser\nfo\nnoitroporP\n)OGNAT(\nsnoiger\nenorp−noitagergga\nB\n25 *** *** *** * * ns ns *** *** *** *** *** ***\n20\n15\n10\n5\n0\n15−20 20−25 25−30 30−35 35−40 40−45 45−50 50−55 55−60 60−65 65−70 70−75 75−80\nGenome GC content (%)\nhtgnel\nretsulc\nACH\n)sAA\nfo\nrebmun(\nC\n*** *** *** *** *** *** *** *** *** *** *** *** ***\n1000\n500\n0\n15−20 20−25 25−30 30−35 35−40 40−45 45−50 50−55 55−60 60−65 65−70 70−75 75−80\nGenome GC content (%)\nhtgneL\n)snodoc\nfo\nrebmun(\nD\nCDSs iORFs\nL443: Add in brackets what this species actually is.\nThe sentence has been edited for greater clarity:\n“E. histolytica and Rotaria sp. have only one de novo gene despite homology profiles similar to A. thaliana\n(99 de novo genes), whereas C. javanica harbors 63 de novo genes despite a poor homology detectability.”\nL460 & L478: Definitions of “effort zone,” “minimal effort,” and “GC space” seem unnecessary and may not\naid in explaining the results.\nWe agree with the Reviewer’s comment and have revised the sentence accordingly.\n“These results enable us to highlight two zones where GC contents are either extreme (< 35% GC and > 50%\nGC) or intermediate (35% GC – 50% GC) (Figure 5).”\nL472: Strange reference formatting—check reference No. 58.\nAs part of our efforts to reduce the length of the manuscript, this statement has been removed from the\nrevised manuscript, where we focus instead on the percentage of hydrophobic residues (31.7%) and GAN\nfrequencies (11.6%) calculated for the CDSs across the 3,379 genomes.\n“Interestingly, achieving the 31.7% hydrophobic residue content observed in CDSs across the 3,379 genomes\nrequires ~44% GC by chance, and although the CDS GAN frequency (~11.6%) cannot be reached randomly\nat any GC level, its maximum occurs near 50% GC (Figure S16).”\nL496: What was at odds with what in these five publications?\nWe thank the Reviewer for this comment. Our initial wording may have unintentionally suggested that the\ncited studies were contradictory, which was not our intention. We have revised the sentence to avoid this\nimplication and to better reflect that our work builds upon and integrates previous observations.\n“Altogether, this study offers a computational framework that unifies previous observations on de novo\ngene emergence6,23,37,41,86 by analyzing 3,379 eukaryotic genomes across broad GC-content and\nphylogenetic spectra.”\nReferences\n1. Kumar, S. et al. TimeTree 5: An Expanded Resource for Species Divergence Times. Mol. Biol. Evol. 39,\nmsac174 (2022).\n2. Romero, M., Cerritos, R. & Ximenez, C. Horizontal Gene Transfers from Bacteria to Entamoeba Complex:\nA Strategy for Da6ng Events along Species Divergence. J. Parasitol. Res. 2016, 3241027 (2016).\n49\n3. Cissé, O. H. et al. Genomic insights into the host specific adapta6on of the Pneumocys6s genus.\nCommun. Biol. 4, 305 (2021).\n4. Gladyshev, E. A., Meselson, M. & Arkhipova, I. R. Massive horizontal gene transfer in bdelloid ro6fers.\nScience 320, 1210–1213 (2008).\n5. Eyres, I. et al. Horizontal gene transfer in bdelloid ro6fers is ancient, ongoing and more frequent in\nspecies from desicca6ng habitats. BMC Biol. 13, 90 (2015).\n6. Craig, R. J. et al. The Chlamydomonas Genome Project, version 6: Reference assemblies for ma6ng-type\nplus and minus strains reveal extensive structural muta6on in the laboratory. Plant Cell 35, 644–672\n(2023).\n7. López-Cortegano, E., Craig, R. J., Chebib, J., Balogun, E. J. & Keightley, P. D. Rates and spectra of de novo\nstructural muta6ons in Chlamydomonas reinhard6i. Genome Res. 33, 45–60 (2023).\n8. Jąkalski, M. et al. Compara6ve genomic analysis of retrogene repertoire in two green algae Volvox\ncarteri and Chlamydomonas reinhard6i. Biol. Direct 11, 35 (2016).\n9. Roginski, P., Grandchamp, A., Quignot, C. & Lopes, A. De Novo Emerged Gene Search in Eukaryotes with\nDENSE. Genome Biol. Evol. 16, evae159 (2024).\n10. Letcher, B., Maciuca, S. & Iqbal, Z. Role for gene conversion in the evolu6on of cell-surface an6gens of\nthe malaria parasite Plasmodium falciparum. PLOS Biol. 22, e3002507 (2024).\n11. Claessens, A. et al. Genera6on of An6genic Diversity in Plasmodium falciparum by Structured\nRearrangement of Var Genes During Mitosis. PLOS Genet. 10, e1004812 (2014).\n12. Freitas-Junior, L. H. et al. Frequent ectopic recombina6on of virulence factor genes in telomeric\nchromosome clusters of P. falciparum. Nature 407, 1018–1022 (2000).\n13. Lo(cid:143)us, B. et al. The genome of the pro6st parasite Entamoeba histoly6ca. Nature 433, 865–868 (2005).\n14. Lorenzi, H. et al. Genome wide survey, discovery and evolu6on of repe66ve elements in three\nEntamoeba species. BMC Genomics 9, 595 (2008).\n15. Yue, J.-X. et al. Contras6ng evolu6onary genome dynamics between domes6cated and wild yeasts. Nat.\nGenet. 49, 913–924 (2017).\n50\n16. Scannell, D. R. et al. The Awesome Power of Yeast Evolu6onary Gene6cs: New Genome Sequences and\nStrain Resources for the Saccharomyces sensu stricto Genus. G3 Bethesda Md 1, 11–25 (2011).\n17. Vakirlis, N. et al. De novo emergence of adap6ve membrane proteins from thymine-rich genomic\nsequences. Nat. Commun. 11, 781 (2020).\n18. Shen, X.-X. et al. Tempo and Mode of Genome Evolu6on in the Budding Yeast Subphylum. Cell 175,\n1533-1545.e20 (2018).\n19. Wacholder, A. et al. A vast evolu6onarily transient translatome contributes to phenotype and fitness.\nCell Syst. 14, 363-381.e8 (2023).\n20. Müller, A., MacCallum, R. M. & Sternberg, M. J. E. Structural Characteriza6on of the Human Proteome.\nGenome Res. 12, 1625–1641 (2002).\n21. Lindeløv, J. mcp: An R Package for Regression With Mul6ple Change Points. Preprint at\nhCps://doi.org/10.31219/osf.io/fzqxv (2020).\n22. Vehtari, A., Gelman, A. & Gabry, J. Prac6cal Bayesian model evalua6on using leave-one-out cross-\nvalida6on and WAIC. Stat. Comput. 27, 1413–1432 (2017).\n23. Dosztányi, Z. Predic6on of protein disorder based on IUPred. Protein Sci. Publ. Protein Soc. 27, 331–340\n(2018).\n24. Linding, R., Schymkowitz, J., Rousseau, F., Diella, F. & Serrano, L. A Compara6ve Study of the Rela6onship\nBetween Protein Structure and β-Aggrega6on in Globular and Intrinsically Disordered Proteins. J. Mol.\nBiol. 342, 345–353 (2004).\n25. Heames, B. et al. Experimental characteriza6on of de novo proteins and their unevolved random-\nsequence counterparts. Nat. Ecol. Evol. 7, 570–580 (2023).\n26. Slavoff, S. A. et al. Pep6domic discovery of short open reading frame–encoded pep6des in human cells.\nNat. Chem. Biol. 9, 59–64 (2013).\n27. Ruiz Cuevas, M. V. et al. Most non-canonical proteins uniquely populate the proteome or\nimmunopep6dome. Cell Rep. 34, 108815 (2021).\n28. Ingolia, N. T., Lareau, L. F. & Weissman, J. S. Ribosome Profiling of Mouse Embryonic Stem Cells Reveals\nthe Complexity and Dynamics of Mammalian Proteomes. Cell 147, 789–802 (2011).\n29. Na, C. H. et al. Discovery of noncanonical transla6on ini6a6on sites through mass spectrometric analysis\nof protein N termini. Genome Res. 28, 25–36 (2018).\n30. Papadopoulos, C. et al. The ribosome profiling landscape of yeast reveals a high diversity in pervasive\ntransla6on. Genome Biol. 25, 268 (2024).\n31. Zheng, E. B. & Zhao, L. Protein evidence of unannotated ORFs in Drosophila reveals diversity in the\nevolu6on and proper6es of young proteins. eLife 11, e78772 (2022).\n32. Barrera-Redondo, J., Lotharukpong, J. S., Drost, H.-G. & Coelho, S. M. Uncovering gene-family founder\nevents during major evolu6onary transi6ons in animals, plants and fungi using GenEra. Genome Biol.\n24, 54 (2023).\n33. Buchfink, B., Reuter, K. & Drost, H.-G. Sensi6ve protein alignments at tree-of-life scale using DIAMOND.\nNat. Methods 18, 366–368 (2021).\n34. Mistry, J., Finn, R. D., Eddy, S. R., Bateman, A. & Punta, M. Challenges in homology search: HMMER3\nand convergent evolu6on of coiled-coil regions. Nucleic Acids Res. 41, e121 (2013).\n35. Lin, Z. et al. Evolu6onary-scale predic6on of atomic-level protein structure with a language model.\nScience 379, 1123–1130 (2023).\n36. Barrio-Hernandez, I. et al. Clustering predicted structures at the scale of the known protein universe.\nNature 622, 637–645 (2023).\n37. van Kempen, M. et al. Fast and accurate protein structure search with Foldseek. Nat. Biotechnol. 42,\n243–246 (2024).\n38. Zhang, Y. & Skolnick, J. TM-align: a protein structure alignment algorithm based on the TM-score.\nNucleic Acids Res. 33, 2302–2309 (2005).\n39. Löytynoja, A. Phylogeny-aware alignment with PRANK. in MulGple Sequence Alignment Methods (ed.\nRussell, D. J.) 155–170 (Humana Press, Totowa, NJ, 2014). doi:10.1007/978-1-62703-646-7_10.\n40. Hubisz, M. J., Pollard, K. S. & Siepel, A. PHAST and RPHAST: phylogene6c analysis with space/6me\nmodels. Brief. Bioinform. 12, 41–51 (2011).\n41. Vakirlis, N., Acar, O., Cherupally, V. & Carvunis, A.-R. Ancestral Sequence Reconstruc6on as a Tool to\nDetect and Study De Novo Gene Emergence. Genome Biol. Evol. 16, evae151 (2024).\n52\n42. Lange, A. et al. Structural and func6onal characteriza6on of a puta6ve de novo gene in Drosophila. Nat.\nCommun. 12, 1667 (2021).\n43. Bungard, D. et al. Foldability of a Natural De Novo Evolved Protein. Struct. Lond. Engl. 1993 25, 1687-\n1696.e4 (2017).\n44. Matsuo, T. et al. Secondary Structure of Human De Novo Evolved Gene Product NCYM Analyzed by\nVacuum-Ultraviolet Circular Dichroism. Front. Oncol. 11, 688852 (2021).\n45. Aubel, M. et al. High-throughput Selec6on of Human de novo-emerged sORFs with High Folding\nPoten6al. Genome Biol. Evol. 16, evae069 (2024).\n46. Middendorf, L., Ravi Iyengar, B. & Eicholt, L. A. Sequence, Structure, and Func6onal Space of Drosophila\nDe Novo Proteins. Genome Biol. Evol. 16, evae176 (2024).\n47. Peng, J. & Zhao, L. The origin and structural evolu6on of de novo genes in Drosophila. Nat. Commun.\n15, 810 (2024).\n48. Papadopoulos, C. et al. Intergenic ORFs as elementary structural modules of de novo gene birth and\nprotein evolu6on. Genome Res. 31, 2303–2315 (2021).\n49. Chen, J. et al. The Rapid Evolu6on of De Novo Proteins in Structure and Complex. Genome Biol. Evol.\n16, evae107 (2024).\n50. Ingolia, N. T., Ghaemmaghami, S., Newman, J. R. S. & Weissman, J. S. Genome-wide analysis in vivo of\ntransla6on with nucleo6de resolu6on using ribosome profiling. Science 324, 218–223 (2009).\n51. Ingolia, N. T. et al. Ribosome Profiling Reveals Pervasive Transla6on Outside of Annotated Protein-\nCoding Genes. Cell Rep. 8, 1365–1379 (2014).\n52. Chen, J. et al. Pervasive func6onal transla6on of noncanonical human open reading frames. Science\n367, 1140–1146 (2020).\n53. Carvunis, A.-R. et al. Proto-genes and de novo gene birth. Nature 487, 370–374 (2012).\n54. Iyengar, B. R. & Bornberg-Bauer, E. Neutral Models of De Novo Gene Emergence Suggest that Gene\nEvolu6on has a Preferred Trajectory. Mol. Biol. Evol. 40, msad079 (2023).\n55. Grandchamp, A., Czuppon, P. & Bornberg-Bauer, E. Quan6fica6on and modeling of turnover dynamics\nof de novo transcripts in Drosophila melanogaster. Nucleic Acids Res. 52, 274–287 (2024).\n56. Grandchamp, A. et al. Popula6on genomics reveals mechanisms and dynamics of de novo expressed\nopen reading frame emergence in Drosophila melanogaster. Genome Res. 33, 872–890 (2023).\n53\n57. Tretyachenko, V. et al. Random protein sequences can form defined secondary structures and are well-\ntolerated in vivo. Sci. Rep. 7, 15449 (2017).\n58. Béni6ère, F., Lefébure, T. & Duret, L. Varia6on in the fitness impact of transla6onally op6mal codons\namong animals. Genome Res. 35, 446–458 (2025).\n59. Dill, K. A. Theory for the folding and stability of globular proteins. Biochemistry 24, 1501–1509 (1985).\n60. Shen, M., Davis, F. P. & Sali, A. The op6mal size of a globular protein domain: A simple sphere-packing\nmodel. Chem. Phys. Leu. 405, 224–228 (2005).\n61. Schavemaker, P. E., Śmigiel, W. M. & Poolman, B. Ribosome surface proper6es may impose limits on the\nnature of the cytoplasmic proteome. eLife 6, e30084 (2017).\n62. Faure, G. & Callebaut, I. Comprehensive Repertoire of Foldable Regions within Whole Genomes. PLOS\nComput. Biol. 9, e1003280 (2013).\n63. Bitard-Feildel, T. & Callebaut, I. Exploring the dark foldable proteome by considering hydrophobic amino\nacids topology. Sci. Rep. 7, 41425 (2017).\n64. Bitard-Feildel, T., Lamiable, A., Mornon, J.-P. & Callebaut, I. Order in Disorder as Observed by the\n‘Hydrophobic Cluster Analysis’ of Protein Sequences. Proteomics 18, e1800054 (2018).\nImpact of GC content on de novo gene birth\nRoginski et al.\nReviewer #1 (Remarks to the Author):\nThe authors have thoroughly and convincingly addressed all of my comments and concerns. I\ncommend them for the clarity and depth of their responses and I am happy to recommend this\nmanuscript for publication.\nThe only thing I would consider changing is the word \"reversal\" in the new Figure 5B legend. (reversal\nimplies going back to a past state, which I don't think fits here) Perhaps \"transition\" ?\nWe thank the Reviewer for their positive comments. We agree that “reversal” does not accurately\nreflect our intended meaning and have replaced it with “transition” throughout the manuscript.”\nReviewer #2 (Remarks to the Author):\nThe revised manuscript by Roginski et al. has been significantly improved. The scientific rigor, clarity,\nand figure quality have all been enhanced. I appreciate that the authors thoughtfully addressed my\ncomments as well as those of the other reviewers.\nOne minor point for the abstract: “Depending on GC content, these microproteins exhibit distinct\nstructural properties, suggesting different cellular impacts of pervasive expression.” This sentence\nmay need to be revised, as the study does not directly examine pervasive expression, and not all\nspecies exhibit it. A possible alternative could be: “...suggesting different cellular impacts if non-genic\nregions are pervasively expressed.”\nWe thank the Reviewer for their positive comments. We acknowledge that our previous formulation\ncould imply that we investigated translation in the studied genomes, which is not the case. We\ntherefore edited the sentence following the Reviewer’s suggestion.\n“Depending on GC content, these microproteins exhibit distinct structural properties, suggesting\ndifferent cellular impacts if non-genic regions are pervasively expressed.”\nReviewer #3 (Remarks to the Author):\nThe authors have implemented most of our comments, and the manuscript now reads more\nsmoothly and engages more thoroughly with the relevant literature. Technical aspects of the study\nhave also improved. We only have a few minor, outstanding comments:\nOverall, we thank the Reviewer for their positive feedback. We address all the points raised below.\nAuthors: “Genomes can be broadly subdivided into coding and noncoding regions, conceptually\nclassified as distinct and impermeable entities.”\n1\n- The term “subdivided” seems unnecessary; “divided” would suffice. The use of “impermeable” is\nproblematic, since the study addresses the transitions between coding and noncoding regions.\nI suggest:\n“Genomic regions can be divided into coding and noncoding regions, conceptually considered as\ndistinct entities.”\nThis sentence should also be supported with a citation. While this concerns functionality beyond\ncoding/noncoding annotation, see Graur et al. 2013, GBE regarding as how distinct annotations can\nbe.\nOur intent in this opening sentence was to introduce the historical and conceptual view of genomes\nas composed of two distinct domains (coding and noncoding), before highlighting that recent OMICS\nstudies have revealed the expression of presumed noncoding regions and, ultimately, a transition\nacross these boundaries.\nWe revised the sentence for clarity as follows:\n“Genomes can be broadly divided into coding and noncoding regions, conceptually considered as\ndistinct entities.”\nBecause this opening statement refers to the traditional way genomes are conceptualized rather\nthan to debates on genomic functionality (e.g., ENCODE), we prefer not to cite Graur et al. (2013),\nwhich focuses on functional interpretations beyond the intent of this sentence and might risk\nconfusing the reader.\nAuthors: “Using computational methods, we showed that yeast intergenic ORFs encode a wide\ndiversity of structural potentials, including aggregation-prone peptides, and a majority of ORFs with\nfolding potential. The latter are not necessarily expected to form well-defined 3D structures but\ncontain elementary structural blocks that may promote folding in specific contexts (e.g., binding with\na partner, chaperone assistance) or the evolution of more stable structures.”\n- Instead of “structural potentials”, I recommend “structural elements”. Similarly", "reviewer_rephrase": "They\nare not necessarily expected to form well-defined 3D structures on their own.", "reviewer_follows": "n“On the other hand, using computational methods, we showed that yeast intergenic ORFs encode a wide\ndiversity of structural elements, including aggregation-prone peptides, and a majority of ORFs with\nfolding potential23. The latter are not necessarily expected to form well-defined 3D structures on their\nown but contain elementary structural blocks that may promote folding in specific contexts (e.g., upon\nbinding to a partner or membrane, or through chaperone assistance39) and could support the evolution\nof more stable structures.”\nAuthors (L31–33): “We show that de novo genes preferentially emerge from GC-rich ORFs with folding\npotential, revealing that GC content and foldability, coordinated by the genetic code structure,\ntogether shape the emergence of novel genes.”\n- This sentence is redundant. “Coordinated by the genetic code structure” does not add clarity.\nI suggest:\n“We show that de novo genes preferentially emerge from GC-rich ORFs with folding potential.”\nWe disagree that the mention of “coordination by the genetic code structure” is redundant and prefer\nto retain it as this concept is central to our work. Our results not only show that de novo genes tend to\nemerge from GC-rich ORFs with folding potential, but also reveal that the relationship between these\ntwo properties (which stems from the structure of the genetic code) influences the emergence of novel\ngenes. To clarify this point, we have slightly revised the sentence as follows:\n“We show that de novo genes preferentially emerge from GC-rich ORFs with folding potential,\nrevealing that the interplay between GC content and foldability - rooted in the structure of the genetic\ncode - shapes the emergence of novel genes.”\nAuthors (L43): “primarily influenced by nucleotide-level constraints”\n- Suggest rephrasing as:\n“primarily constrained at the nucleotide level.”\nUsing both influenced and constrained is redundant.\nDone.\n(L78): Replace “combinations” with “patterns.”\nDone.\n(L90): The word “Alternatively” is intended in relation to the preadaptation theory, but could\ncurrently be misread as referring to the preceding experimental results. Please state this explicitly.\nWe thank the Reviewer for pointing this out. We have revised the sentence as follows:\n“In contrast to the preadaptation model, the proto-gene model proposes a gradual transition, where de\nnovo genes evolve from transitory proto-genes arising from pervasive expression of nongenic\nsequences6.”\n(L445–462): The phrasing “increased/increasing GC content” could confuse readers, implying an\nevolutionary increase over time. Please rephrase more precisely, e.g.:\n“This higher GC content in reconstructed ancestors…”\nWe thank the Reviewer for this subtle remark. We replaced all occurrences of “increase in GC\ncontent/increased GC-content etc” with “higher GC content (in reconstructed ancestors/of ancestral\nORFs).”\nThe correct structural study of the de novo AFGP is:\nApologies, Her C, Yeh Y, Krishnan VV. The ensemble of conformations of antifreeze glycoproteins\n(AFGP8): A study using nuclear magnetic resonance spectroscopy. Biomol. Ther. 2019; 9(6): 235. — not\nBaalsrud et al. 2018., apologies.\nWe thank the Reviewer for bringing this study to our attention; we now cite it alongside the other\nexperimental characterization of de novo proteins.\n“While computational studies on the structural properties of de novo proteins have accumulated23,26,42,43,\nexperimental analyses remain scarce, with only a few de novo proteins structurally characterized44–48.“\n4"}, "author_response": ""}]

Decision letter