Download the carrot genome annotation file and the total sequence file in the carrot genome database (https://plants.ensembl.org/Daucus_carota/Info/Index, genome assembly: GCA001625215.1). The Arabidopsis TCP gene family protein sequence was downloaded from the TAIR website (https://www.arabidopsis.org/). The protein sequences of rice TCP gene family were downloaded from RGAP website (https://rice.uga.edu/). The protein sequences of celery TCP gene family were downloaded from celery website (http://celerydb.bio2db.com). The hidden Markov model file of PF03634 was downloaded from the Pfam website (http://pfam.xfam.org/) as a reference, and the HMMER v3.3.2 search was used to search the carrot genome to obtain the carrot TCP protein sequence. Using Arabidopsis thaliana and rice as reference species, Blastp was used to obtain homologous sequences in carrot (E-value < 10-10), then redundant and incomplete sequences were removed, and finally obtained accurate carrot TCP candidate genes. The physicochemical properties of TCP protein, including amino acid number, molecular weight, isoelectric point, instability index, aliphatic index and hydrophilicity, were analyzed by Protein parameter Calc of TBtools. The subcellular localization of TCP protein was predicted by WoLF-PSORT website (https://wolfpsort.hgc.jp/).

The chromosome location information was obtained from the carrot genome annotation file. The chromosome distribution map of TCP gene family members was drawn by TBtools, and the carrot TCP gene family (DcTCPs) were named according to its position on the chromosome. The homology analysis and collinearity analysis of carrot TCP gene family (DcTCPs) were carried out by using Advanced Circos in TBtools software (Chen et al., 2020). The ratio of non-synonymous substitution rate (Ka) to synonymous substitution rate (Ks) was calculated to evaluate the selection pressure of duplicated gene pairs.

MEGA 11 software was used for phylogenetic analysis. MUSCLE was used for sequence alignment of TCP proteins from carrot, Arabidopsis and rice. Neighbor-Joining algorithm was used to construct phylogenetic tree (Bootstrap = 1000), and the Poisson model for amino acid substitution. The resulting phylogenetic tree was visualized through the online tool iTOL (https://itol.embl.de). The carrot TCP protein sequence was submitted to the MEME online tool (https://meme-suite.org/meme/tools/meme) for conservative motif analysis (motif number = 10, motif width = 10–300 aa). The annotation information of carrot TCP gene was submitted to TBtools software for intron and exon analysis, and the gene structure was visualized.

The TCP gene promoter sequence was submitted to the PlantCare website (https://bioinformatics.psb.ugent.be/webtools/plantcare/html/) for cis-acting element analysis using the 2000 bp upstream of the gene initiation site ATG as the promoter region of the gene. GO functional enrichment analysis of TCP genes was performed using the EggNOG5.0 database (http://eggnog5.embl.de/#/app/home) and visualized using TBtools.

The interaction miRNAs of carrot TCP family members were predicted by psRNATarget online software (https://www.zhaolab.org/psRNATarget/). The Expect value was set to 5, and Arabidopsis miRNAs and some carrot miRNAs were used as reference. The String online website (https://www.string-db.org/, version 11.5) was used to predict the protein interaction network of carrot TCP protein with a confidence score threshold of 0.7, and Arabidopsis thaliana was set as the reference species.

Transcript abundance analysis during carrot root development Based on the transcriptome data (SRR 2177455), RPKM (Reads Per Kilobase per Million mapped reads) was used as an indicator to measure the transcript or gene expression level. The heat map was drawn by Multiexperiment Viewer software (https://mev.tm4.org/).

The expression levels of DcTCP genes in carrot root at different growth stages (30d, 60d, and 90d) were determined using Hieff qPCR SYBR Green Master Mix (Yeason, Shanghai, China). Design specific primers using Primer Premier 6.0 software (Supplementary Table S1), DcActin was used as the reference gene, and its expression stability was verified using geNorm and NormFinder algorithms (Wang et al., 2016). The total reaction volume for RT-qPCR was 20 µL, which included 10 µL of SYBR Premix Ex Taq, 2 µL of cDNA template, 7.2 µL of ddH2O, and 0.4 µL of each forward and reverse primer. The thermal cycling program was as follows: initial denaturation at 95°C for 5 min, followed by 40 cycles of denaturation at 95°C for 10 s and annealing at 60°C for 30 s. The relative expression levels of the genes were calculated using the 2^−ΔΔCT method (Pfaffl, 2001). Each sample was performed for three biological replicates. Statistical analysis was performed using SPSS 25.0 software, and significant differences were determined by Duncan’s multiple range test at p < 0.05.

A total of 50 TCP gene family members were identified in the carrot genome and designated as DcTCP1 to DcTCP50 in sequence according to their distribution on chromosomes (Table 1; Supplementary Table S2). The 50 identified DcTCP family genes encode proteins with amino acid lengths ranging from 55 to 496 residues. Their molecular weights span 6056 to 53,512.5 Da, with an average of 28,857.83 Da. The theoretical isoelectric points (pI) of these proteins ranged from 4.73 to 9.8, and their instability indices varied from 33.81 to 65.9—with 90% of the proteins classified as unstable. The aliphatic indices of the DcTCP proteins fall between 47.38 and 83.43, and all proteins exhibited a grand average of hydropathicity (GRAVY) value less than 0, indicating they were hydrophilic. Subcellular localization prediction revealed that all DcTCP family members were localized in the nucleus.

To further clarify the chromosomal distribution of DcTCP genes in carrot, a chromosomal localization map of the DcTCP gene family was constructed using the bioinformatics tool TBtools (Figure 1A). Analysis results indicated that the 50 DcTCP genes were unevenly distributed across all 9 carrot chromosomes, and these genes were designated as DcTCP1 to DcTCP50 based on their sequential positions on the chromosomes. Among the 9 chromosomes, chromosome 1 harbored the largest number of DcTCP genes (16 in total), whereas chromosomes 5 and 9 contained the fewest, with only 2 DcTCP genes each. For duplication event classification, genes located on the same chromosome with an intergenic distance of less than 200 kb were defined as tandem duplicates; all other duplicated gene pairs were classified as segmental duplicates (Cheung et al., 2003). To compare the evolution of TCPs among plants, we conducted for TCP comparisons of 9 species and constructed schematics of plant evolution (Figure 1B). Higher plants have far more members of the transcription factor family than lower plants. The number of transcription factor families varies greatly among different species, the most abundant being PCF, followed by CIN. The carrot TCP family has 50 members, and the PCF subfamily has 36 members, accounting for 72% of the family.

The results revealed that the DcTCP gene family contained 22 pairs of tandem duplicated genes and 8 pairs of segmental duplicated genes (Figure 2A). Notably, DcTCP29 formed duplicated pairs with two distinct members, namely DcTCP30 and DcTCP36; all other duplicated genes in the family exhibited a one-to-one pairing pattern. The Ka/Ks of 8 pairs of segmental duplicated genes are less than 1, which were purify selection (Supplementary Table S3). Whole-genome duplication (WGD) was the predominant type in carrot (32.0%), These results demonstrate that WGD played an important role in TCP gene expansion in carrot, which is supported by the previous suggestion that they underwent two WGD events since their divergence from lettuce (Song et al., 2020).

Additionally, to clarify the evolutionary relationships of TCP genes, the synteny analysis were performed between carrot and two model plants, Arabidopsis thaliana (Figure 2B) and rice (Figure 2C). Results showed that carrot shares 29 pairs of syntenic genes with A. thaliana, whereas only 4 pairs of syntenic genes were identified between carrot and rice. This observation indicated that carrot exhibits higher TCP gene homology with A. thaliana than with rice.

To investigate the phylogenetic relationships of TCP genes between carrot and rice, a phylogenetic tree was constructed using the maximum likelihood (ML) method. This tree included 50 DcTCP genes from carrot, 27 AgTCP genes from celery and 24 AtTCP in Arabidopsis thaliana (Figure 3). Phylogenetic analysis revealed that TCP genes were clustered into three distinct subfamilies—PCF, CIN, and CYC/TB1—consistent with the subfamily classification of TCP genes in other plant species. Among these subfamilies, the PCF subfamily contained the largest number of members, with 36 DcTCP genes, 16 AgTCP genes and 13 AtTCP genes. The CYC/TB1 subfamily had the fewest members, with a total of 17 genes. We identified three TCP genes—namely, DcTCP34 and AgTCP2, —that clustered together with AtTCP3, suggesting that they are also related to cotyledon fusion (Koyama et al., 2010). Overall, the TCP family members from carrot, Arabidopsis thaliana and celery showed similar evolutionary patterns, suggesting that the plant TCP gene family may have originated from a common ancestral gene and retained conserved evolutionary characteristics during speciation.

The diversity of exon/intron structures and protein domain architectures played a crucial role in the evolution of gene families. For the 50 DcTCP family members identified in carrot, most (47 out of 50) contained only one exon; only DcTCP45 has two exons, and DcTCP11 has three exons. None of the DcTCP genes harbor introns, indicating a relatively simple gene structure. Further analysis of the protein domains encoded by the DcTCP family revealed that all 50 DcTCP proteins contained the canonical TCP domain—a signature feature of the TCP gene family (Supplementary Figure S1), confirming their classification as TCP family members.

To explore conserved motifs among DcTCP proteins, the MEME (Multiple Em for Motif Elicitation) tool was used to identify 10 conserved motifs (Figure 4). Most DcTCP proteins (42 out of 50) contained at least four of these motifs. Among the 10 motifs, Motif 1 was the most widely distributed, presented in all DcTCP proteins except DcTCP9 and DcTCP13. Motif 2 was detected in 31 DcTCP proteins, while Motif 6 was the least abundant, present in only 5 DcTCP proteins. Additionally, spatial analysis showed that most of these conserved motifs were localized within the first 200 amino acid (aa) residues of the DcTCP proteins.

PlantCARE analysis revealed the presence of abundant responsive regulatory elements within the promoter region of the DcTCP genes (Figure 5). The identified cis-acting elements were categorized into four distinct groups based on their functional annotations. The first group comprised hormone-responsive elements, including but not limited to the CGTCA-motif, P-box, TCA-element, ABRE, TGA-element, and AuxRR-core. The second group consisted of photoresponse-associated elements, such as the I-box, G-box, ACE, and MRE. The third group encompasses growth and development regulatory elements, which included the CAT-box, O2-site, CCAAT box, and circadian element. The fourth group was defined as stress-responsive elements, featuring the ATBP-1, LTR, TC-rich repeats, MBS, and MBS I. Further comparative analysis demonstrated that each DcTCP promoter contained conserved elements, including MYB, P-box, CAAT-Box, and LTR. Notably, promoters of members belonging to the carrot PCF subfamily harbored all the predicted cis-acting elements. In contrast, elements such as GCN4-motif, CCAAT-box, and ACE were not detected in promoters of the CYC/TB1 subfamily. These findings collectively suggested that PCF subfamily members in carrot exhibit more extensive regulatory functions.

To more comprehensively delineate the potential functions of TCP genes in carrot, detailed functional annotation and classification were performed using the EggNOG 5.0 tool (Figure 6). Functional enrichment analysis revealed that DcTCPs were involved in a broad spectrum of biological processes, including regulation of morphogenesis of a branching structure, Regulation of secondary shoot formation, and inflorescence development. In terms of molecular functions, DcTCPs primarily exhibited binding activity, and transcription factor activity, among others. For cellular component classification, nucleus was identified as the major categories.

As a crucial family of transcription factors, the TCP gene family exerted its core functions primarily through protein-protein interactions. As illustrated in Figure 7, protein-protein interaction (PPI) network prediction was performed using Arabidopsis thaliana as a reference organism. The results demonstrated that all 28 member proteins of the DcTCP family were capable of interacting with other proteins, leading to the identification of a total of 214 protein pairs. Notably, DcTCP4, DcTCP3, DcTCP14, and DcTCP10 exhibited the highest number of interacting proteins. cTCP33 and DcTCP37 maintained high expression levels during the entire root development stage (30 d to 95 d), which corresponds to the key period of carrot root hypertrophy and nutrient accumulation; this suggests that these two genes may regulate cell expansion and storage substance accumulation in fleshy roots. Collectively, these findings suggested that these specific DcTCP proteins may participated in diverse biological regulatory mechanisms by recruiting or interacting with a larger repertoire of target proteins.

MicroRNAs (miRNAs) are key regulatory molecules in plants, primarily mediating post-transcriptional gene expression regulation. To identify miRNAs targeting the coding regions of DcTCP genes, prediction analysis was performed using the psRNA Target online tool, with detailed results summarized in the Table 2 (Supplementary Table S5). A total of 12 miRNAs were identified to be involved in the regulation of DcTCP genes. Among these, miR319 exhibited the highest number of target DcTCP genes (7 targets), followed by miR172 (6 targets). In contrast, miR5072 had the fewest targets, regulating only the DcTCP34 gene. Regarding the mechanism of action, the predicted regulatory mode of most identified miRNAs on their target DcTCP genes was direct cleavage. Notably, miR402 was the only miRNA predicted to exert its regulatory effect through translational inhibition.

By analyzing and calculating the transcriptome data, the transcriptional abundance of the DcTCPs gene in carrot during the development process was obtained. Thirty-three DcTCPs gene expressions were discovered. As shown in Figure 8, DcTCPs were widely expressed at different developmental stages. DcTCP33, DcTCP37 and DcTCP47 showed relatively high expression levels throughout the entire developmental stage of carrot, while DcTCP39 and DcTCP48 exhibited relatively low expression levels throughout the developmental stage. DcTCP4, DcTCP7, DcTCP10, DcTCP11, DcTCP15, DcTCP16, DcTCP25, DcTCP27, DcTCP28, DcTCP36, DcTCP42 and DcTCP45 gradually increased with the stage of development. We selected six genes for RT-qPCR verification, and the results were consistent with the transcriptome data (Figure 8B). It indicated that these genes may be involved in the formation and development of carrot root.