Identification and functional analysis of five genes that encode distinct isoforms of protein phosphatase 1 in Nilaparvata lugens

Ten distinct cDNAs encoding five different protein phosphatases 1 (PPP1) were cloned from Nilaparvata lugens. NlPPP1α and NlPPP1β are highly conserved whereas NlPPP1-Y, NlPPP1-Y1 and NlPPP1-Y2 are lowly conserved among insects. NlPPP1α and NlPPP1β exhibited a ubiquitous expression, while NlPPP1-Y, NlPPP1-Y1, and NlPPP1-Y2 were obviously detected from the 4th instar nymph to imago developmental stages in males, especially detected in internal reproductive organ and fat bodies of the male. Injection nymphs with dsRNA of NlPPP1α or NlPPP1β was able to reduce the target gene expression in a range of 71.5–91.0%, inducing a maximum mortality rate of 95.2% or 97.2% at 10th day after injection and eclosion ratio down by 65.5–100.0%. Injection with dsNlPPP1Ys targeted to NlPPP1-Y, NlPPP1-Y1and NlPPP1-Y2 was able to induce a maximum mortality rate of 95.5% at 10th day after injection, eclosion ratio down by 86.4%. Knock-down one of the male-biased NlPPP1 genes has no effect on survival and eclosion ratio. Injection of 4th instar nymph with dsNlPPP1Ys led to reduced oviposition amount and hatchability, down by 44.7% and 19.6% respectively. Knock-down of NlPPP1-Y1 or NlPPP1-Y2 gene did not significantly affect oviposition amount but significantly affected hatchability. The results indicate that the male-biased NlPPP1 genes have overlapping functions in N. lugens development, and NlPPP1-Y1 and NlPPP1-Y2 may play important roles in spermatogenesis and fertilization. The dsNlPPP1β and dsNlPPP1Ys in this study could be the preferred sequence in RNAi and low-conserved male-biased NlPPP1 genes could be potential target for N. lugens control.

RNA isolation, sequence amplification and analysis. Total RNA was isolated from N. lugens at different developmental stages and from different tissues using an RNeasy Mini Kit (Qiagen, Hilden, Germany) according to the manufacturer's instructions. The RNA was quantified and the quality verified by NanoDrop 2000 spectrophotometer (Thermo Fischer Scientific, Bremen, Germany). A total of 500 ng RNA was used for reverse transcription in a 10 μL reaction with the ReverTra Ace qPCR RT Master Mix with gDNA Remover Kit (ToYoBo, Osaka, Japan). Synthesized cDNA was diluted tenfold and used as template for quantitative PCR.
The PPP1 gene was amplified with primer pairs which based on our transcriptome database from whole bodies of N. lugens and designed using National Center for Biotechnology Information (NCBI) primer design tool (www.ncbi.nlm. nih.gov/tools/primer-blast). All the primer used in this study were synthesized by Invitrogen Co., Ltd Shanghai China and listed in Table 1. The polymerase chain reaction (PCR) procedure was as follows: 94 °C for 5 min followed by 40 cycles of 94 °C for 30 s, 58 °C for 45 s, and 72 °C for 60 s. The samples were then incubated for 10 min at 72 °C. The PCR products were gel-purified and cloned into the PCR 2.1 Topo vector (Invitrogen, China) and then the plasmids from positive colonies were sequenced with the M13 primer pair on ABI Prism 3100 DNA sequencer (Invitrogen Co., Ltd Shanghai China).The cDNA sequences were analysis with BLAST against N. lugens genome (Nilaparvata lugens (taxid: 108931)).The open reading frame (ORF) was determined using ORF Finder (https ://www.ncbi.nlm.nih.gov/gorf/gorf.html). The translated amino acid sequence was used as a query to identify homologous proteins and compared with other PPP1 deposited in GenBank using the BLASTp tool (https ://blast .ncbi.nlm.nih.gov/Blast .cgi). The molecular weight (Mw) and isoelectric point (pI) of NlPPP1 were calculated by the Compute pI/Mw tool (https ://web.expas y.org/compu te_pi/). The phylogenetic tree of PPP1 was constructed using the maximum likelihood method with MEGA 5.0(https :// megas oftwa re.net/). All the PPP1 sequences from N. lugens were aligned in a multiple sequences alignment using CLUSTAL X and edited with GeneDoc software. The phosphorylation sites were predicted using Netphos3.0 server (https ://www.cbs.dtu.dk). www.nature.com/scientificreports/ Expression analysis by real-time quantitative PCR (RT-qPCR). NlPPP1 transcript levels were quantified in tissues and development stages with specific primers for NlPPP1α, NlPPP1β, and NlPPP1-Ys which were designed based on the cDNA sequence obtained. The RT-qPCR (20 μL per reaction) used 3.0 μL cDNA template, 0.4 μL of each primer (10 mM) and 10 μL SYBR Premix (Toyobo, Japan). RT-qPCRs were carried out using an ABI 7,500 Real-time PCR system (Applied Biosystems, Carlsbad, CA, USA) in a two-step reaction (3 min denaturation at 95 °C, 40 cycles 10 s denaturation at 95 °C, 30 s annealing/extension at 60 °C) followed by a melt curve analysis at the end of the run. Each experiment consisted of 3 separate biological replicates, each of which was comprises 3 technical replicates. Relative expression levels were calculated using the 2 −∆∆CT method 29 . The N. lugens housekeeping genes for β-actin (FJ948574) and 18S rRNA (JN662398) were used as the reference genes 30 . Fold induction values of target genes were calculated with the ΔΔCt equation and normalized to the mRNA level of target genes in control which were defined as 1.0.
RNAi interference and bioassay. The dsRNA synthesis, microinjection experiment and bioassay was conducted based on our previously described method 31 . A PCR method using the plasmid with NlPPP1 insert as template was used to generate the templates for the dsRNA synthesis. The GFP gene (ACY56286) was used as a control. Approximately 70 ng of dsRNA was injected into each newly moulted third-instar nymphs or fourth-instar nymphs (3-day old) N. lugens. The survived nymphs in each treatment were selected and reared on 60-to70 day-old rice variety TN1 in one cage. A total of 175 nymphs (5 replicates, 35 individuals in each replicate) for each treatment were used for each dsRNA injection. RNAi efficiency for dsNlPPP1α and dsNlPPP1β was assessed at 4 days post injection by RT-qPCR. Total RNA was extracted from 5 nymphs sampled from each treatment and each replicate. RNAi efficiency for male-biased gene was assessed at 4 days after emergency by RT-qPCR, internal reproductive organs dissected from approximately 5 males were used in each sample. The survival rates of the injected 3rd-instar nymphs were observed at 24 h intervals with duration of 10 days. Once the 4th-instar nymphs after injection emerged, each female was matched with one male and each pair was allowed to reproduce separately. For the male-biased gene, the injected male was matched with untreated female. In total, 15 single pairs per gene were successfully mated. The number of newly hatched nymphs was recorded every other day until no more nymphs were observed for two successive days. The number of unhatched eggs was also recoded under a light microscope. Eggs were scraped from the leaf sheaths and blades using a pin. All analyses were performed with the data procession system (DPS) of Tang and Feng 32 . Duncan's tests were used to determine differences between the treatment and control. Values of P < 0.05 were considered significantly, all values were expressed as mean ± SEM.  NlPPP1α  CTC ACT TCG TTC GTT CGC ATT  AGA CCT ACT CCA GGT AGC CTT   NlPPP1β  CGT AGA CGT CGG TCT GTG TG  ACC TAA TGC TTT GCT GTT CCTT   NlPPP1-Y  ATT TAC GCG GGT CTT GTG AAC  CGT TCT CGG TCC TTC TTC CTCTT   NlPPP1-Y1  GTT GTA TGC TGG ACA ACA CTG  GAA AAA GCA TCT CTT AAA CCCG   NlPPP1-Y2  TCA ATC TTG AAC CCT TTG TGT GAA  GCG CTT TTT ATC TCG TTC  www.nature.com/scientificreports/ Enzyme-linked immune sorbent assay (ELISA) analysis. Two days after emergence, a total of 40 male internal reproductive organs (4 replicates, 10 individuals in each replicate) were randomly selected and dissected from the treatment group injected with dsNlPPP1Ys and the control group injected with dsGFP. A total of 20 female internal reproductive organs (4 replicates, 5 individuals in each replicate) were randomly selected and dissected from females mated with injected with dsNlPPP1Ys and females mated with injected with dsGFP. The male and female internal reproductive organs were homogenized with a glass tissue grinder in 200 μL 0.8% NaCl buffer. Homogenates were centrifuged at 10,000×g for 15 min at 4 °C, and the supernatants were used for total protein (TP) determination with the Insect TP Elisa Kit (Boshen Biotech Co ltd, Nanjing Jiangsu China). The optical density was read at 450 nm on a Sunrise ELISA reader (Tecan, Mannedorf, Switzerland).

Results
Sequence analysis of NlPPP1. Based on the assembled transcriptome which constructed in our laboratory, five different pair of PCR primers were designed and used to clone the PPP1 gene from N. lugens. Totally 10 cDNA clones were isolated. Based on the comparative study of their nucleotide and deduced amino acid sequences with those reported, these cDNA clones were named NlPPP1α, NlPPP1β, NlPPP1-Y1, NlPPP1-Y2 and NlPPP1-Y. There are 6 transcript variants were identified for NlPPP1-Y named as NlPPP1-Y-X1-6 (listed in Table 2). Sequence analysis showed the difference among the six transcripts was mainly caused by diverse inserts pattern (Fig. 1). The cDNA sequences were blasted against NCBI and N. lugens genomic data. Ten cDNA sequences have one or more blast hits with ≥ 90% query cover and E-value 0.0. Three blast hits (XM_022340879, XM_022340880 and XM_022340881) for NlPPP1-Y transcript variants were retained, which share common 5-UTR (582 bp in length) and ORF (1,005 bp in length) and three different lengths 3-UTR (233, 66 and 32 bp respectively). The 5-UTR is the same with NlPPP1-Y-X6. Analysis of the genomic position and structure showed that NlPPP1α, NlPPP1β, NlPPP1-Y, NlPPP1-Y1 and NlPPP1-Y2 contains 7, 3 and 4 exon and is located at scaffold 2,314, 31,341, 1,294, 10 and 1,537 in the N. lugens genome respectively. No intron was found in cDNA sequences of NlPPP1-Y1 and NlPPP1-Y2. The blast information, deduced amino acid lengths, molecular weights and isoelectric points of the NlPPP1 cDNA clones are illustrated in Table 2.
Phylogenetic and sequence alignment of PPP1. Blastp searches against the NCBI database revealed orthologues of NlPPP1 from other insects. The primary structures of the deduced amino acid sequences PPP1 from N. lugens were compared. Three signature motifs GDxHG, GDxVDRG, and GNHE (G, glycine; D, aspartic acid; x, any amino acid; H, histidine; V, valine; R, arginine; N, asparagine; E, glutamic acid) of PPP family within the catalytic domain were found. The conservative catalytic domain starting from the second α-helix (α1) and ending with the last β-strand (β14) as reported by Goldberg 32 were showed in Fig. 2A. All PPP1s contain a Thr-Pro-Pro-Arg (TPPR) amino acid sequence segment at their carboxyl terminal, which is a consensus sequence for phosphorylation by cyclin-dependent kinases (Cdks) demonstrated in somatic cells [33][34][35] . In Drosophila and N. lugens, this TPPR segment is retained in DmPP1α (PP1α-96A), NlPPP1α, and NlPPP1β, but absent in DmPP1α1 (DmPP1α-87B), DmPP1α2 (DmPP1α-13C), and male-biased NlPPP1-Ys. As a characteristic of most Ser / Thr / Tyr protein kinases, a large number of potential phosphorylation site were also identified in NlPPP1, among them, there are 9 to 17 Ser sites, 5 to 9 Thr sites, and 4 to 6 Tyr sites (Fig. 2B). The number and composition of phosphorylation sites for each gene are also different. The minimum number of potential phosphorylation sites is 19 in NlPPP1-α and the maximum is 30 in NlPPP1-Y. At the nucleic acids level, the identities between DmPP1α and NlPPP1α, NlPPP1-Y1, NlPPP1-Y2 and NlPPP1-Y are 70.1%, 63.1%, 62.5% and 52.3%, respectively.   Table 3). Phylogenetic tree was constructed using the MJ method to evaluate the molecular evolution relationships of the five PPP1s of N. lugens and other PPP1s from representative insect species. Five NlPPP1s were clustered into three classes (α, β and male-baised) (Fig. 3).
Expression characteristic of NlPPP1 in N. lugens. The developmental expression profile of five NlPPP1 in N. lugens was determined using RT-qPCR. NlPPP1α and NlPPP1β were expressed in all developmental stages and both sexes (Fig. 4A, B). In contrast, the transcription of NlPPP1-Y, NlPPP1-Y1 and NlPPP1-Y2 were restricted from 4th instar nymph to adult males whereas in females were nearly undetectable. The trace transcription of NlPPP1-Y and NlPPP1-Y1 were observed from 1 to 3 instar nymph (Fig. 4C-E). www.nature.com/scientificreports/ We then investigated the expression pattern in various tissues dissected from adults, including salivary glands(SG), fat bodies(FB), guts(GT), legs(LG), male internal reproductive organ(MIRO), and female internal reproductive organ (FIRO). The RT-qPCR results demonstrated that NlPPP1α and NlPPP1β showed significant higher expression in guts than in other tissues (Fig. 5A, B). To investigate the tissue-specific expression of malebiased NlPPP1s, total RNA was isolated from male tissues including SG, LG, GT, FB, and MIRO, for RT-qPCR analysis. The transcripts of NlPPP1-Y, NlPPP1-Y1 and NlPPP1-Y2 showed exclusive expression in MIRO, with relative expression levels in MIRO 13,7, and 18 -fold higher than in the MFB respectively (Fig. 5C-E). Their transcripts were also detected in tissues SG, GT and LG at very low level. The trace expression may be caused by fat bodies contamination in dissection. dsRNA sequence analysis. Six dsRNAs targeted to different NlPPP1 genes or different region were synthesized. The dsRNA name, length, targeted position and the largest length of 100% similarity stretch between the dsRNA and NlPPP1s were listed in Table 4. The specific dsRNAs were designed based on 3′-UTR of NlPPP1β, 5′-UTR of NlPPP1-Y and 5′-UTR of NlPPP1-Y1 respectively. The sequence of dsNlPPP1α showed 100% identity with 34 bp stretch in NlPPP1β. dsNlPPP1-Y2 with 78 bp stretch in NlPPP1-Y1, dsNlPPP1Ys showed over 20 bp stretch in NlPPP1-Y1 and NlPPP1-Y2. Moreover, dsNlPPP1Ys also showed over 50 bp stretch with more than 95% and 85% identity in NlPPP1-Y1 and NlPPP1-Y2 respectively.
Influence of injection dsNlPPP1α and dsNlPPP1β on survival rate and fecundity. Injection of dsNlPPP1α and dsNlPPP1β caused a significant decrease in the survivorship of N. lugens. The survival rate at 3rd day after injection was significantly lower in nymphs injected with dsNlPPP1α (83.1 ± 1.8%) and dsNlPPP1β (55.2 ± 3.5%) than with the dsGFP (100.0%). Ten days after injection, the survival rate of nymphs decreased to 4.8 ± 1.6% (dsNlPPP1α) and 2.8 ± 1.2% (dsNlPPP1β) (Fig. 6A). Almost no nymphs injected with dsNlPPP1α did reach to the adult stage. In the treatment of dsNlPPP1β, eclosion ratio was significantly reduced to 21.1 ± 5.8%, when compared to the dsGFP control (61.2 ± 13.2%) (Fig. 6B).
Fourth-instar nymphs were treated with dsRNA to assess the effects on reproduction. Because of high mortality, among 15 pairs, only 3 females injected with dsNlPPP1α and 4 females injected with dsNlPPP1β remained alive and available for oviposition. The mean number of eggs produced by dsNlPPP1α and dsNlPPP1β treated parent pairs was 196 and 133.3, respectively, significantly lower than that of control dsGFP (304 eggs/female).
NlPPP1α and NlPPP1β transcript levels in nymphs after the injection of dsNlPPP1α at 4th day were significantly down-regulated (Fig. 6C). dsNlPPP1α injection treatments led to reduced NlPPP1α and NlPPP1β expression approximately by 91.0% and 71.5% respectively, compared to the dsGFP control. The transcript levels of NlPPP1β decreased by 73.0% after the injection of dsNlPPP1β. However, no significant differences on transcript level of NlPPP1α in specimen injected dsNlPPP1β were observed.

Discussion
One of the most widespread mechanisms of post-translational regulation of proteins is the addition of phosphate by protein kinase, this phosphorylation is antagonized by protein phosphatases. The antagonistic actions of protein kinases and protein phosphatases are of equal importance in determining the degree of phosphorylation of each substrate protein 8 . Five different isforms of PPP1 defined by three signature motifs GDxHG, GDxVDRG, and GNHE within the conserved 30 kDa catalytic domain were identified from N. lugens. The constitutively expressed NlPPP1α and NlPPP1β were highly conserved which has higher to 86% identity with DmPPP1α and DmPPP1β respectively at the amino acid level. Down-regulation of NlPPP1α and NlPPP1β transcription resulted in 90% mortality in ten days and no nymph emergence. The dsNlPPP1α sequence from ORF region of NlPPP1α has 34 bp stretch of 100% identity with NlPPP1β, which contribute to the transcription reduction of NlPPP1β in nymphs injected with dsNlPPP1α. Only down-regulation the NlPPP1β with dsNlPPP1β sequence from 3′-UTR Figure 5. The expression pattern of NlPPP1s in various tissues from N. lugens. NlPPP1s expression value of MIRO was converted to 1. Duncan's tests were used to determine differences between tissues. The histogram bars (mean ± SEM) labeled with the same letter are not significantly different at p < 0.05. www.nature.com/scientificreports/ region of NlPPP1β resulted in 90% mortality in ten days and 80% reduction of eclosion ratio. Our results showed that silencing NlPPP1α and NlPPP1β are semilethal, therefore they are essential genes in N. lugens. The male-biased expressed NlPPP1Y, NlPPP1-Y1 and NlPPP1-Y2 were more divergent than non-sex biased PPP1 gene. Increasing evidences suggest that genes related to sex and reproduction change much faster between species than those limited to survival 36 . As demonstrated in Drosophila, male-biased genes evolve faster than unbiased genes in both coding sequence and expression level [37][38][39][40] . Specific silencing either NlPPP1-Y1 or NlPPP1-Y2 gene resulted in no or only slight mortality. Specific silencing of NlPPP1-Y with dsNlPPP1-Y resulted in 40% mortality, whereas silencing of 3rd instar nymphs using dsNlPPP1Ys designed against ORF of NlPPP1-Y resulted in 90% mortality and 80% reduction of eclosion rate, in which three male-biased NlPPP1-Y, NlPPP1-Y1 and NlPPP1-Y2 were silenced. This result suggests that this group of phosphatases has overlapping function allowing the compensation for the lack of one or the other member of the gene family. In D. melanogaster 41 , Heliothis virescens 42 and N. lugens 43 , male accessory gland proteins transferred to adult females via mating can regulate egg maturation and stimulate oogenesis, ovulation and oviposition. Our result showed that females mated with NlPPP1-Y, NlPPP1-Y1 and NlPPP1-Y2 silenced males led to a reduction in eggs amount and hatchability, whereas mated with NlPPP1-Y1 or NlPPP1-Y2 silenced males only led to a reduction in hatchability. Silenced NlPPP1-Y had no significant effect on both the oviposition amount and hatching rate. This result suggested that NlPPP1-Y1 and NlPPP1-Y2 may play more important roles in spermatogenesis and fertilization, and NlPPP1-Y mainly involved in development of male N. lugens. The physiological role of the male-specific phosphatases is still elusive, although the location and timing of their transcription as well as the conservation of their male-biased expression www.nature.com/scientificreports/ hint a specific role in reproduction. So the fine function of NlPPP1 in male development, sperm development and fertility needs to be clarified. RNAi silences gene expression through the production of small interfering RNAs (siRNAs). In Caenorhabditis elegans, the pairs having high degree of sequence similarity with the RNAi clones (100% over 25 bp, ≥ 94% over 50 bp, ≥ 89% over 100 bp, ≥ 84% over 200 bp and ≥ 81% over 300 bp) were predicted to exhibit off-target crossreaction 44 . The question of over how much length and how much similarity is necessary to observe off-target cross-reaction remains open in N. lugens. In our study, all target NlPPP1 genes were successfully silenced by using dsRNAs containing over 25 bp with 100% similarity. The efficient RNAi effect was also observed between NlPPP1-Y1 with dsNlPPP1Ys and NlPPP1-Y2 with dsNlPPP1Ys. NlPPP1-Y1 has 21 bp with 100% similarity or 54 bp with 95% similarity with the dsNlPPP1Ys. NlPPP1-Y2 has 20 bp with 100% similarity or 59 bp with 85% similarity with the dsNlPPP1Ys. RT-qPCR showed the expression of NlPPP1-Y1 and NlPPP1-Y2 reduced by 44.1% and 54.8% respectively in males injected with dsNlPPP1Ys when compared with males injected with dsGFP.
The overuse of conventional synthetic insecticides caused not only the serious detrimental effect on the environment but also the emergence of pest insect resistance to insecticides 2 . RNAi-based pest control strategies are emerging as environment friendly and species-specific alternatives for the use of conventional pesticides 45 . As the critical molecular switch in cell, protein phosphatase is a preferably considered target when designing RNAi-based pest control strategy as it affects numerous proteins dephosphorylation. NlPPP1-Y, NlPPP1-Y1, and NlPPP1-Y2 have low homology to known PPP1. The dsNlPPP1Ys sequence presents 100% similarity stretch with other organism is short than 10 bp. This means that the sequence of dsNlPPP1Ys varies greatly among insect species and the possibility of off-target effects is tiny. The dsNlPPP1Ys also showed ability to inducing high mortality rate, low eclosion rate and fecundity by interfering three male-biased NlPPP1 genes. In the application of RNAi to conserved genes at the cDNA level, the 3′-UTR is a good candidate sequence 46 . The dsNlPPP1β designed based on 3′-UTR also showed ability to inducing high mortality rate, low eclosion ratio and fecundity by silencing NlPPP1β genes. Therefore, NlPPP1-Ys would be a high efficient potential target gene used for N. lugens control. The selected dsNlPPP1β and dsNlPPP1Ys can be the preferred target sequence used for N. lugens control by means of RNAi.  Table 5. Mating with dsNlPPP1Ys treated males led to reduced protein content in IRO. The data in the table are means ± SEM (N = 4). A total of 40 male internal reproductive organs (4 replicates, 10 individuals in each replicate) and 20 female internal reproductive organs (4 replicates, 5 individuals in each replicate) were dissected from males or females 2 days after emergence, respectively. Different letters indicate significant difference between dsGFP treated control group and dsNlPPP1Ys treated group at p < 0.05.