For both species, we obtained abdominal tissues of third-instar larvae at 5 days after burrowing into the soil for pupation (i.e., prepupae, about 2 days before pupation) and the genital parts of pupae at 2 days after pupation. Development of genital discs was assumed to have started during the prepupal stage, and the copulatory piece and vaginal appendix were assumed to be formed from 2 days after pupation. The tissues of abdominal segments A9 to A11 of larvae and genital tissues of pupae were dissected, fixed in RNAlater, and stored in a freezer until RNA extraction. The sex of a larva was determined by the length difference of the polymerase chain reaction products of dsx gene due to the sex-specific isoforms. The sex of a pupa was determined by the external morphology of the abdominal tip. For each species and sex, three individuals were obtained (total, 24 individuals; table S2). Total RNA was extracted using a Qiagen RNeasy kit. For each sample, paired-end libraries with 170-bp insertions were constructed and sequenced on the Illumina HiSeq 2000 platform at Beijing Genomics Institute (Shenzhen, China).

To identify DEGs between species in each sex and stage, we performed a differential expression analysis (see Supplementary Methods for details). To detect the predominant functional categories of the DEGs, we performed an enrichment analysis for each sex and stage using DAVID Bioinformatics Resources (52) with the D. melanogaster data as the background.

We checked whether the candidate genes in the species-monophyly regions detected by ARGweaver analysis were DEGs. We also listed genes involved in the basic development of insect genitalia and appendage-patterning genes tested in the development of various appendages in beetles. We determined whether these genes were differentially expressed (see table S6 for details).