Coordinated circRNA/lncRNA-Mediated Post-Transcriptional Regulation Balances Anthocyanin and cis-Zeatin Metabolism in Lycium ruthenicum

1 Introduction

Lycium ruthenicum Murray (black goji berry), a perennial shrub that belongs to the Solanaceae family, is an economically important medicinal plant endemic to arid and semi-arid regions of Central Asia [1]. Its fruits are renowned for their exceptionally high content of bioactive compounds, particularly anthocyanins, flavonoids, alkaloids, and amino acids. Petanin—(petunidin-3-O-rutinoside-(p-coumaroyl)-5-O-glucoside) is the primary anthocyanin component, accounting for over 80% of the total pigment content, alongside derivatives such as cyanidin, pelargonidin, and kaempferol aglycones [2,3]. These structurally unique compounds endow the berries with numerous biological activities. Anthocyanins are water-soluble pigments that demonstrate potent antioxidant properties and exhibit protective effects in preclinical models of Alzheimer’s disease, inflammatory bowel disease, and cardiovascular disorders [4,5,6]. Therefore, elucidating the molecular mechanisms underlying anthocyanin accumulation in black goji is necessary for harnessing the plant’s pharmacological potential and implementing molecular breeding strategies to enhance fruit quality. Furthermore, this work offers a unique perspective for understanding the adaptive evolutionary mechanisms of secondary metabolism regulation in plants under abiotic stress conditions.

Anthocyanin biosynthesis is a branch of the flavonoid metabolic pathway. The synthesis is catalyzed by a coordinated series of enzymes encoded by core structural genes, including phenylalanine ammonia-lyase (PAL), chalcone synthase (CHS), and chalcone isomerase (CHI) [7,8,9]. At the transcriptional level, this process is primarily regulated by the conserved MYB-bHLH-WD40 (MBW) protein complex [10]. Zong et al. [11] reported that the content of the MYB transcription factor LrAN2 was 35.66-fold higher in black goji than in Lycium barbarum. They also found that the functional diversity and high expression of LrAN2 likely contribute to anthocyanin accumulation. In tomato (Solanum lycopersicum), a related species, the MYB transcription factor SlAN2-like promotes anthocyanin accumulation in the peel and flesh of ‘Indigo Rose’ fruit exposed to light by activating the anthocyanin biosynthetic pathway and enhancing the expression of the bHLH protein SlAN1 [12]. Beyond transcriptional regulation, phytohormones serve as crucial signaling molecules. They are key regulators of anthocyanin synthesis, with such hormonal control exhibiting differences according to the species and tissue. For instance, auxin levels markedly influence PAP1 expression, in transgenic tobacco callus overexpressing the Arabidopsis PAP1 (At1g56650) gene, modulating the transcription factor MYB75 and ultimately altering anthocyanin content [13]. Abscisic acid (ABA) and light stimulate the expression of the transcription factor LcMYB1 in the pericarp of litchi (Litchi chinensis) by interacting with cis-elements in its promoter region [14]. This, in turn, activates the downstream structural gene LcUFGT to drive anthocyanin synthesis. By contrast, the synthetic cytokinin forchlorfenuron (CPPU) suppresses both LcMYB1 expression and anthocyanin accumulation [14].

However, the regulation of plant secondary metabolism is a multi-layered and networked process. The precise balance of secondary metabolism relies not only on the control by protein-coding genes and hormonal signals at the transcriptional level but also on the coordinated action of post-transcriptional regulatory networks [15]. In recent years, non-coding RNAs (ncRNAs) have emerged as key players in transcriptional as well as post-transcriptional regulation and have been demonstrated to participate in anthocyanin synthesis in various plant species [16,17,18]. Nevertheless, this regulatory dimension remains unexplored in black goji. Unlike protein-coding mRNAs, ncRNAs lack direct coding capacity but can modulate gene expression through diverse mechanisms. Among these, circular RNAs (circRNAs) are a class of single-stranded RNA molecules characterized by a covalently closed circular structure. The absence of free ends renders them more stable than their linear counterparts [19]. Their primary regulatory modes include functioning as miRNA sponges and interacting with RNA-binding proteins (RBPs), enabling precise post-transcriptional gene regulation [20]. Long non-coding RNAs (lncRNAs) are non-coding RNA molecules longer than 200 nucleotides. They primarily participate in the regulation of metabolic pathways during plant growth and development and stress responses by influencing chromatin states, transcriptional efficiency, or mRNA processing via cis- or trans-acting mechanisms [21,22]. Both circRNAs and lncRNAs help regulate flavonoid and anthocyanin biosynthesis in various plant species. For instance, specific circRNAs/lncRNAs have been shown to modulate anthocyanin accumulation in pepper and strawberry by targeting key biosynthetic genes [23,24]. Furthermore, the expression levels of certain circRNAs/lncRNAs exhibit correlations with anthocyanin content in Quercus mongolica and potato, suggesting potential regulatory roles [25,26]. However, whether these mechanisms can be applied to black goji remains to be determined.

No anthocyanin-associated circRNAs or lncRNAs have been systematically identified in black goji, nor has the crosstalk between these ncRNAs and hormonal signals been elucidated. Therefore, the transcriptional and post-transcriptional regulatory mechanisms responsible for its efficient anthocyanin accumulation remain unclear. Here, we used mature fruits of anthocyanin-rich black goji (BS4) and wild white goji (WS4) for the analyses. By integrating whole-transcriptome sequencing, targeted phytohormone metabolomics, and molecular docking, we identified 72 differentially expressed circRNAs and 5064 differentially expressed lncRNAs. Our analysis revealed that the circRNA plantcirc_001165 and the lncRNA MSTRG.28703 co-target and activate Hydroxycinnamoyl-CoA shikimate/quinate hydroxycinnamoyltransferase (LrHCT), a key gene in anthocyanin synthesis. Furthermore, MSTRG.28703 trans-regulates Zeatin O-xylosyltransferase 1 (LrZOX1), promoting the accumulation of cis-zeatin (cZ). Molecular docking confirmed the high-affinity interactions between the corresponding proteins and their ligands. Hence, this work elucidates an ncRNA-mediated mechanism involving the coordinated regulation of LrHCT and LrZOX1 to balance anthocyanin and cZ metabolism. In addition, petanin, the predominant anthocyanin in black goji, possesses potent antioxidant and anti-inflammatory activities [27]. Thus, given that oxidative stress and inflammatory responses are often pivotal factors in the pathogenesis of chronic diseases, we employed network pharmacology approaches to predict the potential multi-target mechanisms by which petanin may intervene in Alzheimer’s disease (AD) and inflammatory bowel disease (IBD). Therefore, through integrated multi-omics analysis, this study also aims to uncover the ncRNA-mediated regulatory mechanism governing the balance of anthocyanin and cZ metabolism, and provide a theoretical foundation for the functional application of black goji anthocyanins.

2 Materials and Methods

2.1 Whole-Transcriptome Sequencing

2.1.1 Plant Materials

Black and white fruit samples were collected from Nuomuhong Farm, Zongjia Town, Qaidam Basin, Qinghai Province, China (longitude: 96°27.2856′ E, latitude: 36°26.9625′ N; altitude: 2856.6 m). Sampling was conducted in September 2021, during which fresh fruits at the fully mature stage (S4) were harvested from different plants grown under identical cultivation conditions. Upon collection, the samples were immediately snap-frozen in liquid nitrogen and stored at −80°C for subsequent RNA extraction and sequencing analysis. All samples were preserved at the School of Life Sciences, Qinghai Normal University.

2.1.2 Total RNA Extraction

Total RNA was extracted from both black goji and white goji using the RNAprep Pure Plant Kit (a polysaccharides & polyphenolics-rich spin-column) (Tiangen, Beijing, China). The quality of the extracted total RNA was verified using a NanoDrop spectrophotometer and agarose gel electrophoresis. All procedures were performed using three biological replicates.

2.1.3 Library Construction

Total RNA from white and black goji fruits that met the quality criteria was sent to Shanghai Panomix Biomedical Technology Co., Ltd. for sequencing. Strand-specific libraries were constructed using a ribosomal RNA depletion method [28]. The prepared libraries were then subjected to paired-end 150 bp (PE150) sequencing on an Illumina platform. Raw sequencing reads were processed to remove low-quality bases and low-complexity sequences from both ends, resulting in high-quality clean reads.

2.1.4 Data Quality Control and Sequence Alignment

Adapter sequences at the 3′ ends and reads with an average quality score below Q20 were removed using fastp (v0.22.0). Subsequently, the reference genome index was built using HISAT2 (v2.1.0), and the paired-end clean reads were aligned to the black goji reference genome (GCA_041430385.1, downloaded from the NCBI database) using HISAT2 [29].

2.1.5 Bioinformatics Analysis of ncRNAs

ncRNA Identification

Twenty bp sequences from both ends of the reads that failed to align in the initial HISAT2 mapping were extracted as anchoring sequences for circRNA detection. These anchors were realigned to the genome using Bowtie2 (v2.5.1). Following alignment of the anchor sequences for each sample, circRNAs were identified using find_circ (v1.0) [30]. High confidence circRNAs were filtered based on the following criteria: only one clear breakpoint; anchor sequence overlap ≤ 2 bp; mismatch ≤ 2 bp; supported by ≥2 reads covering more than half of the samples; at least one anchoring sequence per read pair having an alignment score > 35 to the genome; circRNA length < 100 kb.

For lncRNA identification, transcripts were first assembled from the HISAT2 alignment results using StringTie (v2.2.1) [31]. Transcripts with ambiguous strand orientation were removed. Subsequently, transcripts were filtered to exclude those with uncertain strand direction, length < 200 bp, or minimum read coverage < 3. Furthermore, only transcripts belonging to the “x”, “u” or “i” categories were retained. (“x” indicates exonic overlap with a reference transcript on the opposite strand; “u” indicates unknown, intergenic transcripts; “i” indicates transcripts falling entirely within an intron of a reference transcript). Finally, the protein-coding potential was assessed using three tools in parallel: PLEK (v1.2), CNCI (v1.0), and PfamScan (v1.6.4). Novel transcripts predicted to lack protein-coding potential by all three tools were defined as high confidence lncRNAs.

ncRNA Expression Quantification and Differential Expression Analysis

Reads counts for circRNAs from the find_circ results were tallied, and circRNA expression levels were normalized using as transcripts per million (TPM). For lncRNAs, reads counts were obtained at the transcript level using StringTie, and expression was normalized as fragments per kilobase of transcript per million mapped reads (FPKM). Differential expression analysis between the two comparison groups was performed using DESeq2 (v1.38.3) [32]. The differential expression of ncRNAs was assessed separately, with the screening criteria set as |log2(fold change)| > 1 and p < 0.05.

Enrichment Analysis of Source Genes for Differentially Expressed ncRNAs

GO enrichment analysis for the source genes of differentially expressed ncRNAs was conducted using topGO (v2.50.0) [33]. KEGG pathway enrichment analysis for these source genes was performed using the clusterProfiler software package (v4.6.0) [34].

Target Gene Prediction for Differentially Expressed lncRNAs

Cis-acting regulation is the effect of a lncRNA on nearby target genes. Protein-coding genes located within 100 kb upstream or downstream of a lncRNA were screened as potential cis-targets. Trans-acting regulation indicates the ability for recognition based on expression level correlations. Putative trans-acting lncRNA–mRNA pairs were identified by screening for correlations with an absolute correlation coefficient |r| > 0.95 and p < 0.05.

2.2 Candidate Gene Screening and Bioinformatics Analysis

The open reading frame (ORF) of Lr11g005079 was predicted using the online tool ORF Finder (https://www.ncbi.nlm.nih.gov/orffinder/), and its amino acid sequence was deduced (parameters: default codon table, minimum ORF length ≥ 100 bp). A phylogenetic tree of homologous proteins was constructed using the neighbor-joining method in MEGA11.0 [35].

2.3 Targeted Metabolomic Profiling of Phytohormones

This study employed targeted mass spectrometry to quantitatively analyze 21 phytohormones, ABA, 1-aminocyclopropanecarboxylic acid (ACC), brassinolide (BL), castasterone (CS), dihydrozeatin (DHZ), gibberellin A1 (GA1), gibberellin A3 (GA3), gibberellin A4 (GA4), gibberellin A7 (GA7), indole-3-acetic acid (IAA), isopentenyl adenine (iP), isopentenyl adenosine (iPR), jasmonic acid (JA), jasmonoyl-isoleucine (JA-Ile), cis-12-oxophytodienoic acid (12-OPDA), salicylic acid (SA), typhasterol (TY), cZ, trans-zeatin (tZ), cis-zeatin riboside (cZR), and trans-zeatin riboside (tZR).

2.3.1 Metabolite Extraction

Six samples (three biological replicates each of BS4 and WS4) were placed on ice to thaw (0–4°C). Immediately after thawing, approximately 100 mg of each sample was weighed using an AL104 analytical balance (Mettler Toledo, Greifensee, Switzerland) and combined with mixed isotope-labeled internal standards and pre-cooled extraction solvent (acetonitrile/water/formic acid, 80:19:1, v/v/v). The mixture was vortexed using a QT-1 vortex mixer (Qite Analytical Instrument, Shanghai, China) and then ultrasonicated in a JP-100 ultrasonic bath (Jiemeng Cleaning Equipment, Shenzhen, China) for 25 min at 4°C. After centrifugation at 14,000 rcf for 20 min at 4°C, using a 5430R centrifuge (Eppendorf, Hamburg, Germany), the supernatant was transferred to an Ostro 25 mg 96-well plate (Waters, Milford, MA, USA) and filtered using a positive pressure manifold. The filtrate was collected into 2 mL centrifuge tubes and stored at −80°C.

2.3.2 Chromatography-Mass Spectrometry Analysis

Sample separation was performed using an Agilent 1290 Infinity LC ultra-high-performance liquid chromatography (UHPLC) system (Agilent Technologies, Santa Clara, CA, USA). Samples were maintained at 4°C in the autosampler. Chromatographic separation was achieved on an ACQUITY UPLC BEH C18 column (2.1 mm × 100 mm, 1.7 μm; Waters, Milford, MA, USA) with the temperature set at 45°C. Mobile phase A consisted of 0.05% formic acid in water, and mobile phase B was 0.05% formic acid in acetonitrile. The flow rate was 400 μL/min with an injection volume of 4 μL. The LC gradient was programmed as follows: 0–1 min, from 2% to 10% B; 1–10 min, from 10% to 70% B; 10–11 min, from 70% to 95% B; 11–11.1 min, from 95% to 2% B; 11.1–13 min, 2% B.

Mass spectrometric analysis was conducted on a 6500+ QTRAP mass spectrometer (SCIEX, Framingham, MA, USA) operating in both positive and negative ion modes. The following electrospray ionization (ESI) source parameters for positive ion mode were as used: source temperature, 550°C; ion source gas 1 (Gas1), 55 psi; ion source gas 2 (Gas2), 50 psi; curtain gas (CUR), 30 psi; IonSpray voltage floating (ISVF), 4500 V. For the negative ion mode, the parameters were: source temperature, 550°C; Gas1, 55 psi; Gas2, 50 psi; CUR, 30 psi; ISVF, −4500 V. Analytes were detected in the multiple reaction monitoring (MRM) mode.

2.4 Molecular Docking

The three-dimensional structures of DHZ, tZ, tZR, cZ, cZR, and petanin were obtained from the PubChem database (https://pubchem.ncbi.nlm.nih.gov/) using the 439631, 449093, 6440982, 688597, 13935024, and 6540681 compound IDs (CIDs), respectively. The protein structures of LrHCT and LrZOX1 were retrieved from the UniProt website (https://www.uniprot.org/). Blind docking was performed using the CB-DOCK2 online platform [36] with default parameters to analyze the optimal binding sites and binding affinity between the proteins and compounds. A binding energy lower than or equal to −5.0 kcal·mol−1 (1 kcal = 4.184 J) was considered indicative of a favorable binding capability between the compound and the target protein. The protein-protein interaction structure of LrHCT and LrZOX1 was predicted using AlphaFold3 (https://golgi.sandbox.google.com/) with default parameters and visualized using PyMOL (v3.1).

2.5 Network Pharmacology Analysis

2.5.1 Prediction and Screening of Potential Petanin Targets

The SMILES notation of petanin C[C@H]1[C@@H]([C@H]([C@H](C(O1)OC[C@@H]2[C@H]([C@@H]([C@H](C(O2)OC3=C([O+]=C4C=C(C=C(C4=C3)OC5[C@@H]([C@H]([C@@H]([C@H](O5)CO)O)O)O)O)C6=CC(=C(C(=C6)OC)O)O)O)O)O)O)O)OC(=O)/C=C/C7=CC=C(C=C7)O was obtained from the PubChem database (https://pubchem.ncbi.nlm.nih.gov/) by searching for the keyword ‘petanin’. This SMILES string was then submitted to the SwissTargetPrediction (https://www.swisstargetprediction.ch/), SEA (https://sea.bkslab.org/), and CTD (https://ctdbase.org/) databases to predict potential targets of petanin.

2.5.2 Screening of AD and IBD Targets and Identification of Intersection Genes

Target genes associated with Alzheimer’s disease (AD) and inflammatory bowel disease (IBD) were retrieved from the GeneCards (https://www.genecards.org/, relevance score > 4) and OMIM (https://www.omim.org/) databases using the keywords “Alzheimer’s disease” and “inflammatory bowel disease”, respectively. The obtained AD and IBD target gene sets were then intersected with the potential petanin targets. The intersection was visualized using InteractiVenn (https://www.interactivenn.net/index2.html) [37]. The resulting intersection genes were considered as potential target genes for petanin intervention in AD and IBD.

2.5.3 Construction of Protein-Protein Interaction (PPI) Networks for Petanin-IBD Intersection Targets

The intersection target genes identified above were input into the STRING database (https://cn.string-db.org/), with the organism set to Mus and a confidence threshold of 0.4. The filtered interaction data were exported and imported into Cytoscape (v3.9.1). The “network analyzer” function was used to calculate the degree centrality of the interacting proteins and visualize the PPI network.

2.5.4 GO and KEGG Enrichment Analysis

GO and KEGG enrichment analyses of the key target genes were performed using the Bioinformatics online platform (https://www.bioinformatics.com.cn/).

2.6 Statistical Analysis

Statistical analysis was performed using GraphPad Prism (v10.1). Data are presented as the mean ± standard error of the mean (SEM). Differences between the two groups were assessed using a t-test for independent samples. Statistical significance was set at p < 0.05.

3 Results

3.1 Quality Assessment of Transcriptome Sequencing Data

We performed sequencing on three biological replicates of mature fruits from each of WS4 (white goji) and BS4 (black goji) to assess the reliability of the transcriptome data (Fig. 1A,B). For WS4, an average of 104,999,689 raw reads was obtained, yielding 103,505,081 clean reads after filtering. For BS4, an average of 105,972,212 raw reads yielded 104,497,323 clean reads. For both groups, the Q20 scores exceeded 98%, the Q30 scores exceeded 96%, and the guanine-cytosine (GC) content was above 40% (Table 1), indicating high base accuracy consistent with plant genomic characteristics. The clean reads were then aligned to the reference genome GCA_041430385.1 (Table 2). The overall alignment rate was greater than 92.01%, with a unique alignment rate exceeding 88.11% and a multi-alignment rate below 11.89%, ensuring the specificity of sequence matching. In summary, the sequencing data were of high quality and strong usability, and provided a solid foundation for subsequent ncRNA identification and differential expression analysis.

3.2 Differential Expression and Functional Pathway Analysis of circRNAs

We identified 647 high-confidence circRNAs, using the find_circ software tool (detailed in Section 2.1.5). The analysis of circRNA expression abundance revealed a uniform TPM density distribution for both BS4 and WS4, with a low proportion of lowly and highly expressed circRNAs. This observation indicates a generally stable circRNA expression level across the two sample groups (Fig. 1C). Furthermore, principal component analysis (PCA) demonstrated high clustering of biological replicates within each group and clear separation between groups. This result confirms that the observed differences in circRNA expression primarily originate from genotypic variation and that the experimental reproducibility is reliable (Fig. 1D). Notably, the overall expression level of circRNAs in BS4 was substantially higher than in WS4, suggesting a potential role for circRNAs in the regulation of the distinctive anthocyanin synthesis metabolism in black goji.

Applying the criteria of |log2(FC)| > 1 and p < 0.05, we identified 72 differentially expressed circRNAs, of which 31 were upregulated and 41 were downregulated in BS4 relative to WS4 (Fig. 1E). Gene Ontology (GO) and Kyoto Encyclopedia of Genes and Genomes (KEGG) pathway enrichment analyses were performed on the 111 corresponding source genes. In the GO functional classification, these genes were primarily enriched in three major categories: “cellular component (CC)” (e.g., potassium ion-transporting ATPase complex), “molecular function (MF)” (e.g., glucan exo-1,3-beta-glucosidase activity), and “biological process (BP)” (e.g., stabilization of membrane potential) (Fig. S1A). The KEGG pathway enrichment analysis identified 23 substantially enriched pathways. Among these, the “photosynthesis–antenna proteins” and “basal transcription factors” pathways exhibited the highest enrichment scores. The “photosynthesis–antenna proteins” pathway is involved in the absorption and conversion of luminous energy. Its enrichment suggests that the differentially expressed circRNAs may contribute to enhanced photosynthetic efficiency in BS4, potentially providing more energy for anthocyanin synthesis (Fig. S1B).

Then, we performed cluster analysis on the upregulated circRNAs in BS4 to further pinpoint key circRNAs associated with anthocyanin synthesis. This analysis revealed that plantcirc_001165 was expressed only in BS4, with no detectable signal in WS4 (Fig. 1F). Its source gene, Lr11g005079, was annotated as a key gene in the core anthocyanin synthesis pathways: flavonoid biosynthesis (ko00941) and phenylpropanoid biosynthesis (ko00940). Quantitative experimental validation showed that the relative expression levels of both plantcirc_001165 and Lr11g005079 were significantly higher in BS4 than in WS4 (where both were undetectable, p < 0.05), and their expression trends were highly consistent (p < 0.05). This result suggests that plantcirc_001165 may directly participate in the flavonoid biosynthetic pathway in black goji by regulating its source gene Lr11g005079.

The above results delineate a circRNA-mediated regulatory pathway for anthocyanin synthesis. However, the transcriptional and post-transcriptional regulation of plant secondary metabolism typically relies on the coordinated actions of multiple classes of ncRNAs. Therefore, our subsequent analysis focuses on the differential expression profiles of lncRNAs and investigates their potential relationships with Lr11g005079 and phytohormone metabolism, aiming to systematically elucidate the regulatory mechanisms governing the balance between anthocyanin and hormone metabolism.

Figure 1: Phenotype and transcriptomic analysis of black goji fruits. (A) White goji (WS4). (B) Black goji (BS4). (C) Density distribution of the circRNA expression levels (TPM). (D) Principal component analysis (PCA) of the circRNA expression profiles. (E) Volcano plot of differentially expressed circRNAs. (F) Heatmap of upregulated circRNAs (The asterisk next to plantcirc_001165 highlights it as the key circRNA that is specifically expressed in the black goji (BS4) group).

3.3 Differential Expression and Functional Pathway Analysis of lncRNAs

Although the circRNA plantcirc_001165 may be involved in anthocyanin synthesis by regulating the source gene (Lr11g005079), a single class of ncRNA cannot fully explain the metabolic regulatory network. Thus, to decipher the role of lncRNAs in this process, we employed three software tools—PLEK, CNCI, and PfamScan (detailed in Section 2.1.5)—to screen for high-confidence lncRNAs, and identified a total of 7882 (Fig. S2A). The analysis of lncRNA expression abundance revealed uniform FPKM density curves for both BS4 and WS4, with a low proportion of lowly and highly expressed lncRNAs. This result indicates generally stable lncRNA expression levels across the two sample groups, with no extreme expression biases (Fig. S2B). This outcome, consistent with the circRNA results, further confirms the robust quality of our transcriptomic data.

A total of 5064 differentially expressed lncRNAs were identified using the criteria of |log2(FC)| > 1 and p < 0.05. Among these, 3147 were upregulated and 1917 were downregulated in BS4 compared to WS4, accounting for 62.1% of the differential lncRNAs (Fig. 2A). Notably, the number of differentially expressed lncRNAs was substantially greater than that of circRNAs (72), and the proportion of upregulated lncRNAs was substantially higher. This finding suggests that lncRNAs may play a more central role than circRNAs in transcriptional and post-transcriptional regulation underlying the divergence in anthocyanin metabolic pathways.

Cis-target genes of the differentially expressed lncRNAs were predicted, yielding 40,702 potential target genes (detailed in Section 2.1.5). Then, GO and KEGG enrichment analyses were performed to elucidate the functional associations of these target genes. In the GO functional classification, they were enriched in terms such as “nuclear membrane” (CC), “serine-type endopeptidase inhibitor activity” (MF), and “calcium ion transport into cytosol” (BP), among others (Fig. S2C). KEGG analysis identified 127 enriched pathways. Notably, three secondary metabolic pathways were significantly (p < 0.05) enriched: zeatin biosynthesis (ko00908, 34 mRNAs), flavonoid biosynthesis (ko00941, 23 mRNAs), and phenylpropanoid biosynthesis (ko00940, 74 mRNAs) (Fig. 2B). The significant enrichment of the zeatin biosynthesis pathway represents a unique feature of lncRNA regulation (not observed in the circRNA enrichment analysis). This observation suggests that lncRNAs may simultaneously link anthocyanin synthesis and hormone metabolism, providing a clue for subsequent investigation into the mechanisms balancing these two metabolic processes.

We focused on the key anthocyanin synthesis gene Lr11g005079 to validate the synergistic regulatory roles of ncRNAs. Consistent with the circRNA analysis results, this gene was also enriched in the flavonoid and phenylpropanoid biosynthesis pathways. This finding confirms that both the circRNA and lncRNA regulatory networks target this core gene and converge on the anthocyanin synthesis pathway. By further screening lncRNAs associated with the phenylpropanoid biosynthesis pathway and generating an expression heatmap, we identified the cis-acting lncRNA MSTRG.28703, which targets Lr11g005079. MSTRG.28703 exhibited significantly higher (p < 0.05) expression in BS4, and its expression trend was synchronized with that of Lr11g005079 (Fig. 2C). This result suggests that MSTRG.28703 may activate Lr11g005079 through cis-regulation, promoting efficient anthocyanin accumulation.

These findings establish the regulatory role of lncRNAs, particularly MSTRG.28703, in anthocyanin synthesis. MSTRG.28703 not only synergistically targets Lr11g005079 along with circRNAs, but may also be linked to hormone metabolism through the zeatin biosynthesis pathway. However, the evolutionary conservation of Lr11g005079 and the specific function of its encoded protein remain unclear. Therefore, the following section describes phylogenetic and functional annotation analyses performed to clarify the evolutionary position and enzymatic function of this gene, providing a molecular foundation for the in-depth exploration of the metabolic mechanisms.

Figure 2: Analysis of differentially expressed lncRNAs. (A) Volcano plot of differentially expressed lncRNAs. (B) KEGG pathway enrichment analysis of the cis-target genes of lncRNAs. (C) Heatmap showing the expression of lncRNAs associated with the phenylpropanoid biosynthesis pathway (The asterisk next to MSTRG.28703 denotes this crucial cis-acting lncRNA, which targets Lr11g005079 and exhibits significantly higher expression in BS4.).

3.4 Phylogenetic Analysis and Functional Annotation of the Key Gene LrHCT

We conducted phylogenetic analysis and functional annotation to elucidate the evolutionary position and function of Lr11g005079. The full-length genomic sequence of this gene, comprising 1498 bp, was retrieved from the black goji reference genome. One high-confidence open reading frame, designated ORF39, was identified using the online tool NCBI ORFfinder. Phylogenetic analysis revealed that the amino acid sequence encoded by Lr11g005079 clusters within a clade containing homologous sequences from species belonging to the genera Lycium (e.g., L. ferocissimum and L. barbarum), Anisodus (e.g., A. acutangulus), and Populus (e.g., P. alba) (Fig. 3A). All sequences within this clade are annotated as members of the BAHD acyltransferase family. The core function of the BAHD acyltransferase family is to catalyze the acylation of secondary metabolites, such as flavonoids and phenylpropanoids. This modification involves the transfer of acyl groups to metabolite molecules, regulating their stability, water solubility, and bioactivity. These results suggest that the protein encoded by Lr11g005079 has a conserved function and possesses acyltransferase activity. This protein is likely involved in anthocyanin synthesis in black goji by catalyzing the acylation of anthocyanin precursors.

KEGG pathway mapping was performed on homologous genes to further delineate the specific role of Lr11g005079 within metabolic pathways. The analysis revealed the involvement of this gene in three core secondary metabolic pathways: (1) the biosynthesis of secondary metabolites (ko01110); (2) flavonoid biosynthesis (ko00941); (3) phenylpropanoid biosynthesis (ko00940) (Fig. 3B). Further annotation of key enzymes within these pathways showed that five out of eight homologous genes could be mapped to the critical shikimate enzyme O-hydroxycinnamoyltransferase (HCT, EC: 2.3.1.133). HCT catalyzes the conjugation of p-coumaroyl-CoA with shikimate/quinate to produce key intermediates such as caffeoyl-CoA. Caffeoyl-CoA serves as an essential substrate for the synthesis of anthocyanin precursors, such as cyanidin and petunidin aglycones. This result indicates that Lr11g005079 is a key gene encoding HCT. This gene directly participates in anthocyanin synthesis in black goji by regulating the metabolic flux through the phenylpropanoid–flavonoid pathway.

These findings clarify the function of Lr11g005079 (designated LrHCT) and provide downstream target support for the cis-regulation by its upstream lncRNA, MSTRG.28703. Given the substantial enrichment of lncRNA target genes in the zeatin biosynthesis pathway (ko00908), we hypothesize that lncRNAs may establish a metabolic balance network between anthocyanins and zeatin by dually regulating LrHCT and the associated zeatin metabolism. Therefore, subsequent analysis integrates hormone metabolome data to elucidate the lncRNA-mediated co-regulatory module for anthocyanin–zeatin metabolism.

Figure 3: Phylogenetic and functional annotation analysis of Lr11g005079. (A) Phylogenetic tree of the protein encoded by Lr11g005079 (The red star marks ORF39, which represents the high-confidence amino acid sequence encoded by our key gene Lr11g005079.). (B) KEGG mapping of Lr11g005079 and its homologous genes.

3.5 lncRNA Regulatory Network and Hormone Quantification

Having established that Lr11g005079 encodes the key anthocyanin synthesis enzyme LrHCT, and given the indication from Section 3.3 that lncRNA MSTRG.28703 may cis-regulate LrHCT, we investigated whether this lncRNA is also linked to zeatin metabolism. We predicted the related trans-target genes of MSTRG.28703, and combined this information with hormone quantification for validation. Trans-target prediction yielded 323 co-expressed mRNAs (detailed in Section 2.1.5). KEGG pathway annotation of these mRNAs revealed that Lr01g000261 was annotated as LrZOX1. Although this gene is annotated within the zeatin biosynthesis pathway (ko00908), its biological function is to catalyze the O-xylosylation of zeatin, representing a critical node in hormone metabolism and homeostasis regulation. Furthermore, both LrZOX1 and LrHCT (Lr11g005079) belong to the overarching secondary metabolites biosynthesis pathway (ko01110). Moreover, as the trans- and cis-target gene of lncRNA MSTRG.28703 respectively, they exhibited significant synergy in their expression patterns, suggesting a potential metabolic connection between them.

We performed an intersection analysis between the 34 mRNAs from the zeatin biosynthesis pathway and the 40,702 cis-target genes of the differentially expressed lncRNAs to further identify lncRNAs that cis-regulate LrZOX1. This approach allowed us to identify two key lncRNAs: MSTRG.66 and MSTRG.64. Both cis-regulate LrZOX1 and are concurrently enriched in the phenylpropanoid biosynthesis pathway (ko00940). This forms a potential co-regulatory module targeting both anthocyanin and zeatin metabolism, with a functional division comprising three main parts: (1) The core lncRNA MSTRG.28703, which cis-regulates the key gene for anthocyanin synthesis, LrHCT, which is highly expressed in BS4. (2) The cooperative lncRNAs MSTRG.66/MSTRG.64, which cis-regulate the zeatin modification gene LrZOX1 and are highly expressed in WS4. (3) LrHCT and LrZOX1, which regulate anthocyanin synthesis and zeatin activity, respectively. Expression validation confirmed that all members within this module showed significant differential expression between BS4 and WS4 (p < 0.05, Fig. 4A). Furthermore, MSTRG.28703 exhibited a significant negative correlation with MSTRG.66/MSTRG.64 (p < 0.05, Fig. 4B). This result suggests that the module may regulate the expression balance of LrHCT and LrZOX1 simultaneously through an antagonistic mechanism where MSTRG.28703 suppresses MSTRG.66/MSTRG.64.

LrZOX1 catalyzes the conjugation of active zeatin with UDP-xylose, forming inactive storage forms, such as zeatin-O-xyloside [38,39,40]. Among the active zeatin isomers, cZ is particularly noteworthy. Although it exhibits lower bioactivity than tZ in stimulating cell division and promoting growth, it plays a pivotal role in maintaining the stability of the photosynthetic apparatus and mediating stress responses [41]. This physiological function aligns closely with the stress-resistance characteristics of anthocyanins. Targeted quantitative analysis of 21 phytohormones (e.g., auxins and cytokinins) was performed on BS4 and WS4 fruits to verify whether the differential expression of LrZOX1 affected cZ content (detailed in Section 2.3). The results revealed that the differences in free-form zeatin-type compounds—namely DHZ, cZ, and tZR—were the most pronounced among hormones (p < 0.05, Fig. 4C). The content of DHZ was 12.1-fold higher, cZ was 2.5-fold higher, and tZR was 1.2-fold higher in BS4 than in WS4. These findings confirm the specific accumulation of cZ derivatives in BS4. Furthermore, a significant negative correlation (p < 0.05) was observed between the levels of these active zeatin compounds and LrZOX1 expression, validating the hypothesis reduced LrZOX1 expression in BS4 leads to greater accumulation of cZ. These results suggest that LrZOX1 modulates the active pool of cZ and indicate that cZ may function as a signaling molecule involved in anthocyanin accumulation in black goji.

These results elucidate the mechanism by which the lncRNA module co-regulates anthocyanin and cZ metabolism at both the transcriptional and metabolic levels. The next section describes the use of molecular docking to further validate this metabolic balance model from a structural perspective of protein–ligand interactions.

Figure 4: Analysis of the lncRNA regulatory module and associated hormone levels. (A) Expression levels (bar plot) of MSTRG.28703, MSTRG.66, MSTRG.64, and Lr01g000261. (B) Correlation heatmap of all lncRNAs associated with the phenylpropanoid biosynthesis pathway. (C) Bar plot showing the quantified levels of selected phytohormones. Statistical significance is indicated as follows: ****p < 0.0001, ***p < 0.001, **p < 0.01, *p < 0.05.

3.6 Molecular Docking Validation of Key Protein–Ligand Interactions

Previous findings have confirmed that the lncRNA regulatory module modulates the expression of LrHCT and LrZOX1 through an antagonistic effect, concomitant with a notable increase in multiple free-form zeatins in BS4. However, the molecular underpinnings of this regulatory module—specifically, whether zeatin family small molecules directly interact with LrHCT/LrZOX1, whether LrHCT and LrZOX1 directly interact with each other, and whether petanin serves as a substrate for LrHCT—still lack structural evidence.

Therefore, we used molecular docking to investigate the interactions between key proteins and ligands (zeatin-type small molecules, petanin) and potential protein–protein interactions. The three-dimensional structures of LrHCT and LrZOX1 were obtained from the UniProt and PubChem databases, along with the structures of DHZ, tZ, tZR, cZ, cZR. Blind docking was performed using the CB-DOCK2 online platform. The results showed that all five zeatin molecules formed stable complexes with both LrHCT and LrZOX1, with binding energies below −5.0 kcal·mol−1 (details in Section 2.4). Among them, cZR exhibited the strongest binding affinity to both proteins (LrHCT: −8.5 kcal·mol−1; LrZOX1: −8.7 kcal·mol−1) (Fig. 5A). This outcome verifies the substrate specificity of zeatin for LrZOX1, with cZR as the preferred substrate, and suggests that cZ may bind to LrHCT as a signaling molecule, potentially promoting anthocyanin synthesis by activating LrHCT enzyme activity.

To investigate whether LrHCT and LrZOX1 cooperate in metabolic regulation through direct interaction, we predicted their complex structure using AlphaFold3 and visualized it using PyMOL. The results indicate that LrHCT and LrZOX1 can form a stable heterodimer. Their interaction is mediated by two key hydrogen bonds: the lysine 201 (LYS201) of LrZOX1 forms hydrogen bonds with the arginine 99 (ARG99) and glutamate 101 (GLU101) of LrHCT, with bond lengths of 3.1 Å and 3.5 Å, respectively (Fig. 5B). The prediction confidence parameters (plDDT > 70, pTM + ipTM = 0.86) suggest a reliable structural basis for direct protein–protein interaction between LrHCT and LrZOX1. Considering their negatively correlated expression patterns, we hypothesize that LrHCT and LrZOX1 may cooperatively regulate the balance of anthocyanin and cZ metabolism through an interaction-coupled antagonistic mode.

Furthermore, given that LrHCT is a key acyltransferase in anthocyanin synthesis and petanin is an acylated anthocyanin, we verified whether petanin serves as a catalytic product or a regulatory factor of LrHCT. The three-dimensional structure of petanin was retrieved from the PubChem database and docked with LrHCT. The results showed a binding energy of −8.9 kcal·mol−1 between them, along with the formation of four key hydrogen bonds. Specifically, the methionine 1 (MET1), glutamate 8 (GLU8), glutamine 10 (GLN10), and glycine 13 (GLY13) of LrHCT formed hydrogen bonds with the glycosyl or acyl groups of petanin (bond lengths > 2.5 Å, Fig. 5C). Considering the acyltransferase function of LrHCT, these results suggest that petanin may be the end product catalyzed by LrHCT, or it may enhance LrHCT enzymatic activity through binding. This result provides a clear theoretical basis for subsequent in vitro enzyme activity assays.

These findings provide structural evidence for the binding potential between LrHCT/LrZOX1 and zeatin-type molecules as well as petanin. Furthermore, they offer theoretical support for the metabolic balance mediated by the lncRNA regulatory module and suggest that its catalytic end product, petanin, may possess broader biological activities. Thus, to further investigate the potential application of petanin, we employed a network pharmacology approach to systematically predict its potential therapeutic targets and mechanisms of action in AD and IBD.

Figure 5: Molecular docking analysis of key interactions. (A) Heatmap of binding affinities from molecular docking between LrHCT/LrZOX1 and zeatin-type small molecules. (B) Schematic representation of the docking between LrHCT and LrZOX1. (C) Schematic representation of the docking between petanin and LrHCT.

3.7 Network Pharmacology Analysis of Petanin Intervention in Alzheimer’s Disease and Inflammatory Bowel Disease

To assess the potential bioactivity of petanin, we first predicted a total of 75 petanin-related targets. Subsequently, we screened 10,321 AD-related target genes and 9530 IBD-related target genes. The intersection between petanin targets and disease targets yielded 52 overlapping genes for AD (Fig. 6A) and 50 for IBD (Fig. 6B). These overlapping target genes were then uploaded to the STRING database to construct protein–protein interaction (PPI) networks, which were visualized using Cytoscape (Fig. 6C and Fig. S3A). Notably, the top five key target genes ranked by interaction degree were identical between the two diseases, namely tumor necrosis factor (TNF), interleukin-6 (IL6), epidermal growth factor receptor (EGFR), hypoxia-inducible factor 1α (HIF1A), and prostaglandin-endoperoxide synthase 2 (PTGS2). These shared targets are closely associated with inflammatory response and oxidative stress, suggesting that modulation of this core network may represent a common molecular basis for petanin intervention in these disorders.

We performed GO and KEGG enrichment analyses to systematically decipher the biological functions of these shared targets. The GO analysis revealed a high degree of overlap in the biological processes associated with the potential effect of petanin on AD and IBD. Enriched terms common to both included “response to oxidative stress” and “regulation of inflammatory response.” These commonly enriched terms align closely with the antioxidant activity of petanin. Thus, we speculate that petanin may mitigate oxidative stress and inflammatory damage by modulating these processes to scavenge excess reactive oxygen species (ROS) in the brain (AD) and intestinal tract (IBD). This could represent one of the core mechanisms underlying its potential intervention in both diseases. Despite the substantial number of shared terms, petanin also exhibited disease-specific regulatory potential. GO terms unique to AD were enriched in processes such as “neuron projection terminus” (Fig. 6D). By contrast, terms unique to IBD were concentrated in categories including “regulation of peptidyl-tyrosine phosphorylation” (Fig. S3B). KEGG enrichment analysis further corroborated the high consistency in the mechanistic pathways of these target genes. The 10 most-enriched pathways were identical between the two disease analyses. These included the “AGE-RAGE signaling pathway in diabetic complications” and “HIF-1 signaling pathway” (p < 0.05) (Fig. S4A,B). These results suggest that petanin may exert its potential therapeutic effects on AD and IBD through a multitarget, multipathway approach that coordinately modulates the inflammation–oxidative stress axis.

Figure 6: Network pharmacology analysis of petanin for Alzheimer’s Disease and Inflammatory Bowel Disease. (A) Venn diagram showing the intersection of predicted petanin targets with Alzheimer’s Disease (AD) targets. (B) Venn diagram showing the intersection of predicted petanin targets with Inflammatory Bowel Disease (IBD) targets. (C) Protein-protein interaction (PPI) network of the overlapping targets between petanin and AD, highlighting the key hub genes. (D) Gene ontology (GO) enrichment analysis of the petanin-AD intersection targets.

4 Discussion

4.1 Coordinated Targeting of LrHCT by circRNAs and lncRNAs in Black Goji

Many previous studies have established that circRNAs and lncRNAs regulate gene expression indirectly by constructing competing endogenous RNA (ceRNA) networks [42,43,44]. By contrast, in this study, the circRNA plantcirc_001165, the lncRNA MSTRG.28703, and their target gene LrHCT (Lr11g005079) exhibited a strong co-expression pattern in BS4. All three were highly expressed only in the anthocyanin-rich BS4 and were nearly undetectable in WS4. Hence, we propose the existence of a coordinated circRNA/lncRNA/target gene regulatory mode in black goji. Distinct from the indirect mechanism of ceRNA networks, this mode may involve the simultaneous activation of the core enzyme gene by both classes of ncRNAs, potentially enabling the rapid and efficient initiation of anthocyanin synthesis.

Zong et al. [11] demonstrated that LrAN2 is a key transcriptional regulator for anthocyanin accumulation in black goji. The LrHCT identified in this study is a core enzymatic gene in the flavonoid synthesis pathway [45,46,47]. Notably, HCT can catalyze a reverse reaction [48], enabling the plant to convert chlorogenic acid to caffeoyl-CoA. Hence, we propose that, LrAN2 likely activates LrHCT transcription in black goji by binding to its promoter region, while plantcirc_001165 and MSTRG.28703 facilitate efficient anthocyanin accumulation by stabilizing LrHCT mRNA or enhancing its translational efficiency at the post-transcriptional level. The proposed regulatory network integrates knowledge on anthocyanin regulation in black goji, offering a more comprehensive perspective for understanding its high-anthocyanin trait.

4.2 The Regulatory Network Constructed by lncRNAs Mediates Cross-Talk between Anthocyanin and Zeatin Metabolic Pathways

We posit that MSTRG.28703, which is highly expressed in black goji, regulates metabolic balance through a dual mechanism. On one hand, it cis-activates LrHCT to promote anthocyanin synthesis. On the other hand, it trans-represses the expression of MSTRG.66 and MSTRG.64. We hypothesize that this repression may occur through competitive binding to the promoter region of LrZOX1 or by sequestering miRNAs targeting MSTRG.66/MSTRG.64, indirectly suppressing LrZOX1 transcription. This mechanism would reduce the inactivation of active zeatin, leading to the accumulation of free-form zeatins, such as cZ.

CZ can specifically accumulate under drought stress, contributing to the maintenance of cellular physiological activity and regulating the balance between plant development and stress responses [41]. Meanwhile, anthocyanins can scavenge the excess ROS generated by drought stress [27]. This aligns closely with the evolutionary adaptation of black goji to arid environments. Notably, although cZ acts as a signaling molecule to effectively activate cytokinin-induced genes [49], the regulation of anthocyanin biosynthesis by cytokinins exhibits significant functional duality across different species. In banana [50], cytokinins primarily exert an inhibitory effect, hindering anthocyanin accumulation by upregulating B-type response regulators that repress the expression of genes such as dihydroflavonol 4-reductase (MaDFR). However, our results indicate that black goji fits a positive synergistic regulation model. Consistent with the mechanism elucidated in Arabidopsis thaliana [51], cytokinins act as essential enhancers of blue light signaling, significantly amplifying anthocyanin accumulation via the cryptochrome 1 (CRY1)-dependent pathway. Given the long-term exposure of black goji to the high-radiation environment of the Qinghai-Xizang Plateau, we speculate that the accumulated cZ functions not as an inhibitor, but as a key stress response signaling molecule, producing a synergistic effect with high-intensity environmental light signals to maximize photoprotection.

Cytokinin signaling modulates anthocyanin biosynthesis by directly upregulating specific structural genes in some plants, such as CHS and DFR in Arabidopsis [52] and PAL in Eucalyptus [53]. In contrast, cZ in black goji may tend to enhance the acylation process mediated by LrHCT, thereby improving the stability and bioactivity of anthocyanins. This difference likely stems from variations in the steps of the anthocyanin biosynthetic pathway across species. This study links ncRNA regulation, hormone metabolism, and stress adaptation for the first time, providing a novel case for understanding the ecological and evolutionary significance of plant secondary metabolism.

4.3 Molecular Docking Reveals the High-Affinity Binding of LrHCT and LrZOX1 to Ligands and Their Potential Cooperative Interaction Mechanism

Protein-protein docking analysis revealed that the assembly of a stable LrHCT-LrZOX1 heterodimer core is dependent on three critical residues: LYS201, ARG99, and GLU101. Notably, LYS201 within LrZOX1 establishes a hydrogen-bonding network with ARG99 and GLU101 in LrHCT. The electrostatic interaction between the positively charged LYS201 and the negatively charged GLU101 is key to maintaining dimer stability, while ARG99 further reinforces the complex structure through hydrogen bonding. Notably, GLU101 is located near the predicted acyltransferase active site of LrHCT (based on homology modeling). We speculate that dimer formation may regulate the binding efficiency of LrHCT for p-coumaroyl-CoA by altering the spatial conformation of its active site, achieving dynamic coupling between anthocyanin synthesis and cZ metabolism.

In the ligand binding analysis, cZR exhibited the lowest binding energy with both LrHCT and LrZOX1 (from −8.3 to −8.7 kcal·mol−1). The additional ribose group in its structure can form an extra hydrogen bond with the SER35 of LrHCT, representing the key structural basis for its higher binding affinity than those of other zeatin derivatives. This also clarifies the preferred substrate type for LrZOX1.

The glycosyl and acyl groups of petanin form four hydrogen bonds with the MET1, GLU8, GLN10, and GLY13 of LrHCT. The binding site precisely overlays the predicted acyl-donor binding pocket of LrHCT, based on the structural conservation of the BAHD acyltransferase family [54]. This indicates that petanin is more likely the catalytic end product of LrHCT. The enzyme likely synthesizes petanin by transferring a p-coumaroyl group onto the 5-O-glucoside of an anthocyanin precursor. This inference aligns with experimental results identifying petanin as the major anthocyanin component in black goji [2]. Furthermore, the strong binding between petanin and LrHCT (binding energy −8.9 kcal·mol−1) may also regulate the enzymatic activity of LrHCT through a feedback inhibition mechanism, preventing the metabolic burden caused by excessive anthocyanin accumulation.

4.4 Network Pharmacology Predicts the Multitarget Intervention Potential of Petanin in Alzheimer’s Disease and Inflammatory Bowel Disease

GO analysis revealed that the common targets of petanin for AD and IBD intervention were enriched in biological processes such as ROS metabolism, response to oxidative stress, and regulation of inflammatory response. These findings are consistent with previous experimental evidence [5,55]. In vitro cellular experiments have confirmed that petanin protects PC12 cells from oxidative stress by mitigating lipid peroxidation and modulating endogenous antioxidant enzymes [56]. In vivo studies have demonstrated the ability of petanin to alleviate intestinal inflammation in a DSS-induced mouse colitis model by reducing the expression of IL-6, TNF-α, and other markers [57]. The network pharmacology predictions in the present study further clarify the core targets and pathways underlying these phenotypic effects, thus contributing to a more complete chain of evidence.

Regarding disease-specific mechanisms, petanin targets for AD were enriched in the biological process “neuron projection terminus”. This process is directly associated with synaptic dysfunction in AD pathology [58]. The activation of the core target EGFR can promote neurite outgrowth and synaptogenesis [59]. Thus, we speculate that petanin may alleviate Aβ deposition-induced synaptic damage by inhibiting the overactivation of EGFR [60]. For IBD, the enriched regulatory term “peptidyl–tyrosine phosphorylation” is related to intestinal epithelial barrier repair. The core target HIF-1α can promote intestinal epithelial cell proliferation and the expression of tight junction proteins (e.g., ZO-1) by regulating tyrosine kinase-associated signaling pathways [61,62]. This observation aligns with the potential role of petanin in mitigating intestinal mucosal damage.

Hierarchical analysis of the core pathways indicates that the AGE–RAGE signaling pathway serves as the central upstream node for petanin’s intervention in both diseases. The overactivation of this pathway upregulates downstream targets such as IL6, TNF, and PTGS2 via NF-κB/MAPK signaling [63]. Petanin may directly bind to key receptors within the AGE-RAGE pathway (e.g., RAGE), blocking ligand binding and inhibiting downstream inflammation. The HIF-1 signaling pathway functions as a secondary node, regulating neuronal hypoxia tolerance in AD and intestinal barrier repair in IBD [64,65]. Together, these pathways form a coordinated multi-pathway intervention network for petanin.

5 Conclusions

By integrating multi-omics analysis and computational biology, this study proposes a working model in which circRNAs and lncRNAs coordinately regulate LrHCT and LrZOX1 to mediate the metabolic balance between anthocyanins and cZ. This study provides a new transcriptional and post-transcriptional regulatory perspective for understanding fruit coloration in black goji. Furthermore, the network pharmacology analysis of petanin offers a theoretical prediction for its potential therapeutic efficacy in neurodegenerative diseases and intestinal inflammation. The key findings include: (1) The circRNA plantcirc_001165 and the lncRNA MSTRG.28703 co-target and activate LrHCT, a key gene encoding an enzyme for anthocyanin synthesis, thereby promoting anthocyanin accumulation; (2) MSTRG.28703, in coordination with MSTRG.66/MSTRG.64, regulates the expression of LrZOX1, leading to the accumulation of cZ; (3) Molecular docking confirms the high-affinity binding of LrHCT and LrZOX1 with zeatin-type molecules and petanin; (4) Network pharmacology predicts that petanin may exert potential therapeutic effects in AD and IBD by modulating oxidative stress and inflammation-related pathways, such as the AGE-RAGE and HIF-1 signaling pathways. Future studies will prioritize functional validation and ncRNA perturbation experiments. Furthermore, we aim to explore the potential of utilizing genome editing or ncRNA-based regulatory strategies to enhance the fruit quality of black goji.