LLPSRGs (n=3,633) were retrieved from the Data Resource of LLPS (DrLLPS; http://llps.biocuckoo.cn/). Transcriptomics profiles and clinical records from an OSCC cohort (n=302), whose primary lesion sites were the tongue, base of the tongue, buccal mucosa, oral cavity, floor of the mouth and gingiva, were obtained from The Cancer Genome Atlas (TCGA; http://www.cancer.gov/ccg/research/genome-sequencing/tcga). The GSE213862 dataset from the Gene Expression Omnibus (GEO; http://www.ncbi.nlm.nih.gov/geo/) was also acquired, comprising transcriptomics data from 47 pairs of OSCC tissues and normal mucosal tissues (15).

Consensus clustering analysis was conducted using the R package ‘ConsensusClusterPlus’ (16) to stratify TCGA-OSCC samples into distinct molecular clusters (C1 and C2) based on the expression profiles of LLPSRGs. Cluster stability was evaluated by analyzing the cumulative distribution function (CDF) of the consensus matrix across multiple k-values. Differentially expressed genes (DEGs) between C1 and C2 were identified using the R package ‘DESeq2’ (17), based on the criteria of |log2 fold change (FC)|>1 and P<0.05. Survival analysis between C1 and C2 were conducted using the R package ‘survival’ (https://github.com/therneau/survival).

TCGA-OSCC cohort was randomly divided into training and validation subsets at a ratio of 7:3. Univariate Cox regression and LASSO regression analyses (implemented via the R packages ‘survival’ and ‘glmnet’ (https://glmnet.stanford.edu) were applied to the training cohort to screen seven LLPSRGs significantly associated with OSCC prognosis. A risk score model was established using LASSO-derived coefficients, and patients in both the training and validation cohorts were stratified into high- or low-risk groups based on the median risk score. Model performance was evaluated through Kaplan-Meier curve analysis with Wilcoxon test, time-dependent receiver operating characteristic (ROC) curves, principal component analysis (PCA) and t-distributed stochastic neighbor embedding (t-SNE) using ‘survival’, ‘timeROC’ (https://github.com/cran/timeROC), ‘Rtsne’ (https://github.com/jkrijthe/Rtsne) and the related visualization tool ‘ggplot2’ (https://github.com/tidyverse/ggplot2). The prognostic nomogram was constructed by integrating the continuous risk score with established clinical variables: Tumor-Node-Metastasis (TNM) stage, histological grade and age, using R package ‘rms’ (https://github.com/harrelfe/rms).

Differential expression analysis between the high- and low-risk groups was conducted using the R package ‘DESeq2’ (16) and ‘limma’ (https://bioinf.wehi.edu.au/limma), based on the criteria of P<0.05 and |log2FC|>1. The Minus-vs.-Add (MA) plot visualizing differential expression patterns between risk groups was generated using the R package ‘limma’. Gene Ontology (GO), Kyoto Encyclopedia of Genes and Genomes (KEGG) enrichment analyses and Hallmark GSEA based on the Molecular Signatures Database (MsigDB; http://www.gsea-msigdb.org/gsea/msigdb) of DEGs was conducted using R packages ‘enrichplot’ (https://github.com/YuLab-SMU/enrichplot) and ‘clusterProfiler’ (https://github.com/YuLab-SMU/clusterProfiler).

The scRNA-seq data analyzed in the present study were accessible via the GEO accession number, GSE234933 (18). Based on clinical metadata from GSE234933, samples were categorized into the primary lesion (sample count, 19; cell count, 123,094) and metastatic lesion (sample count, 8; cell count, 38,987) groups. The R package ‘Seurat’ (https://github.com/satijalab/seurat; version 4.3.0) was used for initial quality control removing cells with <200 detected genes. Validated cell data were integrated into a Seurat object for downstream analysis. Distinctly expressed genes were detected with single-cell resolution using the FindAllMarkers tool implemented in Seurat.

The GEPIA2 website (http://gepia2.cancer-pku.cn/) (19) was used to compare the relative expression of model genes in TCGA-OSCC. The IDRs of model genes were identified using IUPred2, while their disordered binding regions were characterized via ANCHOR2 analysis through the IUPred2A platform (https://iupred2a.elte.hu/) (20). LLPS propensity scoring of the model genes was performed using PhaSePred (http://predict.phasep.pro/) (21), incorporating the following computational algorithms: catGRANULE, PdPS-10fea, SaPS-10fea, PdPS-8fea, SaPS-8fea, DeepPhase, Phos, Charged residue, SEG, Hydropathy, ESpritz-DisProt, PScore and PLAAC.

PHLDB2 protein solutions (cat. no. Ag27331; Proteintech Group, Inc.) was conjugated with Alexa Fluor 488 (cat. no. ab236553; Abcam). Fluorescently labeled recombinant PHLDB2 protein was purified and reconstituted in NaCl solutions at physiologically relevant concentration gradients, ranging from 50 to 500 mM. Subsequently, the preparations were immediately combined with a phase separation buffer containing 10% glycerol, 50 mM Tris-HCl (pH 7.4), 1 mM DTT and 10% PEG-8000, serving as a molecular crowding agent. After incubation at room temperature for 5 min, the mixture was pipetted onto glass slides prepared with coverslips secured by double-sided adhesive spacers. Imaging was performed using a confocal microscope (LSM980; Zeiss AG) equipped with a 63X objective lens.

The present study employed two human OSCC lines: SCC1 and HSC3, obtained from the University of Michigan (Ann Arbor, MI, USA). HSC3 cells were cultured in DMEM (cat. no. C11995500BT; Gibco; Thermo Fisher Scientific, Inc.) and SCC1 cells were cultured in DMEM/F-12 medium (cat. no. C11330500BT; Gibco; Thermo Fisher Scientific, Inc.). Both media were supplemented with 10% FBS (cat. no. 10099-141C; Gibco; Thermo Fisher Scientific, Inc.) and 1% penicillin-streptomycin (cat. no. 15140122; Gibco; Thermo Fisher Scientific, Inc.). Cells were maintained at 37°C in a humidified incubator containing 5% CO₂.

The oligonucleotide sequences for siRNA and the PHLDB2-EGFP plasmid are provided in Table SI. SCC1 and HSC3 OSCC cells at 60–80% confluence were transfected with PHLDB2-targeting siRNAs or control siRNAs (50 nM) and the PHLDB2-EGFP plasmid (1 µg/ml) using Lipofectamine^® RNAiMAX (cat. no. 13778150; Invitrogen; Thermo Fisher Scientific, Inc.) according to the manufacturer's protocol. Transfection was performed at 37°C for 48 h and subsequent experiments were conducted 24 h post-transfection.

FRAP analysis was performed using an LSM980 confocal microscope (Zeiss AG). SCC1 cells expressing EGFP-tagged PHLDB2 were visualized with a 63X oil-immersion objective. Phase-separated condensates were selectively bleached using a 488-nm laser at 98% intensity. Post-bleaching recovery was tracked by capturing sequential images at 1-sec intervals for 100 sec. Droplet fluorescence intensities were analyzed using ZEN Microscopy Software (version 3.0; Carl Zeiss AG), and recovery kinetics were plotted with GraphPad Prism 9.0 (Dotmatics).

Total RNA isolation was carried out using TRIzol^® reagent (cat. no. 15596018; Invitrogen; Thermo Fisher Scientific, Inc.), followed by cDNA synthesis using the PrimeScript RT Reagent Kit (cat. no. RR047A; Takara Bio, Inc.) according to the manufacturer's protocol. qPCR was conducted on a Roche LightCycler 96 platform (Roche Diagnostics) employing SYBR Green PCR Master Mix (cat. no. 11201ES08; Shanghai Yeasen Biotechnology, Co., Ltd.). The qPCR amplification protocol consisted of an initial denaturation step at 95°C for 5 min, followed by 40 cycles of denaturation at 95°C for 10 sec, annealing at 60°C for 20 sec and extension at 72°C for 20 sec. GAPDH was used as the reference gene for normalization. Relative gene expression was calculated using the 2^−ΔΔCq method (22). Primer sequences, designed against GenBank (https://www.ncbi.nlm.nih.gov/genbank) reference data, are compiled in Table SII.

Protein extraction was performed with ice-cold RIPA lysis buffer (cat. no. CW2333S; Jiangsu CoWin Biotech Co., Ltd.) containing protease inhibitor cocktail. Protein concentrations were determined using a BCA quantification kit (cat. no. CW0014S; Jiangsu CoWin Biotech Co., Ltd.) according to standardized protocols. The lysates (control-si 5.26 mg/ml; PHLDB2si-1 4.94 mg/ml; PHLDB2si-2 5.04 mg/ml) were electrophoretically separated by SDS-PAGE on 10% gels, followed by semi-dry transfer to PVDF membranes (cat. no. ISEQ00010; MilliporeSigma). The membranes were then subjected to blocking with 5% milk in TBS with 0.05% Tween at room temperature for 1 h, and incubated with primary antibodies at 4°C for 16 h. After washing and incubation with HRP-conjugated secondary antibodies at 25°C for 1 h, immunoreactive bands were detected using ECL substrate (cat. no. WBKLS0500; MilliporeSigma) and semi-quantified with a GeneGnome XRQ imaging system (Syngene; software version, GeneSys 2023). The following primary antibodies were used: Rabbit polyclonal anti-PHLDB2 (cat. no. ab234885; 1:1,000; Abcam), E-cadherin (cat. no. ab231303; 1:2,000; Abcam), N-cadherin (cat. no. ab76011; 1:5,000; Abcam), vimentin (cat. no. ab92547; 1:1,000; Abcam) and mouse monoclonal anti-GAPDH (cat. no. ab8245; 1:5,000; Abcam). The following secondary antibodies were used: Goat Anti-Rabbit IgG H&L (HRP) (cat. no. ab205718; 1:10,000; Abcam) and Goat Anti-Mouse IgG H&L (HRP) (cat. no. ab205719; 1:10,000; Abcam).

Cell proliferation was assessed using the CCK-8 assay. HSC3 and SCC1 cells (2×10³ cells/well) were plated and treated with 10% CCK-8 reagent (Sigma-Aldrich; Merck KGaA) in DMEM. Absorbance was measured after a 2 h-incubation to determine cell viability. For colony formation analysis, after 1 week of culturing cells seeded at an initial density of 1,000 cells per well in 6-well plates, colonies were fixed with 4% paraformaldehyde for 20 min and stained with 0.5% crystal violet (cat. no. C3886; Sigma-Aldrich; Merck KGaA) for 20 min at room temperature, and images were then captured. Colonies defined as cell aggregates containing ≥50 cells within a 0.3–1.0 mm diameter were quantified using ImageJ software (version 1.54p; National Institutes of Health).

Cells were plated in ultra-low attachment 96-well plates at graded densities (1×10³, 5×10², 2×10², 1×10² and 5×10¹ per well) using serum-free DMEM/F-12 containing 2% B27 (cat. no. 17504044; Gibco; Thermo Fisher Scientific, Inc.), 20 ng/ml hEGF (cat. no. 105-04; PrimeGene) and 10 ng/ml hFGF (cat. no. 104-02; PrimeGene). Cultures were maintained under 5% CO₂ at 37°C with quintuplicate wells/density. After a 10-day incubation, images of spheroids were captured using an inverted light microscope (Zeiss Axio Observer; Carl Zeiss AG), and quantified using ImageJ software (version 1.54p; National Institutes of Health). Cancer stem cell frequency was determined using ELDA computational analysis (http://bioinf.wehi.edu.au/software/elda/).

For wound healing evaluation, OSCC cells were cultured to 90% confluence in 6-well plates, and a linear wound was introduced to monolayers using a 200-µl pipette tip. Migration progression was documented at 0, 12, 16 and 24 h intervals through phase-contrast light microscope (Zeiss Axio Observer; Carl Zeiss AG), with wound closure quantified using ImageJ software (version 1.54p; National Institutes of Health). For Transwell invasion assays, the upper chambers were precoated with Matrigel (cat. no. 354234; Corning, Inc.) at 4°C and incubated at 37°C for 30 min. Subsequently, 2×10⁴ cells/well were loaded into Matrigel-precoated upper chambers (cat. no. 354234; Corning, Inc.) in serum-free medium. The lower chambers contained 10% FBS as a chemoattractant. After a 48-h incubation at 37°C, cells on the membrane were fixed with 4% paraformaldehyde at room temperature for 20 min and stained with 0.1% crystal violet for 20 min at room temperature. Invasive cells were counted under bright-field microscopy.

RNA was extracted from SCC1 cells, with RNA purity quantified by Nanodrop (Thermo Fisher Scientific, Inc.) and integrity assessed using an Agilent 4200 TapeStation (Agilent Technologies, Inc.). Libraries were constructed from 1 µg RNA using the TruSeq Stranded mRNA Library Prep Kit (cat. no. 20020594; Illumina, Inc.) according to the manufacturer's protocol. Briefly, mRNA underwent fragmentation using divalent cations, followed by first-strand cDNA synthesis with reverse transcriptase and second-strand synthesis using DNA polymerase I/dNTPs. Double-stranded cDNA was processed through end repair, A-tailing and adapter ligation (200–300 bp size selection). Libraries were size-selected (200–300 bp) using AMPure XP beads. After 12 cycles of PCR amplification, final libraries were purified and quality-controlled using an Agilent 2100 Bioanalyzer High Sensitivity DNA Assay (cat. no. 5067-4626; Agilent Technologies, Inc.). The loading concentration of the final library was 14 nM for control-si group and 15 nM for PHLDB2-si group. Quality-controlled libraries were pooled and diluted to 1.8 nM in 10 mM Tris-HCl (pH 8.5) based on Qubit dsDNA HS Assay measurements (cat. no. Q32854; Thermo Fisher Scientific, Inc.). Paired-end sequencing (2×150 bp) was performed on an Illumina NovaSeq 6000 system using the NovaSeq 6000 SP Reagent Kit v1.5 (100 cycles; cat. no. 20027464) (both from Illumina, Inc.). Bioinformatics analysis identified DEGs between the control group and the PHLDB2-knockdown group using DESeq2 (|log2FC|>1; P<0.05), followed by functional enrichment analysis to characterize biological pathways using the R package ‘clusterProfiler’ (https://github.com/YuLab-SMU/clusterProfiler). Core PI3K-Akt signaling components were retrieved from the KEGG PI3K-Akt signaling pathway gene set. Differential expression of these components between control and PHLDB2 knockdown groups was analyzed with the R package ‘limma’, followed by visualization as a heatmap using the R package ‘ggplot2’.

The Sun Yat-sen University (SYSU)-OSCC cohort retrospectively analyzed clinicopathological data from 51 treatment-naïve patients treated between June 2016 and June 2018. The cohort comprised 35 men and 16 women, with a mean age of 52 years (range: 26–77 years). These data were used for validation in the present study. All participants underwent surgical resection of primary tumors with cervical lymph node dissection at the Hospital of Stomatology, SYSU (Guangzhou, China) with no presurgical adjuvant therapy. Clinicopathological parameters and survival data were retrieved through clinical documentation and structured telephone follow-ups. Tumor staging followed the Union for International Cancer Control 8th edition TNM guidelines (23). The present study received approval from the ethics committee of the Hospital of Stomatology, SYSU (approval no. KQEC202013) and all participants provided written informed consent.

IHC staining was performed on 4-µm-thick paraffin-embedded tissue sections. Sections were blocked with 10% goat serum (cat. no. SL038; Beijing Solarbio Science & Technology Co., Ltd.) in PBS for 30 min at room temperature. Primary antibodies were diluted as follows: Anti-PHLDB2 (1:50; cat. no. ab234885; Abcam) and IgG isotype control antibody (1:100; cat. no. ab172730; Abcam), and incubated overnight at 4°C. Following primary antibody incubation, sections were incubated with the appropriate secondary antibody (1:100 dilution) for 30 min at room temperature. IHC staining was then conducted with anti-rabbit/mouse HRP/DAB (ABC) detection kits (cat. nos. ab64261 and ab64259; Abcam) according to the manufacturer's protocol. Secondary antibodies were diluted as follows: Goat Anti-Rabbit IgG H&L (1:5,000; cat. no. ab205718; Abcam) and Goat Anti-Mouse IgG H&L (1:5,000; cat. no. ab205719; Abcam). Tissue section images were acquired using an Aperio AT2 scanner (Leica Biosystems) with its integrated microscope and analyzed with ImageJ software (version 1.54p; National Institutes of Health). IHC scoring followed a four-tier system based on the staining intensity: 0, Negative; 1, weakly positive; 2, positive; and 3, strongly positive.

Data were processed and analyzed with R (v4.2.2; http://www.r-project.org) and GraphPad Prism 9.0. Statistical significance was evaluated via two-tailed tests. Comparisons between two groups were performed using unpaired Student's t-test for normally distributed data or Mann-Whitney U tests for non-normally distributed data. Paired comparisons between two groups employed paired Student's t-test for normally distributed data or Wilcoxon signed-rank tests for non-normally distributed data. Multi-group comparisons were analyzed by one-way ANOVA with Tukey's post-hoc testing for normally distributed data with equal variance or Kruskal-Wallis test with Dunn's post-hoc correction for non-normally distributed data. Normality was assessed by the Shapiro-Wilk test and homogeneity of variance was verified using Levene's test. Continuous variables were examined by Pearson's correlation, and categorical variables were evaluated by Spearman's correlation. Analysis of paired tumor-normal tissues in TCGA-OSCC was conducted using the R package ‘TCGAplot’ (https://github.com/tjhwangxiong/TCGAplot). P<0.05 was considered to indicate a statistically significant difference.

The workflow of the present study is illustrated in Fig. 1. Through consensus clustering analysis of 3,633 LLPSRGs from the DrLLPS database of transcriptomics data from 302 OSCC cases in TCGA, the samples were successfully stratified into two biologically distinct molecular clusters, C1 and C2 (Fig. 2A). Evaluation of the CDF curve confirmed optimal consistency and cluster stability at k=2 (Fig. S1A). Further analysis revealed that LLPSRGs were significantly upregulated in C1 compared with in C2 (Fig. 2B), and C1 patients exhibited significantly shorter overall survival (OS; Fig. 2C). GO and KEGG enrichment analyses demonstrated that DEGs in the C1 cluster were predominantly enriched in biological processes associated with tumor progression, including ‘cell junction assembly’, ‘PI3K-Akt signaling pathway’ and ‘Wnt signaling pathway’ (Fig. S1B and C). These findings suggested that dysregulated LLPSRGs may influence OSCC prognosis by modulating key oncogenic pathways.

Based on these insights, seven LLPSRGs significantly associated with OSCC prognosis were identified using Cox proportional hazards regression (Fig. 2D) and their predictive efficacy was further validated via LASSO regression (Fig. 2E). To establish a robust model, TCGA-OSCC cohort was randomly divided into training (n=211) and validation (n=91) groups at a 7:3 ratio. A novel prognostic prediction model was developed through LASSO regression analysis using the training set data, ultimately yielding the following mathematical formula: RiskScore=Exp(SYT1) × 0.126671488410991 + Exp(MYO1B) × 0.116232649674581 + Exp(PHLDB2) × 0.111345681109828 + Exp(CORO2B) × 0.0992669799054353 + Exp(CDH13) × 0.0701283734822018 + Exp(DDX21) × 0.0342177920472229 + Exp(PKLR) × −0.202003545600354. After stratifying TCGA-OSCC patients into high- and low-risk groups, the gene expression heatmap reveals distinct expression patterns of model genes demonstrating clear distinctions between the two subgroups (Fig. 2F).

To validate the predictive performance of the developed model, risk scores were calculated for all patients, with samples stratified into high- and low-risk groups based on median values (Fig. 2G and K). PCA and t-SNE revealed distinct spatial separation between the two groups in gene expression profiles (Fig. S2A-D), confirming the biological discriminative power of the model. Survival analysis demonstrated markedly higher mortality in the high-risk group (Fig. 2H and L) and significantly reduced OS rates compared with those in the low-risk group (Fig. 2I and M). Time-dependent ROC curves showed area under the curve (AUC) values of 0.727, 0.737 and 0.604 for 1-, 3- and 5-year survival predictions in the training set, and 0.672, 0.663 and 0.712 in the validation set, indicating high sensitivity and specificity (Fig. 2J and N).

To enhance clinical applicability, the risk score was integrated with traditional clinicopathological parameters (TNM stage, histological grade and age) to construct a nomogram model (Fig. 2O). Calibration curves demonstrated strong concordance between predicted and observed survival probabilities at 1, 3 and 5 years (Fig. 2P). The nomogram model exhibited significantly improved AUC values compared with the risk score alone (Fig. 2Q), highlighting its superior accuracy in individualized prognosis prediction. These results underscore the potential of the LLPSRG-based prognostic model to refine clinical risk stratification and inform targeted therapeutic strategies for OSCC.

To investigate the mechanistic basis of survival disparities between the high- and low-risk cohorts, comparative expression profiling of LLPSRGs across both groups was conducted. The results revealed 401 significantly upregulated genes and 60 downregulated genes in the high-risk group (P<0.05; |log2FC|>1; Fig. 3A). The MA plot demonstrated no significant association between gene expression levels and fold changes (Fig. 3B), ruling out technical bias as a cause of false-positive results for highly expressed genes. According to the core characteristics of LLPS, proteins capable of phase separation require local critical concentrations and rely on IDRs to promote droplet-like condensate formation. Therefore, the functional properties of upregulated genes were focused on.

Hallmark GSEA based on the MsigDB confirmed that EMT was significantly enriched [normalized enrichment score (NES), 2.45; Fig. 3C]. KEGG pathway analysis further demonstrated significant enrichment of these genes in the ‘PI3K-Akt signaling pathway’, ‘ECM-receptor interaction’ and ‘focal adhesion’ (Fig. 3D), suggesting their potential role in promoting cell migration and invasion through ECM remodeling.

Consistent with the aforementioned findings, GO analysis revealed significant enrichment of upregulated genes in biological processes involving ECM remodeling and cell adhesion regulation (Fig. 3E). Notably, GO enrichment analysis of downregulated genes revealed pronounced enrichment in pathways governing epithelial differentiation and keratinization processes (Fig. 3F), underscoring the suppression of epithelial lineage commitment, which corroborates the previously observed EMT pathway activation.

To investigate the impact of the seven LLPSRGs in the prognostic model on the malignant progression of OSCC at a single cell level, systematic analyses were conducted using the public scRNA-seq dataset, GSE234933. Primary tumor samples (n=19, total cells=123,094) and metastatic tumor samples (n=8, total cells=38,987) were selected based on clinical information. Risk scores for all cells were calculated using the established prognostic model formula. Results revealed significantly higher risk scores in metastatic tumor cells compared with primary tumor cells (Fig. 4A). When stratifying cells into high- and low-risk groups based on the median risk score, the metastatic tumor group exhibited a significantly higher proportion of high-risk cells (Fig. 4B; 61.9 vs. 44.2% in primary tumors).

Given the aggressive phenotype and elevated risk scores of metastatic cells, subsequent analyses focused on the metastatic cell population. Unsupervised clustering and canonical marker genes were used to annotate major cell types (Fig. 4C), with marker gene expression heatmaps (such as fibroblast marker COL1A1 and cytotoxic T cell marker CD8B) confirming accurate classification (Fig. 4D). Risk profiling across all cell clusters through computing risk model algorithms revealed pronounced risk score enrichment within Uniform Manifold Approximation and Projection (UMAP)-defined tumor cell niches (Fig. 4E), guiding subsequent investigations to focus on tumor cell subpopulations.

Tumor cells were subdivided into the following five functional subtypes: i) Cell cycle; ii) epithelial differentiation; iii) hypoxia; iv) immune suppression; and v) partial (p)EMT, based on canonical markers (Fig. 4F). Heatmaps of subtype marker genes (such as CDC20 marking the cell cycle subtype and KRT5 marking epithelial differentiation) validated classification accuracy (Fig. 4G). Risk scores were calculated across tumor cell subpopulations (Fig. 4H), revealing a pronounced risk elevation in the pEMT subpopulation (Fig. 4I). Consistently, when tumor cells were stratified by median risk score, high-risk tumor cells exhibited a higher proportion of the pEMT subtype (Fig. 4J; 16.2 vs. 10.3% in low-risk), a subtype strongly associated with tumor invasion and metastatic potential (19). Mechanistically aligning with prior pathway enrichment findings, specific activation of the EMT signature pathway was also observed in the pEMT subpopulation (Fig. 4K; NES, 1.39), indicating that the EMT activation detected in bulk transcriptomics analyses predominantly originates from the pEMT subpopulation within tumor cells.

In summary, single-cell analyses demonstrated that high-risk cell populations foster a pro-metastatic environment in OSCC through activation of pEMT tumor cell subpopulations.

To delineate core oncogenic drivers among the seven prognostic model genes that drive the malignant progression of OSCC, the GEPIA2 website (20) was used to compare the relative expression of the seven LLPSRGs in TCGA-OSCC (Fig. S3A), prioritizing four highly expressed candidates, MYO1B, DDX21, CDH13 and PHLDB2. RT-qPCR of the seven genes in the OSCC cell lines verified the aforementioned results (Fig. S3B and C). Subsequently, integrative prognostic mapping for TCGA-OSCC by GEPIA2 (Fig. S3D) further refined the key regulator candidates to MYO1B, CDH13 and PHLDB2.

To evaluate their phase separation potential, the amino acid sequences of these proteins were obtained from the UniProt database (https://www.uniprot.org/) and uploaded to the IUPred2A platform (17) to analyze their IDR propensity. The prediction revealed markedly higher phase separation potential for PHLDB2 compared with MYO1B and CDH13 (Figs. 5A, S3E and S3F). In addition, LLPS scoring for PHLDB2, MYO1B and CDH13 by PhaSePred (18) supported the aforementioned results (Figs. 5B, S3G and S3H), prompting further in-depth investigation of PHLDB2.

To validate whether PHLDB2 undergoes LLPS, its ability to form condensates in vitro was first examined. In LLPS buffer containing 50 mM NaCl, Alexa Fluor 488-labeled PHLDB2 protein formed phase-separated droplets ranging 2.5–15 µm in diameter, with droplet size positively associated with protein concentration (2.5–15 µM; Fig. 5C) and inversely associated with NaCl concentration (50–500 mM; Fig. 5D). Notably, live-cell imaging of the SCC1 OSCC cell line transfected with PHLDB2-EGFP plasmids demonstrated the formation of cytoplasmic droplets (Fig. 5E). FRAP assays further confirmed the liquid-like properties of these condensates, with fluorescence recovering to 31.4% at 5 sec, 52.6% at 30 sec and 62.7% at 60 sec post-bleaching (Fig. 5F and G). These results conclusively demonstrated that PHLDB2 forms LLPS condensates in OSCC cells.

To investigate the regulatory role of PHLDB2 in the biological behavior of OSCC cells, two sequence-specific siRNAs were designed to knockdown PHLDB2 expression in SCC1 and HSC3 cell lines via transient transfection. RT-qPCR and western blotting results demonstrated that both siRNAs significantly reduced PHLDB2 mRNA and protein levels (Fig. 6A and B), confirming the high efficiency and specificity of the knockdown system.

CCK-8 assays revealed a significant decline in proliferation rates over time in the PHLDB2 knockdown groups (Fig. 6C and D). Tumor-sphere formation assays further showed that the number of tumor-spheres and their average diameter were significantly reduced in the knockdown groups compared with those in the control groups (Fig. 6E-H). Limiting dilution assays confirmed impaired self-renewal capacity of tumor stem cells following PHLDB2 knockdown (Fig. 6I and J). Additionally, colony formation assays indicated a significant decrease in colony numbers in the knockdown groups (Fig. 6K and L), suggesting the critical role of PHLDB2 in maintaining tumor stemness.

For tumor invasion, Transwell assays showed a significant reduction in matrix gel-penetrating cells in the PHLDB2 knockdown groups (Fig. 6M and N). Wound healing assays further revealed decreased wound closure rates in knockdown cells (Fig. 6O-R). These results systematically demonstrated that PHLDB2 promotes OSCC malignant progression by enhancing cell proliferation, stemness maintenance and invasive potential.

To elucidate the molecular mechanisms by which PHLDB2 promotes malignant phenotypes in OSCC cells, a multi-dimensional data cross-validation strategy was employed that systematically compared RNA-seq data from PHLDB2 knockdown experiments with transcriptomics data from TCGA-OSCC cohort. For TCGA-OSCC public data, patients were stratified into high- and low-expression groups based on PHLDB2 levels, followed by differential gene analysis identifying 3,147 distinctively upregulated genes and 15,398 distinctively downregulated genes in the PHLDB2 low-expression group (Fig. 7A). Subsequent Hallmark and GO GSEA based on MsigDB were performed on these DEGs, collectively named the PHLDB2-TCGA pathway enrichment results (Fig. 7C). SCC1 cells were transiently transfected with PHLDB2-si1 siRNA to knockdown PHLDB2 expression, and three biological replicates were conducted for both the experimental group (PHLDB2-si) and the control group (control-si), followed by RNA-seq. Differential expression analysis (thresholds, P<0.05 and |log2FC|>1) confirmed reliable knockdown efficiency of PHLDB2-si, identifying 80 distinctively upregulated genes and 360 distinctively downregulated genes in the PHLDB2 knockdown group (Fig. 7B). Subsequent Hallmark and GO GSEA based on MsigDB were performed on these DEGs, collectively named the PHLDB2-RNAseq pathway enrichment results (Fig. 7C).

To delineate core pathway convergence in PHLDB2-deficient molecular contexts, a cross-dataset intersection of the PHLDB2-RNAseq and PHLDB2-TCGA pathway enrichment results was performed (Fig. 7C). Consensus pathway analysis of the intersecting 152 molecular signatures demonstrated significant suppression of the Hallmark EMT gene set in PHLDB2-depleted samples (NES, −2.833; false discovery rate, 0.023; Fig. 7D), consistent with the Hallmark enrichment results from prognostic model genes (Fig. 3C). Western blot validation of core EMT markers confirmed these findings (Fig. S4), showing decreased N-cadherin and vimentin protein levels accompanied by increased E-cadherin expression in the PHLDB2-si groups compared with those in the control-si groups. GO enrichment analysis further demonstrated that pathways related to ECM remodeling and the cell-substrate junction including ‘extracellular structure organization’, ‘extracellular matrix organization’, ‘cell-substrate junction assembly’ and ‘cell-substrate junction organization’, were markedly downregulated (Fig. 7E), mechanistically implicating that PHLDB2 depletion impaired EMT activation dynamics in OSCC.

To elucidate the molecular mechanisms underlying PHLDB2-mediated EMT regulation, KEGG GSEA was conducted for downregulated genes in PHLDB2-depleted samples, which revealed significant enrichment of the ‘PI3K-Akt signaling pathway’ (Fig. 7F), recapitulating the pathway convergence pattern observed in prior prognostic model gene signatures. Subsequent expression profiling of core PI3K-Akt signaling components in PHLDB2-knockdown models demonstrated marked downregulation of PIK3CA (Fig. 7G). Further validation through RT-qPCR in PHLDB2-knockdown OSCC cell lines confirmed the transcriptional repression of PIK3CA (Fig. 7H and I). Previous research has identified PIK3CA as a critical molecular regulator of EMT activation (24,25). Integrated analysis of the multidimensional datasets in the present study revealed that PHLDB2-mediated upregulation of PIK3CA may constitute a novel mechanistic axis driving EMT progression in OSCC pathogenesis.

To clarify the expression characteristics and clinical significance of PHLDB2 in OSCC clinical samples, multicenter data were systemically integrated. Analysis of TCGA database revealed that PHLDB2 mRNA expression was significantly higher in OSCC tumor tissues (n=302) compared with that in adjacent tissues (n=67; Fig. 8A), with consistent differences observed in paired tumor-adjacent samples (Fig. 8B). This trend was further validated in the independent cohort GSE213862 (n=43), where PHLDB2 expression in OSCC tissues was significantly elevated compared with that in normal mucosal tissues from healthy individuals (Fig. 8C). Notably, analysis of TCGA-OSCC cohort demonstrated that PHLDB2 expression was significantly higher in advanced-stage tumors (T3/T4) compared with that in early-stage tumors (T1/T2) (Fig. 8D), and patients with high PHLDB2 expression had significantly shorter OS compared with that in the patients with low PHLDB2 expression (Fig. 8E).

To further validate protein-level expression characteristics, IHC analysis was performed on an OSCC cohort (n=51) from the Hospital of Stomatology, SYSU. Negative control experiments using non-specific IgG demonstrated assay specificity (Fig. S5). Results indicated significantly higher PHLDB2 protein expression in tumor tissues compared with those in adjacent normal tissues (Fig. 8F-H). Additionally, high PHLDB2 expression was significantly associated with lymph node metastasis (Fig. 8I) and poor tumor differentiation (G2/3 vs. G1; Fig. 8J). These data collectively suggested that elevated PHLDB2 expression may be closely linked to malignant progression and poor prognosis in OSCC.