Light and Dark Conditions Reveal Distinct Components of the Phosphate Starvation Response in Soybean Roots

1. Introduction

Phosphorus (P) is an essential macronutrient for all living organisms. As a component of nucleic acids, phospholipids, and ATP, P supports the structural and energetic foundations of life and plays key roles in energy transfer, membrane stability, and signaling in plant cells [1]. Most P is taken up as inorganic phosphate (P_i), but in many soils P_i is tightly bound to minerals or organic matter, limiting its availability for uptake [1]. To sustain crop productivity, large amounts of P fertilizers are applied, most of it mined from non-renewable geological deposits formed over millions of years. However, intensive mining and fertilizer use have led to declining global reserves and to environmental pollution through runoff and soil erosion [2,3,4].

Improving plant P acquisition and utilization is therefore critical for sustainable agriculture. Soybean (Glycine max), a major source of plant-based protein and oil [5], is particularly relevant because of its global importance and its strong dependence on P to sustain nitrogen fixation and yield. In soybean, progress toward improving P efficiency has come from quantitative trait loci (QTL) analyses [6] and molecular and transcriptomic studies [7,8,9]. These studies have revealed diverse responses to P deficiency, including changes in root architecture and exudation [10], altered leaf P-allocation, and improved P-use efficiency [8].

Under P deficiency, root system architecture changes, typically showing shorter primary roots and more lateral roots and root hairs, which increase the root surface area exposed to soil P_i [10,11]. Roots also secrete organic acids that dissolve P_i and release phosphatases and nucleases to free P_i from organic compounds [1,10,12]. At the cellular level, high-affinity P_i transporters are upregulated to improve uptake efficiency [13,14,15], while membrane lipid remodeling conserves P by shifting from phospholipids to P_i-free membrane lipids [16]. Central carbon metabolism also becomes more P_i-efficient by favoring PP_i-dependent reactions that bypass ATP [17,18]. Although P deficiency is defined at the whole-plant level, its transcriptional and metabolic responses are mediated by changes in cellular P_i availability and P_i-dependent signaling pathways.

These metabolic responses are integrated with signaling networks that coordinate P_i perception and redistribution across the plant. SPX-domain proteins and PHR transcription factors form a core regulatory module [19]. Under P_i-sufficient conditions, SPX proteins bind inositol pyrophosphates and repress PHR activity; under P_i deficiency, this repression is relieved, allowing PHRs to activate P_i-responsive genes [20]. Soybean homologs such as GmPHR1, GmSPX1, and GmSPX3 modulate P_i uptake and homeostasis [21], while PHO1 family members export P_i into the xylem and contribute to whole-plant P_i distribution [22].

Long-distance communication between shoots and roots also plays a central role in coordinating P acquisition. Sucrose, the main product of photosynthesis, functions both as a carbon source and as a signal integrating plant carbon status with nutrient responses [23,24,25]. Under P_i deficiency, increased sucrose flow from leaves to roots stimulates the expression of P_i acquisition genes, including high-affinity transporters and acid phosphatases [26,27]. In soybean, sucrose can trigger rapid transcriptomic shifts within minutes, activating hormone, calcium, and ROS signaling before canonical P_i-starvation genes respond [28]. Other systemic signals, including auxin, ethylene, strigolactones, and microRNAs such as miR399 and miR827, further integrate carbon status with P_i allocation and root development [29,30,31,32]. However, the extent to which P_i-responsive transcription depends on photosynthetically derived sucrose, as opposed to sucrose-independent local P_i sensing, remains unclear.

In addition to systemic signaling, root tips and epidermal regions are also able to sense P_i status locally [33,34]. Because sucrose transport from shoots to roots depends on photosynthesis, comparing P_i-deficient plants grown in light versus darkness offers a practical way to probe the influence of carbon availability on P_i-responsive transcription. Darkness strongly reduces photosynthetic carbon production and sucrose transport from shoot to root but also induces broader metabolic and physiological changes. Thus, light–dark comparisons provide an imperfect but useful framework to explore how carbon availability shapes different components of the P_i-starvation response.

In this study, we used Oxford Nanopore cDNA sequencing to analyze root transcriptional responses to 30 hours of P_i deprivation under light and dark conditions in hydroponically grown soybean. Because sucrose transport from shoot to root depends on photosynthesis, comparing P_i-deficient plants grown in light versus darkness offers a practical way to probe the influence of carbon availability on P_i-responsive transcription. Consistent with this rationale, darkness restricts photosynthetic sucrose production and has been shown to reduce sucrose levels in roots [35], although additional metabolic and physiological effects of darkness are expected. Specifically, we sought to i) test sucrose’s regulatory role in coordinating the phosphate starvation response (PSR) in soybean and ii) identify sucrose-independent components of the P_i deficiency transcriptome.

2. Results

2.1. Overview of Experimental Design and Data Quality

Soybean plants were grown hydroponically for two weeks and then exposed for 30 h to light or darkness under phosphate-sufficient (+P) or phosphate-deficient (−P) conditions (Figure 1a). Darkness prevents photosynthesis and limits shoot-to-root sucrose transport. Each treatment was performed with four biological replicates to ensure reproducibility across experiments.

Oxford Nanopore cDNA sequencing generated 4.8–6.9 × 10⁶ passed reads per sample, with mean quality scores of Q10–11 and N50 values of ~800–1,000 bp (Table 1, see Table S1 for additional sequencing and mapping metrics). Across all samples, sequencing yielded ~90 million passed reads (mean 5.6 million per sample) with >90% mapping efficiency. After filtering out low-abundance genes (≥10 raw reads in ≥3 replicates), 25,348 genes were retained for downstream analyses. These represent ~54% of the 47,068 protein-coding genes annotated in the Glycine max v4.0 genome. Together, these results indicate that the sequencing depth and replicate consistency were adequate for capturing major transcriptional differences across light, dark, +P, and −P treatments.

2.2. RNA-Seq Reveals Major Effects of Light and Minor Effects of P_i Deficiency

Principal component analysis (PCA; Figure 1b) showed that light condition was the primary source of transcriptomic variation: PC1 explained 76% of the variance and separated all light-grown samples from dark-grown samples. Phosphate status contributed a smaller component of variation, captured along PC2 (7%). Under light, +P and −P samples separated clearly, whereas these groups overlapped more extensively in darkness, suggesting that transcriptomic responses to P_i deficiency were reduced or less distinct when plants were grown without light.

Treatment-specific differences were examined using UpSet plots (Figure 1c–d). Consistent with the PCA, the number of differentially expressed genes (DEGs) was substantially higher between light and dark conditions than between +P and −P conditions within either light or dark. However, despite the strong light–dark effect, a distinct set of >300 genes was upregulated under −P in light-grown roots. In contrast, only four genes passed significance thresholds for −P upregulation in darkness, with 13 emerging under relaxed criteria (see Section 3.2).

This sharp contrast indicates that P_i-deficiency transcriptional responses are markedly attenuated in darkness. Because darkness restricts photosynthetic carbon supply and also alters broader metabolic and hormonal processes, the reduced −P response likely reflects a combination of limited carbon availability and additional physiological changes associated with dark-grown roots.

2.3. Light Enables Genome-Wide Activation of the Phosphate Starvation Program

To visualize transcriptional responses to P_i deficiency under each light condition, we compared −P and +P samples using MA and volcano plots (Figure 2). MA plots show overall fold-change patterns across expression levels, while volcano plots emphasize genes with larger expression shifts and stronger statistical support. An exploratory threshold of adjusted p < 0.1 was used in MA plots to highlight subtle changes, whereas volcano plots apply a more stringent adjusted p < 0.05 cutoff.

Under light, −P treatment induced extensive transcriptional changes, with hundreds of genes significantly up- or downregulated across a wide dynamic range. This pattern is consistent with activation of the canonical phosphate starvation response (PSR) in roots under light conditions. In darkness, however, only a small number of genes passed significance thresholds, suggesting that transcriptomic responses to −P are markedly reduced under dark conditions.

Because darkness restricts photosynthetic sucrose production and alters broader aspects of root physiology, the attenuated −P response likely reflects a combination of reduced carbon supply and other dark-associated changes. Despite these factors, the strong contrast between light and dark conditions is clear: light-grown roots show robust −P-induced transcriptional reprogramming, whereas dark-grown roots display a much smaller response.

2.4. Canonical Phosphate Starvation Response Genes Are Induced Only in Light

To examine how transcriptional responses to P_i deficiency differ between light and dark conditions, we compared the genes most strongly affected by −P in each treatment (Figure 3 and Figure 4). Under light, roots showed clear activation of well-established phosphate-starvation markers. Among the top 50 upregulated genes were PHT1;14 (high-affinity Pi transporter) and two PTEN2α (phosphatidylinositol 3,4,5-trisphosphate 3-phosphatase) genes. Additional canonical markers were present within the broader set of ~300 upregulated −P genes, including PHO1, SAP (secreted acid phosphatase), and PIP5K (phosphatidylinositol 4-phosphate 5-kinase). Together, these patterns indicate a coordinated transcriptional response to P_i deficiency under light.

To assess whether the lack of −P induction in darkness reflected biological suppression or low expression, we examined normalized counts for select canonical PSR genes (Table 2). Canonical PSR markers such as PHT1, PHO1, and PTEN2α showed robust induction under −P under light, but only small or no changes under −P in darkness. Importantly, all canonical PSR genes examined in Table 2 exceeded the ≥10-read threshold in at least 3 replicates, including dark-grown roots. For example, PHO1 showed 24 normalized counts in +P dark roots and 32 under −P dark, while PHT1 showed 944 and 807 normalized counts, respectively. These values demonstrate that canonical PSR genes are expressed at substantial levels even in darkness.

Only four genes passed significance thresholds for upregulation in −P in the dark, including an L10-interacting MYB domain-containing transcription factor, while numerous metabolic and stress-associated genes were downregulated. This reduced response is consistent with a general attenuation of P_i-deficiency transcriptional changes in dark-grown roots.

2.5. Sensitivity Analysis Reveals a Small Set of P_i-Responsive Genes in Darkness

To ensure that the minimal P_i response observed in darkness was not an artifact of stringent statistical filtering, we reanalyzed the dark dataset using relaxed criteria (|log₂FC| > 0.5, adjusted p < 0.1). This analysis identified 13 upregulated genes (including the original four), among them PHR11, two L10-interacting MYB domain–containing transcription factors, and several genes related to sugar metabolism (Table 3). These subtle changes support the interpretation that local or basal sensing mechanisms remain active in darkness, whereas systemic activation of the phosphate-starvation program requires light-derived carbon.

2.6. Functional Enrichment Highlights Light-Dependent Transcriptional Reprogramming During P_i Deficiency

To characterize the biological processes underlying the observed transcriptional patterns, we performed functional enrichment analyses using KEGG annotations (Figure 5). As expected from the large number of −P-responsive genes under light conditions, multiple pathways were significantly enriched among genes upregulated under −P in the light (Figure 5a). These included signaling and regulatory categories such as MAPK and plant hormone signaling, as well as several lipid- and inositol-related processes. In contrast, enrichment among −P-upregulated genes in darkness was highly limited: only glycolysis/gluconeogenesis passed the threshold for significance (Figure 5b).

Among downregulated genes, −P in light showed enrichment in categories such as oxidative phosphorylation, ubiquitin-mediated proteolysis, and ATP-dependent chromatin remodeling (Figure 5c), consistent with downregulation of energetically demanding or housekeeping pathways during P_i stress. −P in dark displayed a different pattern, with downregulated genes enriched in flavonoid biosynthesis, plant hormone signal transduction, and plant–pathogen interaction (Figure 5d).

Analyses of the light–dark comparisons (Supplementary Figures S1–S2) revealed similar global trends: photosynthesis- and energy-related pathways were enriched under light but reduced under darkness, reinforcing the dominant effect of light on root transcriptome structure.

2.7. Protein–Protein Interaction Network Identifies Core Response Modules

To examine potential functional relationships among genes upregulated under −P in the light, we constructed a protein–protein interaction network using STRING (Figure 5e). The network showed significant enrichment for known interactions (p < 1 × 10⁻⁴ to p < 1 × 10⁻¹⁶), indicating that these genes form coordinated functional modules rather than independent responses.

Clustering using the Markov Cluster Algorithm (MCL) resolved several groups organized by annotated biological functions. The largest clusters were associated with chromatin remodeling and transcriptional regulation, containing multiple transcription factors and chromatin-associated proteins, features consistent with broad regulatory control of the P_i response. A second major cluster corresponded to carbon and nitrogen metabolism, encompassing enzymes from glycolytic and amino acid pathways, in line with the metabolic reprogramming revealed by KEGG enrichment. Additional clusters were linked to signaling processes, including salicylic acid, ethylene, and phosphatidylinositol signaling, as well as calmodulin binding, highlighting crosstalk between P_i signaling and hormonal or calcium-dependent regulation. Finally, a distinct P_i-metabolism cluster contained phosphate transporters and phosphatases, representing direct components of P_i uptake and turnover.

Overall, the protein–protein interaction network highlights that the −P-responsive genes induced under light form several coordinated functional groups spanning regulation, signaling, and metabolism.

2.8. qRT-PCR Assessment of Selected P_i-Responsive Genes

Expression patterns of selected genes were evaluated by qRT-PCR (quantitative Reverse-Transcription PCR), including one gene upregulated under −P in light-grown roots (PHO1), two genes showing −P-associated induction in darkness in the RNA-seq dataset (L10-interacting MYB and SWEET5), and PHR1, which showed a slight but non-significant increase under −P in darkness. qRT-PCR confirmed that PHO1 was significantly upregulated under −P in the light, whereas expression in darkness remained low, consistent with the RNA-seq results. The L10-interacting MYB gene also showed detectable induction under −P in darkness, but not in light.

For PHR1, qRT-PCR reproduced the small increase observed in the RNA-seq dataset, but this change was not statistically significant in either assay. Such low-amplitude shifts should be interpreted cautiously, as they may fall near the limit of detection for both methods. SWEET5 expression in darkness could not be reliably quantified (Ct > 32), preventing validation; this is consistent with its very low expression in the RNA-seq data and illustrates the challenge of validating genes with extremely low transcript abundance. Overall, qRT-PCR results were broadly consistent with the major expression patterns identified by RNA-seq.

Figure 6. qRT-PCR assessment of selected P_i-responsive genes. Relative transcript levels of PHO1, L10-interacting MYB, PHR1, and SWEET5 were measured under phosphate-sufficient (+P) and phosphate-deficient (−P) conditions in light (L) and dark (D). Bars represent the mean ± SE of three technical replicates of a representative biological replicate. Expression values were normalized to a methyltransferase reference gene that showed stable expression across treatments in the RNA-seq dataset. Asterisks indicate significant differences between +P and −P within each light condition (p < 0.01). SWEET5 expression in darkness was near the detection limit (Ct > 32), and corresponding qRT-PCR values should be interpreted with caution.

Figure 6. qRT-PCR assessment of selected P_i-responsive genes. Relative transcript levels of PHO1, L10-interacting MYB, PHR1, and SWEET5 were measured under phosphate-sufficient (+P) and phosphate-deficient (−P) conditions in light (L) and dark (D). Bars represent the mean ± SE of three technical replicates of a representative biological replicate. Expression values were normalized to a methyltransferase reference gene that showed stable expression across treatments in the RNA-seq dataset. Asterisks indicate significant differences between +P and −P within each light condition (p < 0.01). SWEET5 expression in darkness was near the detection limit (Ct > 32), and corresponding qRT-PCR values should be interpreted with caution.

3. Discussion

3.1. Systemic P_i Signaling Enabled by Light

P_i deficiency triggered a broad transcriptional response in soybean roots under light, consistent with activation of the canonical Phosphate Starvation Response (PSR). Genes encoding P_i transporters such as PHT1 and PHO1 were strongly induced, in agreement with their established roles in high-affinity uptake and xylem loading, and with the known relief of PHO2-mediated repression under P_i limitation [36]. The induction of PTEN2α, phospholipases, as well as KEGG enrichment for lipid metabolism, is consistent with the shift from phospholipids to P-free membrane lipids also characteristic of the PSR [21,37]. In total, 364 genes were upregulated under −P in light, spanning transport, signaling, transcriptional regulation, and metabolism. This broad transcriptional activation under light corresponds to the “canonical” component of the P_i starvation response illustrated in Figure 7.

While sucrose-mediated regulation of individual P_i-responsive genes has been well documented in multiple plant species, most extensively in Arabidopsis [23], our results highlight the extent to which sucrose availability enables genome-wide transcriptional reprogramming during P_i deficiency. Light-dependent sucrose production is likely a major contributor to this effect, as photosynthesis supports both carbon availability and long-distance sucrose signaling from shoot to root.

Dark treatment has been widely used as an experimental approach to reduce endogenous sucrose levels in studies investigating sucrose-mediated signaling in plants [24,35,38,39]. Consistent with this rationale, 24 hour dark treatment has been shown to significantly reduce sucrose levels in Arabidopsis roots [35]. Together with these prior observations, our results support the interpretation that reduced sucrose availability limits activation of the canonical PSR in the dark.

In Arabidopsis, sucrose accumulation has been shown to be required not only for P_i-deficiency signaling but also for activation of Fe-deficiency responses, with darkness suppressing FRO2 and IRT1 induction by limiting sucrose delivery to the root [35]. These studies established sucrose as a central systemic signal linking shoot carbon status to root nutrient responses. In soybean, sucrose-feeding experiments similarly demonstrate that exogenous sucrose enhances root growth and induces metabolic signatures resembling P_i starvation [40]. Together with our observations, these findings suggest that abundant shoot-derived sucrose under light enables engagement of the canonical PSR, supporting broad transcriptional activation across multiple functional categories in soybean roots.

3.2. P_i Signaling in Darkness: A Small Sucrose-Independent Response That Includes a MYB Transcription Factor

In contrast to the broad response seen under light, P_i-deficient roots in darkness showed very limited transcriptional activation, with only four genes passing standard significance thresholds. Notably, among this small set, an L10-interacting MYB-domain transcription factor met stringent differential expression criteria (|log₂FC| > 1; padj < 0.05), and its upregulation was independently confirmed by qRT-PCR.

When relaxing statistical stringency (|log₂FC| > 1; padj < 0.1), a total of 13 upregulated genes emerged, including two additional MYB-domain transcription factors: PHR11, a homolog of the central PSR regulator PHR1 [41], and another L10-interacting MYBs. L10-interacting MYB factors have been linked to nutrient and ribosomal stress [42], and their activation in darkness suggests they may participate in a sucrose-independent P_i deficiency response.

Conservative statistical thresholds using adjusted p-values (padj) can reduce the number of detected DEGs, but even relaxed criteria revealed only a small set of −P responsive genes in darkness. This indicates that the minimal PSR observed in dark-grown roots likely reflects true biological attenuation rather than an artifact of stringent statistical filtering.

Local P_i sensing is known to occur in the root apex, involving LPR1/LPR2-dependent Fe redox chemistry, auxin transport, ROS gradients, and ethylene or strigolactone signaling [36,43]. These mechanisms typically involve only small cell populations that would be diluted in bulk RNA-seq and thus may not have been detected in our study.

Additional insight comes from 11 genes consistently downregulated under −P in both light and dark. Because sucrose-dependent systemic signaling is largely inactive in darkness, these genes likely represent a carbon-independent component of the −P response. These include proteins associated with metabolism (e.g., RuBisCO-associated protein, phosphoglycerate mutase–like protein), signaling (TAXIMIN2, a VQ-motif protein), and secondary metabolism (dirigent protein 22).

Taken together, our data support a two-component view of the soybean P_i response. Light enables a broad, carbon-supported PSR, whereas darkness reveals a much smaller, sucrose-independent transcriptional subset (Figure 7).

3.3. Limitations and Future Directions

Although sucrose availability was not directly quantified, our experimental design relied on well-established effects of darkness on photosynthetic sucrose export. Dark treatment has been widely used to reduce endogenous sucrose levels in studies of sucrose-mediated signaling, and has been shown to significantly decrease sucrose concentrations in roots after extended periods of darkness [35]. Using darkness to suppress sucrose delivery therefore provides a tractable way to probe carbon dependencies, but it also introduces additional metabolic and hormonal changes. Future experiments that supplement sucrose under darkness, or selectively disrupt shoot–root carbon transport (e.g., by girdling), will help clarify the extent to which sucrose regulates the transcriptional response to P_i deficiency.

While Nanopore cDNA sequencing provided sufficient depth to capture robust genome-wide responses under light, additional weakly responsive genes may remain below detection thresholds, particularly if expressed at low levels or restricted to specific cell types. Approaches with increased sensitivity or spatial resolution, such as cell type-resolved or spatial transcriptomic methods, will be valuable to determine whether the sucrose-independent gene set marks discrete P_i-sensing domains within the root apex.

Despite these limitations, our study reveals a clear contrast between the extensive P_i-responsive transcriptome activated under light and the minimal sucrose-independent subset detectable in darkness. Understanding how carbon availability modulates P_i uptake, membrane remodeling, and stress acclimation has practical implications for improving P_i-use efficiency, particularly under field conditions where shading or environmental stress may restrict carbon allocation to the root system.

4. Materials and Methods

4.1. Overview

This study investigated transcriptomic responses of soybean (G. max) roots to phosphate (P) availability under light and dark conditions using Oxford Nanopore long-read sequencing. Plants were grown hydroponically under controlled conditions and exposed to P deficiency for 30 hours under light or dark treatments. Root RNA was extracted for full-length cDNA library preparation, sequencing was performed on MinION flow cells, and reads were processed through a reproducible bioinformatics workflow.

4.2. Plant Growth and Treatment

Glycine max cv. Williams 82 seeds (a friendly gift from John Harada and Julie Marie Pelletier, University of California, Davis) were germinated for 3–4 days and transferred to hydroponic containers containing 8 L of half-strength Hoagland solution. The nutrient solution was replenished every 3–4 days. Plants were grown in a temperature-controlled chamber maintained at 27 °C with a 12 h light/12 h dark photoperiod.

After 10 days, cotyledons were removed to deplete stored P_i reserves. At 14 days, plants were subjected to one of four conditions: phosphate-sufficient (+P) light, phosphate-deficient (−P) light, +P dark, and −P dark, each maintained for 30 h. Each treatment was independently regrown four times to generate biological replicates, with each replicate representing a separate hydroponic experiment performed on different days. At harvest, ~5 cm segments of the main root axes, including the root tip region, were excised, immediately flash-frozen in liquid nitrogen, and stored at –80 °C until RNA extraction.

4.3. RNA Isolation and Quality Assessment

Total RNA was isolated from frozen root tissue using the RNeasy Plant Mini Kit (Qiagen, Hilden, Germany) according to the manufacturer’s protocol for “Purification of Total RNA from Plant Cells and Tissues, and Filamentous Fungi.” RNA concentration and integrity were evaluated using a Qubit 4 Fluorometer (Thermo Fisher Scientific, Waltham, MA, USA) with the Qubit RNA High Sensitivity and RNA IQ assays. The RNA IQ assay measures the proportion of large or structured RNA relative to small RNA to assess sample integrity. Only samples with RNA IQ ≥ 7 were used for sequencing; most samples had RNA IQ > 8.

4.4. cDNA Library Preparation and Quantification

For each sample, 500 ng of total RNA was used to synthesize cDNA using the cDNA-PCR Sequencing Kit V14 with Barcoding (Oxford Nanopore Technologies, Oxford, UK), following the manufacturer’s instructions. cDNA library concentration and quality were verified with an Agilent TapeStation (Agilent Technologies, Santa Clara, CA, USA). For each biological replicate, four barcoded cDNA libraries were pooled to a final amount of 50 fmol.

4.5. Nanopore Sequencing

Sequencing was performed using the Oxford Nanopore Technologies cDNA-PCR Sequencing V14 – Barcoding (SQK-PCB114.24) protocol with R10.4.1 flow cells (FLO-MIN114). Flow cells were primed and loaded with pooled cDNA according to the manufacturer’s instructions. Each MinION sequencing run was conducted for 72 h.

4.6. Nanopore Sequencing Analysis

Raw reads were processed in MinKNOW, retaining only passed reads (Q > 7) and demultiplexed during sequencing. Dorado basecalling was tested for comparison, but MinKNOW produced higher read counts and was used for all downstream analyses. Sequencing data were transferred to the Expanse high-performance computing cluster (San Diego Supercomputer Center, USA).

Read orientation and trimming: Full-length cDNA reads were identified and oriented using Pychopper v2.7 (RRID:SCR_018966; https://github.com/epi2me-labs/pychopper) with the -m edlib setting.

Read alignment: Reads were aligned to the G. max reference genome (GCF_000004515.6 v4.0; downloaded 1 Oct 2025) using Minimap2 v2.30-r1287 [44]. Alignments were sorted and indexed with Samtools v1.20 [45].

Gene quantification: Gene counts were obtained using featureCounts v2.0.6 [46] in long-read mode (-L) with single-end, stranded settings (-s 1). The resulting count matrix and summary tables were exported for downstream analysis in R.

Quality control: Mapping statistics were evaluated using samtools flagstat and featureCounts summary reports.

4.7. Differential Expression Analysis

Differential expression analyses were performed in R (v4.4.0) within RStudio. Raw count matrices from featureCounts were imported into the DESeq2 (v1.42.0) package [47] for normalization and differential expression analysis. Before differential expression analysis, genes were prefiltered by requiring at least 10 raw counts in at least 3 biological replicates within any one of the four treatment groups. This removed low-abundance genes prior to DESeq2 dispersion estimation. Normalized counts were obtained using the variance-stabilizing transformation (VST) in DESeq2. Genes with adjusted p < 0.05 and |log₂ fold change| ≥ 1 were considered significantly differentially expressed. Genes used to exemplify read depth (Table 2) were selected based on their well-established functions in P_i uptake, transport, or transcriptional signaling, ensuring that representative and biologically relevant PSR components were included.

Principal component analysis (PCA) was performed on variance-stabilizing transformed (VST) counts using DESeq2. The PCA was calculated using the top 500 most variable genes, and the first two principal components were visualized with ggplot2. Replicate identifiers were added to the plot to assess clustering within and across treatments.

UpSet plots: Overlaps among DEG sets were visualized using UpSetR v1.4.0 [48] with separate plots for up- and downregulated genes.

Heatmaps: Expression patterns were visualized using ComplexHeatmap v2.20.0 [49], with per-gene z-score scaling of VST counts.

GO enrichment: Gene Ontology (GO) enrichment was performed with g:Profiler (organism = “gmax”) on up- and downregulated DEG sets (FDR < 0.05). The background set consisted of all expressed genes. Results were summarized by Biological Process terms with gentle term collapsing and flattened list columns before export.

KEGG enrichment: KEGG pathway analysis was performed using clusterProfiler v4.10.0 [50] mapping Glyma IDs to KEGG Orthology (KO) terms via KEGGREST. Enrichment was tested using DEGs as the query set and all expressed genes as background.

STRING network analysis: Protein–protein interaction (PPI) networks were analyzed using STRING v12.0 [51] with G. max as the reference organism. Upregulated genes under −P light were uploaded to identify functional clusters using the Markov Cluster Algorithm (MCL, default inflation). Significant network enrichment was defined by STRING-calculated p-values ranging from < 1 × 10⁻⁴ to < 1 × 10⁻¹⁶.

4.8. Validation by qRT-PCR

For each sample, 2000 ng of total RNA was treated with RNase-free DNase to remove genomic DNA, then reverse-transcribed using the iScript™ gDNA Clear cDNA Synthesis Kit (Bio-Rad, Hercules, CA, USA). A methyltransferase (LOC100819876) gene was selected as reference based on stable expression across all conditions. Primers were designed in Primer-BLAST, preferentially spanning exon junctions or framing introns, with amplicon sizes of 70–220 bp and minimal 3′ complementarity to prevent primer–dimer formation.

Primer sequences:

Methyltransferase 17 (LOC100819876): F TTGAGCCCTCACGAGGATTC;

R CTTATGTACACCCGGCGCA;

PHO1: F ACAGAGCCAGCTTATATGGACAT;

R CAATGCTGTGCTTGGTGAGA;

L10-interacting MYB-domain protein (LOC100820117): F TCTCAACTGAGCGGAAGTGC; R TCATTTTGCAACCTCCTGCG;

PHR11: F CCAGAATCAGATGACTGAAAGCTTG;

R GTTCAGCTTTTGCAAGTTCCATA;

SWEET5: F AGGACTACTACATCGCTCTACCA;

R TCCTCAATAGTTATCTCAGTCACAG

qPCR was performed using SYBR Green Supermix (Bio-Rad) on a CFX96 instrument under the following program: 95 °C for 30 s; 40 cycles of 95 °C for 15 s and 60 °C for 30 s. A melt curve analysis followed from 65 °C to 95 °C with 0.5 °C increments every 5 s. CFX Maestro Software (Bio-Rad) was used to obtain Cq values and calculate relative expression using the ΔΔCq method.

4.9. Statistical Analysis

All statistical analyses were conducted in R (v4.4.0). For RNA-seq data, differential expression was evaluated with the Wald test in DESeq2, using Benjamini–Hochberg FDR correction (padj < 0.05, |log₂FC| ≥ 1). Multivariate relationships among samples were examined by PCA. For enrichment analyses, overrepresentation of GO terms and KEGG pathways was assessed in clusterProfiler with adjusted p-value < 0.05 as the significance threshold.

5. Conclusions

Soybean roots require photosynthetically derived sucrose to mount the canonical phosphate starvation response (PSR). Under P deficiency in light, abundant sucrose transport from shoot to root supports coordinated induction of P_i transporters, signaling components, lipid-remodeling enzymes, and metabolic pathways. In contrast, darkness, which greatly reduces photosynthetic sucrose transport, suppresses these responses and reveals only a small, sucrose-independent subset that includes a MYB transcription factor and genes linked to sugar metabolism. Together, these findings highlight the importance of light in shaping early P_i-deficiency responses in soybean roots.