- Open Access
Dnmt1 is essential for egg production and embryo viability in the large milkweed bug, Oncopeltus fasciatus
© The Author(s) 2019
- Received: 20 October 2018
- Accepted: 18 December 2018
- Published: 7 January 2019
The function of cytosine (DNA) methylation in insects remains inconclusive due to a lack of mutant and/or genetic studies.
Here, we provide evidence for the functional role of the maintenance DNA methyltransferase 1 (Dnmt1) in an insect using experimental manipulation. Through RNA interference (RNAi), we successfully posttranscriptionally knocked down Dnmt1 in ovarian tissue of the hemipteran Oncopeltus fasciatus (the large milkweed bug). Individuals depleted for Dnmt1, and subsequently DNA methylation, failed to reproduce. Eggs were inviable and declined in number, and nuclei structure of follicular epithelium was aberrant. Erasure of DNA methylation from gene or transposon element bodies did not reveal a direct causal link to steady-state mRNA levels in somatic cells. These results reveal an important function of Dnmt1 seemingly not contingent on directly controlling gene expression.
This study provides direct experimental evidence for a functional role of Dnmt1 in egg production and embryo viability and uncovers a trivial role, if any, for DNA methylation in control of gene expression in O. fasciatus.
- DNA methylation
DNA methylation in insects has been hypothesized to play numerous functional roles including polyphenism, diapause, longevity, and social behavior and caste differentiation [1–21]. However, these hypotheses are based primarily on correlational studies; for example, there is evidence that the dependence of social behavior and caste differentiation on DNA methylation is not absolute [11–13]. Comparative epigenomic studies have provided valuable insights into patterns of DNA methylation across most insects, but have offered limited insights into its functional significance.
In mammals, DNA methylation typically occurs at CG sites, which is established by DNA methyltransferase 3 (Dnmt3) and subsequently maintained by DNA methyltransferase 1 (Dnmt1). Expression of Dnmt1 and Dnmt3 is spatially and temporally dynamic in Mus musculus (house mouse) and more notably during development where DNA methylation undergoes two rounds of genome-wide erasure and reestablishment . Proper regulation of Dnmt1 during development is necessary to maintain DNA methylation at imprinted genes . Furthermore, mutations causing partial or severe loss of function of Dnmt1 cause abnormal development and embryonic lethality . As in mammals, many insects possess Dnmt1 and Dnmt3, the expression of which are spatially and temporally dynamic with elevated transcript levels in both male and female gonads and throughout the course of development [12, 13, 25–27]. Dnmt1 might be necessary for development and embryonic viability in insects, but the functional contribution of DNA methylation is unknown .
In many insects, DNA methylation is localized to moderately and constitutively expressed genes that are highly conserved between species. This has led some to the hypothesis that DNA methylation functions in transcriptional regulation [6, 15, 28, 29]. However, functional tests of DNA methyltransferases and the consequences to DNA methylation in insects are limited to the posttranscriptional knockdown of Dnmt3 in Apis mellifera, which is associated with modest reductions in DNA methylation, differential gene expression, and alternative splicing . Given the limited genetic and mutant studies, the function of insect DNA methyltransferases, and DNA methylation itself, is obscure.
Here, we show that DNA methylation maintained by Dnmt1 is essential for egg production and embryo viability and is not associated with transcription of the cells profiled in the hemipteran Oncopeltus fasciatus (the large milkweed bug) [31, 32]. Oncopeltus fasciatus has been an established laboratory model organism in fields ranging from physiology and development to evolutionary ecology for > 60 years . Its use as a model in functional genomics is apparent: RNA interference (RNAi) is highly effective across different life history stages (reviewed in ), and the recent release of the genome allows for “-omic” level exploration and hypothesis testing . Posttranscriptional knockdown of Dnmt1 in O. fasciatus led to reduced egg viability, fecundity, and aberrant follicular epithelium, and thus failure to produce a successive generation. Despite finding levels of methylated CG (mCG) within coding regions reduced by 83.55%, we found no evidence for DNA methylation directly affecting transcription in the cells examined. Our results suggest Dnmt1 is required for reproduction in O. fasciatus that is likely mediated by a gene-regulatory independent function. Oncopeltus fasciatus represents a fruitful model species for functional studies of DNA methylation, and continuation of studies in this system will unravel the insect epigenome and its functional consequences.
Dnmt1 is necessary for reproduction
Posttranscriptional knockdown of Dnmt1 affects egg development and viability. In the early stages following injection, females injected with ds-dnmt1-1 did not lay fewer eggs (F = 2.91, df = 1; p = 1.1e−01) (Additional file 1: Fig. S4). Although eggs laid by females injected with ds-dnmt1-1 within the first 8 days post-injection looked typical, they were significantly less likely to develop than the eggs laid by control females (χ2 = 8.470; df = 1; p = 4.0e−03). Furthermore, a mean of 93% of the eggs laid by control females initiated development whereas only a mean of 4% of the eggs laid by ds-dnmt1-1 females initiated development (Fig. 1b). Although eggs laid by the control females that initiated development were viable and hatched, the few eggs laid by ds-dnmt1-1-injected females that initiated development were not viable and failed to hatch.
Posttranscriptional knockdown of Dnmt1 affects egg production following the first 10 days post-injection. By 10–12 days, females injected with ds-dnmt1-1 have mainly stopped laying eggs. While ovaries dissected from control females at this stage have intact eggs in the oviduct, the oviducts of ds-dnmt1-1-injected females are either empty or filled with a mass of what appears to be yolk (Fig. 1c), indicating a fault in production of a functional chorion. Analysis of ovarian cell structure indicates that knockdown of Dnmt1 transcripts affects the cell structure of the follicular epithelium, which is responsible for production of the chorion and vitelline envelope. Nuclei of the follicular epithelium in Dnmt1 posttranscriptional knock-downed females are aberrant and fewer in number than that of the control females (Fig. 1c). Therefore, both through disruption of embryonic development from eggs produced within a week of injection and cessation of egg production, reproduction is compromised following knockdown of Dnmt1 transcripts, preventing a successive generation.
Posttranscriptional knockdown of Dnmt1 successfully and severely reduces mCG in ovaries
High coverage single-base resolution DNA methylomes from ds-dnmt1-1 and control ovaries were generated to understand the impact of the loss of Dnmt1. Greater than 75% reductions of mCG were observed across the genome, consistent with the decrease observed in the low coverage experiments (Additional file 1: Fig. S3). The presence of symmetrical mCG—DNA methylation occurring on both DNA strands at a CpG site—is indicative of the presence of a functional Dnmt1. In control individuals, methylation at CpGs is highly symmetrical and highly methylated. However, even though knockdown of Dnmt1 significantly reduces methylation in ovaries, a minority of CpG sites remain symmetrically methylated (Fig. 2c). This suggests that there are some cells within ovaries that have wild-type methylomes. CG methylation is passively lost in the absence of Dnmt1 upon DNA replication. The lack of a complete loss of mCG is expected due to the injection of fully developed individuals with some cells not undergoing further replication and hence the presence of fully methylated genomes that existed prior to the injection of dsRNA.
Decreased levels of mCG are not associated with transcriptome-wide changes in gene or transposable element expression in ovaries
A maximum of 264 differentially expressed genes (DEG) were observed between O. fasciatus ds-dnmt1-1 and control ovaries (Additional file 1: Table S4), whereas no DEGs are observed between O. fasciatus ds-dnmt1-1 and control for gut, head, and thorax (Additional file 1: Table S8). Of genes with differences in DNA methylation between ds-dnmt1-1 and control ovaries (N = 6590), the majority (N = 6484; 98.39%) had no changes to gene expression (Fig. 3d, Additional file 1: Tables S3 and S4). Furthermore, mCG was always reduced in ds-dnmt1-1 compared to control ovaries in the 6484 gene set. Despite genes being unmethylated (UM) in ds-dnmt1-1 ovaries, changes to gene expression occurred in both directions for the 106 differentially DNA methylated and expressed genes (Additional file 1: Fig. S7). Even for genes that by definition are unmethylated (N = 5982) or similarly CG-methylated (N = 21) in ds-dnmt1-1 and control, a similar proportion (N = 6003; 97.53%) was observed to have no difference in gene expression. The lack of an association between DNA methylation and gene expression is further supported when applying more stringent thresholds to the definition of CG-methylated and unmethylated genes (Additional file 1: Fig. S8). Our observations support a trivial role of DNA methylation in gene expression if any at all.
Comparative epigenomic analysis of the O. fasciatus methylome with other insects
DNA methylation has been implicated to play numerous roles in insects [1–21]; however, its exact function remains uncertain. In this study, we observed that posttranscriptional knockdown of Dnmt1 negatively impacts reproduction through a decrease in egg development and egg production in O. fasciatus.
Posttranscriptional knockdown of Dnmt1 transcripts impacted the maternal somatic gonad and egg maturation. The follicular epithelial cells undergo multiple rounds of mitosis during oogenesis . Oocytes that were in the early stages of maturation would be the most affected, as their follicular epithelium would passively lose DNA methylation with each successive round of mitosis. Loss was not compensated by the de novo DNA methyltransferase Dnmt3 as expression was low in ovaries. Dnmt3 expression during development and across tissues has not been as well explored in hemimetabolous insects compared to holometabolous insects . However, limited Dnmt3 expression in ovaries, which contain developing oocytes and embryos, possibly reflects differences in epigenomic reprogramming and imprinting during oogenesis and early embryo development compared to mammals [37–40]. The loss of function of follicular epithelial cells, which includes production of the chorion and vitelline membrane, could explain the loss of oocyte integrity.
The knockdown of Dnmt1 transcripts in O. fasciatus ovaries subsequently did not lead to a causal link between DNA methylation of genes and their mRNA levels in the cells we examined. No correlation between mCG and gene expression in O. fasciatus ovaries and other insects was observed (Fig. 3, Additional file 1: Fig. S6) [6, 13]. Additionally, knockdown of Dnmt3 in Apis mellifera, although resulting in minor reductions to genic DNA methylation, resulted in no global changes to gene expression (Additional file 1: Fig. S6). A similar lack of relationship between loss of genic CG methylation levels and gene expression has been reported in the plants Arabidopsis thaliana through met1/dnmt1 mutants and Eutrema salsugineum and Conringia planisiliqua through natural, independent losses of gene body DNA methylation [41, 42]. Furthermore, in a recent study, maternal transcriptional knockdown of Dnmt1 in the red flour beetle (Tribolium castaneum)—an insect with Dnmt1, but without DNA methylation—resulted in developmental arrest in offspring embryos . However, the effect on gene expression following dnmt1 knockdown in T. castaneum is unknown . Hence, the importance of mCG in O. fasciatus might lie within the orchestration of genome structure during DNA replication, gametogenesis, or cell types not examined in our study and/or Dnmt1 has an alternative function unrelated to the maintenance of CG methylation. We contribute to a growing body of evidence for a potentially alternative, non-catalytic function of Dnmt1 [43, 44].
DNA methylation might be important for the transcriptional regulation of 106 genes that are differentially expressed when their DNA methylation is reduced (Fig. 3). Around 70% of O. fasciatus genes have DNA methylation above the non-conversion rate (Table S3) with the corresponding 106 genes representing < 1% of genes with DNA methylation. Thus, it seems unlikely that the function of DNA methylation is to regulate gene expression in O. fasciatus and possibly other insects. One of these genes (OFAS014839) has lower levels of DNA methylation and expression in ds-dnmt1-1 compared to control and is orthologous to Drosophila melanogaster secretory 23 (Sec23) (NP_730978.1). In D. melanogaster, P-element insertion in the 5′ UTR of sec23 eliminates or severely attenuates gene function and is lethal . Thus, it is possible that one of these 106 genes is responsible for the phenotype, although no obvious candidates stood out on this list, based on the developmental literature around insect oogenesis and early development. Furthermore, the similar phenotype arising in T. castaneum , an insect without DNA methylation, suggests that the phenotype is likely not due to the loss of DNA methylation.
DNA methylation might be required for downstream functions such as proper mitosis, and the segregation of sister chromatids into their respective daughter cells, as there does appear to be aberrant nuclei structure in ds-dnmt1-1 compared to control ovarian cells (Fig. 1c). This phenotype could be the result of epigenomic defects in scaffold/matrix attachment regions (S/MAR), as some correlate with origins of replication . This phenotype could be exacerbated by multiple mitotic divisions of follicular epithelial cells during oogenesis, and/or the presence of holocentric chromosomes . Future work describing the epigenomic contributions to cellular defects are now possible to study in this newly emerging and tractable model species, O. fasciatus.
We suggest that DNA methylation in O. fasciatus and likely other insects is more important for genome structure, integrity or other cellular processes than it is for the regulation of somatic gene expression. Although it is possible that DNA methylation is important for expression control in a rare cell type that was not examined in our study, we instead propose that regulation of expression is not likely the primary role in insects. Investigating the importance of DNA methylation prior to the onset of meiosis and throughout embryo development is a fundamental next step. Oncopeltus fasciatus represents a fruitful model species for functional studies of DNA methylation, and continuation of studies in this system will unravel the insect epigenome and its functional consequences.
A subset of DNA methyltransferase (Dnmt) 1, 2, and 3 sequences were obtained from  for phylogenetic analysis. The subset included only insect species with available MethylC-seq data, and species representatives from the dipteran suborders Brachycera and Nematocera: Acyrthosiphon pisum, Aedes aegypti, Aedes albopictus, Anopheles gambiae, Apis mellifera, Bombyx mori, Camponotus floridanus, Copidosoma floridanum, Culex pipiens quinquefasciatus, Drosophila melanogaster, Harpegnathos saltator, Microplitis demolitor, Nasonia vitripennis, Nicrophorus vespilloides, O. fasciatus, Ooceraea (Cerapachys) biroi, Polistes canadensis, Polistes dominula, Solenopsis invicta, Tribolium castaneum, and Zootermopsis nevadensis. DNA methyltransferases were reassessed in O. fasciatus by using InterProScan v5.23-62.0  to identify annotated proteins with a C-5 cytosine-specific DNA methylase domain (PF00145). Sequence identifiers are located in Additional file 1: Fig. S1. Full-length protein sequences were aligned using PASTA v1.6.4 , and manually trimmed of divergent, non-homologous sequence in Mesquite v3.2 . Full, aligned, and aligned and trimmed sequence alignments, and parenthetical phylogenetic tree are located in Additional file 1: Data S1. Phylogenetic relationship among Dnmt sequences was estimated using BEAST v2.3.2  with a Blosum62 + Γ model of amino acid substitution. A Markov chain Monte Carlo (MCMC) was run until stationarity and convergence was reached (10,000,000 iterations), and a burnin of 1,000,000 was used prior to summarizing the posterior distribution of tree topologies. A consensus tree was generated using TreeAnnotator v2.3.2, visualized in FigTree v1.4.2 (http://tree.bio.ed.ac.uk/software/figtree/) and exported for stylization in Affinity Designer v1.5.1 (https://affinity.serif.com/en-us/).
PCR confirmation for the presence of a single Dnmt1 ortholog in O. fasciatus
To determine whether OFAS015351 and OFAS018396 were two parts of a single Dnmt1 ortholog, we designed one sense primer at the 3′ end of OFAS015351 (Of_DMNT1-1_3603S; Additional file 1: Table S1) and two antisense primers at the 3′ end of OFAS018396 (Of_DNMT1-2_424A and Of_DNMT1-2_465A) (Additional file 1: Table S1). A fourth sense primer (Of_DNMT1-2_1S) (Additional file 1: Table S1) was designed at the 5′ end of OFAS018396 to confirm the size and sequence of this possibly truncated gene annotation.
Polymerase chain reaction (PCR) with primer combinations Of_DMNT1-1_3603S–Of_DNMT1-2_424A, Of_DMNT1-1_3603S–Of_DNMT1-2_465A, and Of_DNMT1-2_1S–Of_DNMT1-2_424A was performed using Q5 Polymerase (New England BioLabs, Ipswich, MA) per manufactures instructions. Thermacycler conditions were 98 °C for 15 s (s) (denaturing), 60 °C for 30 s (annealing), and 72 °C for 30 s (extension) and repeated for 40 cycles. The PCR products were then purified using QIAquick PCR Purification Kit (Qiagen, Venlo, The Netherlands) and sequenced at the Georgia Genomics and Bioinformatics Core (Athens, GA).
Oncopeltus fasciatus cultures were originally purchased from Carolina Biologicals (Burlington, NC). Mass colonies were maintained in incubators under a 12 h/12 h light/dark cycle at 27 °C. Colonies and individual experimental animals were fed organic raw sunflower seeds and provided with ad libitum deionized water. Late instar nymphs were separated from the mass colonies and housed under the same conditions. Nymph colonies were checked daily for newly emerged adults. Adults were separated by sex and kept with food and water for 7–10 days until females reached sexual maturity.
Template for the in vitro transcription of reactions was prepared from a PCR reaction in which T7 phage promoter sequences were added to the gene-specific Dnmt1 primers . For our control sequence, we used the red fluorescent protein (Red) sequence used in previous parental RNAi experiments in O. fasciatus  or buffer. Primer sequences can be found in Table S2. Sense and antisense RNA was synthesized in a single reaction using the Ambion MEGAscript kit (ThermoFisher Sci, Waltham, MA). After purification, the double-stranded RNA (dsRNA) concentration was adjusted to 2 μg/μL in injection buffer (5 mM KCl, 0.1 mM NaH2PO4) . Females were injected with 5 μL of dsRNA between the abdominal sternites using an insulin syringe. Following injection, females were paired with an uninjected male to stimulate oviposition and fertilize eggs. Individual females with their mate were housed in petri dishes with sunflower seeds, water and cotton wool as an oviposition site. The parental RNAi protocol has been reported to result in 100% penetrance by the third clutch of eggs , and this was also our experience.
To test for any potential off-target effect of our RNAi, we replicated the experiment using a second Dnmt1-specific dsRNA. The original dsRNA was targeted to the cytosine-specific DNA methyltransferase replication foci domain (RFD) of Dnmt1 (ds-dnmt1-1). We developed a second dsRNA that was homologous to the DNA methylase domain (AdoMet) of Dnmt1 (ds-dnmt1-2).
Reproductive phenotype screening and analysis
Eggs were collected between days 4–10 post-injection and assessed for development. Oncopeltus fasciatus embryos change from a creamy white color to orange as they develop, which indicates viability. Thus, color change is a useful tool for assessing healthy development. We examined the number of developing eggs at 5 days post-oviposition, at which point viable eggs are clearly distinguishable from inviable eggs, as well as hatching rate of eggs from females injected with double-stranded Dnmt1 (ds-dnmt1-1) and double-stranded Red (ds-red). The data for the number of eggs laid were normally distributed, so differences among RNAi-treated groups were tested using analysis of variance (ANOVA). The development data, however, were not normally distributed and consisted of binary states (developed and not developed), and so differences among treatment groups were analyzed with a generalized linear model (GLM) using a Poisson distribution.
A second set of O. fasciatus females were injected in the same manner as described in “Parental RNAi” to assess posttranscriptional knockdown of Dnmt1 on ovarian structure. Oncopeltus fasciatus females were dissected 10 days after injection. By 10 days post-injection, O. fasciatus females are beginning to stop laying recognizable eggs. Ovaries were removed from O. fasciatus females and placed in 1 × PBS. Whole ovaries were imaged with a Leica DFC295 stereomicroscope using Leica Application Suite morphometric software (LAS V4.1; Leica, Wetzlar, Germany). Dissected ovaries were fixed within 15 min of dissection in 4% formaldehyde in 1 × PBS for 25 min. Fixed ovaries were stained with Hoechst 33342 (Sigma-Aldrich) at 0.5 µg/ml. The stained ovarioles were imaged using a Zeiss LSM 710 Confocal Microscope (Zeiss) at the University of Georgia Biomedical Microscopy Core.
To assess the effectiveness of posttranscriptional knockdown of Dnmt1, females were dissected 11 days post-injection. Ovaries were removed from each female, flash frozen in liquid nitrogen, and stored at − 80 °C until processing. Total RNA (and DNA) was extracted from a single ovary per female using a Qiagen Allprep DNA/RNA Mini Kit (Qiagen, Venlo, The Netherlands) per manufacturer’s instructions. Complementary DNA (cDNA) was synthesized from 500 ng RNA with qScript cDNA SuperMix (Quanta Biosciences, Gaithersburg, MD).
Expression level of Dnmt1 was quantified by quantitative real-time PCR (qRT-PCR). Primers were designed for Dnmt1 using the O. fasciatus genome as a Ref. . Actin and GAPDH were used as endogenous reference genes. Primer sequences can be found in Additional file 1: Table S2. We used Roche LightCycler 480 SYBR Green Master Mix with a Roche LightCycler 480 (Roche Applied Science, Indianapolis, IN) for qRT-PCR. All samples were run with 3 technical replicates using 10 μL reactions using the manufacturer’s recommended protocol. Primer efficiency calculations, genomic contamination testing, and endogenous control gene selection were performed as described by . We used the ΔΔCT method  to examine differences in expression between control and ds-dnmt1-1 injected females.
Whole-genome bisulfite sequencing (WGBS) and analysis of cytosine (DNA) methylation
MethylC-seq libraries for an O. fasciatus ds-dnmt1-1 and control individual were prepared according to the protocol described in  using genomic DNA extracted from ovaries (see “Materials and methods” section “Quantitative RT-PCR”). Libraries were single-end 75 bp sequenced on an Illumina NextSeq 500 machine [56, 57]. Unmethylated lambda phage DNA was used to as a control for sodium bisulfite conversion, and an error rate of ~ 0.05% was estimated. Oncopeltus fasciatus ds-dnmt1-1 and control were sequenced to a depth of ~ 18 × and ~ 21 × , which corresponded to an actual mapped coverage of ~ 9 × and ~ 11 × , respectively. Additionally, low pass (< 1 ×) WGBS from gut, head, ovary, and thorax was performed for three control and ds-dnmt1-1 biological replicates. Blattella germanica was additionally sequenced to generate equal numbers of hemi- and holometabolous insects investigated in this study. However, DNA was extracted from whole-body minus gastrointestinal tract. MethylC-seq libraries were prepared and sequenced identically to O. fasciatus. An error rate of ~ 0.14% was estimated from unmethylated lambda phage DNA. Blattella germanica was sequenced to a depth of ~ 8 × , which corresponded to an actual mapped coverage of ~ 5 × . WGBS data for O. fasciatus and B. germanica can be found on Gene Expression Omnibus (GEO) under accession GSE109199. Previously published WGBS data for A. mellifera [30, 58], Bo. mori , Daphnia pulex , N. vespilloides , and Z. nevadensis  were downloaded from the Short Read Archive (SRA) using accessions SRR1270128, SRR1270129, SRR445803-4, SRR027157-9, SRR1552830, SRR2017555, and SRR3139749. Thus, DNA methylation was investigated for six insects from six different orders spread evenly across developmental groups, and a crustacean outgroup. Additionally, DNA methylation (and gene expression, see “Materials and methods” section “RNA-seq and differential expression analysis”) changes following posttranscriptional knockdown of Dnmt3 in A. mellifera were investigated . WGBS data were aligned to each species respective genome assembly using the methylpy pipeline . In brief, reads were trimmed of sequencing adapters using Cutadapt v1.9  and then mapped to both a converted forward strand (cytosines to thymines) and converted reverse strand (guanines to adenines) using Bowtie v1.1.1 . Reads that mapped to multiple locations, and clonal reads were removed.
Weighted DNA methylation was calculated for CG sites by dividing the total number of aligned methylated reads by the total number of methylated plus unmethylated reads . For genic metaplots, the gene body (start to stop codon), 1000 base pairs (bp) upstream, and 1000 bp downstream were divided into 20 windows proportional windows based on sequence length (bp). Weighted DNA methylation was calculated for each window and then plotted in R v3.2.4 (https://www.r-project.org/). CG sequence context enrichment for each gene was determined through a binomial test followed by Benjamini–Hochberg false discovery rate [64, 65]. A background mCG level was determined from all coding sequences, which was used as a threshold in determining significance with a false discovery rate (FDR) correction. Genes were classified as CG-methylated if they had reads mapping to at least 20 reads mapping to 20 CG sites and a q < 0.05. Using a binomial test can lead to false negatives—highly CG-methylated genes that are classified as unmethylated (UM)—due to a low number of statistically CG-methylated sites (Additional file 1: Table S4). Genes classified as unmethylated, but had a mCG level greater than the lowest CG-methylated gene were dropped from future analyses.
Putative ortholog identification
Best BLASTp hit (arguments: –max_hsps 1 –max_target_seqs 1 evalue 1e−03) was used to identify orthologs between O. fasciatus and other insect species investigated.
RNA-seq and differential expression analysis
RNA-seq libraries for RNA extracted from ovaries of three biological O. fasciatus ds-dnmt1-1 and control replicates at 11 days post-injection were constructed using Illumina TruSeq Stranded RNA LT Kit (Illumina, San Diego, CA) following the manufacturer’s instructions with limited modifications. RNA from ovaries of an additional three biological O. fasciatus ds-dnmt1-1 and control replicates, and three biological O. fasciatus ds-dnmt1-1 and control replicates from gut, head, and thorax were extracted. The starting quantity of total RNA was adjusted to 1.3 µg, and all volumes were reduced to a third of the described quantity. Libraries were single-end 75 bp sequenced on an Illumina NextSeq 500 machine. RNA-seq data for O. fasciatus ds-dnmt1-1 and control can be found on GEO under accession GSE109199. Previously published RNA-seq data for A. mellifera  and Z. nevadensis  were downloaded from the SRA using accessions SRR893029, SRR893030, SRR893035, SRR893036, SRR2954345, and SRR3139740.
Raw RNA-seq FASTQ reads were trimmed for adapters and preprocessed to remove low-quality reads using Trimmomatic v0.33 (arguments: LEADING:10 TRAILING:10 MINLEN:30)  prior to mapping to the O. fasciatus v1.1 reference genome assembly. Reads were mapped using TopHat v2.1.1  supplied with a reference General Features File (GFF) to the O. fasciatus v1.1 reference genome assembly , and with the following arguments: -I 20000 –library-type fr-firststrand –b2-very-sensitive.
Differentially expressed genes (DEGs) between ds-dnmt1-1 and control libraries were determined using edgeR v3.20.1  implemented in R v3.2.4 (https://www.r-project.org/). Genes were retained for DEG analysis if they possessed a counts per million (CPM) ≥ 1 in at least ≥ 2 libraries. Significance was determined using the glmQLFTest function, which uses empirical Bayes quasi-likelihood F tests. Parameter settings were determined following best practices for DEG analysis as described by . Gene expression metrics for A. mellifera, Z. nevadensis, and O. fasciatus ds-dnmt1-1 and control are located in Additional file 1: Tables S4–S8, respectively.
Gene ontology (GO) annotation and enrichment
GO terms were assigned to O. fasciatus v1.1 gene set  through combining annotations from Blast2Go PRO v4.1.9 , InterProScan v5.23-62.0 (arguments: -goterms -iprlookup -appl CDD,Pfam) , and through sequence homology to D. melanogaster using BLASTp (arguments: -evalue 1.0e−03 -max_target_seqs 1 -max_hsps 1). 11,105/19,615 gene models were associated with at least one GO term, and a total of 8190 distinct GO identifiers were mapped. GO terms are found in Table S12. Enriched GO terms in gene groups were evaluated using topGO v2.30.0  implemented in R v3.2.4 (https://www.r-project.org/), and significance (p < 0.05) of terms was assessed using Fisher’s exact test with a weighted algorithm. Gene groups were contrasted to all O. fasciatus genes associated with GO terms (Additional file 1: Table S13).
Transposable element (TE) annotation and expression
TEs were identified using RepeatMasker v4.0.5 (http://www.repeatmasker.org) provided with the invertebrate repeat library from Repbase (http://www.girinst.org/repbase/) (arguments: -lib < Repbase invertebrate library > -no_is -engine wublast -a -inv -x -gff. Following RepeatMasker, neighboring TEs of the same type were collapsed into a single locus within the outputted GFF. The unmodified GFF is located in Additional file 1: Table S10.
To quantify expression from TEs, RNA-seq libraries from O. fasciatus ds-dnmt1-1 and control were independently combined and mapped to the O. fasciatus v1.1 reference genome assembly  using Bowtie2 v2.2.9  with the following arguments: –sensitive. Mapped reads overlapping with the top ten most abundant TEs of ≥ 500 bp in length were identified using the intersect command in BEDTools suite v 2.26.0 . TE expression is quantified as reads per kilobase per million mapped reads (RPKM) for each intersected TE type by counting the number reads and dividing by the mapped library read number in millions. Significance in expression of TEs between ds-dnmt1-1 and control tissues was assessed using the Mann–Whitney test with the alternative hypothesis set to “greater” in R v3.2.4 (https://www.r-project.org/).
AJB, AJM, PJM, and RJS designed the study. PJM, ZS, and ECM maintained animals and performed RNAi experiments and microscopy analyses. WGBS and RNA-seq data were generated by RJS. AJB analyzed and interpreted the WGBS and RNA-seq data and wrote the manuscript. All authors read and approved the final manuscript.
We thank Sam Arsenault (University of Georgia [UGA]), and Mariana Monteiro and Stefan Götz (Blast2GO Team) for assistance with Blast2GO analyses, Nick Rohr and Tina Ethridge for MethylC-Seq and RNA-seq library preparation, Brigitte Hofmeister (UGA) for bioinformatics advice, Mary Goll (UGA) for comments, Kelly Dawe (UGA) for feedback on cellular phenotypes, Muthugapatti Kandasamy of the Biomedical Microscopy Core (UGA) for technical assistance with confocal microscopy, and Cassandra Extavour (Harvard) for providing technical advice and Red plasmids for RNAi. We also thank the Georgia Advanced Computing Resource Center (GACRC), Georgia Genomics and Bioinformatics Core (GGBC) and the Office of Research at the University of Georgia.
The authors declare that they have no competing interests.
Availability of data and materials
The datasets generated and/or analyzed during the current study are available in the Gene Expression Omnibus (GEO) repository under accession GSE109199.
Consent for publication
Ethics approval and consent to participate
This work was funded by the National Science Foundation (NSF) IOS-1354358 (AJM) and UGA CAES Undergraduate Research Award (ZS). RJS is a Pew Scholar in the Biomedical Sciences, supported by The Pew Charitable Trusts.
Springer Nature remains neutral with regard to jurisdictional claims in published maps and institutional affiliations.
Open AccessThis article is distributed under the terms of the Creative Commons Attribution 4.0 International License (http://creativecommons.org/licenses/by/4.0/), which permits unrestricted use, distribution, and reproduction in any medium, provided you give appropriate credit to the original author(s) and the source, provide a link to the Creative Commons license, and indicate if changes were made. The Creative Commons Public Domain Dedication waiver (http://creativecommons.org/publicdomain/zero/1.0/) applies to the data made available in this article, unless otherwise stated.
- Falckenhayn C, Boerjan B, Raddatz G, Frohme M, Schoofs L, Lyko F. Characterization of genome methylation patterns in the desert locust Schistocerca gregaria. J Exp Biol. 2013;216:1423–9.View ArticleGoogle Scholar
- Cullen DA, Cease AJ, Latchininsky AV, Ayali A, Berry K, Buhl J, De Keyser R, Foquet B, Hadrich JC, Matheson T, et al. From molecules to management: mechanisms and consequences of locust phase polyphenism. Adv Insect Phys. 2017;53:167–285.View ArticleGoogle Scholar
- Pegoraro M, Bafna A, Davies NJ, Shuker DM, Tauber E. DNA methylation changes induced by long and short photoperiods in Nasonia. Genome Res. 2016;26:203–10.View ArticleGoogle Scholar
- Reynolds JA. Epigenetic influences on diapause. Adv Insect Phys. 2017;53:115–44.View ArticleGoogle Scholar
- Kozeretska IA, Serga SV, Koliada AK, Vaiserman AM. Epigenetic regulation of longevity in insects. Adv Insect Phys. 2017;53:87–114.View ArticleGoogle Scholar
- Bonasio R, Li Q, Lian J, Mutti NS, Jin L, Zhao H, Zhang P, Wen P, Xiang H, Ding Y, et al. Genome-wide and caste-specific DNA methylomes of the ants Camponotus floridanus and Harpegnathos saltator. Curr Biol. 2012;22:1755–64.View ArticleGoogle Scholar
- Foret S, Kucharski R, Pellegrini M, Feng S, Jacobsen SE, Robinson GE, Maleszka R. DNA methylation dynamics, metabolic fluxes, gene splicing, and alternative phenotypes in honey bees. Proc Natl Acad Sci USA. 2012;109:4968–73.View ArticleGoogle Scholar
- Lockett GA, Kucharski R, Maleszka R. DNA methylation changes elicited by social stimuli in the brains of worker honey bees. Genes Brain Behav. 2012;11:235–42.View ArticleGoogle Scholar
- Ford D. Honeybees and cell lines as models of DNA methylation and aging in response to diet. Exp Gerontol. 2013;48:614–9.View ArticleGoogle Scholar
- He XJ, Zhou LB, Pan QZ, Barron AB, Yan WY, Zeng ZJ. Making a queen: an epigenetic analysis of the robustness of the honeybee (Apis mellifera) queen developmental pathway. Mol Ecol. 2017;26:1598–607.View ArticleGoogle Scholar
- Libbrecht R, Oxley PR, Keller L, Kronauer DJ. Robust DNA methylation in the clonal raider ant brain. Curr Biol. 2016;26:391–5.View ArticleGoogle Scholar
- Bewick AJ, Vogel KJ, Moore AJ, Schmitz RJ. Evolution of DNA methylation across insects. Mol Biol Evol. 2016;34:654–65.PubMed CentralGoogle Scholar
- Glastad KM, Arsenault SV, Vertacnik KL, Geib SM, Kay S, Danforth BN, Rehan SM, Linnen CR, Kocher SD, Hunt BG. Variation in DNA methylation is not consistently reflected by sociality in Hymenoptera. Genome Biol Evol. 2017;9:1687–98.View ArticleGoogle Scholar
- Provataris P, Meusemann K, Niehuis O, Grath S, Misof B. Signatures of DNA methylation across insects suggest reduced DNA methylation levels in Holometabola. Genome Biol Evol. 2018;10:1185–97.View ArticleGoogle Scholar
- Glastad KM, Gokhale K, Liebig J, Goodisman MA. The caste- and sex-specific DNA methylome of the termite Zootermopsis nevadensis. Sci Rep. 2016;6:37110.View ArticleGoogle Scholar
- Lyko F, Foret S, Kucharski R, Wolf S, Falckenhayn C, Maleszka R. The honey bee epigenomes: differential methylation of brain DNA in queens and workers. PLoS Biol. 2010;8:e1000506.View ArticleGoogle Scholar
- Xiang H, Zhu J, Chen Q, Dai F, Li X, Li M, Zhang H, Zhang G, Li D, Dong Y, et al. Single base-resolution methylome of the silkworm reveals a sparse epigenomic map. Nat Biotechnol. 2010;28:516–20.View ArticleGoogle Scholar
- Wang X, Wheeler D, Avery A, Rago A, Choi JH, Colbourne JK, Clark AG, Werren JH. Function and evolution of DNA methylation in Nasonia vitripennis. PLoS Genet. 2013;9:e1003872.View ArticleGoogle Scholar
- Cunningham CB, Ji L, Wiberg RA, Shelton J, McKinney EC, Parker DJ, Meagher RB, Benowitz KM, Roy-Zokan EM, Ritchie MG, et al. The genome and methylome of a beetle with complex social behavior, Nicrophorus vespilloides (Coleoptera: Silphidae). Genome Biol Evol. 2015;7:3383–96.View ArticleGoogle Scholar
- Patalano S, Vlasova A, Wyatt C, Ewels P, Camara F, Ferreira PG, Asher CL, Jurkowski TP, Segonds-Pichon A, Bachman M, et al. Molecular signatures of plastic phenotypes in two eusocial insect species with simple societies. Proc Natl Acad Sci USA. 2015;112:13970–5.View ArticleGoogle Scholar
- Rehan SM, Glastad KM, Lawson SP, Hunt BG. The genome and methylome of a subsocial small carpenter bee, Ceratina calcarata. Genome Biol Evol. 2016;8:1401–10.View ArticleGoogle Scholar
- Feng S, Jacobsen SE, Reik W. Epigenetic reprogramming in plant and animal development. Science. 2010;330:622–7.View ArticleGoogle Scholar
- Hirasawa R, Chiba H, Kaneda M, Tajima S, Li E, Jaenisch R, Sasaki H. Maternal and zygotic Dnmt1 are necessary and sufficient for the maintenance of DNA methylation imprints during preimplantation development. Genes Dev. 2008;22:1607–16.View ArticleGoogle Scholar
- Li E, Bestor TH, Jaenisch R. Targeted mutation of the DNA methyltransferase gene results in embryonic lethality. Cell. 1992;69:915–26.View ArticleGoogle Scholar
- Zwier MV, Verhulst EC, Zwahlen RD, Beukeboom LW, van de Zande L. DNA methylation plays a crucial role during early Nasonia development. Insect Mol Biol. 2012;21:129–38.View ArticleGoogle Scholar
- Kay S, Skowronski D, Hunt BG. Developmental DNA methyltransferase expression in the fire ant Solenopsis invicta. Insect Sci. 2018;25:57–65.View ArticleGoogle Scholar
- Robinson KL, Tohidi-Esfahani D, Ponton F, Simpson SJ, Sword GA, Lo N. Alternative migratory locust phenotypes are associated with differences in the expression of genes encoding the methylation machinery. Insect Mol Biol. 2016;25:105–15.View ArticleGoogle Scholar
- Sarda S, Zeng J, Hunt BG, Yi SV. The evolution of invertebrate gene body methylation. Mol Biol Evol. 2012;29:1907–16.View ArticleGoogle Scholar
- Hunt BG, Glastad KM, Yi SV, Goodisman MA. Patterning and regulatory associations of DNA methylation are mirrored by histone modifications in insects. Genome Biol Evol. 2013;5:591–8.View ArticleGoogle Scholar
- Li-Byarlay H, Li Y, Stroud H, Feng S, Newman TC, Kaneda M, Hou KK, Worley KC, Elsik CG, Wickline SA, et al. RNA interference knockdown of DNA methyl-transferase 3 affects gene alternative splicing in the honey bee. Proc Natl Acad Sci USA. 2013;110:12750–5.View ArticleGoogle Scholar
- Chipman AD. Oncopeltus fasciatus as an evo-devo research organism. Genesis. 2017;55:e23020.View ArticleGoogle Scholar
- Panfilio KA, Vargas Jentzsch IM, Benoit JB, Erezyilmaz D, Suzuki Y, Colella S, Robertson HM, Poelchau MF, Waterhouse RM, Ioannidis P, et al. Molecular evolutionary trends and feeding ecology diversification in the Hemiptera, anchored by the milkweed bug genome; 2017. Preprint. https://doi.org/10.1101/201731. Cited 26 June 2018.
- Lawrence PA. Mitosis and the cell cycle in the metamorphic moult of the milkweed bug Oncopeltus fasciatus; a radioautographic study. J Cell Sci. 1968;3:391–404.Google Scholar
- Liu PZ, Kaufman TC. Hunchback is required for suppression of abdominal identity, and for proper germband growth and segmentation in the intermediate germband insect Oncopeltus fasciatus. Development. 2004;131:1515–27.View ArticleGoogle Scholar
- Bewick AJ, Hofmeister BT, Lee K, Zhang X, Hall DW, Schmitz RJ. FASTmC: a suite of predictive models for nonreference-based estimations of DNA methylation. G3 (Bethesda). 2015;6:447–752.View ArticleGoogle Scholar
- Bonhag PF, Wick JR. The functional anatomy of the male and female reproductive systems of the milkweed bug, Oncopeltus fasciatus (Dallas) (Heteroptera: Lygaeidae). J Morphol. 1953;93:177–283.View ArticleGoogle Scholar
- Obata Y, Kono T. Maternal primary imprinting is established at a specific time for each gene throughout oocyte growth. J Biol Chem. 2002;277:5285–9.View ArticleGoogle Scholar
- Kaneda M, Okano M, Hata K, Sado T, Tsujimoto N, Li E, Sasaki H. Essential role for de novo DNA methyltransferase Dnmt3a in paternal and maternal imprinting. Nature. 2004;429:900–3.View ArticleGoogle Scholar
- Lucifero D, Mann MRW, Bartolomei MS, Trasler JM. Gene-specific timing and epigenetic memory in oocyte imprinting. Hum Mol Genet. 2004;13:839–49.View ArticleGoogle Scholar
- Skvortsova K, Iovino N, Bogdanović O. Functions and mechanisms of epigenetic inheritance in animals. Nat Rev Mol Cell Biol. 2018. https://doi.org/10.1038/s41580-018-0074-2.View ArticlePubMedGoogle Scholar
- Zhang X, Yazaki J, Sundaresan A, Cokus S, Chan SW, Chen H, Henderson IR, Shinn P, Pellegrini M, Jacobsen SE, et al. Genome-wide high-resolution mapping and functional analysis of DNA methylation in Arabidopsis. Cell. 2006;126:1189–201.View ArticleGoogle Scholar
- Bewick AJ, Ji L, Niederhuth CE, Willing EM, Hofmeister BT, Shi X, Wang L, Lu Z, Rohr NA, Hartwig B, et al. On the origin and evolutionary consequences of gene body DNA methylation. Proc Natl Acad Sci U S A. 2016;113:9111–6.View ArticleGoogle Scholar
- Schulz NKE, Wagner CI, Ebeling J, Raddatz G, Buhr MFD, Lyko F, Kurtz J. Dnmt1 has an essential function despite the absence of CpG DNA methylation in the red flour beetle Tribolium castaneum. Sci Rep. 2018;8:16462.View ArticleGoogle Scholar
- Dunican DS, Ruzov A, Hackett JA, Meehan RR. xDnmt1 regulates transcriptional silencing in pre-MBT Xenopus embryos independently of its catalytic function. Development. 2008;135:1295–302.View ArticleGoogle Scholar
- González-Morales N, Mendoza-Ortíz MÁ, Blowes LM, Missirlis F, Riesgo-Escovar JR. Ferritin is required in multiple tissues during Drosophila melanogaster development. PLoS ONE. 2015;10:e0133499.View ArticleGoogle Scholar
- Vaughn JP, Dijkwel PA, Mullenders LH, Hamlin JL. Replication forks are associated with the nuclear matrix. Nucleic Acids Res. 1990;18:1965–9.View ArticleGoogle Scholar
- Drinnenberg IA, de Young D, Henikoff S, Malik HS. Recurrent loss of CenH3 is associated with independent transitions to holocentricity in insects. eLife. 2014. https://doi.org/10.7554/elife.03676.View ArticlePubMedPubMed CentralGoogle Scholar
- Jones P, Binns D, Chang HY, Fraser M, Li W, McAnulla C, McWilliam H, Maslen J, Mitchell A, Nuka G, et al. InterProScan 5: genome-scale protein function classification. Bioinformatics. 2014;30:1236–40.View ArticleGoogle Scholar
- Mirarab S, Nguyen N, Guo S, Wang LS, Kim J, Warnow T. PASTA: ultra-large multiple sequence alignment for nucleotide and amino-acid sequences. J Comput Biol. 2015;22:377–86.View ArticleGoogle Scholar
- Maddison WP, Maddison DR. Mesquite: a modular system for evolutionary analysis. Version 3.31. 2017. http://mesquiteproject.org.
- Bouckaert R, Heled J, Kühnert D, Vaughan T, Wu CH, Xie D, Suchard MA, Rambaut A, Drummond AJ. BEAST 2: a software platform for Bayesian evolutionary analysis. PLoS Comput Biol. 2014;10:e1003537.View ArticleGoogle Scholar
- Ewen-Campen B, Jones TE, Extavour CG. Evidence against a germplasm in the milkweed bug Oncopeltus fasciatus, a hemimetabolous insect. Biol Open. 2013;2:556–68.View ArticleGoogle Scholar
- Cunningham CB, Douthit MK, Moore AJ. Octopaminergic gene expression and flexible social behaviour in the subsocial burying beetle Nicrophorus vespilloides. Insect Mol Biol. 2014;23:391–404.PubMedPubMed CentralGoogle Scholar
- Livak KJ, Schmittgen TD. Analysis of relative gene expression data using real-time quantitative PCR and the 2(-Delta Delta C(T)) Method. Methods. 2001;25:402–8.View ArticleGoogle Scholar
- Urich MA, Nery JR, Lister R, Schmitz RJ, Ecker JR. MethylC-seq library preparation for base-resolution whole-genome bisulfite sequencing. Nat Protoc. 2015;10:475–83.View ArticleGoogle Scholar
- Cokus SJ, Feng S, Zhang X, Chen Z, Merriman B, Haudenschild CD, Pradhan S, Nelson SF, Pellegrini M, Jacobsen SE. Shotgun bisulphite sequencing of the Arabidopsis genome reveals DNA methylation patterning. Nature. 2008;452:215–9.View ArticleGoogle Scholar
- Lister R, O’Malley RC, Tonti-Filippini J, Gregory BD, Berry CC, Millar AH, Ecker JR. Highly integrated single-base resolution maps of the epigenome in Arabidopsis. Cell. 2008;133:523–36.View ArticleGoogle Scholar
- Herb BR, Wolschin F, Hansen KD, Aryee MJ, Langmead B, Irizarry R, Amdam GV, Feinberg AP. Reversible switching between epigenetic states in honeybee behavioral subcastes. Nat Neurosci. 2012;15:1371–3.View ArticleGoogle Scholar
- Asselman J, De Coninck DI, Pfrender ME, De Schamphelaere KA. Gene body methylation patterns in Daphnia are associated with gene family size. Genome Biol Evol. 2016;8:1185–96.View ArticleGoogle Scholar
- Schultz MD, He Y, Whitaker JW, Hariharan M, Mukamel EA, Leung D, Rajagopal N, Nery JR, Urich MA, Chen H, et al. Human body epigenome maps reveal noncanonical DNA methylation variation. Nature. 2015;523:212–6.View ArticleGoogle Scholar
- Martin M, Marcel M. Cutadapt removes adapter sequences from high-throughput sequencing reads. EMBnet J. 2011;17:10–2.View ArticleGoogle Scholar
- Langmead B, Trapnell C, Salzberg SL. Ultrafast and memory-efficient alignment of short DNA sequences to the human genome. Genome Biol. 2009;10:R25.View ArticleGoogle Scholar
- Schultz MD, Schmitz RJ, Ecker JR. ‘Leveling’ the playing field for analyses of single-base resolution DNA methylomes. Trends Genet. 2012;28:583–5.View ArticleGoogle Scholar
- Takuno S, Gaut BS. Body-methylated genes in Arabidopsis thaliana are functionally important and evolve slowly. Mol Biol Evol. 2012;1:219–27.View ArticleGoogle Scholar
- Niederhuth CE, Bewick AJ, Ji L, Alabady MS, Kim KD, Li Q, Rohr NA, Rambani A, Burke JM, Udall JA, et al. Widespread natural variation of DNA methylation within angiosperms. Genome Biol. 2016;17:194.View ArticleGoogle Scholar
- Vleurinck C, Raub S, Sturgill D, Oliver B, Beye M. Linking genes and brain development of honeybee workers: a whole-transcriptome approach. PLoS ONE. 2016;11:e0157980.View ArticleGoogle Scholar
- Bolger AM, Lohse M, Usadel B. Trimmomatic: a flexible trimmer for Illumina sequence data. Bioinformatics. 2014;30:2114–20.View ArticleGoogle Scholar
- Kim D, Pertea G, Trapnell C, Pimentel H, Kelley R, Salzberg SL. TopHat2: accurate alignment of transcriptomes in the presence of insertions, deletions and gene fusions. Genome Biol. 2013;14:R36.View ArticleGoogle Scholar
- Robinson MD, McCarthy DJ, Smyth GK. edgeR: a bioconductor package for differential expression analysis of digital gene expression data. Bioinformatics. 2010;26:139–40.View ArticleGoogle Scholar
- Chen Y, Lun ATL, Smyth GK. From reads to genes to pathways: differential expression analysis of RNA-Seq experiments using Rsubread and the edgeR quasi-likelihood pipeline. F1000Res. 2016;5:1438.PubMedPubMed CentralGoogle Scholar
- Götz S, García-Gómez JM, Terol J, Williams TD, Nagaraj SH, Nueda MJ, Robles M, Talón M, Dopazo J, Conesa A. High-throughput functional annotation and data mining with the Blast2GO suite. Nucleic Acids Res. 2008;36:3420–35.View ArticleGoogle Scholar
- Alexa A, Rahnenfuhrer J. topGO: enrichment analysis for gene ontology. R package version 2.30.0. 2016.Google Scholar
- Langmead B, Salzberg S. Fast gapped-read alignment with Bowtie 2. Nat Methods. 2012;9:357–9.View ArticleGoogle Scholar
- Quinlan AR, Hall IM. BEDTools: a flexible suite of utilities for comparing genomic features. Bioinformatics. 2010;26:841–2.View ArticleGoogle Scholar