Acetyltransferase Enok regulates transposon silencing and piRNA cluster transcription.

The piRNA pathway is a highly conserved mechanism to repress transposon activation in the germline in Drosophila and mammals. This pathway starts from transcribing piRNA clusters to generate long piRNA precursors. The majority of piRNA clusters lack conventional promoters, and utilize heterochromatin- and HP1D/Rhino-dependent noncanonical mechanisms for transcription. However, information regarding the transcriptional regulation of piRNA clusters is limited. Here, we report that the Drosophila acetyltransferase Enok, which can activate transcription by acetylating H3K23, is critical for piRNA production from 54% of piRNA clusters including 42AB, the major piRNA source. Surprisingly, we found that Enok not only promotes rhino expression by acetylating H3K23, but also directly enhances transcription of piRNA clusters by facilitating Rhino recruitment. Taken together, our study provides novel insights into the regulation of noncanonical transcription at piRNA clusters and transposon silencing.


Introduction
In a wide range of organisms, repressing the activation of transposon insertions is essential for maintenance of genome stability [1]. Small  important roles in silencing transposons in eukaryotic genomes [2]. Also, mammals and Drosophila utilize the PIWI-interacting RNA (piRNA) pathway to achieve transcriptional and post-transcriptional silencing of transposons in the germline [3,4]. In Drosophila, the piRNA pathway starts from transcription of the 142 piRNA clusters, usually ranging from 50 to a few hundred kilobases and containing multiple copies of truncated or full-length transposons, which produces long piRNA precursors [5][6][7]. The long RNA precursors would then be processed through slicer-and Zucchini (Zuc)-dependent mechanisms into mature 23-29 nt piR-NAs that get loaded to the Piwi protein [8]. In addition, another two PIWI-clade Argonaute proteins, Ago3 and Aubergine (Aub), function in the ping-pong cycle to specifically amplify piRNAs against active transposons [9]. Guided by complementary piRNAs, Ago3 and Aub can mediate degradation of transposon transcripts, and the Piwi-piRNA complex can also direct heterochromatin formation at the loci of transposons and piRNA clusters by recruiting epigenetic factors [10][11][12][13], resulting in effective repression of transposons both transcriptionally and post-transcriptionally.
The Drosophila piRNA clusters can be divided into two classes: uni-strand, which produces piRNAs mainly from one genomic strand, and dual-strand, which produces piRNAs from both genomic strands. Transcription of the uni-strand clusters is proposed to be similar to the canonical transcription of protein-coding genes, as they contain clear promoter structures with enriched H3K4me2 and peaks of RNA polymerase II (Pol II) [14]. In contrast, dualstrand clusters lack clear Pol II promoter regions and are enriched for the heterochromatic H3K9me3 mark. Therefore, these clusters undergo noncanonical transcription that utilizes multiple internal initiation sites via heterochromatin-and Rhino (Rhi)-dependent mechanisms [6,15].
Rhi is the germline-specific heterochromatin protein 1D (HP1D), and it associates with Deadlock (Del) and Cutoff (Cuff) to form the RDC complex [14,16]. The RDC complex is recruited to dual-strand clusters by the interaction between H3K9me3 and the chromodomain of Rhi. At dual-strand clusters, the RDC complex licenses and promotes their transcription through four mechanisms. First, Del interacts with the germline-specific paralog of transcription factor IIA (TFIIA)-L, Moonshiner (Moon), and in turn recruits TFIIA and the TATAbox binding protein (TBP)-related factor TRF2 for transcription initiation [6]. Second, the RDC complex has been shown to suppress the splicing of piRNA cluster transcripts, which is proposed to facilitate piRNA production [16]. Third, Cuff recruits the transcription-export (TREX) complex to nascent transcripts to promote efficient transcription at piRNA clusters [17]. Fourth, Cuff interferes with recruitment of the cleavage and polyadenylation specificity factor (CPSF) complex, and therefore prevents premature termination during transcription of piRNA precursors [18]. While the positive roles of the RDC complex in noncanonical transcription of piRNA clusters were studied extensively, further transcriptional regulation upstream to the recruitment of this complex to piRNA clusters is still unclear.
The KAT6 acetyltransferases are highly conserved from budding yeast to mammals, and preferentially acetylate histone H3 among the four core histones [19]. The fly KAT6, Enok, has been shown to function as the major acetyltransferase for establishing the H3K23ac mark, which plays activating roles in transcription of genes [20,21]. H3K23ac has been suggested to destabilize the interaction between H3K27me3 and the chromodomain of Polycomb [22], and therefore may contribute to transcription activation. In the ovarian germline, Enok is important for the maintenance of germline stem cells [23], and is required for proper polarization of oocytes by promoting expression of the actin nucleator spir [20]. In this study, we report a novel role for Enok in the piRNA pathway. Mutating or knocking-down enok in the ovarian germline led to derepression of transposons and reduction in levels of piRNAs produced from a subset of piRNA clusters including the major piRNA source 42AB. We show that Enok binds to and acetylates H3K23 in the 5' region of rhi, and is required for its normal expression levels in the ovary. We further show that Enok is also required for proper Rhi recruitment to a subset of piRNA clusters to promote their transcription. Therefore, Enok contributes to proper transposon silencing in the germline by promoting transcription of rhi and piRNA clusters.

Transcript levels of 7 transposons are increased upon loss of functional Enok in the germline
The corresponding author has previously utilized the FLP/FRT/ovoD system to generate enok mutant germline clones to investigate the roles of Enok in the germline [20]. In this system, only egg chambers containing 16 homozygous mutant germline cells can develop beyond stage 4 of oogenesis [24]. Two enok mutant alleles, enok 1 and enok 2 , were used in the germline clone analysis, and both alleles encode mutant Enok proteins with abolished acetyltransferase function [20]. The enok 1 allele contains a W765-to-amber mutation, and produces a truncated Enok protein lacking more than 50% of the histone acetyltransferase (HAT) domain and the entire C-terminal region. The enok 2 allele contains a C746Y mutation inside the HAT domain that mutates a conserved cysteine residue, and therefore produces a full-length Enok protein with compromised HAT function. It was found that the enok mutant germline clone egg chambers showed defects in oocyte polarization. Therefore, RNA-seq analysis was performed using egg chambers of stages 5-14 from enok mutant or wild-type (WT) germline clone ovaries to identify genes involved in oocyte polarization that are transcriptionally regulated by Enok. Surprisingly, in addition to genes up-or down-regulated in enok mutants, we also found that several transposon families were activated in the enok 1 and enok 2 germline clone ovaries as compared with the WT control (Fig 1). Given the high degree of correlation between each independent replicate in our next generation sequencing (NGS) data sets (S1, S2, and S3 Figs), results from one representative replicate are shown in all figures demonstrating a Genome Browser view. The transposon families activated upon loss of functional Enok in the germline include both the long terminal repeat (LTR) and the non-LTR transposons. For example, as shown in Fig 1A, the LTR transposon Burdock (left panel) and the non-LTR transposon HeT-A (right panel) were activated in both the enok 1 and enok 2 mutants compared with WT ovaries. Furthermore, among the 126 transposon families, 7 families were significantly activated in enok mutant ovaries as compared with the WT control ( Fig 1B and 1C, and S1 Table). Taking these results together, we concluded that Enok is important for proper silencing of 7 transposon families in the germline.
Since the piRNA pathway is the major mechanism to suppress transposon activation in the germline, we examined the piRNA levels in WT and the two enok mutant germline clone ovaries by small RNA-seq analysis (S2 Table). While the levels of piRNAs mapping to genic regions were similar in WT and enok mutants (Fig 2A; left panel), the levels of piRNAs that mapped to a subset of piRNA clusters were decreased in enok mutant ovaries as compared with the WT control (Fig 2A; middle and right panels). We further examined how many piRNA clusters were down-or up-regulated in enok mutants using DESeq2 (S2 Table). Because the levels of piRNAs uniquely mapping to cluster 70, 107, 110 and 128 were zero in WT samples, these four clusters were omitted, resulting in total 138 piRNA clusters in this analysis. Among these 138 piRNA clusters, 63 (46%) and 67 (49%) clusters showed significantly decreased levels of piRNAs in both the enok 1 and enok 2 mutants compared with WT ovaries in the sense and antisense strands, respectively (S4A Fig). Also, 56 clusters (41%) showed reduced piRNA levels in both the sense and antisense strands in both enok mutants (S4A Fig; top right panel). In contrast, only 11 (8%) clusters had increased piRNA levels in both the sense and antisense strands in both enok mutants (S4A Fig; bottom right panel). Therefore, our data suggest that Enok mainly plays a positive role in the piRNA biosynthesis.
The major piRNA source 42AB and cluster 29 are among the most down-regulated clusters (Figs 2A and S4B). As shown in S5A Fig, levels of sense and antisense piRNAs mapping to 42AB were both reduced in enok mutants, and the sense piRNAs were affected to a larger extent than the antisense species. The ping-pong amplification of piRNAs mapping to 42AB was also compromised in enok mutant germline clone ovaries as compared with the WT control. We further examined the piRNAs mapping to transposons. While the enok mutants only showed mild reductions in antisense piRNAs mapping to all transposons (S5B Fig), antisense piRNAs mapping to the highly overexpressed Burdock transposon were more severely reduced (S5C Fig). This result suggests that our small RNA-seq data are in line with the results from RNA-seq analysis ( Fig 1C) for a couple of transposons including Burdock.
After the piRNA clusters were first annotated [5], further studies have identified more regions in the fly genome that give rise to piRNAs [14]. Therefore, all the genomic regions giving rise to piRNAs were divided into 1 kb bins and defined as piRNA source loci (piRNA-SL). These piRNA-SL includes >91% of 1 kb bins defined by previously annotated piRNA clusters and extends the size of piRNA source regions more than 3-fold [14]. We further examined how mutations of enok affect levels of piRNAs mapping to the 3 types of piRNA-SL defined in the previous study [14]: Rhi-dependent (RD-SL), Rhi-independent (RI-SL) and soma-specific (SO-SL). A population of heterochromatic 1 kb bins not giving rise to piRNAs (het non-SL) was also included as the control group [14]. As shown in Fig 2B, levels of piRNAs mapping to RD-SL were significantly decreased in the two enok mutants as compared with the WT control. In contrast, piRNAs mapping to RI-SL, SO-SL or het non-SL were not down-regulated in enok mutants ( Fig 2B). Therefore, Enok may play an important role in the production of piRNAs that are enriched in germline cells and dependent on Rhi.
Next, we determined the respective Enok dependency (piRNA fold change upon enok mutation) of each RD-SL (S6A Fig) [14]. Based on this dependency, the majority (87%) of RD-SL were split into two populations. The log2 values of piRNA fold changes from the first population (n = 1074) were less than -1.5 in both the enok 1 and enok 2 mutants as compared with the WT control (S6A and S6B Fig), and therefore piRNAs from this population were dependent on Enok (Enok-dependent source loci; ED-SL). The second population (n = 4501) produces piRNAs that were not or only weakly dependent on Enok (Enok-independent source loci; EI-SL) (S6A Fig). Notably, only 404 RD-SL showed log2 values of piRNA fold changes more than 1.5 in both enok mutants compared with the WT control (S6C Fig), again suggesting that Enok mainly plays a positive role in piRNA production. As shown in Fig 2C, piRNAs from ED-SL were severely reduced in both enok mutants, while EI-SL piRNAs were unaffected. Thus, we concluded that Enok may be important for piRNA production from~20% of the Rhi-dependent loci (S3 Table). Enok is important for piRNA production from a subset of piRNA clusters. (A) Loss of functional Enok in the germline reduced levels of piRNAs from a subset of piRNA clusters but not from genic regions. The levels of piRNAs that uniquely mapped to genic regions (left panel) and piRNA clusters (middle and right panels) in enok mutant germline clone ovaries versus the WT control are plotted in a log scale. Selected piRNA clusters that were severely down-regulated in enok mutants are indicated by red dots. (B-C) A box-plot displaying the distribution of piRNA levels for all 1 kb bins belonging to the indicated groups. Center line, median; box limits, upper and lower quartiles; whiskers, the 5th and 95th percentile (outliers not shown). In (A-C), genotypes of females are as described in Fig 1. In (B-C), P-values were calculated using Wilcoxon Rank sum test. https://doi.org/10.1371/journal.pgen.1009349.g002

Enok regulates the expression levels of three genes involved in piRNA biosynthesis
We sought to investigate the mechanisms by which Enok contributes to the Rhi-dependent piRNA production. First, we examined the expression levels of genes known to be involved in the piRNA pathway. A previous transcriptome-wide RNAi screen has identified top 100 genes that are important for the germline piRNA pathway [25]. Among them, 9 genes showed a RPKM value less than 5 in our RNA-seq data obtained from WT ovaries, and therefore were excluded from the gene list. In the remaining 91 genes, expression levels of 5 genes, rhi, Brother of Yb (BoYb), maelstrom (mael), shutdown (shu) and CG12721/moon, were significantly downregulated in the enok mutant germline clone ovaries as compared with the WT control (Figs 3A and S7A). In contrast, levels of other known key regulators of the piRNA pathway, including piwi, were largely unaffected upon loss of Enok in the germline (S7A and S7B Fig, and S4 Table). Immunostaining in the enok mutant germline clone ovaries also revealed reductions in the protein levels of Rhi and Mael (S7C and S7D Fig) similar to those in their mRNA levels (S7A Fig). However, the enok mutant ovaries still expressed more than 47% of the wild-type levels of BoYb, mael, shu and moon mRNA, which is similar to or better than a haploid condition. On the other hand, the expression levels of rhi showed a~75% reduction in enok mutants as compared with WT ( Fig 3A and S4 Table). Since heterozygotes of most piRNA pathway components show normal transposon silencing in the germline [26], our results suggest that Enok may contribute to proper piRNA biosynthesis partly by promoting the expression of rhi. Nevertheless, we cannot exclude the possibility that the combination of those mild reductions in BoYb, mael, shu and moon levels may also have effects on piRNA production.
To examine whether enok and rhi mutants show similar effects on piRNA production, we first analyzed the abundance of antisense piRNAs encoded by each transposon family (S8A Fig). Among the 7 transposon families activated in enok mutants (Fig 1C), levels of antisense piRNAs encoded by Burduck, HeT-A and Bari1 were reduced in enok mutants as compared with the WT control, while the other 4 transposon families only showed subtle decreases in antisense piRNAs (S8A Fig). The decreases in levels of primary piRNAs in enok mutants may be masked by ping-pong mechanism with participation of maternal transposon-specific piR-NAs, and the increased RNA levels of blood, TAHRE, GATE and HMS-Beagle in enok mutants ( Fig 1C) may be determined by another mechanism. We further analyzed the fold change in transposon family expression in enok mutants (S8B Fig; y axis, enok1/WT or enok2/WT) and the fold change in antisense piRNAs encoded by the same transposon family (S8B Fig; x axis, enok1/WT or enok2/WT) as performed in the previous study for the rhi mutant [15]. The results for enok mutants (S8B Fig) show some similarity to those for the rhi mutant [15]: some of the highly overexpressed transposons (Burdock and HeT-A) also showed large reductions in antisense piRNAs, while expression of blood increased 10-12 fold but the total antisense piRNA pool was reduced by only 12%-25%. This similarity between enok and rhi mutants further supports the hypothesis that Enok may promote piRNA biosynthesis by regulating rhi expression.
To further examine whether Enok binds to and directly regulates transcription at these 5 genes involved in piRNA biosynthesis, Enok ChIP-seq analysis was performed in wild-type ovaries using two independently raised α-Enok antibodies (S9A Fig). ChIP-seq results showed that Enok was localized to the 5' regions of rhi, mael and shu (Figs 3B and S9B). This result is also in line with the previous report that Enok was enriched at the 5' end of mael and depleting Enok in S2 cells reduced Pol II occupancy at mael 5' [20]. Furthermore, ChIP-qPCR analysis showed that Enok levels at the 5' end of rhi, mael and shu were reduced upon germline-specific knockdown of enok ( Fig 3C). Therefore, our results suggest that Enok may contribute to It was previously demonstrated that Enok is the major enzyme for establishing the H3K23ac mark in Drosophila, and the Enok-mediated H3K23ac mark plays positive roles in transcriptional regulation [20]. Therefore, we asked whether H3K23ac is enriched at rhi and mael. To answer this question, ChIP-qPCR analysis of H3K23ac was carried out in the ovary. As shown in Fig 3D and 3E, the H3K23ac mark was enriched 2-3 fold in the 5' regions of rhi and mael relative to the intergenic region. Moreover, knocking down enok in ovaries reduced the H3K23ac level in the 5' region of rhi (Fig 3E), supporting that Enok acetylates H3K23 at rhi 5'. On the other hand, since the Enok level at mael 5' only reduced 43% upon enok knockdown (Fig 3C), the H3K23ac level at mael 5' did not reduce in enok knockdown ovaries ( Fig 3E). However, it was previously shown that depleting Enok in S2 cells reduced H3K23ac levels in the 5' and 3' regions of mael [20]. Taken together, we concluded that Enok may promote the expression of rhi and mael by acetylating H3K23.
Rhi is the germline-specific HP1D. It forms the RDC complex with Del and Cuff, and is recruited to the H3K9me3-enriched dual-strand piRNA clusters to license noncanonical transcription and to suppress piRNA precursor splicing [6,16]. Since rhi is severely down-regulated in enok mutant ovaries, we asked whether levels of RNAs derived from the Rhi-dependent 42AB and the Rhi-independent cluster 2 are regulated by Enok. To answer this question and to verify our RNA-seq results obtained from enok mutant germline clone ovaries, RNAi against enok in the ovarian germline was carried out using two different UAS-shRNA-enok transgenes driven by MTD-Gal4 (P{w[+mC] = otu-GAL4::VP16.R}1; P{w[+mC] = GAL4-nos. NGT}40; P{w[+mC] = GAL4::VP16-nos.UTR}CG6325 MVD1 ). Then, the RNA levels in the enok knockdown and control ovaries were examined by RT-qPCR. Consistent with the RNA-seq results in enok mutants, knocking down enok in the germline resulted in activation of the Burdock transposon as compared with the EGFP RNAi-1 control (S10A Fig). Also, the rhi mRNA levels were significantly reduced in both the enok RNAi-1 and enok RNAi-2 ovaries when RNAi against enok was performed at 29˚C for 3 days (S10A Fig), suggesting that Enok indeed regulates the expression of rhi. The levels of RNAs derived from 42AB and cluster 2 were further examined by strand-specific RT-qPCR. Similar to the rhi mutant [15], depletion of Enok expression plotted against the expression levels for the enok 1 (left panel) or the enok 2 (right panel) mutant ovaries versus the WT control. RPKM, reads per kilobase per million. Ninety-one genes (RPKM > 5 in WT) involved in the piRNA pathway are marked with black circles. Genotypes of females are as described in Fig 1. (B) Genome Browser view of Enok ChIP-seq data at rhi, BoYb, mael, shu and CG12721/moon. ChIP-seq experiments were performed in ovaries from Oregon-R (OreR) in two independent replicates, and results from one representative replicate are shown. Input (IN) and Immunoprecipitation (IP) tracks are shown. Enrichment of Enok was observed in the 5' region of rhi, mael and shu. The FlyBase Genes track is shown to represent the genes, showing the exon (wider bars) and intron (blue line) structure of each annotated coding transcript, with direction of transcription indicated. in the germline led to decreased levels of RNAs derived from 42AB but increased levels of RNAs from cluster 2 (S10B Fig). This result suggests that Enok is important for transcription of the Rhi-dependent 42AB cluster.
We further analyzed the RNA-seq reads uniquely mapping to piRNA clusters to examine whether levels of RNAs transcribed from each cluster were decreased in the enok mutant germline clone ovaries. Six piRNA clusters (42AB, cluster 5, 29, 78, 131 and 137) were significantly down-regulated (adjusted P-value < 0.05) in the two enok mutant ovaries (S5 Table). Also, as shown in S11 Fig, the four severely down-regulated clusters in RNA-seq analysis (42AB, cluster 5, 29 and 78; red dots) match the severely down-regulated clusters identified by small RNA-seq (Fig 2A; red dots). On the other hand, the RNA-seq reads uniquely mapping to three severely down-regulated clusters labeled in Fig 2A (clusters 11, 33 and 42) were not significantly affected. Cluster 11 contains very low levels of unique RNA-seq reads, and this may be why the adjusted p-value is not significant (>0.05; S5 Table). Clusters 33 and 42 both contain an expressed gene within them, and therefore decreases in the relatively low levels of piRNA precursor transcripts may be masked by the mRNA levels from those two genes. Taken together, these results suggest that the defective transcription of a subset of piRNA clusters in enok mutants indeed resulted in reduced piRNA production from those clusters.

Enok has roles in promoting transcription of piRNA clusters that are independent of its role in regulating rhi expression
We noticed that the rhi mRNA level in the enok RNAi-1 ovaries was reduced only~30% (S10A Fig), suggesting that the effects on RNA levels transcribed from 42AB in the enok RNAi-1 ovaries may be independent of the role of Enok in promoting rhi expression. To test this possibility, RNAi against enok was carried out at 25˚C for 1 day and then at 29˚C for 2 days to knockdown enok more mildly, as the activity of Gal4 is higher at 29˚C than 25˚C. Under this condition, we were able to obtain efficient enok knockdown without affecting the mRNA levels of rhi, mael and BoYb in the enok RNAi-1 ovaries (Fig 4A). Supporting our hypothesis, strand-specific RT-qPCR analysis showed that levels of RNAs derived from 42AB, except for cl1-32-minus, were significantly decreased in the enok RNAi-1 ovaries as compared with the Luciferase control ( Fig 4B; left panel). On the other hand, levels of RNAs derived from cluster 2 were not decreased upon depletion of Enok ( Fig 4B; right panel). Since knocking down enok in the germline without affecting rhi levels decreased the levels of RNAs from 42AB, our results suggest that, in addition to its role in regulating rhi levels, Enok may play another role in facilitating transcription of piRNA clusters including 42AB.
We further examined whether knocking down enok in the germline can indeed reduce the levels of RNAs derived from 42AB without affecting the bulk protein levels of Rhi. As shown in S10A Fig, even when RNAi against enok was performed at 29˚C for 3 days, the mRNA levels of rhi and mael in the enok RNAi-1 ovaries were only decreased 28% and 13%, respectively. Consistent with these small changes in the rhi and mael levels, immunostaining showed only mild reductions in staining signals of Mael and Rhi in the enok RNAi-1ovaries as compared with the EGFP RNAi-1 control (S12A and S12B Fig). Since the staining signals of endogenous Rhi were relatively weak and the α-Rhi antibody does not work for western blotting, we further utilized the GFP-rhi transgene expressing GFP-tagged Rhi in the context of its endogenous control region as reported in the previous study [14]. Both western blotting and fluorescence microscopy showed that the global levels of GFP-Rhi were largely unaffected upon germlinespecific knockdown of enok as compared with the control (Luciferase) (S12C-S12E Fig). Also, bulk levels of the GFP-tagged Del, which forms a complex with Rhi [14], were similar in the Luciferase and enok RNAi-1 ovaries (S12C Fig; lanes 4-6). Taken together, we concluded that knocking down enok in the germline can decrease levels of RNAs from 42AB (S10B Fig; enok RNAi-1) without affecting the global Rhi levels (S12B-S12E Fig).

Enok is important for Pol II occupancies at ED-SL
To further test the hypothesis that Enok may regulate transcription at piRNA clusters, ChIPseq analysis of Pol II was performed in ovaries. The Pol II ChIP-seq was carried out in ovaries with the UAS-EGFP RNAi-1 (GFP_KD) or UAS-enok RNAi-1 (enok_KD) transgene driven by MTD-Gal4 to examine whether Enok is important for Pol II occupancies at piRNA clusters. As shown in Fig 5A, while sites of Pol II occupancy at cluster 2 were not affected upon depletion of Enok in the germline (right panel), those at 42AB were decreased in the enok_KD ovaries compared with the GFP_KD control in two independent replicates (left panel; light yellow shades). This result suggests that Enok has a role in promoting transcription of a subset of dual-strand clusters including 42AB.
We further examined the changes in Pol II occupancies at ED-SL upon enok knockdown. Consistent with our small RNA-seq results (Fig 2C), Pol II occupancies at ED-SL were significantly reduced upon enok knockdown, while those at EI-SL were unaffected (Fig 5B). In addition, RNA-seq reads uniquely mapping to ED-SL were also significantly decreased in the two enok mutants (Fig 5C). Since knocking down enok without affecting Rhi levels (S12C-S12E Fig) decreased Pol II occupancies at ED-SL (Fig 5B), we concluded that Enok may have a role, which is independent of its role in promoting rhi expression, in enhancing transcription of these Enok-dependent loci.

Enok is critical for the recruitment of Rhi to ED-SL
Next, we sought to investigate how Enok regulates transcription of piRNA clusters. As Enok ChIP-seq analysis did not detect Enok peaks at 42AB (Fig 5A) or enrichment of Enok at ED-SL relative to EI-SL (S13A Fig), this regulation is unlikely to be mediated by the H3K23ac mark. H3K23 and H3 ChIP-seq analyses in GFP_KD and enok_KD ovaries further confirmed that the H3K23ac levels mapping to piRNA clusters or transposons were not decreased upon enok knockdown (S13B Fig and S6 Table). Therefore, we explored other possibilities. Rhi anchors the RDC complex at piRNA clusters by recognizing the H3K9me3 mark, and has been reported to play multiple key roles in seeding and facilitating productive transcription at all dual-strand piRNA clusters [6,[14][15][16][17][18]. Although the global protein levels of Rhi were largely unaffected upon enok knockdown (S12C-S12E Fig), it is possible that Enok may regulate Rhi recruitment to Enok-dependent piRNA clusters. Indeed, knocking down enok in the germline severely reduced both the transcript levels and Rhi occupancy at 42AB (Figs 6A, S14A and S14B). Since Rhi recognizes H3K9me3, we further examined the effects of enok knockdown on the H3K9me3 levels at 42AB by ChIP analysis. Upon depletion of Enok in the germline, levels of H3K9me3 at the cl1-A and cl1-32 sites in 42AB were reduced only 10%-20% as compared with the control EGFP RNAi-1 ovaries (S14C Fig). It has been shown previously that the H3K9me3 levels were decreased~65% at the cl1-A site when piwi was knocked down during embryogenesis, and the Rhi levels at the same site were also reduced~65% [12]. Therefore, thẽ 85% reduction in Rhi occupancy at 42AB upon enok knockdown (Fig 6A; cl1-A) is unlikely to be mediated by the~15% reduction in H3K9me3 levels at 42AB (S14C Fig; cl1-A). Also, H3K9me3 levels in the 5' regions of Burdock and HeT-A remained similar in the EGFP RNAi-1and enok RNAi-1 ovaries (S14C Fig). These results suggest that depletion of Enok may not have great impacts on H3K9me3 levels at piRNA clusters or transposons.
To further assess the effects of enok knockdown on Rhi recruitment to piRNA clusters/ source loci, Rhi ChIP-seq analysis was performed in the GFP_KD and enok_KD ovaries. Consistent with the Rhi ChIP-qPCR result (Fig 6A), knocking down enok in the germline resulted in decreased Rhi occupancies across 42AB (Fig 6B). Also, Rhi occupancies at ED-SL were significantly reduced upon enok knockdown, while those at EI-SL were unaffected ( Fig 6C). As the global protein levels of Rhi were largely unaffected in enok knockdown ovaries (S12C-S12E Fig), this result suggests that Enok may promote transcription of these Enok-dependent piRNA source loci by facilitating Rhi recruitment. Notably, Rhi occupancies at ED-SL were higher than those at EI-SL in the control knockdown ovaries (Fig 6C; GFP_KD). However, the H3K9me3 levels (GEO accession GSE55824) were similar between ED-SL and EI-SL (Fig 6D). Therefore, Enok may have a role in enhancing Rhi recruitment/binding to the Enok-dependent loci as compared with the Enok-independent loci. Since knocking down enok in the germline led to reduced Rhi occupancy at 42AB, we sought to examine whether overexpression of rhi can rescue this defect in the enok knockdown ovaries. To this end, the UAS-RNAi or UAS-Luciferase construct and an UASp-rhi:GFP transgene [15] were overexpressed in the ovarian germline by the nos-Gal4-VP16 driver (S15A and S15B Fig). Interestingly, while the UASp-rhi:GFP transgene was overexpressed up to~64-fold of the endogenous rhi levels (Figs 7A and S15C), the defective transcription at 42AB in the enok RNAi-1 ovaries could not be rescued by overexpressing GFP-rhi (Figs 7B and S15D). Moreover, ChIP analysis of Rhi showed that overexpression of GFP-rhi was not able to rescue the reduced Rhi occupancy at 42AB upon enok knockdown (Fig 7C). We also asked whether the regulation of Rhi recruitment by Enok is dependent on its HAT activity. To answer this question, we examined whether overexpressing rhi can rescue the defective transcription at 42AB in the enok 2 mutant germline clone ovaries, as the enok 2 allele encodes a full-length Enok protein with compromised HAT function. Interestingly, overexpression of GFP-rhi failed to rescue the decreased levels of RNAs derived from 42AB in the enok 2 mutant (S15E and S15F Fig). These results strongly suggest that Enok plays a critical role, which is dependent on its HAT activity and cannot be bypassed by overexpressing rhi, in Rhi recruitment to some specific Enok-dependent piRNA source loci.
Taken together, our results suggest that Enok contributes to proper transposon silencing in the ovarian germline through three pathways summarized in Fig 7D. First, Enok acetylates H3K23 at genes involved in the piRNA production, mainly rhi, and promotes their expression (left). Second, Enok regulates the recruitment of Rhi to enhance transcription at piRNA clusters, facilitating the synthesis of piRNA precursors (middle). Third, as only 3 out of 7 activated transposon families (Fig 1C) showed decreased antisense piRNA levels in enok mutants (S8A

Discussion
We report herein a novel role for Enok in suppressing the activation of transposons in the germline. Loss of functional Enok in the ovarian germline resulted in activation of 7 transposon families ( Fig 1C). This amount of activated transposon families in enok mutant ovaries is comparable to the 17 families activated in the rhi mutant [15]. Our RNA-seq analysis showed a~75% reduction in the mRNA levels of rhi in enok mutant germline clone ovaries as compared with the WT control (Figs 3A and S7A). Knocking down enok in the ovarian germline using two different UAS-shRNA-enok constructs also reduced the mRNA levels of rhi as compared with two different control fly lines (Figs 4A and S10A). In addition, Enok ChIP-seq analysis revealed that Enok is localized to the 5' region of rhi (Fig 3B), and the Enok-dependent enrichment of H3K23ac at the 5' end of rhi suggests that the Enok-mediated H3K23ac mark promotes rhi expression, contributing to proper piRNA production (Fig 3D and 3E). Indeed, enok mutant germline clone ovaries showed decreased levels of piRNAs that mapped to RD-SL (Fig 2B). However, not all RD-SL showed decreased piRNA levels in enok mutants. About 20% of the 6426 RD-SL showed reduced piRNA levels in enok mutants (Fig 2C and S3 Table). Therefore, the remaining 25% of rhi levels in enok mutant ovaries (Fig 3A) may be sufficient to support transcription of the RD-SL that were not affected by loss of Enok. More strikingly, knocking down enok in the germline, without affecting the global protein levels of Rhi (S12B-S12E Fig), reduced Rhi occupancies at ED-SL but not at EI-SL (Fig 6C). This result suggests that Enok regulates Rhi recruitment specifically at ED-SL. The enok and rhi mutants show similar effects on the fold changes in transposon family expression and in antisense piRNAs (S8B Fig) [15]. However, among the top 24 most highly overexpressed families in rhi [15], loss of Enok in the germline specifically activates 7 families. This specificity suggests that these 7 families may be more sensitive to reductions in Rhi recruitment to a subset of piRNA source loci. Taken together, Enok may contribute to fine-tuning transcription of piRNA clusters by modulating rhi expression and by regulating Rhi recruitment to Enok-dependent piRNA source loci. Three genome-wide RNAi screens have been reported before, but two of them were specifically performed in ovarian somatic cells [27,28]. Knocking down enok in ovarian somatic cells using the tj-Gal4 driver did not activate the soma-dominant transposon, Gtwin (S16A Fig), suggesting that Enok may be dispensable for transposon silencing in the soma. In the genome-wide screen in the germline [25], the enok RNAi construct (KK108400) is a long hairpin RNA. The efficiency of knocking down enok by long hairpin RNAs is lower than by short hairpin RNAs in the germline, even in the presence of additional Dicer-2 (S16B Fig). Czech et al. indeed showed that knocking down enok weakly activated the blood and Burdock transposons (z-scores of -0.5 and -0.74, respectively). However, this activation effect did not reach the threshold (z-score of -1.5 or lower) applied in the screen. We used two different short hairpin RNA constructs against enok to deplete Enok in the germline, and therefore we were able to detect stronger activation of transposons (Figs 4A and S10A), possibly due to better knockdown efficiencies.
Enok is the major enzyme responsible for the abundant H3K23ac mark [20,21]. It was previously demonstrated that Enok is localized to the 5' end of its target genes, spir and mael, and promotes their expression by acetylating H3K23 [20]. Here we further report rhi and a subset of RD-SL (defined as ED-SL) as novel targets that are transcriptionally regulated by Enok. Intriguingly, while the 5' region of rhi is enriched with Enok and H3K23ac (Fig 3B and 3D), Enok is not enriched at ED-SL relative to EI-SL (S13A Fig). Also, knocking down enok in ovaries reduced the H3K23ac levels at the 5' end of rhi ( Fig 3E) but not at piRNA clusters (S6 Table). These results suggest that Enok facilitates rhi expression by acetylating H3K23, but regulates the transcription of ED-SL through other mechanisms. Notably, knocking down enok in the ovarian germline severely reduced the Rhi occupancy to sites in 42AB (Fig 6A and 6B), while global protein levels of Rhi and the H3K9me3 levels at 42AB were largely unaffected (S12B-S12E and S14C Figs). Therefore, Enok is likely to promote transcription of ED-SL by regulating Rhi recruitment.
The transcription of dual-strand piRNA clusters utilizes noncanonical heterochromatindependent internal initiation. Transcription initiation at these clusters was proposed to take place by the H3K9me3-bound RDC complex recruiting the germline-specific paralog of TFIIA-L, Moon [6]. In this study, we show that Enok is important for both Rhi and Pol II occupancies at a subset of RD-SL (defined as ED-SL) (Figs 5B and 6C), suggesting that Enok can facilitate transcription of these piRNA source loci. As Rhi is highly enriched across the entire 42AB cluster [14] and no Enok peaks were detected within 42AB (Fig 5A), Enok is unlikely to regulate the Rhi occupancy at 42AB by directly recruiting it. Also, our Co-IP assay failed to detect interaction between Enok and the overexpressed Rhi in ovaries. Interestingly, the HAT activity of Enok is critical for transcription of 42AB even when rhi is overexpressed (S15E and S15F Fig). Therefore, it is possible that Enok may play a role in acetylating some factors that are required for Rhi recruitment, or it may have an indirect role in Rhi recruitment by promoting expression of other genes with yet unidentified functions in the piRNA pathway. Notably, while knocking down enok in the germline decreased the RNA levels transcribed from both genomic strands at cl1-A and from the sense strand at cl1-32, RNA levels transcribed from the antisense strand at cl1-32 was not affected by depletion of Enok (Figs 4B and S10B). Thus, within dual-strand clusters, Enok may regulate the internal initiation in specific regions. Taken together, our study provides novel information regarding noncanonical transcription and transposon silencing in the germline.

Fly strains and culture
Drosophila melanogaster was crossed and grown on standard media at 25˚C unless stated otherwise. To generate germline clone ovaries, larvae were heat shocked twice at 37˚C for 40 min, with 30 min break at 25˚C in between, for two consecutive days starting from the second day after hatching. Females were conditioned with wet yeast for 2-3 days before ovary dissection. MTD-Gal4 driven RNAi and overexpression of transgenes were carried out at 29˚C for 2-3 days. Females were conditioned with wet yeast for 1-2 days before ovary dissection.

RNA purification and RT-qPCR
RNA isolation and cDNA synthesis were performed essentially as described previously with minor modifications [20]. Briefly, total RNA was isolated using TRIzol (Invitrogen) or REzol (Protech), treated with DNase (Invitrogen), and subjected to cDNA synthesis using Super-scriptIII reverse transcriptase (Invitrogen) or M-MuLV reverse transcriptase (Protech). For examination of the mRNA levels of transposons or genes, qPCR was performed in triplicate for each biological replicates, and the data were analyzed using the standard curve method. The mRNA levels of rp49 were used as the normalization control. For the strand-specific RT-qPCR analysis for RNAs derived from piRNA clusters, qPCR data were analyzed by the ΔΔCt method, using the transcript levels of rp49 in each sample as the internal normalization standard, which was set as 1. For each biological replicate, the mean of technical triplicates was used. The p-values were calculated using these means by T-tests

RNA-seq analysis
RNA-seq data were published previously [20], and were retrieved from the ArrayExpress database (www.ebi.ac.uk/arrayexpress; accession number E-MTAB-2521). The RPKM values for each gene were calculated as described in the previous study [20]. The expression levels of transposon families were calculated with piPipes 1.5.0, using the default RNA-seq pipeline with UCSC dm3 as the reference genome (0 mismatches; multiple mappers assigned using Expectation-Maximization algorithm; normalized to genome-mapping reads) [29]. Reads mapping to the transposon sequences from flyBase were normalized to transposon length and library depth (total mapped reads) to obtain the RPKM values for each transposon family. Transposon families activated in enok mutant ovaries were identified by TEtranscripts [30] using the cutoff of adjusted P-value < 0.05, fold change > 1.5 and RPKM > 1 in the enok 1 and enok 2 data sets. For S4 Table, only the unique mappers were used (library depth for RPKMs: uniquely mapped reads), and the RNA-seq data were analyzed by edgeR [31] using UCSC dm3 and dm6, respectively (cutoff: FDR < 0.05). For S5 Table, only the unique mappers were used, and reads mapping to piRNA clusters were analyzed by DESeq2 [32] using UCSC dm3 (cutoff: adjusted P-value < 0.05).

Small RNA-seq analysis
Ovaries were dissected from 3-5 day old flies. Total RNA was isolated from egg chambers of stages 5-14 using TRIzol (Invitrogen) followed by DNase treatment (Invitrogen), and subjected to small RNA-seq analysis.
Small RNA libraries were generated using TruSeq Small RNA Sample Preparation Kit (RS-200-0012; Illumina) with the following modification. Size selection was done by purifying 18-30 nt long small RNAs from a polyacrylamide gel. A 2S blocking DNA oligo for Drosophila rRNA was added (at 65˚C; 5 min) between the step of 3' primer Stop solution incubation and the step of 5' primer incubation. The resulting libraries were further size-selected by excising 135-160 nt long libraries from a polyacrylamide gel before sequencing on a HiSeq2500 (Illumina).
Reads generated were 51-base-pair (bp) single-end, directional using the Illumina protocol. The abundance of piRNAs mapping to each piRNA cluster or genic region shown in S2 Table were calculated with piPipes 1.4.12, using the default small RNA-seq pipeline with UCSC dm3 as the reference genome (0 mismatches; only unique mappers; normalized to miRNAs) [29]. Mapping statistics are shown in S1 Fig. The piRNA clusters with decreased or increased levels of piRNAs mapping to them in enok mutant ovaries were identified by DESeq2 [32] using the cutoff of adjusted P-value < 0.05.
For the 1 kb bin analysis, only the unique mappers were used. The definition of RD-SL, RI-SL, SO-SL and het non-SL loci was obtained from the previous study [14]. After normalized to library depth (uniquely mapped reads), the mapped reads from 3 independent replicates were pooled and normalized to mappability [14], and bins with zero reads were excluded from the analysis. The Wilcoxon Rank sum test with continuity correction was used to calculate P values. As mentioned in the previous study [14], a P value of <2.2 × 10 −16 is the maximum R (wilcox.test) can accurately calculate and therefore the P values less than 2.2 × 10 −16 were shown as "P < 2.2e-16".
For the 1 kb bin analysis, only the unique mappers were used, and the definition of RD-SL, RI-SL, SO-SL and het non-SL loci was obtained from the previous study [14]. After normalized to library depth (uniquely mapped reads), the mapped reads from independent replicates were pooled and normalized to mappability [14], and bins with zero reads were excluded from the analysis. The Wilcoxon Rank sum test with continuity correction was used to calculate P values, and the P values less than 2.2 × 10 −16 were shown as "P < 2.2e-16".

Correlation analysis of the NGS data
The Spearman's correlation coefficients between each replicate/sample in our RNA-seq, small RNA-seq and ChIP-seq data sets were calculated with the plotCorrelation tool from deepTools using multiBamSummary, 10 kilobases by default [34], additional parameters "-skipZeros". The Principal Component Analysis (PCA) plot was generated using R.

Immunofluorescent staining
Immunofluorescent staining of ovaries was carried out as previously described [20]. Ovaries were dissected in Grace's insect medium. The following primary and secondary antibodies were used: α-

Primers
Primers used in strand-specific RT reaction are as follows:  rp49. (B) The same total RNA samples used in (A) were subjected to strand-specific RT-qPCR analysis for RNAs derived from 42AB. The location of amplicon used in the qPCR reaction is as indicated in S10B Fig. (C) The H3K9me3 levels at 42AB (cl1-A and cl1-32), the 5' region of Burdock (Burdock 5') and the promoter of HeT-A (HeT-A pro) in ovaries were analyzed by ChIP-qPCR. The H3K9me3 IP signals were first normalized to the histone H3 IP signals (H3K9me3/H3), and then the H3K9me3/H3 values were normalized to the mean of H3K9me3/H3 values obtained for the rp49 loci in EGFP RNAi-1 ovaries, which was set as 1. The location of amplicons used for 42AB is the same as indicated in S10B Fig. The amplicon for Burdock 5' is located at the position 152-283 of a full-length Burdock insertion, and that for HeT-A pro is located at the position 4094-4240 of the HeT-A{}4795 insertion [37]. In (A-C), data represent the mean of three biological replicates +/-SD. � P < 0.05, �� P < 0.01, ��� P < 0.001 (Student's t-test). Genotypes of the females are as described in S10 Whole cell extracts were prepared from ovaries and subjected to western blotting. β-tubulin was used as the loading control. (B) Ovarioles containing the germarium and egg chambers up to stage 9 were stained with DAPI. Bars: 50μm. (C) RT-qPCR analysis of ovaries was used to examine the expression levels of the indicated genes. The mRNA levels were normalized to the levels of rp49. (D) The same total RNA samples used in (C) were subjected to strand-specific RT-qPCR analysis for RNAs derived from 42AB. The location of amplicon used in the qPCR reaction is as indicated in S10B Fig. (E) RT-qPCR analysis of ovaries was used to examine the expression levels of the indicated gene or transposon. The mRNA levels were normalized to the levels of rp49. (F) The same total RNA samples used in (E) were subjected to strand-specific RT-qPCR analysis for RNAs derived from 42AB. The location of amplicon used in the qPCR reaction is as indicated in S10B Fig. In (C-F), data represent the mean of three biological replicates +/-SD except for enok 2 , which represent two biological replicates. � P < 0.05, �� P < 0.01, ��� P < 0.001 (Student's t-test). n.s.: not significant.  Table. Abundance of piRNAs mapping to piRNA clusters and genic regions. The abundance of piRNAs mapping to the sense or antisense strand of each piRNA clusters or genic regions obtained from the small RNA-seq analysis of triplicate samples of the WT ovary or the enok mutants are shown. The unique reads mapping to each piRNA cluster or genic region were normalized to miRNAs and defined as the abundance of piRNAs. The results of DESeq2 analysis on piRNA reads mapping to piRNA clusters are also shown. Genotypes are as described in S1 Table. (XLSX) S3 Table. Abundance of piRNAs mapping to ED-SL and EI-SL. The abundance of piRNAs mapping to ED-SL or EI-SL obtained from the small RNA-seq analysis of triplicate samples of the WT ovary or the enok mutants are shown. The unique reads mapping to each ED-SL or EI-SL were normalized to library depth and mappability, and defined as the abundance of piR-NAs. The location, mappability and transposon content of each ED-SL or EI-SL are also shown. Genotypes are as described in S1 Table. (XLSX) S4 Table. RPKMs of genes involved in the piRNA pathway. The list of genes involved in the piRNA pathway was obtained from the previous transcriptome-wide RNAi screen [25]. The RPKM values obtained from the RNA-seq analysis of triplicate samples of the WT ovary or the enok mutants are shown (using UCSC dm3 or dm6 as the reference genome). Differential analysis was performed using edgeR. Genotypes are as described in S1 Table. (XLSX) S5 Table. RNA levels at piRNA clusters. Unique mappers from the RNA-seq analysis of triplicate samples of the WT ovary or the enok mutants mapping to each piRNA cluster were analyzed using DESeq2. The six clusters that were significantly down-regulated (adjusted pvalue < 0.05) in the two enok mutants are highlighted in purple. Genotypes are as described in S1 Table. (XLSX) S6 Table. H3K23ac levels at piRNA clusters and transposons. Unique mappers from the H3K23ac and H3 ChIP-seq analyses of two replicates of the control (GFP_KD) or the enok knockdown (enok_KD) ovaries mapping to each piRNA cluster or transposon were analyzed using DESeq2. The H3 ChIP-seq data were used as the normalization control for the H3K23ac ChIP-seq data. The results of differential analysis of H3K23ac levels in enok_KD compared with GFP_KD are shown. Transposons or piRNA clusters with decreased H3K23ac levels upon enok knockdown were defined as log2FoldChange < -1 and adjusted p-value < 0.05.