Experimental materials were the purified offspring S1 (WO108-4) and S2 (WO59-1) of the doubled haploid hybrid population strains selected and bred by the P. haitanensis Germplasm Improvement and Application Research Laboratory of Jimei University, Fujian Province, China. Conchocelis samples without conchosporangia (0.05 g fresh weight) were added to a cell culture bottle containing 250 mL sterilized seawater and then cultured in an incubator set at 29°C with a 9-h light (1,000–1,500 lux)p : 15-h dark cycle to induce conchocelis maturation. The culture medium was supplemented with 1.36 mg L⁻¹ NaH₂PO₄ and 14 mg L⁻¹ NaNO₃ and refreshed every 7 d.

Conchocelis maturation was analyzed using a microscope (Eclipse80i; Nikon, Tokyo, Japan). Briefly, 10 randomly selected algal pellets were examined in each bottle, with 25 fields of view per algal pellet. To simplify the examination, algae were separated using a dissecting needle and placed on a glass slide. Using the microscope (200× magnification), the number of conchosporangia and filamentous conchocelis in each field of view was recorded.

Conchocelis maturation, which was determined on the basis conchosporangium formation, was calculated using the following formula: number of conchosporangia/(total number of conchocelis) × 100%. The conchocelis maturation of S1 in 49 d is 90%, while S2 will not reach 90% maturity until 91 d. The conchosporangia developmental process can be divided into four stages: filamentous conchocelis (S1-0d and S2-0d), early maturation stage (S1-1d and S2-1d), middle maturation stage (S1-7d, S2-7d, and S2-49d), and late maturation stage (S1-49d and S2-91d) (Lin et al. 2021).

S1 and S2 conchocelis samples were collected on days 0, 1, 7, and 49. For S2, conchocelis samples were also collected on day 91. Three samples (0.5 g each) were collected at each time point and immediately frozen in liquid nitrogen (Each sample represents an independent conchocelis culture). For each sample, 100 mg of vacuum-dried powder was accurately weighed and dissolved in 1 mL of pre-chilled 70% methanol solution containing 0.1 mg L⁻¹ lidocaine as an internal standard. The mixture was stored overnight at 4°C and vortexed three times during this period. Subsequently, the samples were centrifuged at 10,000 ×g for 10 min at 4°C. The supernatant was carefully transferred to an EP tube and filtered through a 0.22 μm microporous membrane. From each sample, 35 μL was aliquoted and pooled to create a quality control (QC) sample, while an additional 60 μL was transferred to an injection vial for further analysis.

For chromatographic separation, 2 μL of each sample was injected into a Waters ACQUITY UPLC HSS T3 C18 column (2.1 mm × 100 mm, 1.8 μm; Waters Corp., Milford, MA, USA) maintained at 40°C with a flow rate of 0.4 mL min⁻¹. The mobile phases consisted of acidified water (0.04% acetic acid, phase A) and acidified acetonitrile (0.04% acetic acid, phase B). The gradient elution program was as follows: 95 : 5 (phase A/phase B) at 0 min, transitioning to 5 : 95 at 11.0 min, maintaining 5 : 95 until 12.0 min, and returning to 95 : 5 at 12.1 min, which was held until 15.0 min. The column effluent was directed to an ESI-triple quadrupole-linear ion trap (Q TRAP) mass spectrometer.

Mass spectrometry analysis was performed using a Q TRAP 6500 System (AB Sciex, Marlborough, MA, USA) equipped with an ESI-Turbo Ion-Spray interface operating in positive ion mode. The system was controlled by Analyst 1.6.1 software (AB Sciex). Key operational parameters included an ESI source temperature of 500°C, an ion spray voltage of 5,500 V, and a curtain gas pressure of 25 psi. Collision-activated dissociation was set to the highest level. Multiple reaction monitoring (MRM) scans were conducted with optimized declustering potential and collision energy (CE) for each MRM transition. The mass-to-charge (m/z) range was set between 50 and 1,000 (Chen et al. 2013). Missing values were imputed using half the minimum value of all metabolites across all samples. Feature retention criterion: the feature is detected in ≥50% of samples within each sample group.

Data processing, including peak detection, calibration, and filtering, was performed using Analyst 1.6.1 software (Applied Biosystems Sciex, Foster City, CA, USA). Metabolites were identified by querying internal and public databases such as MassBank, KNApSAcK, HMDB (Wishart et al. 2017), MoTo DB, and METLIN (Zhu et al. 2013). Identification was based on matching m/z values, retention times, and fragmentation patterns with reference standards. Multivariate statistical analyses, including principal component analysis (PCA) and orthogonal partial least squares-discriminant analysis (OPLS-DA), were employed to compare metabolite profiles across samples. The OPLS-DA models were rigorously validated using a 5-fold cross-validation strategy (k = 5) and assessed for statistical significance via cross-validated analysis of variance (CV-ANOVA) (p < 0.05). Model quality was evaluated using the cumulative R²X, R²Y, and Q² values, and 200 permutation tests were conducted to guard against overfitting. Changes in metabolite levels among experimental groups were visualized and analyzed. Significant differentially abundant metabolites were identified using a combination of multivariate (OPLS-DA; Variable Importance in Projection [VIP] ≥ 1) and univariate (t-test; p < 0.05) statistical approaches. Metabolites were further annotated using the Kyoto Encyclopedia of Genes and Genomes (KEGG) database to assign corresponding Compound IDs (Cid).

Algal samples were collected as described above, with four samples collected at each time point (Each sample is a separate culture of conchocelis). For each sample, 100 mg was mixed with 750 μL chloroform : methanol (2 : 1) solution in a 2 mL EP tube, which was vortexed at −20°C for 30 s, after which 100 mg glass beads were added before grinding at 60 Hz for 1 min using a grinder. The grinding was repeated two times and then the tube was placed on ice for 40 min. Next, 190 μL ddH₂O was added and then the solution was vortexed for 30 s and placed on ice for 10 min. The tube was centrifuged at 12,000 rpm for 5 min and the lower organic phase was aspirated, after which 500 μL chloroform : methanol (2 : 1) solution was added to the upper residue. Following a centrifugation at 12,000 rpm for 5 min, the lower organic phase was aspirated. The two aspirated organic phases were combined, concentrated, and dried under vacuum conditions. Prior to the liquid chromatography-mass spectrometry analysis, 200 μL isopropanol was added to the lipid extract, which was vortexed for 30 s. The supernatant was filtered through a 0.22 μm membrane and the collected filtrate was analyzed. The reproducibility of this analysis was assessed using a QC sample comprising 20 μL of each sample (Narváez-Rivas and Zhang 2016).

The chromatographic separation was performed using a Thermo Ultimate 3000 instrument with an ACQUITY UPLC BEH C18 column (1.7 μm, 2.1 × 100 mm) and an autosampler set at 8°C. Briefly, 2 μL sample was added at 50°C at a flow rate of 0.25 mL min⁻¹ and the gradient elution was completed with mobile phases comprising acetonitrile : water (60 : 40) (0.1% formic acid + 10 mM ammonium formate) (C) and isopropanol : acetonitrile (90 : 10) (0.1% formic acid + 10 mM ammonium formate) (D). The gradient elution program was as follows: 0–5 min, 70–57% C; 5–5.1 min, 57–50% C; 5.1–14 min, 50–30% C; 14–14.1 min, 30% C; 14.1–21 min, 30–1% C; 21–24 min, 1% C; 24–24.1 min, 1–70% C; 24.1–28 min, 70% C (Triebl et al. 2017).

ESI-MSⁿ analyses were performed using a Thermo QExactive Focus mass spectrometer (Thermo Fisher Scientific, Waltham, MA, USA), with a spray voltage of 3.5 and 2.5 kV in positive and negative modes, respectively. Sheath and auxiliary gas were set at 30 and 10 arbitrary units, respectively. The capillary temperature was 325°C. The Orbitrap analyzer scanned over a mass range of m/z 150–2,000 for a full scan at a mass resolution of 35,000. Data dependent acquisition tandem mass spectrometry (MS/MS) analyses were completed using HCD scan. The normalized CE was 30 eV. On the basis of dynamic exclusion, some unnecessary information was removed from MS/MS spectra.

Integrate lipidomic data with transcriptomic data (Lin et al. 2021) to screen out metabolites and genes that have significantly changed in lipid metabolic pathways (Data publicly available). The samples for lipidomics and transcriptomics were cultured separately. Biological replicates were paired based on identical culture conditions (i.e., the same strain and the same time point). The data for the 0, 1, 7, and 49 d samples correspond between the two datasets. Conchocelis was defined as fully mature and entering the late maturation time point upon reaching 90% maturity. Since the S2 strain reached full maturity on 91 d in lipidomics and on 84 d in transcriptomics under this criterion, we paired its 91 d lipidomics data with its 84 d transcriptomics data for analysis.

To identify shared metabolic pathways between genes and metabolites, we queried the KEGG database and analyzed the correlation patterns of genes and metabolites within these pathways. Utilizing comprehensive transcriptomic and metabolomic datasets, we performed bidirectional orthogonal projections to latent structures (O2PLS) analyses using the OmicsPLS package (Bylesjö et al. 2007, Bouhaddani et al. 2016). This approach allowed us to integrate transcriptomic and metabolomic data by calculating Pearson correlation coefficients. Gene-metabolite pairs were ranked based on the absolute values of their correlation coefficients in descending order.

For further analysis, the top 50 genes and metabolites with the highest correlation coefficients were selected and visualized using a heatmap generated with the pheatmap package in R (R Foundation for Statistical Computing, Vienna, Austria). Additionally, the top 250 gene-metabolite pairs exhibiting absolute Pearson correlation coefficients greater than 0.5 were used to construct a metabolite-transcript network. This network was analyzed and visualized using the igraph package in R, providing insights into the interactions and regulatory relationships between genes and metabolites within the shared pathways.

To evaluate the quality of metabolomics data, we used positive and negative ion modes for quantitative analysis of metabolites and verified the repetition and stability of the data by integral correction graphs and total ion chromatogram overlap graphs of QC samples (Supplementary Fig. S1). The results show that metabolomic data have high repetition and reliability.

P. haitanensis strains S1 (early-maturing conchocelis) and S2 (late-maturing conchocelis) were examined via MRM, with all samples analyzed qualitatively. A total of 484 metabolites were detected (Supplementary Table S1). To assess the reproducibility and stability of the metabolomic data, all samples were included in a PCA (Fig. 1A), which revealed clustering between samples at different time points, good reproducibility between samples, and close principal components. All metabolites were clustered after the z-score standardization of the data using the R package pheatmap (v1.0.12). The corresponding heatmaps were plotted (Fig. 1B) to visualize metabolite contents in all samples. Although S1-0-2 was apparently an outlier according to the PCA plot and metabolite clustering analysis, the Pearson correlation coefficient plot of all samples (Supplementary Fig. S2) showed that the correlation between S1-0-2 and the samples within the same group was 0.848 (i.e., high correlation). Therefore, the data in this group did not need to be excluded for further analyses. Because the clustering among the duplicate samples was good, this set of data was appropriate for the subsequent analysis.

Differentially abundant metabolites were screened to detect differences in metabolism. Using VIP ≥ 1 and p < 0.05 as thresholds, 123 significant differentially abundant metabolites were identified (Fig. 1C). A comparison of these metabolites between S1 and S2 on days 0, 1, 7, and 49 of the maturation period detected 107 metabolites with differing abundances in the two strains (VIP ≥ 1, p < 0.05), including lipids, amino acids and their derivatives, organic acids, nucleotides and their derivatives, alkaloids (Fig. 2A). Lipids and amino acids (and their derivatives) were the two most common differentially abundant metabolites (48.54 and 19.42% of the differentially abundant metabolites, respectively) (Fig. 2B). The levels of proline and glutamic acid in S1 were significantly higher than those in S2 at different maturation time points (Fig. 2A).

Because lipid metabolites were the main differentially abundant metabolites between the two strains, we conducted a lipid metabolomic analysis of the conchocelis of the two strains at different maturation time points. The OPLS-DA model showed excellent separation between groups, with high explained variance and robust predictive ability. The validity of the model was confirmed by a significant CV-ANOVA p-value and a permutation test (n = 200), which demonstrated that the model was not overfitted (Supplementary Table S2). The QC performance metrics for the lipid metabolites are detailed in Supplementary Table S3. The Pearson correlation coefficient plot of all samples in the lipidomics dataset (Supplementary Fig. S3) showed high correlations among samples within each group. We combined a multivariate statistical analysis and a univariate statistical analysis to screen for significant differentially abundant lipid metabolites between groups. A total of 800 lipid metabolites were identified (Supplementary Table S4), including 301 glycerolipids, 167 sphingolipids, 20 sterol lipids, 184 glycerophospholipids, 8 fatty acyls, 100 glucosylsphingosines, 3 prenol lipids, and 17 other lipids (Supplementary Table S5). Notably, the abundance of 372 lipid metabolites differed significantly among samples (VIP ≥ 1; a full list with p-values and q-values is available in Supplementary Table S6).

The top 50 differential genes with correlation coefficients with differential lipid metabolites are shown in a heat map (Supplementary Fig. S4). By comparing S1 and S2 on days 0, 1, 7, and 49, a total of 221 differentially abundant lipid metabolites were identified. The main KEGG pathways enriched among these metabolites were glycerophospholipid metabolism, alpha-linolenic acid metabolism, arachidonic acid metabolism, linoleic acid metabolism, biosynthesis of secondary metabolites, and glycosylphosphatidylinositol (GPI)-anchor biosynthesis (Supplementary Fig. S5A), of which glycerophospholipid metabolism was the most enriched pathway. Transcriptomic data generated in an earlier study by Lin et al. (2021) were analyzed to screen for genes that were differentially expressed in the conchocelis of the two strains at the same time points. The identified genes were included in a KEGG enrichment analysis, which indicated the main enriched metabolic pathways included endocytosis, PI signaling system, mitogen-activated protein kinase signaling pathway–plant, plant hormone signal transduction, and glycerophospholipid metabolism (Supplementary Fig. S5B). According to an association analysis of the lipid metabolome and transcriptome, the shared enriched KEGG pathways among the identified genes and metabolites included glycerophospholipid metabolism, alpha-linolenic acid metabolism, arachidonic acid metabolism, linoleic acid metabolism, biosynthesis of secondary metabolites, and GPI-anchor biosynthesis (Supplementary Table S7). In addition, there were 108 candidate genes and 48 lipid metabolites associated with the glycerophospholipid metabolic pathway (Supplementary Table S7). These 48 lipid metabolites included 25 PCs, 13 LPCs, and 10 phosphatidylethanolamines (PEs). Correlation analysis was performed using screening criteria of Pearson correlation coefficient absolute value > 0.5 and significance p-value < 0.05, with the results shown in Supplementary Table S8 & Fig. S6. Accordingly, the expression of differentially expressed genes involved in PC synthesis- and metabolism-related pathways was analyzed in the two strains during maturation (Fig. 3A). The analyzed genes included those encoding lysophosphatidylcholine acyltransferase, phosphatidylethanolamine N-methyltransferase, choline kinase, DGK, phosphatidate cytidylyltransferase, phosphatidylserine decarboxylase, phospholipase C (PLC), phosphatidate phosphatase, phospholipase A (PLA), lysophospholipase 1, cytosolic phospholipase A2, phospholipid-N-methyltransferase, and phospholipid : diacylglycerol acyltransferase, diacylglycerol acyltransferase, glycerol-3-phosphate acyltransferase, and lysophosphatidic acid acyltransferase. The expression levels of these genes were up-regulated in S1 on days 1 and 7 (relative to the corresponding levels on day 0), but were down-regulated on day 49. In S2, the expression of these genes was up-regulated on days 1, 7, 49, and 91 (relative to the corresponding levels on day 0) (Supplementary Table S9). Moreover, these genes were more highly expressed in S1 than in S2 on day 1 (Fig. 3B), indicating that S1 activated the PC anabolic pathway earlier than S2 during the maturation process (mean Δlog₂FC = −1.43, 95% confidence interval [CI]: −1.6 to −1.27) (Table 1). PLA-mediated hydrolysis of PC releases free FAs, including arachidonic acid and α-linolenic acid. The arachidonic acid, α-linolenic acid, γ-linolenic acid, and linolenic acid contents were higher in S1 than in S2 on days 0 and 1, but the opposite trends were observed on day 7 (Fig. 3C).

TG lipase (Unigene0042736-SDP1) gene expression was higher in S1 than in S2 on days 7 and 49 (Fig. 3B), implying that S1 may metabolize TG into FAs faster than S2. A comparison of the lipid metabolites in S1 and S2 on days 0, 1, 7, and 49 revealed 257 differentially abundant lipid metabolites, with TG accounting for the largest proportion (Fig. 4). On day 0, the TG content was significantly higher in S1 than in S2 (mean Δlog₂ FC = −1.87, 95% CI: −2.08 to −1.65) (Table 2). However, the TG content changed as the conchocelis maturation process was initiated, with significantly higher levels in S2 than in S1 on day 7 (mean Δlog₂FC = 1.37, 95% CI: 1.20 to 1.54) (Table 2). On day 49, the difference in the TG contents of the two strains narrowed. These findings suggest that TG accumulation and subsequent consumption are closely linked to the initiation of conchocelis maturation. According to the heatmap of the changes in differentially abundant lipid metabolites in the S1 conchocelis during maturation, TG accumulated on day 7, but was subsequently metabolized as the maturation process proceeded (Supplementary Fig. S7A). By contrast, the heatmap of the changes in differentially abundant lipid metabolites in the S2 conchocelis during maturation showed that TG also accumulated on day 7, but was subsequently metabolized relatively slowly (Supplementary Fig. S7B).

Several major lipid metabolites among the differentially abundant lipid metabolites underwent a statistical analysis (Fig. 5), including sulfoquinovosylmonoacylglycerol (SQMG), digalactosyldiacylglycerol (DGDG), lysobisphosphatidic acid, monogalactosyldiglycerol (MGDG), sulfoquinovosyldiacylglycerol, monogalactosylmonoacylglycerol, ceramide, sphingosine, phosphatidylethanol (PEt), coenzyme, diacylglycerol (DG), FA, monoglyceride, PC, TG, PI, and phosphatidylserine. In S1, TG accounted for 0.2, 20.9, 42.9, and 2.2% of the total lipid metabolite content on days 0, 1, 7, and 49, respectively (Fig. 5A–D). In S2, TG accounted for 29.8, 30.1, 38.8, 47.6, and 41.8% of the total lipid metabolite content on days 0, 1, 7, 49, and 91, respectively (Fig. 5E–I). Thus, the S1 conchocelis accumulates TG in the early maturation stage, but the accumulated TG is subsequently metabolized to provide the energy required for conchocelis maturation.

MGDG and DGDG are the main components of the thylakoid membrane in plants. The DGDG/MGDG ratio is an important factor related to the structure and function of photosynthetic organs because it reflects the permeability of the chloroplast membrane and the stability of the membrane bilayer. The sum of the category peak areas of all, DGDG or MGDG detected was counted, and the ratio was then calculated. During the conchocelis maturation process, the DGDG/MGDG ratio of S1 continued to increase (from 0.443 to 4.761). However, the DGDG/MGDG ratio in S2 exhibited an initial non-significant change from 1.502 (day 0) to 1.257 (day 1), peaked at 2.753 on day 7, and then experienced a subsequent drop to 0.334 by day 49 (Fig. 6).