Contrary to dogma, evolutionarily young and dynamic genes can encode essential functions. We find that evolutionarily dynamic ZAD-ZNF genes, which encode the most abundant class of insect transcription factors, are more likely to encode essential functions in Drosophila melanogaster than ancient, conserved ZAD-ZNF genes. We focus on the Nicknack ZAD-ZNF gene, which is evolutionarily young, poorly retained in Drosophila species, and evolves under strong positive selection. Yet we find that it is necessary for larval development in D. melanogaster. We show that Nicknack encodes a heterochromatin-localizing protein like its paralog Oddjob, also an evolutionarily dynamic yet essential ZAD-ZNF gene. We find that the divergent D. simulans Nicknack protein can still localize to D. melanogaster heterochromatin and rescue viability of female but not male Nicknack-null D. melanogaster. Our findings suggest that innovation for rapidly changing heterochromatin functions might generally explain the essentiality of many evolutionarily dynamic ZAD-ZNF genes in insects.
Although organisms display enormous phenotypic diversity, their cellular organization and early development are highly conserved across broad taxonomic ranges (Miklos and Rubin, 1996). Such widespread conservation has led to a commonly held view that an ancient, conserved genetic architecture encodes fundamental biological functions, which appears to be largely borne out by comparative genomics (Miklos and Rubin, 1996; Krylov et al., 2003). However, recent studies demonstrated that 30% of 185 evolutionarily young genes in D. melanogaster (Chen et al., 2010)?have acquired roles in development, cell biology, and reproduction that render them essential for viability or fertility (Lee et al., 2017; Long et al., 2013; Ross et al., 2013; Lee et al., 2019). Sometimes, evolutionary turnover of genes underlying essential cellular processes can occur, as seen in the evolution of kinetochore proteins (Drinnenberg et al., 2016). Even when they are retained over long evolutionary periods, genes encoding essential functions can evolve unexpectedly rapidly across plants and animals (Malik and Henikoff, 2001; Talbert et al., 2004). Thus, at least a subset of essential functions is encoded by rapidly evolving genes or genes that are subject to high genetic turnover. The functional basis of this correlation, which runs counter to dogma, is unclear.
To study this unexpected class of rapidly evolving, essential genes further, we focused on the highly dynamic ZAD-ZNF gene family, which encodes several essential transcription factors. ZAD-ZNF proteins contain a conserved N-terminal ZAD (Zinc-finger-associated domain), a linker, and a C-terminal domain that includes tandem C2H2 zinc fingers (Chung et al., 2002; Lespinet et al., 2002). The ZAD facilitates protein-protein interactions but does not have DNA-binding ability, whereas the C2H2 zinc fingers often mediate sequence-specific DNA binding (Jauch et al., 2003). Unlike the ZAD and ZNF domains, which are homologous between different ZAD-ZNF proteins, the linker domains are highly variable across ZAD-ZNF proteins in both sequence and length and have no discernible structural motifs. ZAD-ZNF genes arose in the ancestor of vertebrates and arthropods, but dramatically expanded within insect lineages (Chung et al., 2007), becoming the most abundant class of TFs in many genomes, including in D. melanogaster (Chung et al., 2002). However, ZAD-ZNF gene repertoires can vary quite extensively across insect lineages via gene gains and losses (Chung et al., 2002; Chung et al., 2007); the causes and consequences of this gene dynamism are poorly understood.
Insect ZAD-ZNF proteins might act analogously to mammalian KRAB-Zinc finger (KZNFs), providing an explanation for their gene dynamism (Chung et al., 2002; Chung et al., 2007). KZNF genes arose in the ancestor of tetrapods and expanded through lineage-specific gene amplifications to become the most abundant class of transcription factors present in mammalian genomes (Emerson and Thomas, 2009; Imbeault et al., 2017). KZNF proteins contain an N-terminal KRAB (Krüppel-associated box) domain and C-terminal arrays of C2H2 zinc-finger domains that define their DNA-binding specificities. Many KZNF genes play a critical role in genome defense through recognition and repression of transposable elements (Wolf and Goff, 2009; Rowe et al., 2010). As a result, KZNF genes involved in genome defense are subject to rapid genetic turnover and positive selection of their DNA-binding domains (Thomas and Schneider, 2011). Similar functions might have driven the diversification of ZAD-ZNF genes in insects (Chung et al., 2002; Chung et al., 2007).
Despite their abundance, only a few ZAD-ZNFs have been functionally characterized, with studies restricted to D. melanogaster. Approximately half of all ZAD-ZNF genes are highly expressed in ovaries and early embryos (Chung et al., 2002), where they might play crucial roles for fertility or development. Biochemical characterization of 21 ZAD-ZNF proteins shows that they bind to unique DNA consensus sequences, with putative targets in the regulatory regions of specific target genes (Krystel and Ayyanathan, 2013). For example, the ZAD-ZNF protein M1BP (Motif 1 Binding Protein) binds core promoters and promotes the expression of numerous housekeeping genes (Baumann et al., 2017; Li and Gilmour, 2013). Similarly, ZAD-ZNF gene Grauzone promotes cortex expression and is necessary for meiotic progression during oogenesis (Chen et al., 2000; Harms et al., 2000; Page and Orr-Weaver, 1996). Finally ZAD-ZNF genes Molting Defective, Ouija Board, and Séance encode proteins that promote the transcription of heterochromatin-embedded genes, Spookier and Neverland, required for larval progression (Uryu et al., 2018). In contrast, other ZAD-ZNF proteins repress rather than drive transcription. For example, the ZAD-ZNF protein Oddjob localizes to pericentromeric heterochromatin and is required for gene silencing (position-effect variegation) (Swenson et al., 2016). Another ZAD-ZNF protein encoded by CG17801 helps repress HetA and Blood transposable elements in the ovary (Czech et al., 2013).
Not all ZAD-ZNF functions are directly related to transcription. For example, the ZAD-ZNF proteins ZIPIC, Zw5, and Pita help organize chromatin architecture (Gaszner et al., 1999; Zolotarev et al., 2016) whereas the ZAD-ZNF protein Trade Embargo mediates the initiation of meiotic recombination during oogenesis (Lake et al., 2011). Finally, some ZAD-ZNFs might not function in the nucleus at all. Even though the Weckle ZAD-ZNF protein possesses zinc-finger domains, it localizes to the plasma membrane instead of nuclear chromatin, where it interacts with the Toll-MyD88 complex to help establish the anterior-posterior axis of the developing embryo (Chen et al., 2006).
The ZAD-ZNF gene family in Drosophila is ideal for studying the relationship between genetic innovation and essentiality because of its involvement with essential cellular processes despite rapid evolutionary dynamics. In order to study this relationship between rapid evolution and essentiality, we leveraged extensive phylogenomic and population genetics datasets in Drosophila species. We also?took advantage of genome-wide screens for phenotype and tools for cytological and genetic analyses in D. melanogaster. Using evolutionary analyses, we identified the D. melanogaster ZAD-ZNF genes that had been either subject to genetic turnover or positive selection. Although only a few ZAD-ZNF genes have undergone positive selection, we found that these genes are more likely to be required for viability or fertility in D. melanogaster than slowly-evolving ZAD-ZNF genes. We focused on the characterization of one of these positively-selected ZAD-ZNF genes: Nicknack (CG17802, Nnk). We show that Nicknack is essential for larval development in D. melanogaster despite being evolutionarily young and differentially retained among Drosophila species. Nicknack belongs to a small cluster of ZAD-ZNF paralogs, the best characterized of which is Oddjob (CG7357, Odj) (Swenson et al., 2016). We show that both Odj and Nnk encode heterochromatin-localizing proteins. Although Odj broadly localizes to heterochromatin, we found that Nnk predominantly localizes to discrete foci within heterochromatin. Surprisingly, despite a strong signature of positive selection, we found that the protein encoded by the divergent D. simulans Nicknack ortholog can still localize to heterochromatin in D. melanogaster cells. Furthermore, D. simulans Nicknack can significantly rescue the viability of D. melanogaster Nicknack-null females but is unable to rescue Nicknack-null males. Based on our functional and cytological analyses, we conclude that rapidly changing requirements for heterochromatin function likely drove essential innovation of ZAD-ZNF genes such as Nicknack and Oddjob in Drosophila.
We searched Flybase to identify all genes in D. melanogaster that encode a ZAD domain (PF07776; Pfam database, Pfam.org). We found 91 ZAD-ZNF genes distributed across Chromosomes 2, 3, and X. 37 ZAD-ZNF genes occur in 13 gene clusters, containing two or more tandemly-arrayed ZAD-ZNF genes. Of these, many ZAD-ZNF genes share intron/exon structures with their neighbors and likely arose via segmental duplication. In contrast, seven ZAD-ZNF genes (CG3032, CG4318, CG9215, CG44002, CG17361, CG7963, CG17359) lack introns found in their closest relatives; we infer that these genes were likely born via retrotransposition.
Although ZAD-ZNF proteins are defined as having both a ZAD and a ZNF domain, further analysis of ZAD-containing proteins using NCBI’s Conserved Domain Database revealed significant variation in the ZNF domains. We found that ZAD-containing proteins have an average of?six C2H2 domains. However, some ZAD-containing proteins have no C2H2 domains at all (dlip, dbr, CG15435, CG31109, CG31457), whereas others contain up to 23 C2H2 domains (CG11902). In addition to the ZAD in the N-terminus, 13 ZAD-ZNF proteins possess another N-terminal domain, such as a Sodium-Calcium exchanger domain (CG12391) or ASF1 histone chaperone domain (CG10321). The functional significance of these additional domains is unknown.
Next, we surveyed ZAD-ZNF genes in all 12 previously sequenced and annotated Drosophila genomes. These species represent a range of evolutionary divergence from D. melanogaster, from just a few million years (e.g., D. simulans), to more than 40 million years (e.g., D. virilis) (Drosophila 12 Genomes Consortium, 2007). Our analysis reveals a surprisingly wide range in number of ZAD-ZNF genes across different Drosophila species genomes (Figure 1). For example, we found that the D. persimilis genome encodes 130 ZAD-ZNF genes whereas the D. willistoni genome encodes only 75 ZAD-ZNF genes (Supplementary file 1). Such analyses are dependent on the state of completion and annotation of individual Drosophila species’ genomes. Thus, this number may be an underestimate. To complement these analyses, we assessed the apparent age of each D. melanogaster ZAD-ZNF gene by examining the presence of orthologs in the 12 annotated genomes of Drosophila species (flybase.org). We found that 74 of 91 D. melanogaster ZAD-ZNF genes arose in or prior to the common ancestor of all Drosophila species, 40 million years ago. Of these 74 ZAD-ZNF genes, 61 have been preserved over Drosophila evolution, whereas 13 have been lost in at least one lineage or species. We estimate that 17 ZAD-ZNF genes arose via gene duplication after the last common ancestor of D. melanogaster and D. virilis. At least 3 ZAD-ZNF genes found in D. melanogaster (CG4318, CG17612, neu2) originated via gene duplication less than 10 million years ago. Our findings complement previous large-scale surveys that identified rapid changes in ZAD-ZNF gene repertoires within insect genomes (Chung et al., 2002; Lespinet et al., 2002; Chung et al., 2007).
These rapid changes in ZAD-ZNF gene repertoires suggested that selection might favor their genetic innovation. We investigated whether evolutionary retention is a predictor of essentiality. We took advantage of the fact that knockdown or knockout phenotypes have been characterized for almost all D. melanogaster genes. Indeed, 85 of 91 D. melanogaster ZAD-ZNF genes have associated phenotypic outcomes (summarized in Supplementary file 1). Of these, knockdown or knockout of 22 ZAD-ZNF genes showed complete lethality or sterility in D. melanogaster, whereas the other 63 were determined not to be essential.
Of the 61 ZAD-ZNF genes with orthologs retained in all 12 annotated Drosophila species’ genomes, we found that 14 are essential for either fertility or viability in D. melanogaster, whereas 42 genes are not (five genes have no phenotypic data available). In comparison, we found that 8 of the 30 genes not universally conserved in Drosophila species are essential, whereas 21 are not (one has no phenotypic data available). Thus, somewhat surprisingly, we find that genes not globally retained over Drosophila evolution are just as likely to encode a necessary function in D. melanogaster as genes that have been strictly maintained over 40 million years of Drosophila evolution (8:21 versus 14:42, p=0.8, two-tailed Fisher’s exact test) (Table 1). These findings further support the idea that gene families subject to rapid evolutionary turnover may become involved in essential functions.
Gene duplication and loss is only one form of genetic innovation. We also analyzed the D. melanogaster ZAD-ZNF genes for signatures of recent positive selection. We performed McDonald-Kreitman tests of ZAD-ZNF orthologs found in both D. melanogaster and the closely-related species, D. simulans. The McDonald-Kreitman test assesses sequence diversity within a species versus divergence?between species by comparing the ratio of non-synonymous (amino-acid altering, or replacement) to synonymous substitutions fixed during the divergence of the two species (Dn: Ds), to that of non-synonymous to synonymous polymorphisms within a species (Pn: Ps) (McDonald and Kreitman, 1991). The Pn: Ps ratio is a proxy for functional constraint acting on a gene within species and is expected to be similar to Dn: Ds between species under the null hypothesis. However, a higher than expected number of fixed non-synonymous changes would indicate the action of adaptive evolution during species divergence (McDonald and Kreitman, 1991).
We took advantage of previous efforts that sequenced the genomes of hundreds of D. melanogaster strains and a reference D. simulans strain (Lack et al., 2016; Langley et al., 2012) to perform the McDonald-Kreitman test using the Popfly server (popfly.uab.cat) (Hervas et al., 2017). Of the 91 D. melanogaster ZAD-ZNF genes, only CG2202 is absent in D. simulans. We found that 12 out of the remaining 90 ZAD-ZNFs show evidence for recent adaptive evolution, i.e., have an excess of fixed non-synonymous changes (Dn) (Table 2). Using D. yakuba as an outgroup species, we also polarized fixed differences between D. melanogaster and D. simulans to assess whether the lineage leading to D. melanogaster showed evidence of positive selection in these 12 ZAD-ZNF genes. We found that 5 out of 12 genes showed evidence of positive selection with this polarized McDonald-Kreitman test (Table 2).
Subsequently, we used a domain-restricted analysis to define which of the three protein domains (ZAD, linker, or ZNF) were subject to positive selection. Although we found evidence of domain-specific positive selection in eight cases (two in the case of the polarized McDonald-Kreitman test), we were unable to define the domain subject to positive selection in four cases, due to small numbers of intra-species polymorphisms (Table 2). One gene (CG7386) showed signatures of adaptive evolution in its ZAD domain, whereas three genes (CG2712, CG7386, CG10321) showed evidence of adaptive evolution in the ZNF domain. Unexpectedly, we found that the linker domain has evolved under positive selection in six ZAD-ZNF genes. The biochemical function of these linker domains is largely unexplored as they lack predicted structural motifs (Chung et al., 2002; Chung et al., 2007), are highly variable in length and sequence even between orthologs, and are predicted to encode intrinsically disordered domains (Ishida and Kinoshita, 2007). Our findings implicate the poorly characterized linker region in mediating the adaptive potential of several ZAD-ZNF genes.
Using this signature of positive selection in a subset of ZAD-ZNF genes, we evaluated whether positively selected genes encode essential functions. Intriguingly, we found that 8 of 12 ZAD-ZNF genes that have evolved under positive selection are essential, whereas only 14 of the remaining 73 genes that have been phenotypically assayed are essential (six have no phenotypic data available). Thus, contrary to the dogma, we find that ZAD-ZNF genes that are subject to positive selection are more likely to be essential for viability or fertility (8:4 versus 14:59, two-tailed Fisher’s exact test, p=0.0016) (Table 1). Our findings not only imply that several essential ZAD-ZNF genes evolve rapidly, but also raise the possibility that rapid evolution of some ZAD-ZNF genes might be critical for organismal viability.
To further investigate the biological basis for the correlation between positive selection and gene essentiality, we decided to focus on one cluster of ZAD-ZNF genes on chromosome 3 of D. melanogaster (Figure 2A). This cluster of five genes contains two of the eight positively-selected, essential genes in D. melanogaster: Oddjob (Odj, CG7357) and CG17802; both also show evidence of positive selection in the polarized McDonald-Kreitman test (Table 2). Odj is the best-characterized gene in this cluster. The Oddjob protein interacts with Heterochromatin Protein 1a (HP1a), dynamically localizes to heterochromatin, and is required for heterochromatin-mediated gene silencing, or position-effect variegation (Swenson et al., 2016). Knockdown or over-expression of Odj is lethal in D. melanogaster, although the functional basis for this lethality is still unknown (Schertel et al., 2015). In keeping with the Oddjob nomenclature theme of ‘James Bond henchmen,’ and the mutant phenotype (described below), we renamed CG17802 as Nicknack, or Nnk. This cluster of ZAD-ZNF genes also includes three genes that do not evolve under positive selection (CG17801, CG17806, and CG17803) in D. melanogaster (Figure 2A). Although previous studies found that the knockdown of CG17801 in the germline led to de-repression of HetA and Blood transposable elements in the ovary, CG17801 knockdown did not significantly impair fertility or viability (Czech et al., 2013). The other two ZAD-ZNF genes in this cluster (CG17803, CG17806) have not been previously investigated.
We evaluated the functional roles of all five ZAD-ZNF genes in the Odj-Nnk cluster by generating knockdowns using an Actin5C-Gal4 driver and RNAi constructs specific to each gene (Figure 2B) to generate ubiquitous knockdown of individual genes in the resulting zygotes. Consistent with previous results (Schertel et al., 2015), we found that the ubiquitous knockdown of Odj or Nnk is lethal. In contrast, we found that knockdowns of CG17801, CG17803 and CG17806 are viable (Figure 2B). Based on these results, we conclude that the two positively-selected members of this ZAD-ZNF cluster – Odj and Nnk – are both essential, confirming the unexpected correlation we previously observed between positive selection and gene essentiality in the ZAD-ZNF gene family (Table 1B).
To gain deeper insight into the evolutionary dynamics of the Odj-Nnk cluster, we identified orthologs of these genes using reciprocal TBLASTN searches with each of the five D. melanogaster genes as queries. We searched both the originally-sequenced, well-annotated 12 Drosophila genomes (Drosophila 12 Genomes Consortium, 2007), as well as eight additional genomes that were subsequently sequenced to sample the melanogaster group within the Sophophora subgenus more densely (Chen et al., 2014). In several cases, we were not able to confidently assign the TBLASTN hits to orthologous groups because they matched closely to more than one D. melanogaster gene, or there were several putative hits within a single genome, or because the hit contained only the ZAD or ZNF domain.
To more rigorously identify orthologs, we conducted phylogenetic analyses (Figure 2—figure supplement 1) using a multiple alignment of the ZNF domain (Supplementary file 2). Our phylogenomic analyses reveal that Odj-Nnk cluster evolution was highly dynamic during the evolution of the Drosophila genus (summarized in Figure 2C). Although Oddjob orthologs are present throughout 40 million years of Drosophila evolution, no other member of the D. melanogaster Odj-Nnk cluster is universally present in Drosophila species. In addition to Odj, CG17801 and CG17806 also date back prior to the origin of Drosophila but unlike Odj, CG17801 and CG17806 have since been lost in some species (Figure 2C). While CG17801 has been lost in the obscura group,CG17806 underwent multiple independent duplication and loss events. CG17803 arose in the ancestor of D. melanogaster and D. ananassae and underwent two independent losses. Finally, Nnk appears to have arisen in the ancestor of D. melanogaster and D. pseudoobscura (~30 mya, Figure 1), and later experienced multiple independent duplications and losses (Figure 2C). We note that our estimates for the age of Nnk are higher than those reported previously (Chen et al., 2010), likely because the prior estimate was based on fewer sequenced species. Based on these analyses, we conclude that, despite its essential function in D. melanogaster, Nnk is not universally conserved in Drosophila species.
Nnk has a dramatic evolutionary history: young evolutionary age, differential retention, positive selection. Yet, it serves an essential function in D. melanogaster based on our and other previous analyses. However, this claim of essentiality has been challenged by two previous findings. First, the screen that originally identified Nnk as an essential gene used the KK RNA-interference (RNAi) collection from the Vienna Drosophila Stock Center (VDRC) (Figure 3A; Chen et al., 2010). Many lines in this ‘KK’ collection were later found to harbor a second-site mutation that caused lethality as a result of ectopic tiptop expression when crossed to GAL4-driver lines (Green et al., 2014; Vissers et al., 2016). As a result of its dependence on the KK lines, the claim of Nnk essentiality remained ambiguous. A second Nnk mutant allele, created by CRISPR-Cas9-mediated mutagenesis (hereafter referred to as CRISPR-null), had a four base pair deletion within the coding sequence of the gene that created a frameshift and a premature stop codon within the linker domain of Nnk (Kondo et al., 2017; Figure 3A). Although this deletion was homozygous lethal, it was unexpectedly viable when paired with a deficiency covering this ZAD-ZNF cluster (Kondo et al., 2017), again challenging the result that Nnk is an essential gene.
Given these ambiguous results and the importance of Nnk for our claims of ZAD-ZNF essentiality despite evolutionary innovation, we re-investigated whether Nnk is essential for viability in D. melanogaster. We found several lines of evidence to support the conclusion that it is indeed essential (Figure 3B–D). First, we began by validating the RNAi line used in the original study (Chen et al., 2010), which first identified Nnk as a young, essential gene. We found that the Nnk VDRC RNAi line does not have an insertion upstream of the tiptop gene, which is associated with lethality in other KK lines (Green et al., 2014; Vissers et al., 2016). Second, we were able to rescue this lethality via complementation with an intact Nnk-mel?rescue transgene. This Nnk-mel?rescue transgene was flanked by regulatory sequences from the endogenous locus (1 kb segments upstream and downstream of endogenous Nnk) and recoded via mutations in synonymous sites to make it resistant to RNAi knockdown (Figure 3B; Supplementary file 4). Endogenous levels of Nnk expression are too low for us to validate expression of the transgene. However, our ability to rescue Nnk knockdown-mediated inviability strongly imply that the recoded transgene is appropriately expressed (Figure 3B).
To rule out any indirect effects arising from RNAi knockdown, we also examined two previously generated mutant alleles of Nnk. The first is a piggyBac insertion two bp upstream of the start codon in the 5’ UTR of Nnk (Schuldiner et al., 2008;?Figure 3C). This piggyBac insertion, which is marked by a DsRed reporter driven by an eye-specific promoter, is flanked by stop codons in all three reading frames, which prevents translation downstream of the insertion. We refer to this insertion as a piggyBac-null allele of Nnk. We found that this insertion allele is homozygous lethal but viability can be fully rescued upon mobilization of the piggyBac element, which repairs the intact 5’ UTR in a ‘scarless’ fashion to restore Nnk function (schematized in Figure 3C). Furthermore, piggyBac-null flies can also be rescued by the Nnk-rescue transgene (Figure 3D; Figure 3—source data 1). Second, we re-examined the previously-generated CRISPR-null Nnk allele (Kondo et al., 2017; Figure 3A). Contrary to previous results, we found that this allele is lethal when paired with a Nnk-spanning deficiency (BL9207, Figure 3D) and also fails to complement the piggyBac-null allele (Figure 3D; Figure 3—source data 1; Supplementary file 4). Moreover, the Nnk-rescue transgene can restore the viability of the CRISPR-null Nnk allele (Figure 3D; Figure 3—source data 1). Based on all these results, we conclude that Nnk is unambiguously an essential gene in D. melanogaster.
Next, we investigated the developmental stage at which Nnk-null progeny die. For this, we crossed Nnk-null heterozygotes to each other (Figure 4A). These flies contain the piggyBac insertion upstream of the Nnk gene (Figure 3C) on one chromosome, along with a balancer chromosome (TM3G) marked with GFP and carrying a wildtype Nnk allele. Progeny homozygous for the balancer chromosome TM3G die as early embryos. Consequently, all larvae lacking GFP expression are Nnk-null homozygotes whereas those that are heterozygous express the GFP encoded on the balancer. We conducted egg-lay experiments from crosses between heterozygote Nnk-null flies and tracked the developmental progression of GFP-negative, Nnk-null progeny relative to their GFP-expressing heterozygote siblings. We found that Nnk-null progeny progress through embryogenesis at the same rate as their heterozygote siblings and are morphologically indistinguishable from heterozygotes until the L1 larval stage. The Nnk-null larvae are able to move toward and consume yeast paste much like their heterozygote siblings. However, when heterozygote siblings molt into the L2 stage 48 hr after egg laying (AEL), Nnk-null larvae do not molt (Figure 4B). Instead, 48 hr AEL, Nnk-null larvae progressively become unable to move or eat, eventually dying by 60 hr AEL.
The Nnk-null larval arrest phenotype is reminiscent of previous findings with the ZAD-ZNF genes Séance, Ouija Board and Molting Defective, which encode proteins necessary for expression of the Spookier and Neverland genes required for ecdysone biosynthesis (Uryu et al., 2018; Neubueser et al., 2005; Komura-Kawa et al., 2015). Defects in any of the Séance, Ouija Board and Molting Defective ZAD-ZNF genes leads to arrest of larval development and death (Uryu et al., 2018; Neubueser et al., 2005; Komura-Kawa et al., 2015). However, this lethality can be rescued by supplementing the diet with ecdysone or by overexpression of ecdysone biosynthetic enzymes, bypassing the requirement for these three ZAD-ZNF genes. We therefore, tested whether dietary supplementation with an ecdysone precursor, cholesterol (7DH, or 7-Dehydrocholesterol) or ecdysone (20E, or 20-Hydroxyecdysone) could bypass or significantly delay the death of Nnk-null larvae (Figure 4C). We found that it could not. In contrast, the same dietary supplementation is able to significantly restore the viability of Npc1a larvae, which lack an essential transporter of cholesterol (Figure 4C). Based on these results, we conclude that Nnk plays an essential role in larval progression in D. melanogaster that is distinct from known steps of the ecdysone biosynthesis pathway.
To further?investigate our developmental findings, we also performed RNA-seq analyses, comparing the transcriptomes of Nnk-null larvae to their age-matched heterozygote siblings. We compared the two genotypes at each of two important time-points. First, we compared transcriptomes 24 hr AEL, when the two genotypes are morphologically indistinguishable L1 larvae. At this timepoint, we find that 249 genes (2.4% of all expressed genes) are differentially expressed between Nnk-null larvae and control larvae, with 116 genes at least two-fold upregulated and 133 genes two-fold downregulated in Nnk-null larvae (Figure 4D). Among genes downregulated in Nnk-null L1 larvae, we find a significant over-representation of functional categories related to proteolysis and sterol transport (both important functions during larval development), as well as dopamine monooxygenase activity (Supplementary file 3). In contrast, we find that genes related to lysosome, cytochrome P450s (also important for larval molts), and eye-related functions are upregulated upon Nnk loss.
Second, we performed comparisons 48 hr AEL when the Nnk-null larvae are significantly smaller and appear developmentally arrested compared to their age-matched controls. At this timepoint, we find that 3027 (28.1% of all expressed genes) genes are differentially expressed, with 1301 genes at least two-fold upregulated and 1726 genes two-fold downregulated in the Nnk-null mutants compared to the age-matched controls (Figure 4E). Thus, there are significantly more genes affected by Nnk loss by 48 hr AEL. Intriguingly, clustering samples by the transcriptional profile of all genes shows that Nnk-null larvae at 48 hr AEL are transcriptionally more similar to control larvae at 24 hr AEL of either genotype than they are to age-matched control larvae (Figure 4F). The transcriptional status of Nnk-null larvae therefore mirrors the phenotypes we observe, displaying a severe developmentally arrested phenotype and transcriptional profile at 48 hr AEL (Figure 4B).
Most of the Odj-Nnk cluster of ZAD-ZNF genes are functionally uncharacterized. Since Odj encodes a protein that is highly enriched in pericentric heterochromatin (Swenson et al., 2016), we speculated that its close paralog, Nnk, might also encode a protein with heterochromatic localization in D. melanogaster cells. To test this possibility, we used transient transfections to introduce epitope-tagged Odj and Nnk genes into D. melanogaster Schneider 2 (S2) cells and induced their expression with a heat-shock promoter (Figure 5A). Upon induction, we confirm that Oddjob has a broad localization pattern within heterochromatin (marked by histone H3 lysine nine methylation, H3K9me3 and outlined with a dashed line,?Figure?5A). We found that the Nnk-encoded protein also localizes to heterochromatin, but its localization is restricted to discrete foci, unlike Oddjob (Figure 5A). Since heat-shock induction can alter chromatin properties in cells, we also employed a complementary transient transfection strategy in which we expressed mCherry epitope-tagged Odj and Nnk genes under the control of a constitutive pCopia promoter in S2 cells. These analyses also revealed a broad heterochromatic localization of Odj contrasting with discrete foci within heterochromatin for Nnk protein (Figure 5B). These foci do not overlap with centromeres (identified by the centromeric histone Cid) or dual-strand piRNA clusters (marked by the piRNA-binding HP1 protein Rhino) (Figure 5—figure supplement 1).
Odj has a broad heterochromatic localization in D. melanogaster cells, which could result from direct interaction with HP1a (Swenson et al., 2016). Indeed, Odj has two potential PxVxL motifs, which are putative interaction sites for HP1a (Smothers and Henikoff, 2000), in the linker and ZNF domains (Figure 5—figure supplement 2A). Investigating Odj orthologs in other Drosophila species revealed that the PxVxL motif in the linker domain is well-conserved, but the one in the ZNF domain is not. We found that mutation of the putative HP1a-interaction site in the linker domain (V164A) converted Odj localization from a broad to a discrete pattern that at least partially overlapped with Nnk (Figure 5—figure supplement 2B). In contrast, mutation of the second putative PxVxL site (V321A) did not significantly affect Odj localization. Our results suggest that the putative PxVxL motif in the Odj linker region is a major contributor to Oddjob’s broad localization to heterochromatin potentially by mediating a direct interaction with HP1a. In the absence of this interaction, Oddjob and Nnk (which lacks a canonical PxVxL motif) localize similarly to discrete foci within heterochromatin (Figure 5B). Our findings suggest that altered protein–protein interactions or DNA-binding specificity via the linker domain may provide a means for functional diversification between closely-related ZAD-ZNF paralogs.
To gain deeper insight into the heterochromatic localization of Nnk and Odj, we performed chromatin immunoprecipitation and sequencing (ChIP-seq) using transient transfection of S2 cells with the epitope-tagged pCopia constructs for Odj and Nnk. Consistent with its cytological localization (Figure 5A and B), we found that Odj was highly enriched throughout the large chromosome regions previously defined as pericentric heterochromatin (Hoskins et al., 2015) whereas Nnk-bound regions are more narrowly distributed within heterochromatin (Figure 5C). Finer-scale analyses suggested that Nnk and Odj signals are enriched close to transcription start sites classified by modENCODE as residing in TSS-proximal chromatin in S2 cells, both in heterochromatin and euchromatin (Figure 5—figure supplement 3;?Ho et al., 2014). However, many of these regions overlap with false positive ‘phantom peaks’ previously identified in ChIP-seq experiments in S2 cells (Jain et al., 2015; Figure 5—figure supplement 3). While some heterochromatic TSSs could still represent true binding, we are more confident that the broader occupancy we observe outside TSSs represents Odj and Nnk’s true localization. Moreover, we did not observe a significant effect of Nnk loss on heterochromatin-embedded gene expression at the larval L1 stage via RNA-seq analyses (Figure 4D, Figure 5—figure supplement 4). Instead, we found that some heterochromatin-embedded repetitive elements are derepressed upon Nnk loss (Figure 5—figure supplement 4). Based on these results, we hypothesize that rather than acting as a transcription regulator for specific heterochromatin-embedded genes like its paralogs Séance, Ouija Board and Molting Defective (Uryu et al., 2018), Nnk may instead repress expression from heterochromatic repeats.
Our analyses revealed Nnk is an evolutionarily young, differentially-retained gene that is essential for larval development in D. melanogaster. Nnk has also evolved under dramatic positive selection in just the 2.5 million-year divergence between D. melanogaster and D. simulans, having undergone 52 fixed non-synonymous differences in the 439 aa protein-coding region (Table 2). Given the intriguing correlation between essentiality and positive selection that we had observed in the ZAD-ZNF genes, we investigated whether positive selection of Nnk has affected its function.
We first assayed the subcellular localization of epitope-tagged Nnk orthologs from D. melanogaster and D. simulans in D. melanogaster S2 cells. We used transient transfections to introduce these genes into S2 cells and used heat shock to drive their expression. We found that both proteins similarly localize to foci within heterochromatin (Figure 6A–B; dashed line marks heterochromatin boundaries). Thus, the rapid evolution of Nnk during D. melanogaster- D. simulans divergence has not dramatically affected its gross subcellular localization.
Next, we examined the consequences of Nnk positive selection on viability in D. melanogaster. We created a D. simulans Nnk ‘rescue’ transgene (Figure 6C). This rescue transgene is similar to the D. melanogaster Nnk-rescue construct (Figure 3B), except that it contains the D. simulans Nnk coding sequence (codon-optimized to D. melanogaster) with 1 kb Nnk-flanking sequences from D. melanogaster (Supplementary file 4). We introduced this D. simulans Nnk-rescue transgene into the same attP site on D. melanogaster X chromosome as the D. melanogaster Nnk-rescue transgene via PhiC31-mediated transgenesis (Figure 6D). The attP-insertion rescue design allowed us to put the transgene in the same genetic location in the Nnk-mel and Nnk-sim rescue crosses, normalizing for variability in expression of transgenes. This allowed a near-isogenic comparison of the D. simulans and D. melanogaster Nnk transgenes’ ability to rescue inviability of Nnk-null D. melanogaster flies, despite their high level of sequence divergence. Unfortunately, low levels of endogenous Nnk expression did not allow us to assess whether the expression levels of both Nnk transgenes were equivalent.
We crossed heterozygous females either containing one or no copy of the (RNAi-resistant) Nnk-rescue transgene and the Nnk-RNAi allele to Act5C-GAL4/CyO-GFP males. We expected?to?find that progeny that did not inherit the rescue transgene would be inviable. In contrast, in the resulting progeny from rescue transgene-bearing females, we expect half of the progeny to inherit the Nnk-rescue transgene (see Materials?and?methods). We found that the Nnk-mel transgene significantly rescued Nnk knockdown compared to no-transgene controls; however, males were recovered at slightly lower levels than females (67 males: 101 females compared to expectation of 1:1 ratio; p=0.08, Fisher’s exact test). We similarly found that the D. simulans Nnk transgene can significantly rescue the lethality caused by knockdown of endogenous D. melanogaster Nnk in females (Figure 6C), although at a slightly lower level than Nnk-mel rescue. In contrast, rescue of male viability by D. simulans Nnk is extremely poor (2 males: 33 females compared to expectation of 1:1 ratio; p<0.0001, Fisher’s exact test) resulting in a severe sex-bias. Thus, Nnk-sim is much worse than Nnk-mel in rescuing male viability (67:101 versus 2:33, p<0.0001, Fisher’s exact test). Based on these findings, we infer that the D. simulans Nnk-rescue transgene specifically fails to rescue Nnk-knockout males, in the presence of the heterochromatin-rich Y chromosome. Our findings suggest that not only is Nnk a positively-selected, essential ZAD-ZNF gene, but also that its positive selection is required for optimal function in the D. melanogaster genome.
In this study, we explored the relationship between genetic innovation and essentiality in the ZAD-ZNF gene family, which encodes the most abundant class of transcription factors in insects. Due to their lineage-specific amplification, protein structure and expression patterns, ZAD-ZNF genes were previously hypothesized to be analogous to the KZNF (KRAB-Zinc Finger) transcription factor-encoding gene family found in vertebrates (Chung et al., 2002; Chung et al., 2007), many of which target transposable element sequences (TEs) inserted in the genome. Just like KZNF genes in mammals, we find strong evidence for evolutionary dynamism in Drosophila ZAD-ZNF genes. For example, we find that only 61 of 91 ZAD-ZNFs found in D. melanogaster are universally retained in most Drosophila species. Furthermore, 12 ZAD-ZNFs have evolved under positive selection during D. melanogaster- D. simulans divergence.
Despite these similarities, however, there are considerable differences between these gene families. First, a direct connection to TEs has only been revealed for one ZAD-ZNF gene in Drosophila; CG17801 has a role in regulating HetA and Blood transposable elements in the female germline (Czech et al., 2013). Second, we find that ZAD-ZNF genes that evolve under positive selection are more likely to encode essential functions in embryonic axial patterning, larval development, and meiosis (Chen et al., 2000; Harms et al., 2000; Page and Orr-Weaver, 1996; Uryu et al., 2018; Lake et al., 2011; Chen et al., 2006; Table 2). In contrast, most KZNF genes that have been shown to be essential for sterility or viability in mammals are slowly evolving (Wolf et al., 2020; Imbeault et al., 2017). Finally, unlike KZNFs (Thomas and Schneider, 2011), positive selection in ZAD-ZNFs is not primarily focused on their DNA-binding C2H2 domains but rather on the poorly-characterized linker domains that connect the ZAD and C2H2 domains (Table 2). We speculate that the linker, which is often comprised of intrinsically disordered domains, may play an important role in chromatin localization of ZAD-ZNF proteins either via direct DNA-binding or via protein–protein interactions (Brodsky et al., 2020; Erijman et al., 2020). For example, the Oddjob linker domain encodes a PxVxL motif that is important for its broad heterochromatic localization, potentially by mediating direct interaction with HP1a (Figure 5, Figure 5—figure supplement 2).
Based on these dissimilarities, we do not favor the possibility that ZAD-ZNF innovation is solely driven by arms-races with TEs. Instead, we favor the hypothesis that the recurrent adaptation of a subset of ZAD-ZNF stems from their roles in pericentromeric heterochromatin organization or regulation of gene expression. Indeed, many genes that encode crucial heterochromatin functions are often critical for viability or fertility, yet are quite variable even among closely-related species (Ross et al., 2013; Klattenhoff et al., 2009; Levine et al., 2016; Parhad et al., 2017; Vermaak and Malik, 2009). Heterochromatin is a gene-poor component of most eukaryotic genomes. Yet, its establishment and maintenance is nevertheless essential for many cellular processes including chromosome condensation and segregation, repression of TEs, and genome stability (Abe et al., 2016; Grézy et al., 2016; Levine et al., 2015; Nambiar and Smith, 2018; Okita et al., 2019; Vernì and Cenci, 2015; Liu et al., 2014; Azzaz et al., 2014; Inoue et al., 2008; Ruiz-Estévez et al., 2014; Verschure et al., 2005; Brennecke et al., 2007; Senti and Brennecke, 2010; Goriaux et al., 2014a; Goriaux et al., 2014b). Thus, the rapid evolution of genes encoding heterochromatin functions might reflect lineage-specific mechanisms to package heterochromatic DNA, or silence TEs. In this context, it is important to note the heterochromatic satellite DNA sequences themselves are among the most rapidly evolving component of Drosophila genomes (Wei et al., 2018; Chakraborty et al., 2020; Jagannathan et al., 2017).
Several characterized ZAD-ZNFs have been found to play key roles at heterochromatin. For example, the Oddjob protein co-immunoprecipitates with HP1a and broadly localizes to heterochromatin (Figure 5; Swenson et al., 2016). Similarly, Séance and Ouija Board control the expression of heterochromatin-embedded genes necessary for larval development (Uryu et al., 2018). Based on these observations, we propose that ZAD-ZNF diversification (marked by gene turnover and positive selection) is driven by the high turnover of sequences embedded within heterochromatin. Although the bulk of heterochromatin is made up of highly repetitive elements such as satellite DNAs and TEs, heterochromatin also harbors many genes that are deeply embedded within heterochromatin (Yasuhara et al., 2005; Yasuhara and Wakimoto, 2006; Eberl et al., 1993; Schulze et al., 2006; Schulze et al., 2005; Devlin et al., 1990). These genes, many of which encode essential functions (Sinclair et al., 2000), require a heterochromatic environment to ensure their correct expression and regulation (Wakimoto and Hearn, 1990). We posit that the constant turnover of flanking and embedded sequence elements such as TEs and satellite DNAs may require constant adaptation of transcription factors required for the proper expression of genes embedded in heterochromatin. ZAD-ZNF adaptation might also be necessary to protect against the inappropriate expression of heterochromatin-embedded elements especially at crucial developmental stages. Because of the high rate of evolutionary turnover of heterochromatic sequences, we hypothesize that ZAD-ZNF genes that are essential in one species could nevertheless be lost in another, either because the target loci of the ZAD-ZNF genes is lost, or because another ZAD-ZNF paralog has acquired this essential regulation function.
By testing the causal link between positive selection and essential function, we also find further support for Nnk’s essential function being related to heterochromatin biology. We find that D. simulans Nnk can significantly rescue Nnk-null inviability in females, but not in males. We speculate that this failure to rescue male viability could be due to the heterochromatin-rich Y chromosome in males. For example, if the D. melanogaster Nnk protein, but not the D. simulans Nnk protein, could appropriately repress the D. melanogaster Y chromosome, this might explain the sex-bias seen in the D. simulans Nnk-rescue progeny. The Y chromosome is itself not essential for viability or sex determination in Drosophila. However, de-repression of the heterochromatin-rich Y chromosome could nevertheless lead to detrimental consequences on larval development. This effect could be direct, leading to inappropriate expression of Y-chromosome-embedded genetic elements that block larval development. Alternatively, this effect could be indirect; de-repression of the Y chromosome could indirectly impact several other chromatin processes genome-wide e.g., due to inappropriate titration of transcription or heterochromatin factors, which could exacerbate an already hypomorphic function of the D. simulans Nnk allele (Branco et al., 1869; Francisco and Lemos, 2014; Wang et al., 2018; Piergentili, 2010; Brown et al., 2020). Finally, different functional optimality of male-specific or female-specific Nnk functions could also drive its rapid evolution (VanKuren and Long, 2018). In any of these scenarios, we hypothesize that rapid co-evolution with heterochromatic sequences might have driven the rapid evolution of Nnk and possibly other heterochromatin-interacting ZAD-ZNF proteins.
Based on our findings, we hypothesize that constant adaptation in ZAD-ZNF genes is driven by rapid alterations in heterochromatin across Drosophila and other insect species. This co-evolutionary arms-race may provide the explanation for the unexpected correlation we find between gene essentiality and innovation in the largest family of transcription factors in insect genomes.
We used the Flybase database (http://flybase.org) to identify all ZAD-containing proteins (Pfam motif PF07776) in 12 sequenced and annotated Drosophila species (Drosophila 12 Genomes Consortium, 2007). Using NCBI’s Conserved Domains search, we identified other domains found in the ZAD-containing proteins (Marchler-Bauer et al., 2017). To estimate the evolutionary age of these ZAD-containing genes, we used OrthoDB to identify orthologs of the 91 ZAD-containing proteins across Drosophila (Zdobnov et al., 2017). We used the divergence between the two most distantly-related species that still encoded and ortholog to calculate the minimum age of each gene using timetree.org (Kumar et al., 2017). If orthologs were identified in basally-branching but not later-branching species by OrthoDB, we defined this as a loss of the ZAD-ZNF gene in the lineage. The repertoire of different ZAD-ZNF genes in different Drosophila species inferred from our analysis is summarized in Supplementary file 1.
We used FlyBase gene summaries (FlyBase.org) and published studies when available to define gene essentiality. Our criteria for essentiality were broad: if there was a lethal allele reported in any assay, we counted the gene as essential for viability in D. melanogaster. Supplementary file 1 summarizes the phenotypic data (and its source) that we used to classify the ZAD-ZNF genes into either essential or non-essential categories.
For the McDonald-Kreitman test, we extracted the gene of interest from D. melanogaster population genetic datasets available through Popfly (www.popfly.com) and removed low frequency (<0.05) variants from the dataset to minimize the effects of false positives and low-frequency variants that may not have been subject to selection (Hervas et al., 2017; Fay et al., 2002). We used a manually trimmed alignment of the D. melanogaster filtered dataset and the reference D. simulans sequence for the McDonald-Kreitman test (http://mkt.uab.es/mkt/; McDonald and Kreitman, 1991; Egea et al., 2008). We also carried out a polarized McDonald-Kreitman test in which we only analyzed non-synonymous and synonymous sites that were fixed along the lineage leading to D. melanogaster, after its divergence from D. simulans, inferred using D. yakuba as an outgroup species, except for the Trem gene where we used D. erecta outgroup instead because the D. yakuba ortholog aligned poorly to D. melanogaster and D. simulans.
We used an Act-GAL4/CyO-GFP driver for ubiquitous knockdown. The RNAi lines used to specifically target Oddjob cluster genes are the VDRC KK or GD lines: Oddjob (27971), CG17803 (38869), CG17801 (29501), CG17806 (40106) and Nicknack (102311). RNAi controls used for the experiment were Cid (43856) and HP1B (26097) for the ubiquitous knockdown. Ubiquitous knockdown of Cid produced no viable progeny. We crossed five virgin females carrying Act-GAL4/CyO-GFP to 3 males of each RNAi line. We allowed the females to lay eggs for 3 days and flipped the flies into fresh vials three times. Each cross was performed in triplicate. Progeny were counted 10–15 days after each cross was set up. A minimum of 20 CyO-GFP males and 20 CyO-GFP females (i.e., control genotype) were required for us to quantify the crosses. Each cross is represented by a point on the graph and shown as a ratio of observed (number of knockdown or rescue progeny counted in each cross) over expected (the number of expected knockdown or rescue progeny from Mendelian segregation inferred by inheritance of the balancer chromosome). If the knockdown has no effect on viability or the rescue is 100% effective, the observed/expected?=?1. If the knockdown of the essential gene is 100% penetrant or the rescue is ineffective, the observed/expected?=?0. Plotting and statistical analyses were conducted using Graphpad Prism 8 software.
Since Oddjob cluster genes experienced numerous independent segmental duplications, it was not possible to determine orthologs by synteny alone. Instead, we used TBLASTN (Gertz et al., 2006) to identify candidate orthologs of Oddjob cluster genes, using the genes in D. melanogaster as queries. We used a reciprocal blast search strategy to identify potential orthologs and further investigated these candidates by making a maximum likelihood phylogenetic tree (LG substitution model in PhyML with 100 bootstrap replicates) of a manually trimmed protein alignment (constructed using Clustal Omega program [Thompson et al., 2002] in the Geneious package, www.geneious.com) (Supplementary file 2, Figure 2—figure supplement 1). We assigned orthologs based on genes that formed a monophyletic clade with the each of the D. melanogaster Oddjob cluster genes. We mapped these Oddjob cluster orthologs back to each genome assembly to determine the composition of the Oddjob locus across Drosophila species. In cases where there were other ZAD-ZNFs present, we blasted them against the D. melanogaster genome to examine if there were any orthologs present in D. melanogaster. If the top hit was not a member of the Oddjob locus, we did not include it in the tree. In the case where there are partial ZAD-ZNFs (containing just the ZAD domain or the zinc-finger domains), we performed a blast search against the D. melanogaster genome. If the top hit was a member of the Oddjob cluster, we assigned orthology by making a phylogenetic tree with all other Oddjob cluster orthologs. Our phylogenomic inferences are dependent on the accuracy and state of completion of genomes from other Drosophila species, which can vary.
We designed a D. melanogaster recoded transgene comprising a 3.3 kb fragment containing the genomic region of Nicknack plus 1 kb upstream and downstream of the start and stop codons, based on the D. melanogaster reference assembly. We recoded the sequence targeted by the VDRC RNAi KK line (103211) by making synonymous changes at each codon. The synonymous changes made the gene resistant to RNAi knockdown but did not alter the amino-acid sequence of Nnk. The recoded region spans 270 base pairs that corresponds to the region targeted by the VDRC RNAi KK line (103211). The resulting sequence was synthesized by GENEWIZ Co. Ltd. (Suzhou, China) and cloned into a plasmid we generated that contains 3xP3-DsRed attP, which produces fluorescent red eyes in the adult to mark the presence of the transgene. To generate the D. simulans transgenic allele, we codon-optimized the D. simulans Nnk coding sequence for the D. melanogaster genome using IDT’s codon-optimization tool. The resulting sequence was synthesized by GENEWIZ Co. Ltd. (Suzhou, China) and swapped for the D. melanogaster Nnk coding sequence in the plasmid described above, using the NEBuilder kit (New England Biolabs). We submitted transgenic constructs to The BestGene Inc (Chino Hills, CA) for injection into the X-chromosome attP18 line (BL 32107) using PhiC31 site-specific integration (Groth et al., 2004).
For the transgene rescue cross, we crossed five virgin female flies bearing one copy of Nnk-rescue transgene (on the X chromosome) and the RNAi allele (on the second chromosome) to 3 Act5C-GAL4/CyO-GFP males. We allowed the females to lay eggs for 3 days and flipped the cross three times. We set up each cross in triplicate and progeny were counted after 10–15 days. The description of the observed/expected calculation can be found in ‘viability studies’ section. Each cross had at least three replicates. All flies were raised at 25°C.
We placed 50–75 flies heterozygous for the Nnk-null allele (Nnk pBac null/TM3G) or for the Nnk CRISPR allele (Nnk CRISPR/TM3G) into a small embryo collection cage containing a grape-juice plate with a thin strip of yeast paste and collected embryos for 3 hr at 25°C. We transferred the larvae to fresh grape-juice plates containing yeast paste daily and scored developmental stage by mouth hook morphology. We used fluorescence to distinguish between heterozygotes (GFP-positive larvae) and homozygotes (GFP-negative larvae). For the trans-heterozygote evaluation, we crossed 30–40 virgin female Nnk CRISPR/TM3G to 10 Nnk pBac null/TM3G males. Crosses were done in triplicate and at least 100 progeny were counted per cross.
We placed 50–75 flies heterozygous for the Nnk-null allele (Nnk pBac/TM3G) into a small embryo collection cage containing a grape-juice plate with a thin strip of yeast paste and collected embryos for 3 hr at 25°C. The first time point was collected 24 hr after egg laying (AEL) and the second 48 AEL. We transferred the larvae to fresh grape-juice plates containing yeast paste daily.
Whole larvae (~30 animals at 24 hr AEL and?~20 animals at 48 hr AEL for each sample; RNA from each time point and genotype was prepared in triplicate) were ground in a 1.5 mL Eppendorf tube containing 50 μL of TRIzol reagent using a DNase, RNase and DNA free 1.5 mL pestle. 450 μL of TRIzol reagent was added after grinding. Immediately, we added 500 μL of chloroform and the tube was inverted gently 2–3 times. We removed the aqueous phase into a fresh tube containing 1 mL of 200 proof EtOH and mixed by inversion. The mixture was then bound to a Zymo-spin column according to the manufacturer's instructions (Zymo Research). We followed the DNase extraction and purification protocol outlined in the RNA Clean and Concentrator kit (Zymo Research). We eluted the RNA in 15 μL of DNase/RNase-free water and immediately placed the samples at ?80°C. We checked the quality of the samples with a 2200 Tapestation (Agilent Technologies) and selected samples that had an RNA integrity number?>9.0 for library preparation. Library construction and Illumina 150 bp paired-end RNA-sequencing were conducted at Novogene Bioinformatics Technology Co., Ltd (Beijing, China).
We used Kallisto (Bray et al., 2016) to quantify abundances of the D. melanogaster reference transcriptome (refMrna.fa for dm6, obtained from UCSC Genome Browser Oct 16th, 2018, which contains 34,114 transcripts). For each transcript, we acquired the gene name using R (org.Dm.eg.db). Kallisto counts were read into R using the tximport package (Soneson et al., 2015) and using the summarizeToGene function, we summarized alternative splice-form counts into a single count per gene. We used DESeq2 to identify differentially expressed genes with adjusted p-value<=0.05 and absolute log2(fold change)>=1 (Anders and Huber, 2010), comparing Nnk-null larvae to controls for each timepoint separately. Before performing each comparison, we excluded genes with low expression (<100 counts total across all samples): this filtering yields 10,574 ‘expressed’ genes in 24 hr AEL larvae, and 10,758 ‘expressed’ genes in 48 hr AEL larvae. The RNA-seq data has been deposited in the SRA database with the accession number PRJNA643855.
We performed Gene Ontology (GO) enrichment analysis on each of four gene lists: 116 over-expressed genes and 133 under-expressed genes in 24 hr AEL larvae, and 1301 over-expressed genes and 1726 under-expressed genes in 48 hr AEL larvae. We used the Bioconductor GOstats package to perform conditional hypergeometric enrichment tests for each of three ontologies (biological process, molecular function and cellular component) (Huber et al., 2015; Falcon and Gentleman, 2007). For the ‘universe’ of all genes examined (i.e., background) we used the corresponding list of ‘expressed’ genes at each developmental timepoint. We report only over-represented categories with p<0.001, and use annotations found in the org.Dm.eg.db Bioconductor package.
To estimate expression levels aggregated across all instances of each repeat type, we took an approach similar to that of Day et al., 2010. Here, rather than mapping reads to the typical reference genome assembly, we constructed a ‘repeat assembly’, where we used RepeatMasker annotations for dm6 (obtained from UCSC) to extract and concatenate all instances of each repeat type, adding 75 bp flanking sequences (half the length of each read) and inserting 150 N bases between each instance. The repeat assembly therefore consists of a single ‘chromosome’ for each repeat type. We then used BWA-mem to map sequences as single reads (not paired-ends) to the repeat assembly (Li, 2013), and samtools to count reads mapping to each repeat type (‘chromosome’). We combined counts for pairs of simple repeats that represent the reverse-complement of one another. We used DESeq2 to perform statistical analyses on raw counts, using the total number of sequenced fragments as size factors (Love et al., 2014). We normalized counts by dividing each by the number of reads sequenced for that sample in millions.
To evaluate the ability of dietary sterols to rescue Nnk-null phenotype, we carried out sterol supplementation as previously outlined (Uryu et al., 2018). Briefly, we mixed together 20 mg of dry yeast in 38?μL of water. To this yeast paste, we added either 2?μL of EtOH (negative control) or 2?μL of EtOH plus cholesterol (Sigma), 7-dehydrocholestrol (Sigma) or 20-hydroxyecdysone. Stocks used for these experiments were balanced over GFP-balancers (Npc1a57/CyO-GFP and Nnk-null/TM3G). The control used for these experiments was the w1118 stock. Eggs were laid on yeasted grape-juice plates for 3 hr at 25°C. 24 hr AEL, GFP-negative larvae were transferred onto grape plates containing fresh yeast paste at 25°C. For 72 hr AEL, larvae were transferred to fresh grape-juice plates containing yeast paste daily and scored for their developmental stage based on the morphology of their mouth hooks. Npc1a57/CyO-GFP stocks were a kind gift from Leo Pallanck (Fluegel et al., 2006).
Oddjob cluster genes were amplified from genomic DNA from D. melanogaster and D. simulans and directionally cloned into pENTER/D-TOPO (ThermoFisher) according to the manufacturer’s instructions. We verified that clones had the appropriate insertions by sequencing. We used LR clonase II (ThermoFisher) to get each Oddjob cluster gene into the Drosophila Gateway Vector destination vector to express the gene of interest with an N-terminal Venus tag under the control of the D. melanogaster Hsp70 promoter (pHVW). For Rhino localization in S2 cells, the Rhino gene was cloned into the Drosophila Gateway Vector destination vector to enable expression of Rhino with an N-terminal 3xFlag tag under the control of the D. melanogaster Hsp70 promoter (pHFW). Nnk and Odj were also cloned into pCopia vectors containing N-terminal and C-terminal mCherry epitope tags, respectively, and HP1a was previously tagged with mCerulean at the N-terminus.
Schneider 2 cells were obtained from the Drosophila Genomics Resource Center (Bloomington, IN, USA) and grown at 25°C in M3+BPYE+10%FCS. For the transfections, one million cells were seeded, and one day later 2 micrograms of plasmid DNA was transfected into cells using Xtremegene HP transfection reagent (Roche) according to the manufacturer’s specifications. For the Rhino localization experiment, cells were co-transfected with 1?μg of GFP-tagged ZAD-ZNF vector and 1?μg of Flag-tagged Rhino vector. Cells were allowed to recover for 24 hr post transfection, heat shocked for 1 hr at 37°C and recovered for 3 hr at 25°C prior to fixation. For the pCopia vectors, transient transfections were conducted on S2 cells using the TransIT‐2020 reagent (Mirus), and live imaging was performed 72 hr later, using an Applied Precision Deltavision microscope and analyzed using SoftworX software.
For the heatshock vector transfections, cells were transferred to coverslips for 30–45 min prior to starting the immunohistochemistry protocol. 0.5% sodium citrate hypotonic solution was added to the coverslip for 10 min to swell cells, which was then spun at 1900 rpm for 1 min in a Cytospin III to remove the cytoplasm from the cells. The sodium citrate was immediately removed and cells were subsequently fixed. For fixation, 4% PFA + PBST (PBS + 0.1% Triton) was added to the cells for 10 min. Coverslips were then washed in PBST and blocked for 30 min in PBST + 3% BSA. Cells were incubated with primary antibody overnight at 4°C in a humid chamber. For immunolocalization, the following dilutions of primary antibodies were used: GFP (Abcam AB13970) 1:1000, H3K9me3 (Abcam AB8898) 1:500 M2 FLAG (Sigma-Aldrich F4042) 1:1000. After washing with PBST three times for 10 min per wash, the following fluorescent secondary antibodies were used at 1:1000 dilution: goat anti-chicken (Invitrogen Alexa Fluor 488, A-11039), goat anti-rabbit (Invitrogen Alexa Fluor 568, A-11011) and goat anti-mouse (Invitrogen Alexa Fluor 568, A-11031). The cells were incubated with 1x DAPI in the final wash and mounted in SlowFade Gold Mounting Medium (ThermoFisher). We imaged cells on a Leica TCS SP5 II confocal microscope with LASAF software and images were processed using ImageJ and were representative of the cell population.
S2 cells transfected with either mCherry-Nnk or Odj-mCherry for 72 hr were fixed with 1% paraformaldehyde for 10 min and sheared with Bioruptor sonicator (Diagenode) to obtain chromatin. Chromatin fragments were confirmed to contain DNA in the 200–500 bp size range using a Bioanalyzer (Agilent). Each immunoprecipitation was performed on chromatin from 2 × 107 S2 cells by overnight incubation with Protein-G Dynabeads and 5?μg of anti-mCherry antibody (Novus). Library construction from immunoprecipitated DNA was conducted using TruSeq sample preparation kits (Illumina). 150 bp paired-end sequences were generated by the Vincent J. Coates Genomics Sequencing Laboratory at UC Berkeley.
We used BWA-mem to map paired reads to version 6 of the D. melanogaster genome assembly (dm6) from which we had removed unplaced scaffolds (Li, 2013). We used deepTools’ bamCompare (with the ‘--binsize 1 --extendReads’ options) to obtain log2 ratios of fragment coverage for matched ChIP and input samples, and visualized those ratios in IGV (Robinson et al., 2011). We further visualized ChIP-seq signal around TSSs using deepTools’ computeMatrix and plotHeatmap tools, using TSS annotations obtained from the TxDb.Dmelanogaster.UCSC.dm6.ensGene BioConductor package, taking the most upstream TSS for genes with alternative start sites. We split TSS annotations according to whether they were within ‘TSS-proximal’ (active) chromatin according to modENCODE’s nine state annotation for S2 cells, obtained from http://intermine.modencode.org and converted from dm3 to dm6 coordinates using UCSC’s liftOver tool (Hinrichs et al., 2006). We further split the ‘active’ TSSs according to whether they are within cytogenetic heterochromatin, using coordinates from Supplementary file 2 of Hoskins et al., 2015 as well as the whole of chromosomes four and Y. The ChIP-seq data has been deposited in the SRA database with the accession number PRJNA644950.
The genetics of a small autosomal region of Drosophila melanogaster containing the structural gene for alcohol dehydrogenase VII Characterization of the region around the snail and cactus lociGenetics 94:679–694.
Human heterochromatin protein 1α promotes nucleosome associations that drive chromatin condensationJournal of Biological Chemistry 289:6850–6861.https://doi.org/10.1074/jbc.M113.512137
GFZF, a Glutathione S -Transferase Protein Implicated in Cell Cycle Regulation and Hybrid Inviability, Is a Transcriptional CoactivatorMolecular and Cellular Biology 38:e00476-17.https://doi.org/10.1128/MCB.00476-17
Sex-specific adaptation and genomic responses to Y chromosome presence in female reproductive and neural tissuesProceedings. Biological Sciences 284:20172062.https://doi.org/10.1098/rspb.2017.2062
The Drosophila Y chromosome affects heterochromatin integrity Genome-WideMolecular Biology and Evolution 20:msaa082.https://doi.org/10.1093/molbev/msaa082
Lineage-specific expansion of the zinc finger associated domain ZADMolecular Biology and Evolution 24:1934–1943.https://doi.org/10.1093/molbev/msm121
The organization and expression of the light gene, a heterochromatic gene of Drosophila melanogasterGenetics 40:129–140.
Evolutionary turnover of kinetochore proteins: a ship of theseus?Trends in Cell Biology 26:498–510.https://doi.org/10.1016/j.tcb.2016.01.005
Control of genetic stability by a new heterochromatin compaction pathway involving the Tip60 histone acetyltransferaseMolecular Biology of the Cell 27:599–607.https://doi.org/10.1091/mbc.E15-05-0316
PopFly: the Drosophila population genomics browserBioinformatics 33:2779–2780.https://doi.org/10.1093/bioinformatics/btx301
Orchestrating high-throughput genomic analysis with bioconductorNature Methods 12:115–121.https://doi.org/10.1038/nmeth.3252
Perturbation of HP1 localization and chromatin binding ability causes defects in sister-chromatid cohesionMutation Research/Genetic Toxicology and Environmental Mutagenesis 657:48–55.https://doi.org/10.1016/j.mrgentox.2008.08.010
PrDOS: prediction of disordered protein regions from amino acid sequenceNucleic Acids Research 35:W460–W464.https://doi.org/10.1093/nar/gkm363
Comparative analysis of satellite DNA in the Drosophila melanogaster Species ComplexG3: Genes, Genomes, Genetics 7:693–704.https://doi.org/10.1534/g3.116.035352
Active promoters give rise to false positive ‘Phantom Peaks’ in ChIP-seq experimentsNucleic Acids Research 43:6959–6968.https://doi.org/10.1093/nar/gkv637
Mapping second chromosome mutations to defined genomic regions in Drosophila melanogasterG3: Genes, Genomes, Genetics 8:9–16.https://doi.org/10.1534/g3.117.300289
New genes often acquire male-specific functions but rarely become essential in DrosophilaGenes & Development 31:1841–1846.https://doi.org/10.1101/gad.303131.117
TimeTree: a resource for timelines, timetrees, and divergence timesMolecular Biology and Evolution 34:1812–1819.https://doi.org/10.1093/molbev/msx116
A thousand fly genomes: an expanded Drosophila genome nexusMolecular Biology and Evolution 3313:3308.https://doi.org/10.1093/molbev/msw195
Rapid evolution of gained essential developmental functions of a young gene via interactions with other essential genesMolecular Biology and Evolution 36:2212–2226.https://doi.org/10.1093/molbev/msz137
Recurrent gene duplication diversifies genome defense repertoire in DrosophilaMolecular Biology and Evolution 33:1641–1653.https://doi.org/10.1093/molbev/msw053
New gene evolution: little did we knowAnnual Review of Genetics 47:307–333.https://doi.org/10.1146/annurev-genet-111212-133301
Small Drosophila zinc finger C2H2 protein with an N-terminal zinc finger-associated domain demonstrates the architecture functionsBiochimica Et Biophysica Acta. Gene Regulatory Mechanisms 2020:194446.https://doi.org/10.1016/j.bbagrm.2019.194446
CDD/SPARCLE: functional classification of proteins via subfamily domain architecturesNucleic Acids Research 45:D200–D203.https://doi.org/10.1093/nar/gkw1129
Piragua encodes a zinc finger protein required for development in DrosophilaMechanisms of Development 144:171–181.https://doi.org/10.1016/j.mod.2016.12.003
Molting defective is required for ecdysone biosynthesisDevelopmental Biology 280:362–372.https://doi.org/10.1016/j.ydbio.2005.01.023
Multiple roles of the Y chromosome in the biology of Drosophila melanogasterThe Scientific World JOURNAL 10:1749–1767.https://doi.org/10.1100/tsw.2010.168
The distribution of and complementation relationships between spontaneous X-linked recessive lethal mutations recovered from crossing long-term laboratory stocks of Drosophila melanogasterMutation Research/Fundamental and Molecular Mechanisms of Mutagenesis 163:115–144.https://doi.org/10.1016/0027-5107(86)90042-4
The piRNA pathway: a fly's perspective on the guardian of the genomeTrends in Genetics 26:499–509.https://doi.org/10.1016/j.tig.2010.08.007
Coevolution of retroelements and tandem zinc finger genesGenome Research 21:1800–1812.https://doi.org/10.1101/gr.121749.111
Multiple sequence alignment using ClustalW and ClustalXCurrent Protocols in Bioinformatics 2:3.https://doi.org/10.1002/0471250953.bi0203s00
Gene duplicates resolving sexual conflict rapidly evolved essential gametogenesis functionsNature Ecology & Evolution 2:705–712.https://doi.org/10.1038/s41559-018-0471-0
Multiple roles for heterochromatin protein 1 genes in DrosophilaAnnual Review of Genetics 43:467–492.https://doi.org/10.1146/annurev-genet-102108-134802
In vivo HP1 targeting causes large-scale chromatin condensation and enhanced histone lysine methylationMolecular and Cellular Biology 25:4552–4564.https://doi.org/10.1128/MCB.25.11.4552-4564.2005
A Drosophila RNAi library modulates hippo pathway-dependent tissue growthNature Communications 7:10368.https://doi.org/10.1038/ncomms10368
Vitellogenesis and fertility in Drosophila females with low ecdysteroid titres; the L(3)3DTS mutationJournal of Insect Physiology 33:137–142.https://doi.org/10.1016/0022-1910(87)90139-9
Variable rates of simple satellite gains across the Drosophila phylogenyMolecular Biology and Evolution 35:925–941.https://doi.org/10.1093/molbev/msy005
Oxymoron no more: the expanding world of heterochromatic genesTrends in Genetics 22:330–338.https://doi.org/10.1016/j.tig.2006.04.008
Patricia J WittkoppSenior and Reviewing Editor; University of Michigan, United States
In the interests of transparency, eLife publishes the most substantive revision requests and the accompanying author responses.
[Editors' note: this paper was reviewed by Review Commons.]
This work is significant because it challenges the widely held belief that essential genes are highly conserved over evolutionary time by demonstrating that a rapidly evolving family of transcription factors has essential functions in Drosophila. It also beautifully demonstrates how an essential gene in D. melanogaster that interacts with heterochromatin has functionally diverged from its ortholog in D. simulans. Overall, the study is impressive for its breadth and depth, and because it provides important new insights into how new gene functions evolve.https://doi.org/10.7554/eLife.63368.sa1
We would like to thank the three reviewers for their time and thoughtful comments on the manuscript. We found their constructive feedback extremely helpful and are grateful for the peer review process in clarifying and strengthening our study. All three reviewers were highly enthusiastic but requested a few additional experiments and clarifications. We have added the following analyses in our revision in response to their comments:
1) Two reviewers commented that polarized McDonald-Kreitman tests would provide better discrimination about whether the positive selection we saw previously in the unpolarized D. melanogaster-D. simulans comparison could be attributed to D. melanogaster alone. We have now added this analysis to our revision. Our revisiting of this analysis also made us aware of an error we had overlooked, which we have since corrected. Correcting this error has made it an even stronger correlation between genes that are evolving under positive selection and those that are required for essentiality. We corrected and revised the tables in our revision.
2) As the reviewers requested, we have also included a new Supplementary file 1, which lists the retention of all D. melanogaster ZAD-ZNF genes across Drosophila species, and the assays done to measure their phenotype.
3) Upon one reviewer’s request, we have changed the presentation of our crosses to show observed/expected ratios rather than numbers of adult flies, which are more accurate and clearer to interpret.
Reviewer #1 (Evidence, reproducibility and clarity):
Kathinathan et al. is an important study with two-fold significance in general. By a thorough search of the ZAD-ZNF genes with their functional importance assayed by knockdown and knockout, an adequately big dataset in Drosophila was created for testing how functional essentiality evolved with the ages of genes. The reported data convincingly show that the ZAD-ZNF genes quickly evolved essential functions in viability and the positive selection is related to essentiality. Further, in a thorough experimental study, it is shown that the rapid heterochromatin-derived sequence evolution is a force to lead to dynamic evolution of the ZAD-ZNF genes with their functional importance. This provides a nice mechanistic model with strong evidence for understanding the evolution of essential functions. Nicknack and Odj are beautiful examples, with detailed analyses that reveal several unexpected insights for its role in the early development of Drosophila.
This work, overall adding significantly novel data and new concepts to the literature on the evolution of gene functions and the evolutionary development, should be published soon in any journal with broad audiences because of the general interests of the topics and findings.
There are several minor issues that can be fixed easily, with one more consideration to interpret the rescue experiment using the D. simulans orthologue.
We thank the reviewer for the positive and constructive feedback and have incorporated their comments to our revision as detailed below.
Introduction: Literature 4, a more relevant one should be cited: Lee et al., 2019. Mol Biol Evol 36: 2212-2226.
We appreciate the reviewer’s suggestion; this citation has been updated.
The recue experiment by the D. simulans orthologue gives the most intriguing result that may need more consideration of interpretation. The two hypotheses of direct and indirect effects are based on the heterochromatin of Y. There can be a third hypothesis based on the sequence difference of the coding sequences between the two species. The MK test revealed a highly excess of Dn (52, Table 2), which indeed is unusually high, given the short divergence time and size of the proteins. The male cannot be rescued while female can, also suggesting a hypothesis of rapid evolution of male functions that lead to the divergence between the two species. The further relevance of this hypothesis is that sexual selection may be a driving force, as was seen from a rate even higher than the rate driven by the adaptive evolution that helps species (both males and females) to adapt to changed environments (for example, VanKuren and Long, 2018). This is not contradictory to the two hypotheses already proposed (Pahe 17-18). After all, the Y chromosome is logically a ground where sexual selection is operating.
We are grateful to the reviewer for making this suggestion and have added this possibility and citation in our Discussion to explain the lack of male rescue.
The current analysis only counted the divergence between D. melanogaster and D. simulans, so it is unknown if the positive selection is acting on D. melanogaster or D. simulans or both. A simple counting by comparing to outgroup, D. yakuba and other related species, using parsimonious principle, the 52 substitutions can be assigned to the D. melanogaster lineage or the D. simulans lineage and then do MK test against the polymorphisms in the two species separately. This testing of the symmetricity of substitutions may help understand the substitution rates in the two species and better understand the rescue experiments.
We agree with this comment and have now added the polarized MK test (Table 2).
The discussion about the heterochromatin emphasizes that the heterochromatin functions are "quite variable even among close-related species". Does this mean the genomic region covered by heterochromatins also evolved quickly in the time scale of close-related species?
We now add a sentence in response to the query that “Heterochromatic DNA sequences themselves are among the most rapidly evolving component of Drosophila genomes” along with appropriate citations.
Materials and methods: The designed techniques in this section is tricky and cleverly designed, but the description is too concise for outsiders. For example, what was the technical purpose to make the synonymous changes? how many codons were involved? A few more words can help the audiences who are not familiar with the Drosophila knockdown and rescue techniques.
We thank the reviewer for this comment and suggestion. We have added additional details to the Materials and methods and Results section to more clearly explain how the experiment was performed and interpreted.
Reviewer #1 (Significance):
This is a breakthrough in the area of evolution of gene functions.
Reviewer #2 (Evidence, reproducibility and clarity):
Authors study the evolutionary dynamics of ZAD-ZNF gene family of transcription factors in Drosophila. These are proteins defined by having both a ZAD (protein–protein interaction) domain and a ZNF (DNA-binding) domain. These is a family of transcription factors that has expanded in insects including flies and there are 91 proteins of this kind in Drosophila. However, although many are conserved and present in every fly genome, some have experienced recent duplications or losses. Some are evolving fast and others not and authors study the correlation between essentiality, mode of evolution and function. The approach is incisive because we still need to understand how often new genes become essential if at all and for what functions. I can also see how in particular some chromatin functions can be essential but evolve fast at the same time. I just would like authors to provide more details in some aspects mainly of their phylogenomics analyses and general functionality inferences.
We thank the reviewer for the excellent summary of the manuscript and kind appraisal of the work.
1) There has been a debate about how well knockdowns recapitulate knockouts phenotype. Authors themselves mentioned some shortcomings of some knockdown approaches. Sometimes knockouts are not consistent with knockdowns. Do authors give the same weight to the knockdown and knockout data when they define essentiality? Do conclusions on essentiality hold even with smaller numbers when only knockout phenotypes are analyzed? Authors might want to provide a supplementary table containing these inferences?
We thank the reviewer for their question. We gave the knockdown and knockout data equal weight in our analysis. Indeed, our detailed analysis of the Nnk gene shows that the initial skepticism about the Nnk knockdown phenotype was unwarranted. Although it is difficult to weigh the relative strengths of available phenotypic data because of the published results, we point out that our results are completely consistent with previously published studies, and we have rigorously addressed the one discrepancy on Nnk essentiality. In response to the reviewer’s suggestion, we have now added a new supplementary file 1 that details the type of assay utilized to determine the essentiality of the ZAD-ZNFs.
2) Authors need to explain how they infer the age of the gene using OrthoDB and how they made sure of the losses in the genes that they classify as not universally conserved in the Drosophila species. Would it be possible to explain or provide a supplementary table containing the details?
We agree and apologize that we were not clear before. We have now added that gene age was determined by identifying the most divergent species that still encodes an ortholog of each ZAD-ZNF gene and using previous estimates of the time of divergence between these species. If orthologs were identified in basally-branching but not later-branching species by OrthoDB, we defined this as a loss of the ZAD-ZNF in that specific lineage.
3) Are the linker protein regions easy to align or do they have indels? Regions that do not align well could lead to the spurious detection of positive selection (Genome Res. 2011. 21(6):863-74). Have authors looked at that?
The linker protein regions between D. melanogaster and D. simulans are generally easy to align. For example, for the nucleotide sequence of Nnk, the linker region was 88.8% identical and had one 21 bp region that is unique to D. simulans; this was removed for the analysis. Linker regions get progressively more difficult to align at larger divergences, however, which is why we focused on pairwise analyses to infer positive selection. We have now appended a related manuscript file showing the high nucleotide identity between D. simulans and D. melanogaster Nnk as an example. We would be happy to append this as an additional supplementary file if the editor desires.
4) How are authors sure that duplicates are not alleles in heterozygous regions or collapsed in some species reference genome? See PLoS Comput Biol 16(7): e1008104 as an example of these effects in genome assemblies. What could be the consequence of those effects? Could this effect change any of the authors conclusions?
The reviewer is correct that our genomic inferences are sensitive to the quality of the sequencing and assembly of contigs. We have now added a sentence explaining this caveat.
5) How do authors envision that this kind of essential gene can be lost in some species?
We envision that this kind of essential gene can be lost in some species either because target(s) of the ZAD-ZNF in question were lost (e.g., regulatory sequences for heterochromatin-embedded genes) or because a paralogous ZAD-ZNF took over essential function. We have added this hypothesis in our revised Discussion.
1) Tables 1 and 2 do not have titles or footnotes. What do authors mean by * or the superscript numbers in the last column of Table 2?
The reviewer may have missed the title and footnotes we had originally attached along with the figure legends. We have now updated these.
Reviewer #2 (Significance):
This work provides a conceptual advance. As I mentioned above, the approach is incisive because we still need to understand how often new genes become essential if at all and for what functions. It is revealing that those functions happen to be sometimes functions that evolve fast at the same time.
Reviewer #3 (Evidence, reproducibility and clarity):
This study analyzes the evolutionary patterns of ZAD-ZNF DNA-binding proteins in Drosophila and performs extensive functional characterization of two members of a particular subfamily, Nicknack and Oddjob. The authors first examine the conservation across Drosophila of the ZAD-ZNF genes found in D. melanogaster and use publicly available phenotypic data to ask whether the genes that show more rapid turnover across the genus or that are evolving under positive selection between D. melanogaster and D. simulans are more likely to have essential functions. Finding marginal evidence of the latter, the authors next focus on a particular cluster of tandemly duplicated ZAD-ZNF genes, which contains both Nnk and Odj. They document extensive gene gain and loss in this cluster across a wider set of species and use whole-organism RNAi to test each gene for effects on viability in D. melanogaster. Finding that knockdown of either Nnk or Odj results in lethality, they responsibly confirm these results using complementary forms of genetic ablation (CRISPR and transposon insertion). They then perform a series of functional experiments that show that Nnk is required for the transition to the L2 larval stage; that Odj and Nnk proteins both localize to heterochromatin, but with different specificities; and that the sequence divergence between D. melanogaster and D. simulans Nnk may specifically affect a part of the protein's function that affects male viability.
Overall, this is a strong study that has well-supported conclusions. The functional analyses of Odj and Nnk clear up some ambiguities from prior papers but also break meaningful new ground, and the overall work is a nice example overall of combining evolutionary genomic and functional genetic analyses. That the authors went to the trouble of rescuing RNAi with re-coding and doing a cross-species complementation experiment is impressive and shows why this research group continues to lead the field of functional molecular evolution. The genomic documentation of ZAD-ZNF genes will be useful for other researchers studying this important gene family, as will the careful analysis of the particular cluster housing Nnk. The claims about the relative functional importance of conserved vs. fast-evolving family members are not as strong because they rely on phenotypic data in FlyBase that is, by necessity, more variable in quality, and because the genome annotations of other Drosophila species are of variable quality. Still, these analyses remain reasonable to include.
We thank the reviewer for the excellent summary of our manuscript and detail the responses to their constructive feedback below.
The story is acceptably complete and comprehensive as is, so I do not feel that any more experiments are necessary. However, I think there is an issue in data analysis that should be addressed. In Figures 2B, 3B and 6D, the authors report their RNAi and genetic rescue results by depicting the raw number of flies that emerge from different replicates of fertility assays. The problem in each assay is that the expected number of flies with the particular genotype of interest depends on the overall productivity of each vial's cross, which can be affected by environmental variables (e.g., food quality, presence of bacterial contamination, health of the female parent, death of any parents during the length of the assay) that vary from vial to vial. Thus, it would be more appropriate to report these data as proportions that take into account the total number of progeny that emerged from the vial and the fraction of progeny expected from Mendelian segregation patterns to be in the target class, assuming no ill effect. Crosses that show the expected number of progeny in the genetic class of interest would have value = 1, while crosses that showed complete non-viability of the genetic class of interest would have a value of 0. I don't think this will affect the major conclusions of any of the experiments in question, and thus of the study, but it would be more accurate and interpretable.
We agree with the reviewer. In all cases, we have a built-in control (siblings that inherited the Balancer chromosome instead of the Gal4 driver or knockdown transgene) that provide the expected number of progeny assuming no knockdown. As suggested by the reviewer, we now report observed/expected ratios instead of raw numbers of flies. This is the most accurate way for us to correct for any idiosyncratic differences between conditions, and the data is more interpretable.
Here are some more minor comments and suggestions:
Some of the genomic data described for all members of the ZAF-ZNF family do not seem to be provided. For example, there is no supplemental table listing the species distribution of the 91 family members present in melanogaster, nor do the authors document which genes are the 28 for which experimental data show lethality or sterility. It's true that all of these data can be found on FlyBase, but it would be more convenient for other groups wishing to study other family members if the authors provided their findings here.
This suggestion was also made by reviewer #1. We now include this information in a new Supplementary file 1 in our revision.
Could the authors include a comment on the strength of the evidence for the essentiality of the 28 genes described in the Results? For example, the manuscript already notes the problem with some KK RNAi lines used by Chen et al., 2010, to study new gene viability, so if the bulk of the data that contribute to this count comes from studies like these, it would diminish the present authors' conclusions.
Reviewer #2 also raised this concern. We now provide details about the evidence for the viability phenotypes that we observed in a new Supplementary file 1. In the vast majority of cases, these data were not dependent on the KK lines (as the fly community has recognized the potential for artifact). Other knockdown sets do not have a similar systemic problem. We note especially that we do not rely on KK lines for any of the positively selected, essential ZAD-ZNF genes.
A couple of sentences that describe evolutionary findings could be worded more carefully. At the end of the first paragraph the authors "raise the possibility that rapid evolution of some ZAD-ZNF genes might be critical for organismal viability." I suppose the literal interpretation of this sentence is possible: if these genes hadn't changed, the ancestral organisms could have died (e.g., due to TE or repetitive region dysregulation). But, it could also be that the genes that eventually underwent positive selection already had essential functions, those functions were refined by the selection, and now modern experiments on them (which probably knocked out/down the whole gene, rather than altering the positively selected sites) have simply revealed their essential functions.
The reviewer is bringing up an important point that we have not explicitly tested whether positive selection per se is required for essentiality. However, we do show that positively selected genes are more likely to be essential. This is contrary to prevailing dogma, which argues that essential functions are more likely to be highly conserved rather than variable across species. We believe that highly unexpected finding justifies our hypothesis that rapid evolution of some ZAD-ZNF genes might be critical for organismal viability. Under the model the reviewer proposes, a previously essential gene would still have to be able to carry out its essential functions in spite of the positive selection it underwent.
Furthermore, our D. simulans Nnk swap experiment directly assesses the impact of positive selection. The failure of D. simulans Nnk to rescue male viability demonstrates the essentiality of the recent adaptive changes for proper Nnk function.
Likewise, the authors write that "62 of 91 ZAD-ZNFs found in D. melanogaster are universally retained." But didn't they find that only 73 of those genes were inferred to be present at the base of the Drosophila phylogeny (and thus were eligible for "universal retention")?
We would like to clarify this concern. While 73 of the genes were inferred to be present at the base of the Drosophila phylogeny, 11 of these genes were inferred to be lost in some species and therefore did not meet our criteria for “universal retention.” We make this clearer in our revision to avoid ambiguity.
There is a small discrepancy in the estimated time of divergence between D. melanogaster and D. pseudoobscura, listed as 35 mya but depicted as ~29 mya in Figure 1, perhaps due to different sources being used for each? The exact time is not critical for the story here, but better to be consistent (or to give a range in the text?).
We thank the reviewer for pointing out this inconsistency. We now use 30 mya as the estimated time of divergence based on the 12 Drosophila genomes paper.
One small weakness is that the authors did not appear to confirm that the genes of the Odj cluster in Figure 1 that did not show viability effects (CG17801 and CG17803) were being effectively knocked down, since a lack of knockdown would lead to the same result.
This is a fair point. However, we cannot address this concern due to the very low levels of expression of these genes and the absence of any tools (antibodies, transgenes) for this.
Also, it is a bit surprising that Nnk expression (and the expression of its transgenes) could not be detected by RT-PCR. I agree with the authors that the other phenotypic data supports appropriate expression of the rescue transgene, but I think they should at least discuss the possibility that differential expression levels of the D. simulans transgene in females vs. males could account for the sex-specific rescue of viability in Figure 6D.
The reviewer raises the important caveat that sex-specific rescue of viability is due to expression levels of the D. simulans transgene. Although this is formally possible, we think this is unlikely because the regulatory regions are shared between the Nnk-mel and Nnk-sim transgenes (and are derived from D. melanogaster); only the protein-coding regions vary. Because the Nnk-mel transgene rescues both males and females, we infer that this regulatory region does not operate in a sex-specific fashion. Thus, if there was a sex-specific expression effect, it would have to be entirely mediated by the protein-coding region.
While I largely trust the Malik Lab to know heterochromatin localization patterns when they see them, Figure 6B does not have any counterstaining to show that the localization of the Nnk transgenes is actually to heterochromatin.
We appreciate the reviewer’s trust! However, to aid other readers, we highlighted the chromocenter, identified using the H3K9me staining, in figure 6B in the panels.
In Figure 3—source data 1, please clarify whether the observed/expected data shown are for one representative replicate experiment, or are combined across all 3 or 9 replicates.
We have updated the legend to clarify that the counts are combined across all replicates.
In Figure 6D, the figure itself shows male progeny as unshaded circles and female progeny as filled circles, while the legend refers to "closed circles" representing male progeny and "open circles" representing female progeny. Please check for consistency and clarity.
Thank you for politely bringing up this inconsistency. The figure and legend have been updated.
Reviewer #3 (Significance): Overall, I am impressed by this study's analytical breadth and depth of functional characterization. I think these findings will appeal to a broad audience, for example: molecular evolutionary biologists interested in the effects of adaptive evolution or gene duplication; biologists who study heterochromatin and its regulation; and, researchers who study zinc finger protein function and evolution across taxa. While previous work has catalogued ZAD-ZNF genes in insects and Drosophila, this study adds the important lens of adaptive evolution to this work. Likewise, while Nnk had undergone cursory, high-throughput functional testing in two prior papers (Chen et al., 2010 and Kondo et al., 2017), this study performs a much more thorough characterization and adds valuable functional comparisons to its paralog, Odj, and between the melanogaster and simulans orthologs.
(For reference, my expertise lies in Drosophila genetic analysis, molecular evolution, and the function and evolution of novel genes. I have passing familiarity with the literature on heterochromatin and techniques like RNAseq and ChIP-seq and think those parts of this paper look reasonable, but their in-depth technical review is best left to those with more expertise.)
We appreciate the reviewer’s rigorous comments to clarify and improve the paper.
I appreciate point #1 raised by reviewer #2; it is a better articulated version of my second minor bullet point. Since both of us hit on this concern, I think it is a worthwhile issue for the authors to address. I also agree with their point #3, regarding potential indel variability in the linker region. I had thought of this as well, but figured that the divergence between mel and sim (the comparison most relevant for the authors' analysis of positive selection) would likely be minimal. However, I agree with reviewer #2 that it makes sense to explicitly ask the authors about this.
We appreciate this concern and have addressed it in response to reviewer #2. Briefly, linker variability is not high enough to affect the pairwise MK analyses but could affect the phylogenomic analyses, which therefore does not rely on linker alignments.
I also like reviewer #1's suggestion to attempt to polarize the MK test to determine if the adaptive evolution occurred on the branch leading to D. melanogaster or D. simulans. If data exist already for D. yakuba to do this analysis, I agree that it would add nicely to the study.
We agree and have added this analysis in our Revision.https://doi.org/10.7554/eLife.63368.sa2
- Bhavatharini Kasinathan
- Bhavatharini Kasinathan
- Bhavatharini Kasinathan
- Gary H Karpen
- Harmit S Malik
- Harmit S Malik
The funders had no role in study design, data collection and interpretation, or the decision to submit the work for publication.
We thank past and present members of the Malik lab for their comments and discussion. In addition, we are especially grateful to Celeste Berg and Barbara Wakimoto for their helpful comments during this project, and to Leo Pallanck for helpful discussions and generous gift of fly strains. We also thank Eric Lai for sharing the CRISPR-null Nnk fly line. We thank Ching-Ho Chang, Michelle Hays, Pravrutha Raman, and Aida de la Cruz for their comments on the manuscript. This work was supported by the following grants: NIH training grant T32 CA009657, NIH F30 CA225077 fellowship, and Julie Tall Achievement Rewards for College Scientists (ARCS) endowment from the Seattle Chapter of the ARCS Foundation (to BK); NIH R01 GM117420 (to GHK), NIH R01 GM074108 and HHMI Investigator grant (to HSM). The funders had no role in study design, data collection and analysis, decision to publish, or preparation of the manuscript. The authors declare that they have no conflicts of interest.
- Patricia J Wittkopp, University of Michigan, United States
- Received: September 23, 2020
- Accepted: October 12, 2020
- Version of Record published: November 10, 2020 (version 1)
? 2020, Kasinathan et al.
This article is distributed under the terms of the Creative Commons Attribution License, which permits unrestricted use and redistribution provided that the original author and source are credited.
Downloads (link to download the article as PDF)
Download citations (links to download the citations from this article in formats compatible with various reference manager tools)
Open citations (links to open the citations from this article in various online reference manager services)
Biomedical and clinical sciences are experiencing a renewed interest in the fact that males and females differ in many anatomic, physiological, and behavioral traits. Sex differences in trait variability, however, are yet to receive similar recognition. In medical science, mammalian females are assumed to have higher trait variability due to estrous cycles (the 'estrus-mediated variability hypothesis'); historically in biomedical research, females have been excluded for this reason. Contrastingly, evolutionary theory and associated data support the 'greater male variability hypothesis'. Here, we test these competing hypotheses in 218 traits measured in >26,900 mice, using meta-analysis methods. Neither hypothesis could universally explain patterns in trait variability. Sex-bias in variability was trait-dependent. While greater male variability was found in morphological traits, females were much more variable in immunological traits. Sex-specific variability has eco-evolutionary ramifications including sex-dependent responses to climate change, as well as statistical implications including power analysis considering sex difference in variance.
Paired fins are a defining feature of the jawed vertebrate body plan, but their evolutionary origin remains unresolved. Gegenbaur proposed that paired fins evolved as gill arch serial homologues, but this hypothesis is now widely discounted, owing largely to the presumed distinct embryonic origins of these structures from mesoderm and neural crest, respectively. Here, we use cell lineage tracing to test the embryonic origin of the pharyngeal and paired fin skeleton in the skate (Leucoraja erinacea). We find that while the jaw and hyoid arch skeleton derive from neural crest, and the pectoral fin skeleton from mesoderm, the gill arches are of dual origin, receiving contributions from both germ layers. We propose that gill arches and paired fins are serially homologous as derivatives of a continuous, dual-origin mesenchyme with common skeletogenic competence, and that this serial homology accounts for their parallel anatomical organization and shared responses to axial patterning signals.