Two weeks after germination, we maintained the pepper (C. annuum cv. CM334) seedlings in a growth room at 24°C ± 1°C with a 16-h light and an 8-h dark photoperiod. At the six-true-leaf stage, we subjected young plants to a temperature of 10°C or 40°C to mimic cold or heat stress. The third or fourth leaves were harvested per replicate at 0, 3, 6, 12, 24, and 72 h after treatment. We collected three biological replicates at each time point per condition. The leaves were flash-frozen in liquid nitrogen and stored in an ultra-low-temperature refrigerator at −80°C for RNA extraction and quantitative real-time polymerase chain reaction (qRT-PCR).

We obtained the genome data of C. annuum cv. CM334 from the pepper website ¹ and downloaded other genome files, such as Arabidopsis thaliana TAIR10, Solanum lycopersicum ITAG4.0, and Solanum tuberosum v6.1, from Phytozome V13. ² The P450 family gene seed file (PF00067) was downloaded from the Pfam online database. ³ We used HMMER 3.0 (Eddy, 2011) software and the local BLASTP program to compare and search sequences containing the P450 protein domain. We uploaded the protein sequences to the online software Pfam, the Simple Modular Architecture Research Tool database (SMART), ⁴ and the Conserved Domain Database (CDD). ⁵ We performed sequence alignment and analysis to remove unannotated genes and redundant sequences. Finally, we obtained all members of the P450 gene family in four species and named the pepper P450 genes according to their positions on the chromosome. The amino acid number, molecular weight (MW), and isoelectric point (PI) were obtained from the ExPASy website, ⁶ and the subcellular localization information was subsequently obtained from the CELLO website. ⁷ We downloaded the Arabidopsis P450 family gene members from the Arabidopsis genome database TAIR website. ⁸

Because of the huge difference in the sequence of P450 proteins in pepper, we selected a protein length longer than 300 amino acids to construct the phylogenetic tree. We used Clustal X to perform multiple sequence alignments of the obtained genes from four species. The phylogenetic tree included multiple plant species and was constructed using MEGA 11 with the NJ method, and the selected model was the Poisson model with 1,000 bootstrap replications (Kumar et al., 2016). We used the ITOL ⁹ web server to beautify the resulting phylogenetic tree. We used the software blastp to reclassify the proteins with a sequence length of less than 300 amino acids according to the constructed phylogenetic tree. We confirmed P450 protein motifs by the Multiple Expectation Maximization for Motif Elicitation (MEME) ¹⁰ program with the motif number set to 10, and all the other parameters were set to default (Bailey et al., 2009). TBtools was used for visualization (Chen et al., 2020).

The chromosome physical location information of the pepper P450 gene family was extracted from the pepper genome annotation file information and visualized with TBtools. We analyzed the tandem duplication events of the P450 gene family using TBtools with MCScan X (Wang et al., 2012; Chen et al., 2020). In addition, we investigated segmental duplication events and the collinearity of gene pairs from different species using TBtools with the MCScanX and BLASTP methods (Wang et al., 2012; Chen et al., 2020). We calculated non-synonymous (Ka) and synonymous (Ks) substitutions between gene pairs using TBtools.

We downloaded the transcriptome data of pepper root, stem, flower, fruit development, and different periods under various viral infections (PRJNA223222), abiotic stresses (PRJNA525913), and hormone treatments (PRJNA634831) from the NCBI Sequence Read Archive (SRA) database ( Supplementary Table 1 ) (Kim et al., 2014; Kang et al., 2020; Lee et al., 2020). We used the genome C. annuum cv. CM334 as the reference genome. We calculated the FPKM value (i.e., fragments per kilobase of the exon model per million mapped reads) using HISAT2 software for alignment and StringTie (v2.1.7) software for the assembly and quantitation with default parameters. We used log (FPKM+1) to calculate the expression level and imported the results into the R package “pheatmap” to draw the expression level of heatmap with the following parameters: scale = “raw,” cluster_rows = TRUE, and cluster_cols = FALSE.

Total RNA was extracted from pepper leaves in various stages after cold and heat treatment of C. annuum cv. CM334 using TIANGEN’s polysaccharides- and polyphenolic-rich RNAprep pure plant plus kit. The online website Primer-Blast ¹¹ was used to design primers for the qRT-PCR experiment ( Supplementary Table 2 ). We used HiScript IIQRT SuperMix from Vazyme to synthesize first-strand cDNA. We performed qRT-PCR with an ABI 7500 real-time PCR system using ChamQ Universal SYBR qPCR Master Mix from Vazyme. We also used the pepper ubiquitin-conjugating protein (CaUbi3) as a control for normalization between samples (Qin et al., 2018). We calculated relative transcript levels using the 2^−ΔΔCT method.

According to previous reports, the CYPome (i.e., the CYP complement of a given species) from angiosperms consists of approximately 300 genes (Nelson and Werck-Reichhart, 2011). In this study, we identified 478 P450 genes from the genome of C. annuum cv. CM334, which we named from CaCYP1 to CaCYP478 according to their location on chromosomes ( Figure 1 ; Supplementary Table 3 ). In the identified proteins encoded by P450 genes, the smallest protein (CaCYP356) and the largest protein (CaCYP154) contained 53 and 1,321 amino acids, with molecular weights of 6.11 and 150.85 kDa, respectively. The theoretical isoelectric point of the coding protein ranged from 4.21 to 11.18. For subcellular localization, 11 CaCYP proteins were predicted to locate in the chloroplast, 203 CaCYP proteins in the cytoplasm, 10 CaCYP proteins in the extracellular membrane, 79 CaCYP proteins in the mitochondria, 32 CaCYP proteins in the nucleus, and 143 CaCYP proteins in the plasma membrane ( Supplementary Table 4 ).

Because of the significant difference in the length of the P450 genes in pepper, the phylogenetic tree could not successfully include all P450 genes. To investigate the evolutionary relationship of the P450 genes, 335 proteins were encoded by P450 genes for which a protein length longer than 300 amino acids was used to construct the phylogenetic tree. Finally, we constructed a phylogenetic tree containing 1,305 P450 genes from four species (A. thaliana, C. annuum, S. lycopersicum, and S. tuberosum) using the MEGA software ( Figure 2 ; Supplementary Figure 1 ). We divided all of the P450 genes into one of the two major clades: A type, which contained the CYP71 clan; and non-A type, which contained eight clans. The remaining 143 CaCYP proteins with a length less than 300 aa were classified through comparison of protein sequence with the P450 protein used for the construction of the phylogenetic tree. Finally, we identified 478 CaCYP genes in the pepper genome. The number in each clan showed significant differences. For example, the CYP71 clan included 305 (63.81% of the total P450 genes in pepper) CaCYP genes, and the CYP51 or CYP711 clan contained only one (0.21%) of the CaCYP genes ( Supplementary Table 3 ).

To better understand their similarity and diversity, we investigated the conserved motif and gene structure using MEME software. We observed a total of 10 putative different motifs in the P450 genes of pepper ( Supplementary Figure 2 ). Almost all of the P450 genes contained at least one motif and the greatest number of motifs was 17 ( Supplementary Figure 3 ). We did not find, however, any motifs in the protein sequence of a few genes, including CaCYP219, CaCYP354, CaCYP468, and CaCYP418. Most of the members in the pepper P450 family contained eight motifs (motif 4, motif 7, motif 2, motif 6, motif 3, motif 5, motif 1, and motif 8) in sequence near the 3′ UTR regions. In addition, the motif composition of each subfamily member was different, but the motif composition of the same subfamily genes was similar. Moreover, compared with other subfamilies, we found that the CYP71 clan contained a specific motif 9. These results showed that the protein structure of each group of CaCYP genes was highly conserved. Except for motif analysis, we also analyzed gene structure and found that variations in the gene structure were also significant. The coding sequence (CDS) length of CaCYP genes ranged from 100 to more than 30,000 bp, and the number of introns or exons ranged from 1 exon and 0 introns to 16 exons and 15 introns, demonstrating that significant variation existed among the gene structure of CaCYP genes.

We identified 417 CaCYP genes distributed on 12 chromosomes, and another 61 CaCYP genes were unmapped scaffolds ( Figure 1 ). Most of the CaCYP genes were concentrated at the terminal of each chromosome, with fewer exceptions of genes located in the intermediate of the chromosome. The number of CaCYP genes on each chromosome was slightly different, and the number ranged from 11 located on chromosome 8 to 72 located on chromosome 3.

According to the duplication event analysis at the genome-wide level, 106 pairs of tandem duplicated genes were discovered to be distributed across the 12 chromosomes ( Figure 1 ). The number of tandem duplicated genes was different on each chromosome and 17, 5, 20, 15, 5, 7, 14, 3, 3, 10, 2, and 5 pairs of tandem duplicated genes were distributed on chromosomes from 1 to 12 ( Supplementary Table 3 ). At the same time, we also observed 23 pairs of segmental duplicated genes on the pepper chromosomes ( Figure 3 ; Supplementary Table 5 ). Moreover, we also evaluated the pressure of a selective constraint on each pair of duplicated CaCYP genes by calculating the Ka and Ks ratios ( Supplementary Table 6 ). The results showed that every gene pair had a Ka/Ks ratio lower than 1, which indicated that the P450 gene family may have undergone strict purification selection during evolution. This result suggested that replicators were evolutionarily conserved and structurally stable and may have had consistent functions.

To further explore the evolutionary relationships of the P450 gene family, we analyzed the collinearity among pepper and other species (tomato, potato, and Arabidopsis) ( Figure 4 ). According to the collinearity analysis between Arabidopsis and pepper, we identified a total of 50 pairs of genes, which indicated that the P450 gene family was significantly amplified before the differentiation of the two species. We observed 160 pairs of collinear genes between tomato and pepper and 178 pairs of collinear genes between potato and pepper ( Supplementary Table 7 ). These results showed that the P450 genes of tomato or potato and pepper were highly evolutionarily conserved. Collinearity analysis showed that CaCYP77 had four pairs of homologous genes and CaCYP325 had three pairs of homologous genes between pepper and tomato. In pepper and potato, CaCYP77 and CaCYP325 each had three pairs of homologous genes. This result indicated that these homologous genes may have played a vital role in the evolutionary process. In addition, some pepper P450 genes matched one gene from potato or tomato, such as CaCYP348, CaCYP325, CaCYP351/Soltu.DM.10G024690.1, CaCYP77, CaCYP178, CaCYP183, and CaCYP411/Solyc04g078270.3.1, which indicated that these genes may have expanded during pepper evolution. Moreover, a large number of genes in three Solanaceae species (e.g., pepper, tomato, and potato) had the same homologous genes, indicating that these plants had great homology.

We investigated the expression profile of CaCYP genes in various tissues and different fruit developmental stages according to published data (Kim et al., 2014). The results showed that in the CYP71 clan, some genes, such as CaCYP57, CaCYP63, CaCYP186, CaCYP113, and CaCYP461, showed high expression at placenta development before the fruit mature stage [at stage 1, stage 2, stage 3, and mature green (MG)], and some genes, such as CaCYP40, CaCYP78, CaCYP175, CaCYP240, CaCYP315, and CaCYP352, showed high expression at pericarp development before the fruit mature stage (at stage 1, stage 2, stage 3, and MG) ( Figure 5 ; Supplementary Figure 4 ). In addition, we also observed high-expression genes both in the placenta and pericarp development after the fruit mature stage (breaker, breaker plus 5 days, and breaker plus 10 days) in the CYP71 clan, which showed that the CYP71 clan may have been significantly involved in fruit development. Except for the CYP71 clan, other clans also were found to be involved in this process. For example, CaCYP56, CaCYP292, CaCYP295, and CaCYP294 in the CYP72 clan were related to placenta development after the fruit mature stage. We also found that some genes were related to pericarp development after the fruit mature stage. For example, CaCYP310, CaCYP303, CaCYP305, and CaCYP309 in clan 72 were related to pericarp development before the fruit mature stage, whereas CaCYP296, CaCYP273, CaCYP291, and CaCYP52 were related to pericarp development after the fruit mature stage. In the CYP74 clan, six CaCYP genes (66.7% of the total CYP74 clan) were related to the pericarp development before the fruit mature stage. In the CYP85 clan, CaCYP456, CaCYP208, CaCYP173, and CaCYP222 were related to placenta development. CaCYP95, CaCYP30, CaCYP275, CaCYP83, CaCYP51, CaCYP81, and CaCYP332 were related to pericarp development. Interestingly, all of the genes in the CYP97 clan were expressed at pericarp development, which showed that this clan may be connate for pericarp development.

To further explore the function of CaCYP genes, we analyzed the RNA-seq data of pepper leaves after cold, heat, sodium chloride (NaCl), and D-mannitol treatment for 0, 3, 6, 12, 24, and 72 h (Kang et al., 2020) to explore the function of CaCYP genes in abiotic stresses. Most of the CaCYP genes were expressed under the four abiotic stresses, indicating that the CaCYP genes were closely related to abiotic stress ( Figure 6 ; Supplementary Figure 5 ). In the CYP71 clan, 70 genes had high-expression levels under cold treatment, of which 11 genes had specific high expression at 72 h cold treatment, and 49 genes were expressed at almost all time points of cold treatment. Additionally, CaCYP15, CaCYP137, CaCYP316, CaCYP219, and CaCYP325 had higher expression at 3 h of heat treatment. CaCYP135, CaCYP32, CaCYP357, and CaCYP154 had higher expression at 12 h of heat treatment. CaCYP114, CaCYP115, CaCYP106, CaCYP422, CaCYP62, and CaCYP120 were highly expressed at 12 h after D-mannitol treatment, whereas CaCYP186, CaCYP130, CaCYP263, CaCYP462, CaCYP421, and CaCYP364 were relatively highly expressed under NaCl treatment. In addition, a total of 16 genes in the CYP72 clan were associated with cold treatment, and most of the remaining genes were expressed under the other three abiotic stresses. Three genes were related to cold treatment and five genes were associated with D-mannitol treatment and NaCl treatment in the CYP74 clan. In the CYP85 clan, eight genes were related to cold treatment, 13 genes were related to heat treatment, and seven genes had high expression levels in D-mannitol treatment and NaCl treatment. In the CYP86 clan, 10 genes were associated with cold stress, whereas CaCYP465, CaCYP128, and CaCYP466 had very high expression levels under NaCl treatment. Finally, we found that some genes had similar expression patterns under heat, D-mannitol, and NaCl treatment, indicating these genes may have been involved in the production of cellular secondary metabolites to regulate cellular osmotic pressure. Moreover, we validated the RNA-seq results by qRT-PCR analysis of 12 CaCYP genes that responded to cold and heat treatment of pepper P450 genes. These selected genes were analyzed for relative transcript abundance ( Figure 7 ). The results showed that the expression trend of these genes was highly consistent with that of RNA-seq analysis. All of these genes were highly expressed under cold or heat treatment. The expression levels of CaCYP328 and CaCYP30 were higher than those of the control after cold stress at 3 h. Heat stress induced the expression of CaCYP225 at 6 h. CaCYP454, CaCYP216, and CaCYP182 showed the highest expression levels after 12 h of cold or heat stress. CaCYP60, CaCYP399, and CaCYP295 transcripts were elevated in response to cold and heat treatment at 24 h. Finally, CaCYP108, CaCYP159, and CaCYP245 showed the highest levels after 72 h of cold or heat stress. These results indicated that the P450 family genes may play an important role in stress tolerance in pepper.

To further explore the function and expression pattern of CaCYP genes, we analyzed the RNA-seq data of pepper leaves treated with four hormones, including methyl jasmonate (MeJA), salicylic acid (SA), ethylene (ET), and abscisic acid (ABA) (Lee et al., 2020) ( Figure 8 ; Supplementary Figure 6 ). High expression levels of CaCYP genes were observed under four hormone treatments. In the CYP71 clan, 11, 5, and 6 CaCYP genes had extremely high expression at 1, 6, and 12 h of MeJA treatment, respectively. We identified 52 genes related to SA treatment, of which 11 genes had higher expression levels during the early stage, 30 genes had higher expression levels during the late stage, and the remaining 11 genes had higher expression levels in every stage after SA treatment. A total of 56 genes were upregulated after ET treatment. Among them, CaCYP244, CaCYP28, and CaCYP64 had higher expression levels at 1 h. CaCYP22, CaCYP325, CaCYP137, CaCYP167, CaCYP101, CaCYP106, and CaCYP422 had higher expression levels at 6 h. CaCYP155, CaCYP42, CaCYP277, CaCYP246, and CaCYP405 had higher expression levels at 12 h. A total of 17 genes were associated with SA in the CYP72 clan, and five and seven genes were associated with MeJA in the CYP74 clan and CYP86 clan, respectively. In addition, eight genes of the CYP85 clan were associated with ET. Moreover, 25 genes also were related to ABA treatment from 3 to 12 h, including CaCYP378, CaCYP384, CaCYP453, and CaCYP69. The analysis results showed that the P450 family genes responded to stress resistance-related hormones, which may be involved in the stress tolerance of pepper.

To explore the potential function of CaCYP genes under biotic stress, we analyzed the RNA-seq data of leaves infected with the tobacco mosaic virus (TMV) P0 strain and pepper mottle virus (PepMoV) at different infection times (Kim et al., 2014) ( Figure 9 ; Supplementary Figure 7 ). The results showed that in the CYP71 clan, 68 genes showed slightly higher expression after the virus infected the leaves for 2 and 3 days compared with the control. Interestingly, these genes showed a strong response to TMV after the virus infected the leaves for 2 and 3 days. Except for the CYP71 clan, other clans also showed higher expression after TMV and PepMoV infection, including a moderate response to PepMoV after incubation; however, these genes showed a strong response to TMV. In the CYP72 clan, approximately 17 CaCYP genes significantly increased in the presence of virus infection, among which the expression level was higher in TMV infection. Eight CaCYP genes were expressed only in the presence of PepMoV infection. Three genes in the CYP85 clan showed higher expression levels after 2 and 3 days of infection. The results showed that the P450 family genes may be closely related to plant disease resistance.