Adults of Pylorgus porrectus (Fig. 1a, b) were collected from Zhongshan Botanical Garden (32°03.27'N, 118°49.85'E), Nanjing, Jiangsu Province, China, 20 April 2022, Cuiqing Gao leg. Adults of P. sordidus (Fig. 1c) were collected from Hongqi Management and Protection Station, Yintiaoling Nature Reserve (31°23.87'N, 109°41.32'E), Wuxi County, Chongqing, China, 1 July 2022, Suyan Cao leg. The specimens were identified based on the morphological characteristics seen under a Zeiss Stereo Discovery V8 Zoom Microscope and deposited in the Insect Collection, College of Forestry, Nanjing Forestry University.

Figure 1. 

Pylorgus species sequenced a, b P. porrectus, dorsal and ventral views c P. sordidus, dorsal view.

The complete genomic DNA was extracted from an adult sample using a Rapid Animal Genomic DNA Isolation Kit (Sangon Biotech, Shanghai, China). Libraries were prepared on an Illumina MiSeq PE300 platform (Sangon Biotech, Shanghai, China). Low-quality and short reads were removed using Fastp v. 0.36 (Chen et al. 2018) to obtain clean reads and ensure rich quality of information analysis.

SPAdes v. 3.15 (Bankevich et al. 2012) was used to assemble the high-quality next generation sequencing data de novo to construct contig and scaffold. After the assembly was completed, we evaluated and quality controlled the assembly results, excluding any contamination that may originate from the host genome in the subsequent analysis, and only retained the scaffolds derived from the genome of the organelle. We used BLASTn to compare the scaffolds with the NCBI library to obtain sequence similarity information, extracted the sequencing depth and coverage information of each scaffold, and manually selected possible target scaffolds after sorting out and comprehensively considering the above information. Then GapFiller v. 1.11 (Boetzer and Pirovano 2012) was adopted to supplement GAP to the contig obtained by splicing, and PrInSeS-G was adopted to carry out sequence correction to correct editing errors and insertion and deletion of small fragments in the splicing process, and finally the complete mitochondrial genome was obtained.

For mitochondrial gene annotation, we used tBLASTn and GeneWise to back-align with near-source reference databases to obtain the coding sequence (CDS) gene boundaries, and MiTFi to obtain the transfer RNA genes (tRNAs) sequence annotation. The non-coding ribosomal RNA genes (rRNAs) were identified by cmsearchrfam alignment and finally summarized into a complete annotation result.

The nucleotide composition and RSCU (relative synonymous codon usage) were calculated using PhyloSuite v. 1.2.2 (Zhang et al. 2020) and MEGA X (Kumar et al. 2018). Strand asymmetry was calculated using the formula: AT-skew = [A−T]/[A+T] and GC-skew = [G−C]/[G+C] (Perna and Kocher 1995). DnaSP v. 5 (Librado and Rozas 2009) was used to calculate the value of Ka (the nonsynonymous substitution rate), Ks (the synonymous substitution rate), and nucleotide diversity.

The mitochondrial genome data of 25 species in Pentatomomorpha were used to reconstruct the phylogenetic relationship of Lygaeoidea, in which 21 species of Lygaeoidea were regarded as ingroup and four species was regarded as outgroup (Table 1). All sequences were standardized and extracted 13 protein-coding genes (PCGs) by PhyloSuite v. 1.2.2 (Zhang et al. 2020). The 13 PCGs of the 25 species were aligned individually using codon-based multiple alignments with MAFFT v. 7.313 software with default settings (Katoh and Standley 2013). Gblocks v. 0.91b software was used to remove the intergenic gaps and ambiguous sites (Talavera and Castresana 2007), and all PCGs sequences were concatenated in PhyloSuite v. 1.2.2. The best partitioning scheme and evolutionary models for constructing Bayesian inference (BI) and maximum-likelihood (ML) trees were selected by PartitionFinder2 (Lanfear et al. 2016), with a greedy algorithm, BIC criterion, and the gene and codon model.

BI phylogenies were inferred using MrBayes v. 3.2.6 (Ronquist et al. 2012) under partition model (2,000,000 generations), in which the initial 25% of sampled data were discarded as burn-in. ML phylogenies were inferred using IQ-TREE (Nguyen et al. 2015) under the Edge-linked partition model for 5000 standard bootstraps with 1000 replicates.

The assembled complete mitogenomes of Pylorgus porrectus and P. sordidus are circular DNA molecules of 15,174 bp and 15,399 bp in length, respectively, which is within the range of the sequenced mitogenomes of Lygaeidae in GenBank (Table 1). These mitogenomes all have a similar typical insect mitogenome structure, closed-circular and double-stranded DNA, containing 13 PCGs, 22 tRNAs, two rRNAs, and a control region (D-loop) (Fig. 2). The sequence of mitochondrial protein-coding genes is the same as that in other Lygaeoidea (Cao et al. 2020). Among the 37 genes, 23 genes (9 PCGs and 14 tRNAs) are on the majority strand (N-strand), while the remaining four PCGs, eight tRNAs, and two rRNA genes are on the minority strand (J-strand).

Figure 2. 

Circular maps of the complete mitogenome of Pylorgus species a P. porrectus b P. sordidus. Different colors indicate different types of genes and regions. Genes in the outer circle are located on the J-strand, and genes in the inner circle are located on the N-strand.

The basic composition of P. porrectus was A = 42.7%, T = 31.8%, G = 9.6%, and C = 15.8%, and of P. sordidus, A = 42.8%, T = 33.1%, G = 9.6%, C = 14.5%. Furthermore, both mitochondrial genome sequences were biased toward A and T. The AT content of P. porrectus was 63.74% and that of P. sordidus was 64.12%. The AT-skew value was greater than 0, whereas the GC skew value was less than 0, indicating that the base composition of P. porrectus and P. sordidus showed a strong A-bias and T-bias (Table 2).

The complete length of the 13 PCGs of P. porrectus and P. sordidus were 10,991 bp and 10,993 bp, respectively. Of these, nine PCGs are located at the N-strand, and the other four PCGs were encoded on the J-strand (Fig. 2). Most PCGs started with ATN except for cox1 and nad4l that began with TTG. Ten PCGs terminated with TAA/TAG, and the remaining three PCGs (cox1, cox2, and cox3) terminated with an incomplete T residue (Tables 3, 4). It has been speculated that these incomplete termination codons can be completed by adding ‘A’ during transcription (Ojala et al. 1981; Lavrov et al. 2002), and do not affect translation.

The RSCU of the two species was calculated (Fig. 3). The codons that were most used TTA-Leu and AGA-Arg. Most of the frequently used codons are composed of A and T, which may be related to the fact that the A-T skewness is higher than the G-C skewness in the PCGs of the two species.

Figure 3. 

RSCU values of Pylorgus species a P. porrectus b P. sordidus. The abscissa represents the type of amino acid translated by the codon, and the ordinate represents the codon bias score calculated for the amino acid. The higher the score, the more the types of codons, and the more active the evolutionary variation of genes in the genome.

The nucleotide diversity (Pi) of the two species based on 13 PCGs was computed (Fig. 4) and ranged from 0.05 to 0.11. Among the PCGs, nad3 (0.11) had the highest Pi values, and nad4l (0.05) had the lowest Pi values, which implies that nad4l is the most conserved gene in Pylorgus.

Figure 4. 

Nucleotide diversity (Pi) of 13 PCGs among two newly sequenced Pylorgus mitogenomes.

The ratios of Ka/Ks for each gene of the 13 PCGs were also computed (Fig. 5). All Ka/Ks values were less than 1 and ranged from 0.01 to 0.13, indicating that the genes have been subjected to purification selection. In particular, the Ka/Ks values were the highest for nad4 and nad5, suggesting that they had the highest evolution speed, and lowest for cox1, indicating the slowest evolution.

Figure 5. 

The ratios of Ka/Ks of 13 PCGs in the mitochondrial genomes of Pylorgus porrectus and P. sordidus.

Eleven gene overlaps were observed in the two mitogenomes, ranging from 1 bp to 25 bp (Tables 3, 4), and nad4l and trnT possessed the longest overlap.

Intergenic spacers were identified in the two mitogenomes, and their lengths ranged from 1 bp to 38 bp (Tables 3, 4). The longest intergenic spacer of 38 in P. sordidus was located between trnH and nad4.

The two mitogenomes both contain the complete set of 22 tRNA genes typical of Lygaeidae mitogenomes, ranging from 60 to 71 bp, which is consistent with previously sequenced mitogenomes of Lygaeidae (Cao et al. 2020; Huang et al. 2021). Fourteen of the 22 tRNAs were on the N-strand, and eight were on the J-strand (Fig. 2).

All tRNA have the typical cloverleaf secondary structure, including the TΨC arm, the amino acid acceptor arm, the anticodon arm, and the dihydrouridine arm. Some of tRNA genes (trnY, trnA, trnS1, trnF, trnH, trnP, and trnV) showed individual base mismatches, which is a common phenomenon in insect mitogenomes (Zhang et al. 2019).

The rrnL genes of the two mitogenomes are located at the intergenic region between trnL and trnV, with lengths that range from 1217 bp to 1221 bp. The rrnS genes are located between trnV and the D-loop, which are both 590 bp in length. Both rRNAs are located on the N-strand.

Phylogenetic relationships within Lygaeoidea were reconstructed based on mitochondrial 13 PCGs using BI and ML methods (Figs 6, 7). A total of 21 Lygaeoidea species were selected as the ingroup and an additional four species from Pyrrhocoroidea, Coreoidea, Rhopalidae, and Alydidae were used as the outgroup. Compared to the ML tree, the BI tree had higher confidence values, and the monophyly of all the studied families was supported except Rhyparochromidae.

Figure 6. 

Phylogenetic tree inferred from ML methods based on 13 PCGs. Nodal support is given as standard bootstrap (%); only values > 70% are shown. The newly sequenced Pylorgus porrectus and P. sordidus mitogenomes are highlighted in dark blue and bold.

Figure 7. 

Phylogenetic tree inferred from BI methods based on 13 PCGs. Nodal support is given as partition model; only values > 0.80 are shown. The newly sequenced Pylorgus porrectus and P. sordidus mitogenomes are highlighted using dark blue and bold typeface.

The clades making up the Lygaeidae had high support values in the BI results and confirmed the monophyly of Lygaeidae (Figs 6, 7). The monophyly of Lygaeidae was also supported in the ML results, but the nodal support is not so high. However, the Lygaeidae clusters as sister to Malcidae in the ML tree, but sister to Geocoridae in the BI tree, implying that the positions of Geocoridae and Malcidae are unstable. The two species of Rhyparochromidae are not clustered together. Neolethaeus assamensis clusters as sister to the Pyrrhocoroidea species, and together they are sister to the remaining ingroups.

No conflict of interest was declared.

No ethical statement was reported.

This research was funded by Major Project of Agricultural Biological Breeding (grant no. 2022ZD0401501), the National Natural Science Foundation of China (grant no. 31402010), and the Highly Educated Talents Foundation in Nanjing Forestry University (grant no. G2014002).

Conceptualization, C.G. and W.D.; methodology, C.G. and W.D.; investigation, C.G. and W.D.; funding acquisition, C.G.; writing—original draft preparation, W.D.; writing—review and editing, C.G. Both authors have read and agreed to the published version of the manuscript.

Cuiqing Gao https://orcid.org/0000-0002-0177-5161

Wen Dong https://orcid.org/0009-0004-6559-808X

All of the data that support the findings of this study are available in the main text or Supplementary Information.