Identification and Characterization of Neuropeptides and Their G Protein-Coupled Receptors (GPCRs) in the Cowpea Aphid Aphis craccivora

Neuropeptides are the most abundant and diverse signal molecules in insects. They act as neurohormones and neuromodulators to regulate the physiology and behavior of insects. The majority of neuropeptides initiate downstream signaling pathways through binding to G protein-coupled receptors (GPCRs) on the cell surface. In this study, RNA-seq technology and bioinformatics were used to search for genes encoding neuropeptides and their GPCRs in the cowpea aphid Aphis craccivora. And the expression of these genes at different developmental stages of A. craccivora was analyzed by quantitative real-time PCR (qRT-PCR). A total of 40 candidate genes encoding neuropeptide precursors were identified from the transcriptome data, which is roughly equivalent to the number of neuropeptide genes that have been reported in other insects. On this basis, software analysis combined with homologous prediction estimated that there could be more than 60 mature neuropeptides with biological activity. In addition, 46 neuropeptide GPCRs were obtained, of which 40 belong to rhodopsin-like receptors (A-family GPCRs), including 21 families of neuropeptide receptors and 7 orphan receptors, and 6 belong to secretin-like receptors (B-family GPCRs), including receptors for diuretic hormone 31, diuretic hormone 44 and pigment-dispersing factor (PDF). Compared with holometabolous insects such as Drosophila melanogaster, the coding genes for sulfakinin, corazonin, arginine vasopressin-like peptide (AVLP), and trissin and the corresponding receptors were not found in A. craccivora. It is speculated that A. craccivora likely lacks the above neuropeptide signaling pathways, which is consistent with Acyrthosiphon pisum and that the loss of these pathways may be a common feature of aphids. In addition, expression profiling revealed neuropeptide genes and their GPCR genes that are differentially expressed at different developmental stages and in different wing morphs. This study will help to deepen our understanding of the neuropeptide signaling systems in aphids, thus laying the foundation for the development of new methods for aphid control targeting these signaling systems.

Neuropeptides are the most abundant and diverse signal molecules in insects. They act as neurohormones and neuromodulators to regulate the physiology and behavior of insects. The majority of neuropeptides initiate downstream signaling pathways through binding to G protein-coupled receptors (GPCRs) on the cell surface. In this study, RNA-seq technology and bioinformatics were used to search for genes encoding neuropeptides and their GPCRs in the cowpea aphid Aphis craccivora. And the expression of these genes at different developmental stages of A. craccivora was analyzed by quantitative real-time PCR (qRT-PCR). A total of 40 candidate genes encoding neuropeptide precursors were identified from the transcriptome data, which is roughly equivalent to the number of neuropeptide genes that have been reported in other insects. On this basis, software analysis combined with homologous prediction estimated that there could be more than 60 mature neuropeptides with biological activity. In addition, 46 neuropeptide GPCRs were obtained, of which 40 belong to rhodopsin-like receptors (A-family GPCRs), including 21 families of neuropeptide receptors and 7 orphan receptors, and 6 belong to secretin-like receptors (B-family GPCRs), including receptors for diuretic hormone 31, diuretic hormone 44 and pigment-dispersing factor (PDF). Compared with holometabolous insects such as Drosophila melanogaster, the coding genes for sulfakinin, corazonin, arginine vasopressin-like peptide (AVLP), and trissin and the corresponding receptors were not found in A. craccivora. It is speculated that A. craccivora likely lacks the above neuropeptide signaling pathways, which is consistent with Acyrthosiphon pisum and that the loss of these pathways may be a common feature of aphids. In addition, expression profiling revealed neuropeptide genes and their GPCR genes that are differentially expressed at different developmental stages and in different wing morphs. This study will help to deepen our understanding of the neuropeptide signaling systems in aphids, thus laying the foundation for the development of new methods for aphid control targeting these signaling systems.

INTRODUCTION
Insect neuropeptides are trace polypeptides produced by neurosecretory glands or neurosecretory cells in insects (1). They are the most abundant and diverse signal molecules and play important regulatory roles in the physiology and behavior in insects such as the growth, development, metabolism, reproduction, and feeding of insects (2). Most of the known neuropeptides are oligopeptides or small protein molecules composed of several to tens of amino acids. Biologically active mature neuropeptides are formed by complex post-translational processing of neuropeptide precursors. Neuropeptides act as signal molecules to exert their physiological functions as neuromodulators or neurohormones by binding to corresponding membrane receptors. Except for a few families of neuropeptides, such as prothoracicotropic hormone (PTTH), eclosion hormone (EH), and neuropeptide-like precursor 1 (NPLP1), the majority of neuropeptide receptors belong to the G protein-coupled receptors (GPCRs) superfamily (3). In recent years, insect GPCRs have been proposed to be used as potential targets to explore new insecticides or behavioral regulators, opening up new ways for the safe control of pests (4) By targeting neuropeptide GPCRs, future insecticides or behavioral regulators can prevent or over-stimulate the normal activity of these proteins, producing a lethal or deleterious effect on pests. Many species of aphids are important agricultural and forestry pests. At the same time, aphids are also model organisms for studying phenotypic plasticity, insect-plant interactions, insect-plant-virus interactions and endosymbionts and thus have attracted much attention from researchers (5).
In recent years, the continuous innovation of next-generation sequencing technology and bioinformatics technology has greatly promoted the prediction and identification of insect neuropeptides and their receptors (6)(7)(8)(9). However, only Acyrthosiphon pisum has so far been comprehensively investigated for neuropeptides and their receptors in aphids (10)(11)(12). Aphis craccivora is a worldwide insect pest that causes harm to many crops such as cowpea, peanut and pea (13,14). In addition to directly sucking plant sap, it can also cause fungal diseases and acts as a vector for various plant viruses such as cowpea mosaic virus (15) and peanut stripe virus (16). At present, A. craccivora is mainly controlled by chemical insecticides, resulting in increasingly severe insecticide resistance and ecological problems (17). So it is imperative to develop alternative new control strategy and identification of neuropeptides and their receptors of A. craccivora will be the first step toward achieving this goal.
In the present study, based on transcriptome sequencing and bioinformatics analysis, we screened the genes encoding neuropeptide precursors and their GPCRs in A. craccivora. We also conducted structural analysis and mature peptide prediction for the identified neuropeptide precursors and classified the neuropeptide receptors into subclasses based on phylogenetic analysis. The expression profiles of these neuropeptides and their GPCRs at different developmental stages were determined by quantitative real-time PCR (qRT-PCR). These identified gene sequences provide the basis for further investigation of the function of A. craccivora neuropeptides and their GPCRs.

Insect Rearing, Sample Collection, and RNA Extraction
Aphis craccivora was collected from peanut fields in Laixi, Shandong Province, China, in 2015. Then the collected aphids were reared for many generations with Vicia faba seedlings in our laboratory under the photoperiod of 14 L: 10 D, at 22 ± 1 • C. For transcriptome sequencing in the next step, the heads of 1,000 wingless parthenogenetic adults were quickly removed into phosphate-buffered saline (PBS) buffer and immediately frozen with liquid nitrogen. To find as many target genes as possible, we also collected 50-60 whole-body aphids from each instar (including wingless and winged morphs) as a mixed sample to generate the second transcriptome. Total RNA extraction was performed from the above two samples using RNAiso Plus reagent (Takara, Dalian, China) following the manufacturer's instructions.

Transcriptome Sequencing, Assembly and Annotation
We used 2 µg total RNA to construct each library with the NEBNext R Ultra TM RNA Library Prep Kit for Illumina R (Illumina, San Diego, CA, USA). Paired-end sequencing was performed using an Illumina HiSeq TM 2,500 by Gene Denovo Biotechnology Co. (Guangzhou, China). After the raw sequencing reads were filtered to remove adaptors, lowquality sequences and reads with N content <5%, de novo assembly were performed using Trinity software with the default parameters (18). Finally, the annotation of unigene sequences were performed using a Blastx search against the NCBI non-redundant protein sequence (Nr), A manually annotated and reviewed protein sequence database (Swiss-Prot), Kyoto Encyclopedia of Genes and Genomes (KEGG), and Clusters of Orthologous Groups (COG) databases with a cut-off E-value of 10 −5 .

Structure and Domain Analysis
For the neuropeptide precursors, SignalP 4.1 (http://www.cbs. dtu.dk/services/SignalP/) was employed to predict the signal peptides (32). The splice sites were predicted using the Known Motif and Insect Models methods of NeuroPred (http:// stagbeetle.animal.uiuc.edu/cgi-bin/neuropred.py) (33) and were corrected based on the processing procedures of known insect neuropeptide precursors. Sulfinator (https://web.expasy. org/sulfinator/) was used to identify the sulfation of tyrosine residues (34). Cyclization of N-terminal glutamine/glutamic acid residue and amidation of the C-terminal glycine residue were predicted based on homology to the reported invertebrate neuropeptides. Multiple alignments of amino acid sequences were performed using Clustal X 2.0 and the results were visualized using GeneDoc software.

Reverse Transcription PCR (RT-PCR) Amplification of Neuropeptide GPCR Genes
Two unigene sequences in the transcriptome data were likely partial fragments from the same gene. RT-PCR was used to amplify and obtain longer sequences. Upstream and downstream primers were designed for each of the fragments using Primer Premier 5.0 software ( Table S1). The cDNA template was synthesized using a PrimeScript TM II 1st Strand cDNA Synthesis Kit (Takara, Dalian, China) The reaction procedure was as follows: pre-denaturation at 94 • C for 3 min; 37 cycles of denaturation at 94 • C for 30 s, annealing at 56 • C for 30 s, and elongation at 72 • C for 1 min; and a final elongation at 72 • C for 5 min. The amplication products were sent to Invitrogen (Shanghai, China) for sequencing. The sequencing results were compared with the corresponding unigene sequences using Blastn to determine sequence identity.

Phylogenetic Analysis
Neuropeptide GPCRs of A. craccivora were compared with the reported homologs from A. pisum, N. lugens, D. melanogaster, B. mori, and T. castaneum. Amino acid sequences with a length shorter than 150 aa were removed from the data set for phylogenetic analysis. Multiple alignments were carried out using MAFFT (37). The maximum-likelihood (ML) phylogenetic trees were constructed using the FastTree software (version 2.1.7) under the default settings (38). Local support values were calculated based on the Shimodaira-Hasegawa (SH) test implemented within FastTree. Metabotropic glutamate receptor (mGluR) of D. melanogaster (encoded by CG11144) served as the outgroup for the phylogenetic tree of A-family GPCRs, and 5-hydroxytryptamine receptor (5-HTR) of D. melanogaster (encoded by CG1056) served as the outgroup for the phylogenetic tree of B-family GPCRs. The phylogenetic trees were finally edited and displayed using FigTree 1.4.3 (http://tree. bio.ed.ac.uk/software/figtree).

Expression Profiling of Neuropeptides and Their GPCRs at Different Developmental Stages
To determine the transcriptional expression levels of neuropeptides and their GPCRs at different life stages of A. craccivora, 20-60 aphids from each instar were collected separately for qRT-PCR analysis. From the third instar, we were able to distinguish winged from wingless aphids by wing buds of the thorax. For this reason, winged and wingless morphs of the aphids were collected separately from the 3rd, the 4th instar nymphs and adults. After total RNA was extracted as described above, cDNA was synthesized by reverse transcribing 1 µg of total RNA using a PrimeScript TM RT reagent Kit with gDNA Eraser (Perfect Real Time) (Takara, Dalian, China) according to the manufacturer's instructions. The gene-specific primers for qRT-PCR were designed using the NCBI online primer design tool (https://www.ncbi.nlm.nih.gov/tools/primer-blast/) ( Table S2). We select 40S ribosomal protein S8 gene (RPS8) as the reference gene to normalize the expression of the target genes (39).
The qRT-PCR reactions were run with SYBR Premix Ex Taq II (Tli RNaseH Plus) (Takara, Dalian, China) on an iQ TM 5 Multicolor Real-Time PCR Detection System (Bio-Rad Hercules, CA, USA). The qRT-PCR program was set as follows: 95 • C for 3 min; 40 cycles of 95 • C for 10 s and 58 • C for 30 s; and followed by a melting cycle (from 60 to 95 • C with an increase of 0.5 • C per cycle). The experiment included 3 biological replicates and 3 technical replicates. The relative expression level of the target genes was calculated using the 2 Ct method (40). Statistical analysis was performed using SPSS 20.0 (Systat Software Inc., London, UK). Statistical significance was evaluated using one-way ANOVA followed by the Tukey HSD test at the 0.05 level.

Transcriptome Sequencing, Assembly, and Annotation
We sequenced the sample from the heads of A. craccivora wingless parthenogenetic adults and 56,634,657 clean reads (8.4 Gbp) were obtained. For the mixed sample, 70,096,246 clean reads (10.3 Gbp) were obtained. After assembly, 60,113 unigenes (N50 = 1,169) were obtained with an average length of 703 bp and a GC content of 38.08%. Among all the unigenes, a total of 45,150 were annotated in the Nr, Swissprot, KOG, and KEGG databases (Figure S1A), and the top 3 species with the most matches in Nr database were A. pisum, Diuraphis noxia, and Diachasma alloeum ( Figure S1B). The transcriptome raw data have been submitted to the Sequence Read Archive (SRA) database at NCBI under the accession number SRP189414.

Identification of Neuropeptide Precursors in A. craccivora
Based on the homology search and Nr-annotation, a total of 40 genes encoding neuropeptide precursors were identified, involving 32 neuropeptide families ( Table 1; Supplementary Material 1) and covering most insect neuropeptide families (Table S3). However, a small number of family members were missing, including CNMamide (CNMa), corazonin (Crz), elevenin, neuroparsin, pigmentdispersing factor (PDF), RYamide (RYa), sulfakinin (SK), sex peptide (SP), trissin, arginine vasopressin-like peptide (AVLP), and adipokinetic hormone/corazonin-related peptide (ACP). Of the 40 neuropeptide genes obtained, 36 had intact ORFs, and the remaining four non-full-length sequences include three insulinlike peptide (ILP)-encoding genes and 1 PTTH -encoding gene ( Table 1). The majority of A. craccivora neuropeptide precursors showed more than 90% similarity with their homologous sequences in other aphids, whereas a few had low similarity with homologous sequences in related species, such as proctolin (Proc), leucokinin (LK), and diuretic hormone 44 (DH44). Besides, orthologs of some A. pisum neuropeptides were not identified from A. craccivora transcriptome, including ion transport peptides (ITPs) and ILPs. However, the transcript encoding ITP-like was found in the A. craccivora transcriptome. Ten different genes coding for ILPs have been found in A. pisum genome but only six of these genes in A. craccivora, i.e., genes encoding ILP1, ILP4, ILP5, ILP7, ILP8, and ILP10. By contrast, the gene coding for NPLP3 possesses orthologs in A. craccivora but not in A. pisum.

Mature Neuropeptide Prediction
During the process of neuropeptide production, a larger inactive neuropeptide precursor is first synthesized. One or several bioactive mature peptides are then formed through a series of enzymatic digestions and post-translational modifications, which are critical to ensure the biological activity and stability of mature neuropeptides (41). The signal peptide cleavage site, sulfation site, amidation site, and other modification sites of A. craccivora neuropeptide precursors were predicted, and the prediction results were revised based on the homology and conserved motif of the known invertebrate neuropeptides. A total of 164 mature peptides (including protein hormones) were predicted from 40 neuropeptide precursors in A. craccivora. Among them, more than 60 mature peptides could be biologically active. Meanwhile, a total of 52 mature peptides were predicted to be amidated at the C-terminus, and 11 sites of 6 mature peptides were predicted as potential sites of sulfation (Figure 1;  Supplementary Material 2).
Most A. craccivora neuropeptides have highly conserved characteristic motifs within the same family. For example, the predicted 8 mature allatostatin A (AstA) peptides contain the Cterminal signature sequence of -Y/FXFGL/Iamide; adipokinetic hormone (AKH) has the sequence motif characterized by pQX 6 Wamide. Some mature neuropeptide sequences are identical in different aphids, such as short neuropeptide F (sNPF) and SIFamide (SIFa) of A. craccivora, A. pisum, and Myzus persicae. Even the sequences of mature peptides of DH31 and DH44 that are up to tens of amino acids are identical in A. craccivora and A. pisum (Supplementary Material 2).
The mature neuropeptides of A. craccivora also exhibit some structural variability, which is mainly reflected by the expansion and truncation of the polypeptide sequences. The former mainly involves Proc, CCHamide (CCHa), and AstB, and the latter occurs in SIFa. The predicted Proc of A. craccivora is identical to that of M. persicae with an N-terminal extension compared to the pentapeptide proctolin (RYLPT) of most other insects (42); the effect of this expansion on Proc biological activity remains unknown. The N-terminus of A. craccivora CCHa1 has an addition of 4 amino acids (KQGA), and the N-terminus of CCHa2 has an additional Gly. There are two AstB peptides with a sequence motif of -WX 7 Wamide in A. craccivora, which has an extension of one amino acid residue in the middle of the polypeptide compared to the classical motif. This subtype of AstB has also been found in A. pisum, M. persicae, A. gossypii, N. lugens, and R. prolixus (6,43). Therefore, it is speculated that -WX 7 Wamide may be a special AstB subtype that is prevalent in hemipteran insects. The predicted A. craccivora SIFa sequence (FRKPPFNGSIFamide) is consistent with that of A. pisum and M. persicae, but lacks an amino acid residue at the N-terminus compared to the classical motif.
In addition, variability is also reflected in the mutation of amino acid residues. The A. craccivora capability (CAPA) gene encodes a CAPA-pvk (periviscerokinin) peptide and a CAPApk (pyrokinin) peptide with the conserved characteristic motif, -PRVamide, and -PRLamide, respectively. This gene also encodes another polypeptide, EGLIAFPRIamide, which is identical to the CAPA-2 peptide of A. pisum and M. persicae (42,44). This subtype of CAPA with -PRIamide is less common in insects. In addition to aphids, it is found only in a few insects such as T. castaneum, A. mellifera, and Camponotus floridanus (29,45,46).

Neuropeptide GPCRs in A. craccivora
The insect neuropeptide GPCRs belong to rhodopsin-like receptors (A-family GPCRs) or secretin-like receptors (B-family GPCRs). A total of 46 putative neuropeptide GPCRs were identified from the A. craccivora transcriptome data, among which 34 contained complete ORFs ( Table 2; Supplementary Material 3). The encoded theoretical proteins had a length of 154-1184 amino acids. Forty of 46 GPCRs belong to the A-family (represented by 7tm_1 or PF0001 in the Pfam database), and six belong to the B-family (represented by 7tm_2 or PF0002 in the Pfam database) ( Table 2). Phylogenetic analysis revealed that all A. craccivora neuropeptide GPCRs had unique orthologous proteins in A. pisum (Figures 2-4).
Leucine-rich repeat-containing GPCRs (LGRs) are a special class of neuropeptide/protein hormone receptors in A-family GPCRs. They act as receptors for glycoprotein hormone or insulin/relaxin-related proteins and play important regulatory roles in development and reproduction of insects (48,49). According to the structure characteristics, LGRs are divided into type A, type B, and type C. The main differences among the 3 types are the number of leucine-rich repeats (LRRs), the presence or absence of a low density lipoprotein receptor-like cysteinerich (LDLa) motif, and type-specific hinge region (48,49). In this study, a total of 5 LGRs (A36-A40) were identified from A. craccivora (Figure 3). Among them, only A40 has a complete ORF. CDD analysis indicated that A40 has 15 LRR motifs and no LDLa motif, showing typical structural characteristics of a type B LGR. A36 has an N-terminus containing 7 LDLa repeats and 8 LRR motifs, consistent with the typical characteristics of a type C LGR. The presence of conserved domains was not detected in other LGRs due to incomplete sequences. Phylogenetic analysis revealed that A. craccivora A37 and A39 were clustered into the branch of GPA2/GPB5 receptor. A40 is the ortholog of D. melanogaster LGR2, which has been confirmed to be a receptor for Bursicon and is also the first deorphanized LGR in invertebrates (50,51). It has been confirmed that dILP8 is a ligand for LGR3 in D. melanogaster (52), and the A38 receptor of A. craccivora is clustered into this branch. Thus, it is speculated that A38 may be an insulin-like receptor. In addition, there is an orphan receptor, A36, clustered into the C2 LGR subfamily (Figure 3).
B-family GPCRs are a small class of receptors that differ structurally and functionally from other families. Their remarkable structural feature is a long N-terminal domain. The B-family can be further subdivided into three subfamilies: subfamily B1, B2, and B3 (53). Neuropeptide GPCRs belong to the B1 subfamily and contain three types of neuropeptide/protein hormone receptors: calcitonin receptor (also known as DH31 receptor), corticotropin releasing factor (CRF)-like diuretic hormone receptor (also known as DH44 receptor), and PDF receptor (19). All these three types of receptors were identified in A. craccivora and involved a total of 6 sequences (B1-B6) (Figure 4). Phylogenetic analysis revealed that A. craccivora B1, B3, and B4 were clustered in the branch of DH31 receptor; B5 and B6 were clustered in the branch of DH44 receptor; and B2 was clustered in the branch of PDF receptor (Figure 4).

Expression Profiles of Neuropeptides and Their GPCRs at Different Developmental Stages of A. craccivora
The expression profiles of A. craccivora neuropeptides and their GPCRs at different developmental stages were determined by qRT-PCR. For most of the measured neuropeptide-encoding genes, the expression abundance showed a trend of decreasing first and then increasing during development, i.e., the expression level was significantly downregulated at the 2 nd and 3 rd instar stages (Figure 5). For the neuropeptide GPCR genes, most were differentially expressed at different developmental stages and in different wing morphs, and moreover, most were more abundant in the adult stage than in other stages (Figure 6). The two TK receptors, A24 and A29, had basically the same developmentalstage expression pattern, while the expression patterns of the two members of NPF receptor family, A31 (NPFR ortholog) and A33 (NPY2R ortholog), were different ( Figure 6). Overall, the gene expression patterns of most neuropeptide-encoding genes and their presumed GPCR-encoding genes were not consistent (Figures 5, 6).

DISCUSSION
Based on RNA-seq technology, we identified a large number of genes encoding neuropeptides and their GPCRs from A. craccivora in this study. The types and numbers of these genes are basically similar to those found in other insects, indicating that RNA-seq can be used as a powerful tool for the excavation of these genes from insects in the absence of whole-genome data. For some neuropeptide signaling systems, A. craccivora seems to lack both receptor and neuropeptide precursor (SK, Crz, AVLP, trissin, and elevenin). For others, the neuropeptide was not found or is lacking, but its presumed GPCR is present (e.g., PDF, RYa, CNMa, etc.). Why were they not found? There may be three main reasons. Firstly, the sequenced samples probably did not cover all stages of the A. craccivora life cycle, and some genes may not have been expressed in these measured samples; in other cases, some genes are expressed at low abundance in certain tissues or at certain developmental stages and the transcripts would be further diluted in the mixed sample. Secondly, due to the lack of strong sequence conservation among some neuropeptide precursors or neuropeptide GPCRs, their clear orthologs could not be found in A. craccivora based on homology searches. Thirdly, the neuropeptides may truly be absent in this species, a phenomenon that is common in insects. For example, sex peptide was reported only in the Drosophila genus; while neuroparsin is missing in D. melanogaster and other Drosophila species (54). Compared to D. melanogaster, B. mori, and T. castaneum, we found neither the genes coding for SK, Crz, AVLP, trissin, and elevenin in A. craccivora nor the genes coding for their receptors. If both neuropeptide and its receptor are truly absent, it is more likely that this signaling system is lacking as a whole based on neuropeptide-receptor co-evolution in this aphid species. Likewise, the former four signaling systems are not found in A. pisum (10,12). These absent pathways are mainly involved in feeding, ecdysis and osmotic homeostasis in other insects (55)(56)(57)(58), suggesting the absense is the result of constant adaptation to the environment during evolution of aphids. SK, Crz, and AVLP neuropeptide signaling systems were found in another hemipteran insect, N. lugens, but the trissin signaling system was absent in this species (6). For the neuropeptide GPCRs, two receptor genes for LK, SIFa and NPF have been identified from A. craccivora, A. pisum, and N. lugens (Figure 2), suggesting that each of these receptor genes duplicated in the evolution. These 3 types of receptors are mainly functionally involved in water balance, sexual behavior, and feeding behavior (59)(60)(61)(62). However, we cannot conclude from these findings that the replication of these genes is a common feature of hemipterans, because they are a highly diverse group of insects. A total of 40 genes enconding neuropeptides were identified from A. craccivora. Compared with another aphid species A. pisum (10), ITP and several ILPs were not detected from A. craccivora. ITP is a member of the crustacean hyperglycaemic hormone/ITP (CHH/ITP) superfamily, which stimulates the ileum to transport Cl − ion and may act as an antidiuretic hormone in insects (63,64). In addition to ITP peptide, the ITP gene also encodes an alternatively spliced variant without C-terminal amidation (known as ITP-like or ITPL) (64). Only the transcript encoding ITPL was found in the A. craccivora transcriptome, whereas the ITP-encoding transcript was not. Generally, invertebrate genomes contain multiple genes encoding ILPs including 10 ILP precursor-encoding genes in the A. pisum genome (10), while only 6 of this family were found in A. craccivora. In contrast to ITPs and ILPs, an NPLP3 ortholog is present in A. craccivora but not in A. pisum. There are 4 NPLP-encoding genes, NPLP1-4, in D. melanogaster (65); NPLP2 and NPLP3 were found in A. mellifera (25). However, only one NPLP-encoding gene, NPLP1, was found in A. pisum (10). We found NPLP1 and NPLP3 in A. craccivora, which is consistent with that previously reported in a hemipteran insect, Diaphorina citri (7), while NPLP1, NPLP3, and NPLP4 are present in another hemipteran insect, N. lugens (6). Interestingly, no genes encoding PTTH orthologs were originally identified in hemipteran insects such as A. pisum and R. prolixus, probably due to insufficient conservation of PTTH sequences (6,10). However, we used the sequences of PTTH precursors reported in other insects as a query to conduct a homology search and found a unigene in the A. craccivora. Moreover, a homologous sequence (ARM65502.1) was also found in A. pisum by the Blast program. PTTH is a polypeptide hormone that stimulates prothoracic glands to synthesize and release ecdysone and it was FIGURE 5 | Relative expression levels of neuropeptide-encoding genes at different developmental stages of A. craccivora. The transcript level was measured via qRT-PCR and normalized against RPS8. Data are presented as the mean ± SD. The data were statistically analyzed by one-way ANOVA followed by Tukey's HSD test. Different lowercase letters indicate significant differences at the 0.05 level. "−1" indicates wingless aphids; "−2" indicates winged aphids.
first cloned and characterized in B. mori (66). Mature B. mori PTTH is a homodimer. Each monomer contains 7 conserved Cys residues, of which 6 Cys form 3 internal disulfide bonds and the 7th Cys forms a disulfide bond between two monomers (66,67). PTTH was later identified in several lepidopterans (68) and D. melanogaster (69). Compared with other neuropeptides, the sequence similarity of PTTH between different insects was relatively lower. However, the position of Cys residues and the folding structure of dimers were highly conservative (69-72). Tanaka et al. (6) used B. mori PTTH as a query to search for the analog in the N. lugens transcriptome and finally found a unigene whose theoretical protein has conserved structural features of PTTH. This is the first PTTH-like gene identified from hemimetabolous insects. We also found a unigene from the A. craccivora transcriptome that lacks the 5 ′ -end of the ORF and may encode a partial sequence of the PTTH precursor. Sequence alignment revealed that the PTTH-like precursor of A. craccivora lacks the Cys residue at position 4 ( Figure S2). The PTTH ortholog in A. pisum (ARM65502.1) also lacks the Cys residue at this position. Notably, all the 7 Cys residues of PTTH monomer are involved in the formation of disulfide bonds, which is required for PTTH to achieve biological activity (70). The effect of the absence of the fourth Cys residue on the formation of the PTTH advanced structure and whether it affects its normal biological function in aphids deserves further investigation.
We also carried out the expression profiling of genes encoding neuropeptides and their receptors in different developmental stages. Unfortunately, the low expression level of some measured genes reduced the accuracy of their qRT-PCR data. Thus, such genes were not included in the results (Figures 5, 6). FIGURE 6 | Relative expression levels of the genes encoding G protein-coupled receptors for the neuropeptides at different developmental stages of A. craccivora. The transcript level was measured via qRT-PCR and normalized against RPS8. Data are presented as the mean ± SD. The data were statistically analyzed by one-way ANOVA followed by Tukey's HSD test. Different lowercase letters indicate significant differences at the 0.05 level. Neuropeptide receptors A1-A40 are abbreviated as A1-A40, respectively; Neuropeptide receptors B1-B6 are abbreviated as B1-B6, respectively. And the receptors are marked with their predicted neuropeptide ligand partners in brackets; orphan receptors are not marked. "−1" indicates wingless aphids; "−2" indicates winged aphids.
The expression and functions of insect neuropeptides may change with development (73). In the present study, the tested neuropeptide genes also showed different levels of expression during development of A. craccivora. However, whether the differential temporal expressions of these genes have a biological role on development needs further investigation. It is worth noting that NPF was expressed in A. craccivora wingless adults at a significantly higher level compared with that in winged adults. NPF is homologous to the vertebrate neuropeptide Y (NPY) and its functional roles are well-known in regulation of feeding behavior in a wide range of insects (60). This neuropeptide was also suggested to play a positive role in the regulation of food intake in A. pisum using RNAi knockdown (74). In view of this, whether the differential expression of NPF between morphs has an impact on feeding behavior of different wing-morph adults or not, worths further research.
Nutritional conditions are key environmental factors that affect the growth and development of insects. The ILPs play a main role in systematically regulating the growth of the body in response to nutritional conditions (75). A prominent feature of insect ILPs is that multiple genes encode different ILPs. For example, 8 ILPs have been identified in D. melanogaster (76,77), while 39 ILPs have been identified in B. mori (78). The synthesis and secretion of different ILPs in the same insect species exhibits different temporal and spatial pattern, and their functions are also different (79). In A. pisum, ILP5 has attracted the attention of researchers due to its abundant expression. It is speculated that this gene may be involved in the rapid development of aphids (10). It has also been reported that ILP5 is involved in embryonic development during wing differentiation in A. pisum (79). From the RPKM values, the expression abundance of ILP5 gene in A. craccivora is much higher than other ILP-encoding genes, which is consistent with that reported in A. pisum (80). The ILP5 transcript was detected from both the head as well as the mixed whole-body samples of A. craccivora. Besides, it was abundantly expressed in the head of wingless adult aphids.
A. pisum ILP5 gene shows the similar expression pattern: it is expressed systemically in A. pisum adults and third instar nymphs with the highest expression level in the head of wingless adults (including the antennae) (81). This expression characteristic of ILP5 in aphids suggests the importance of its biological function and deserves further in-depth investigation.
In this study, the sequence similarity between neuropeptide precursors is very high among different species of aphids. Several neuropeptide genes, such as sNPF, SIFa, AstC, AstCC, and Proc, even encode identical mature active peptides in different species of aphids. For neuropeptide GPCRs, besides high sequence similarity, all of the identified A. craccivora neuropeptide GPCRs have a one-to-one orthologous relationship with the homologs of A. pisum. This result indicates that the aphid neuropeptide GPCRs are also highly conserved reflected by their types and numbers. In summary, the neuropeptide signaling systems in aphids are highly conserved and provide potential insecticide targets for the control of aphids. In addition, most neuropeptide signalings are conserved in both holometabolous and hemimetabolous insects, suggesting that they play key roles in the physiological processes of insects (3). Some neuropeptide signaling systems are lost in specific species, possibly due to functional redundancy that occurs during adaptation and evolution or other neuropeptide signalings taking over their functions (3).

DATA AVAILABILITY STATEMENT
The datasets generated for this study can be found in online repositories. The names of the repository/repositories and accession number(s) can be found in the article/Supplementary Material.

AUTHOR CONTRIBUTIONS
T-XL and M-JQ supervised the whole project. XL and LD performed the experiments and wrote the manuscript with equal contributions. X-JJ, QJ, and C-JQ provided technical support. All authors contributed to the article and approved the submitted version.