Spatiotemporal atlas of internal hemorrhoids in rats elucidated using integrated single-cell RNA sequencing and spatial transcriptomics

Abstract

BACKGROUND

Internal hemorrhoids are a highly prevalent vascular anorectal disease worldwide. However, the pathological mechanisms of disease progression, and particularly early-stage cellular heterogeneity, intercellular communication, and microenvironmental dynamics, remain poorly defined.

AIM

To investigate the spatiotemporal transcriptome profile of internal hemorrhoids development in rats.

METHODS

We established a rat model of internal hemorrhoids and performed integrated single-cell RNA sequencing, spatial transcriptomic, histopathology, and western blot analyses to elucidate the spatiotemporal mechanisms of hemorrhoid pathogenesis.

RESULTS

Macrophages switched from a proinflammatory phenotype to a reparative phenotype; fibroblasts mainly differentiated into growth factor-regulated subsets and profibrotic protomyofibroblasts; and vascular endothelial/smooth muscle cells participated in hemorrhoidal vascular remodeling. Spatial transcriptomic clarified the in situ distribution of key genes and spatiotemporal correlations between immune cell infiltration and extracellular matrix remodeling in early lesions. Multiomics revealed that Cd74⁺ macrophages act as early inflammatory switches, driving Igfbp5⁺ fibroblasts/protomyofibroblasts activation via Col-integrin interactions to trigger the phosphatidylinositol 3-kinase/protein kinase B signaling, synergistically promoting vascular dilation, inflammation, and fibrosis.

CONCLUSION

This study delineates a spatiotemporal transcriptomic landscape of internal hemorrhoid development in rats, revealing a pathogenic cascade from vascular inflammation to extracellular matrix remodeling.

Key Words: Internal hemorrhoids; Single-cell RNA sequencing; Spatial transcriptomics; Cellular heterogeneity; Fibrosis

Core Tip: This study integrates single-cell and spatial transcriptomics to reveal spatiotemporal dynamics in rat internal hemorrhoids. Macrophages shift from pro-inflammatory to pro-reparative phenotypes, while fibroblasts differentiate into growth factor-regulated and profibrotic subsets, driving vascular remodeling. Cd74⁺ macrophages act as early inflammatory switches, activating Igfbp5⁺ fibroblasts/protomyofibroblasts via collagen–integrin interactions and the phosphatidylinositol 3-kinase/protein kinase B signaling, thus promoting vascular dilation, inflammation, and fibrosis. These findings provide novel insights into hemorrhoid progression from inflammation to extracellular matrix remodeling.

INTRODUCTION

Hemorrhoidal disease is among the most prevalent anorectal disorders worldwide, with clinical manifestations and complications severely impairing patients’ quality of life and working capacity. It represents the third most common diagnosis for outpatient gastrointestinal disorders in the United States, with nearly 4 million associated outpatient visits documented annually[1]. In China, this condition has a prevalence of 49.14%, which ranks first among all anorectal disorders[2]. Its detection rate is correlated with body mass index, and females are more frequently affected, likely because of body mass index-related increases in abdominal pressure and sex-specific anatomical features[3]. Internal hemorrhoids are characterized by hematochezia and prolapse. Traditional Chinese medicine attributes hemorrhoid pathogenesis to “meridian relaxation and intestinal irregularity after excessive eating” (Su Wen Sheng Qi Tong Tian Lun)[4], whereas modern theories focus on varicose veins and anal cushion prolapse (abnormal sliding of the Treitz muscle)[5]. In recent decades, an increasing number of studies have explored the molecular pathogenesis of internal hemorrhoids beyond traditional anatomical theories, and current research has focused mainly on three core pathological links: Dysregulated inflammation, disrupted vascular remodeling, and an extracellular matrix (ECM) metabolic imbalance[6]. However, the existing theories lack synergistic explanatory power (e.g., anal cushion prolapse fails to fully explain early hematochezia)[7], and the fundamental pathogenic mechanisms remain incompletely defined. Most studies focus on isolated pathological links or molecular targets and lack a systematic analysis of spatiotemporal cellular/molecular features across disease-associated regions, hindering the development of precision therapies targeting core mechanisms.

Cutting-edge genomics technologies are transformative tools for dissecting complex diseases. Single-cell RNA sequencing (scRNA-seq) resolves cellular heterogeneity and state transitions at unprecedented resolution, while spatial transcriptomics (ST) preserves the tissue architecture to map gene expression in situ[8]. These complementary approaches enable the bridging of molecular changes with anatomical localization, addressing critical gaps in understanding the spatiotemporal progression of the disease. Fibroblasts and macrophages are central to fibrotic disorders. Santacroce and Di Sabatino[9] used scRNA-seq to analyze ileal samples from Crohn’s disease patients and identified NT5E⁺, FAP⁺, CCL11⁺, and FDFR2⁺ fibroblast subsets; moreover, transforming growth factor-β (TGF-β)-induced myofibroblasts can mediate inflammation-driven ECM remodeling[10], epithelial/endothelial cells can transdifferentiate into activated myofibroblasts[11], and chronic inflammation-induced macrophage dysfunction can promote fibrosis[12]. COL10A1⁺ fibroblasts, for example, drive the epithelial-mesenchymal transition and M2 macrophage polarization[13]. Notably, upon exogenous stimulation, hemorrhoidal tissues often exhibit mucosal damage, which may serve as the pathological basis for the recurrent nature of hemorrhoidal disease[14], whereas FBLN5 increases fibroblast activation to modulate hemorrhoid repair[15]. Abnormal collagen metabolism (e.g., a reduced ratio of type I/III collagen) impairs the mechanical stability of the anal cushion[16], and Xiaozhiling injection targets the SphK1-S1P pathway to promote collagen synthesis in fibroblasts[17], highlighting the relevance of fibroblast/macrophage dynamics in hemorrhoid pathogenesis. Gene set enrichment analysis of human hemorrhoid gene sets revealed key shared pathways essential for vascular and intestinal tissue development, including vascular morphogenesis and development, arterial morphogenesis and development, epithelial morphogenesis, and smooth muscle morphogenesis[18]. Here, we integrated scRNA-seq, ST, histopathology, and protein validation to construct a spatiotemporal atlas of hemorrhoidal disease progression in a rat model. Our goals were to identify core regulatory cell subsets and signaling pathways and to provide novel mechanistic insights and technical support for internal hemorrhoid research and clinical translation.

MATERIALS AND METHODS

Experimental animals and model establishment

Eight-week-old SPF-grade male Sprague-Dawley rats (weight: 180-220 g) were purchased from Shanghai Jiesijie Laboratory Animal Co., Ltd. [Shanghai, China; license number: SCXK (Shanghai) 2018-0004]. All the rats were housed in the Experimental Animal Center of Fujian University of Traditional Chinese Medicine, and the animal procedures strictly adhered to the Guiding Opinions on the Humane Treatment of Laboratory Animals[19] and were approved by the Animal Ethics Committee of Fujian Academy of Traditional Chinese Medicine, approval No. FJACMS-PZ-20250016. The rat internal hemorrhoid model established by our team[20] is a multifactor composite modeling method that combines anal dilation + chemical stimulation + fatigue stress and is stable and reproducible. A total of 15 rats were randomly divided into three groups (n = 5 rats per group): The control group (Ctrl), 1 week of modeling (Model-1w group), and 2 weeks of modeling (Model-2w group). The rats in the Ctrl were housed under standard temperature and humidity conditions with free access to food and water. For the model groups, anal dilation was performed first, after which a longitudinal incision with a length of approximately 0.5 cm was made at the 12:00 position while the rats were in the lithotomy position. Subsequently, 0.3 mL of a croton oil mixture (prepared by mixing distilled water, pyridine, ether, and 6% croton oil at a volume ratio of 1:4:5:10) was instilled into the perianal mucosa daily, and standing and swimming experiments were conducted to simulate fatigue. Two modeling cycles (1 week and 2 weeks) were established. The rats in the Model-1w group exhibited prolapse of maternal hemorrhoid regions at 3 o’clock, 7 o’clock, and 11 o’clock, whereas those in the Model-2w group exhibited near-circumferential prolapse. Histopathological observations via hematoxylin and eosin (HE) and Masson’s trichrome staining revealed pathological changes in the hemorrhoidal regions of the model rats, including vascular dilation, capillary proliferation, inflammatory cell infiltration, and submucosal collagen matrix remodeling, which are highly consistent with the pathological characteristics of human internal hemorrhoids.

Sample collection and processing

After modeling, hemorrhoidal tissue and normal rectal tissue were collected from the rats. The samples were divided into four parts and processed as follows: (1) Fixed with 4% paraformaldehyde for paraffin embedding and HE/Masson’s trichrome staining; (2) Embedded in optimal cutting temperature and snap-frozen in liquid nitrogen for the spatial transcriptomic analysis; (3) Fresh tissue was used to prepare single-cell suspensions for scRNA-seq; and (4) Fresh tissue was snap-frozen in liquid nitrogen for western blotting.

HE staining

Fresh hemorrhoidal tissues were fixed with 4% paraformaldehyde, dehydrated with a gradient of ethanol solutions, cleared in xylene, and embedded in paraffin. Continuous 4-μm sections were prepared and stained with HE using a commercial kit (Beyotime Biotechnology Co., Ltd., Shanghai, China) according to the manufacturer’s protocol. An Olympus BX45 optical microscope was used for observation and image collection, with a focus on the degree of inflammatory cell infiltration, vascular dilation, the looseness of the tissue structure, and fiber arrangement in the hemorrhoidal region.

Masson’s trichrome staining

Adjacent serial sections were subjected to Masson’s trichrome staining with a commercial kit (Solarbio Science & Technology Co., Ltd., Beijing, China) to visualize the collagen fibers. Staining was performed according to the manufacturer’s instructions, with collagen fibers stained blue and muscle fibers stained red. The area ratio of collagen fibers in the submucosa was semiquantitatively analyzed using Image-Pro Plus software to evaluate the degree of matrix remodeling.

scRNA-seq

Fresh hemorrhoidal tissues were rinsed with precooled phosphate-buffered saline, minced, and enzymatically digested at 37 °C with constant rotation. The dissociated cells were collected, filtered, centrifuged and resuspended to prepare single-cell suspensions. Cell viability was confirmed to be ≥ 80% by trypan blue staining, and the cell concentration was adjusted to 700-1200 cells/μL. Approximately 1 × 10⁴ cells were loaded to generate gel beads in emulsions, followed by reverse transcription, cDNA purification, and library construction. The qualified libraries were subsequently sequenced on the NovaSeq Xplus/DNBSEQ-T7 platform (PE150 mode). After quality control was performed to select qualified cells from the single-cell sequencing data, standardization, dimensionality reduction, clustering, and cell type annotation were performed with reference to the CellMarker database (http://xteam.xbio.top/CellMarker/) and the published literature. A differential expression analysis and gene set Gene Ontology (GO) and Kyoto Encyclopedia of Genes and Genomes (KEGG) enrichment analyses were performed to compare the transcriptomic characteristics and pathways enriched in differentially expressed genes from different cell subsets between groups.

Spatial transcriptomic sequencing

Fresh tissues were embedded in optimal cutting temperature compound, snap-frozen in liquid nitrogen, and sectioned for the ST analysis. ST libraries were constructed using the 10 × Genomics Visium platform according to a standard protocol and sequenced on the Illumina HiSeq platform. The raw data were processed with FastQC and Space Ranger software for quality control, image alignment, spot gene quantification, and duplicate removal. Spots with ≥ 200 detected genes were retained for the downstream ST analysis, which was performed using the Seurat R package. The expression matrix was normalized via SCTransform, followed by principal component analysis dimensionality reduction, unsupervised clustering, and visualization with uniform manifold approximation and projection. Differentially expressed genes for each cluster were identified with the Find All Markers function, and cell types were annotated with reference to the CellMarker database. GO/KEGG functional enrichment analyses, single-cell subset spatial mapping (SpotLight algorithm), and spatial differential expression analysis (SPARK software) were performed as described in the Statistical Analysis section.

Western blotting

Fresh tissue samples were lysed on ice with radioimmunoprecipitation assay lysis buffer containing protease inhibitors. Total protein was extracted, and the concentration was determined using a bicinchoninic acid protein assay kit. Equal amounts of protein were separated by sodium dodecyl sulfate-polyacrylamide gel electrophoresis and transferred to polyvinylidene difluoride membranes. The membranes were blocked with 5% skim milk for 1.5 hours at room temperature, incubated with primary antibodies at 4 °C overnight, and then incubated with the corresponding secondary antibodies for 1.5 hours at room temperature. The protein bands were visualized with a chemiluminescence system, and the gray values were quantified using ImageJ software, with Tubulin β serving as the internal reference.

Mapping rat internal hemorrhoid datasets to human datasets

We first retrieved immunohistochemical staining profiles of healthy human anal canal tissues from the human protein atlas database to verify the clinical relevance of the key genes identified in our rat model. We further obtained transcriptomic datasets of fibrotic diseases from the Fibrotic Disease-associated RNAome database[21] to assess the generalizability of our identified genes across fibrosis-related pathologies. Moreover, the public human hemorrhoidal disease transcriptome dataset (GSE154650)[18] was downloaded from the National Center for Biotechnology Information Gene Expression Omnibus database, which comprises mRNA sequencing data for hemorrhoid tissues from 20 patients and normal anal cushion tissues from 18 patients with anal fissure as controls. Orthologous gene pairs between humans and rats were mapped using the Ensembl Genes database. The differential expression analysis was conducted using the limma R/Bioconductor package, and weighted gene co-expression network analysis was performed to identify disease-related co-expressed gene modules. The functional enrichment analysis revealed that genes within these modules were enriched mainly in four core biological processes: ECM organization, smooth muscle contraction, epidermal development, and mitochondrial organization. Among them, the M1 module was defined as the core pathogenic module, with significant enrichment in ECM organization and smooth muscle contraction functions, which directly recapitulated the canonical pathological hallmarks of hemorrhoidal disease: Anal cushion descent and connective tissue degeneration. We further delineated the functional roles of multiple high-confidence pathogenic genes, including ANO1 (encoding a calcium-activated chloride channel in intestinal interstitial cells of Cajal), SRPX (encoding an ECM-associated protein), ELN (encoding elastin, the core structural component of anal cushion connective tissue), and MYH11 (encoding smooth muscle myosin heavy chain that mediates smooth muscle contractility). A cross-species conservation analysis between our rat model and human hemorrhoidal specimens further validated these findings.

Statistical analysis

All the statistical analyses were performed using R package, which includes the Seurat, CellChat, SPARK, and ClusterProfiler packages. All animal experimental data were expressed as mean ± SD. For hemorrhoid samples from model rats, low-quality cells were first removed, after which the batch effects were corrected. After principal component analysis dimensionality reduction, Louvain clustering, and uniform manifold approximation and projection visualization of single-cell data, cell subsets were annotated by screening differentially expressed genes (Wilcoxon rank-sum test, log2|FC| > 1, adjusted P value < 0.05). Differences in the proportions of cellular subsets among the three groups were analyzed using the Kruskal-Wallis H test, and P values were corrected by the Benjamini-Hochberg method, with adjusted P value < 0.05 considered to indicate statistical significance. For spatial transcriptomic data, single-cell subset spatial mapping was performed using the SpotLight algorithm, and a spatial differential expression analysis was conducted with SPARK software. Differential gene expression between groups was analyzed by performing the DESeq2 negative binomial distribution test. GO and KEGG functional enrichment analyses were performed using ClusterProfiler software, with a false discovery rate < 0.05 set as the significance threshold. For the cell-cell communication analysis, ligand-receptor interactions were analyzed using CellPhoneDB, with significance determined using a 1000-time permutation test (P < 0.05). Western blots were analyzed using ImageJ software. Differences among groups were analyzed using one-way analysis of variance with Tukey’s post hoc test, with P < 0.05 considered to indicate statistical significance.

Safety statement

No rats died during the experiment.

RESULTS

Single-cell transcriptomic landscape of internal hemorrhoidal tissue

After strict quality control, a total of 37015 high-quality cells were obtained. The clustering analysis identified 33 cell clusters (Figure 1A), which were annotated into 9 cell types (Figure 1B): Fibroblasts (17494 cells, 47.26%), epithelial cells (6467 cells, 17.47%), macrophages (5098 cells, 13.77%), smooth muscle cells (3155 cells, 8.52%), endothelial cells (2311 cells, 6.24%), T cells (1390 cells, 3.76%), dendritic cells (437 cells, 1.18%), plasma cells (371 cells, 1.00%), and basophils (292 cells, 0.79%). An analysis of the dynamic changes in cell proportions revealed significant compositional changes in the major cell types in hemorrhoidal tissues during disease progression (Figure 1C). The proportion of fibroblasts, the most abundant cell type, increased significantly from 41.24% in the Ctrl to 51.07% in the Model-1w group (1.24-fold, adjusted P value < 0.05) and 48.84% in the Model-2w group (1.18-fold, adjusted P value < 0.05). The proportion of macrophages also increased significantly from 12.09% in the Ctrl to 15.44% in the Model-1w group (1.28-fold, adjusted P value < 0.05) and 13.37% in the Model-2w group (1.11-fold, adjusted P value < 0.05). Consistently, the proportion of endothelial cells was significantly increased in both model groups (all adjusted P value < 0.05), whereas the proportions of epithelial cells and vascular smooth muscle cells were significantly decreased (all adjusted P value < 0.05). These results confirmed the initiation of inflammatory activation and matrix remodeling in the early stage of internal hemorrhoid development. Functional clustering of the differentially expressed genes revealed that the highly expressed genes were enriched mainly in functional modules such as collagen synthesis (Col1a1, Col1a2, and Col3a1), immune inflammation (Ccl11 and Tnf), matrix remodeling (Fap, Tgfbi, and Mmp), and smooth muscle regulation (Acta2 and Mgp) (Figure 1D and E).

Figure 1 Single-cell transcriptomic landscape of internal hemorrhoidal. A and B: The X-axis and Y-axis denote dimensionality-reduced components; different cell clusters are distinguished by distinct colors; C: Bar plot showing the relative proportions of 9 major cell types in the control, 1 week of modeling and 2 weeks of modeling groups (n = 5 rats per group). Differences in cell proportions among groups were analyzed using the Kruskal-Wallis H test with the Benjamini-Hochberg correction. Adjust P < 0.05 vs the control group; D: The X-axis indicates normalized gene expression; each dot indicates differential gene expression in each cell; and the Y-axis indicates clusters; E: The X-axis indicates colored clusters; the Y-axis indicates the number of differentially expressed genes per cluster.

Identification and functional analysis of core cell subsets

Further clustering analysis identified heterogeneous subpopulations in fibroblasts, macrophages, vascular smooth muscle cells, and endothelial cells of hemorrhoidal tissues (Figure 2A-L).

Figure 2 Visualization of the subpopulations of internal hemorrhoid cells. A, D, G, and J: The X-axis and Y-axis denote dimensionality-reduced components; different cell clusters are distinguished by distinct colors; B, E, H, and K: Different colors represent distinct cell type identities; C, F, I, and L: The X-axis indicates cell types; the Y-axis indicates normalized gene expression values; and each dot indicates the average expression of the differentially expressed genes in each cell type.

Fibroblasts exhibit high heterogeneity and various differentiation trajectories: Fibroblasts were further clustered into 7 functional subsets: (1) Ccl11-fibroblasts (inflammatory chemotactic subtype) were characterized by the upregulation of Ccl11, Nfkbia, Acta2 (adjusted P value < 0.05), S100a6, Cd9, Gsn, and Mt-atp6 1 week after modeling, synergistically promoting the inflammation-fibrosis transition through the amplification of inflammation, immune cell recruitment, and the energy supply; (2) Igfbp5-fibroblasts (growth factor-regulated subtype) were characterized by the expression of the marker gene Igfbp5, which was upregulated 1 week after modeling, accompanied by the upregulation of S100a4, Pla2g2a, Timp1, and Serpina3a (adjusted P value < 0.05; Figure 3A). This subtype mediates the transition from inflammatory cell infiltration to reparative fibrosis by regulating the balance between ECM synthesis and degradation via growth factors; (3) Protomyofibroblasts (myofibroblast precursors) expressed Pdpn, Aspn, and Fap at high levels and were activated early (1 week after modeling), accompanied by the upregulation of S100a4, Pla2g2a, and Ifi27l2a (adjusted P value < 0.05; Figure 3B). These cells have the potential to differentiate into myofibroblasts; (4) Fascia-fibroblasts (fascia-like structural support subtype) were characterized by high expression of Col3a1, which increased 2 weeks after modeling, accompanied by the upregulation of Nfkbia, Cxcl12, and Jun (adjusted P value < 0.05; Figure 3C); these cells participated in late tissue mechanical support and fibrotic remodeling; (5) Sdk1-fibroblasts (cell adhesion signaling subtype) were characterized by high expression of Sdk1, and are presumably involved in cell adhesion and signal transduction; (6) Cadm2-fibroblasts (spatial localization subtype) were characterized by high expression of Cadm2; these cells are presumably involved in regulating the spatial localization of fibroblasts; and (7) Proliferate-fibroblasts (proliferatively active subtype) expressed high levels of proliferation-related genes, representing a rapidly proliferating fibroblast population under pathological conditions.

Figure 3 Volcano plots of differentially expressed genes. The X-axes indicate the log2 fold change (log2FC) values for gene expression; the Y-axes indicate the adjusted P values (adjust P values). Adjust P < 0.05 was considered to indicate statistical significance. Each dot represents a gene: Red indicates significantly upregulated genes, blue indicates significantly downregulated genes, and gray indicates genes whose expression did not change significantly. Dots farther from the center or with larger values along the Y-axis indicate genes with more significant differences in expression. A: Igfbp5-fibroblast [1 week of modeling (Model-1w) vs control group (Ctrl)]; B: Protomyofibroblast (Model-1w vs Ctrl); C: Fascia-fibroblast [2 weeks of modeling (Model-2w) vs Model-1w]; D: Cd74-expressing macrophage (Model-1w vs Ctrl); E: M2 macrophage (Model-2w vs Model-1w); F: Vascular endothelial cell (Model-1w vs Ctrl); G: Fibroblast-smooth muscle cell (SMC) (Model-1w vs Ctrl); H: Contractile-SMC (Model-1w vs Ctrl); I: Synthetic-SMC (Model-1w vs Ctrl). Ctrl: Control group; Model-1w: 1 week of modeling; Model-2w: 2 weeks of modeling; VEC: Vascular endothelial cell; SMC: Smooth muscle cell.

Macrophage polarization and regulation of inflammation: The macrophages were clustered into 3 subsets: (1) Cd74-expressing macrophages (antigen-presenting proinflammatory subtype) were associated with major histocompatibility complex class II; expressed high levels of Cd74, Cd86, and Tnf; and were activated 1 week after modeling, resulting in the upregulation of S100a4, S100a6, Sod2, and Mt-atp6 expression (adjusted P value < 0.05; Figure 3D). Cd74 is a type II transmembrane glycoprotein, and previous studies have confirmed that it is a marker of M1 macrophage infiltration[22]; thus, Cd74-expressing macrophages serve as a key cell subtype that initiates early inflammation; (2) Monocyte-derived macrophages (monocyte-derived chemotactic subtype) expressed high levels of Clec7a, Tnf, Il6, and other genes, and their recruitment peaked 1 week after modeling to supplement the local inflammatory cell pool; and (3) M2 macrophages (anti-inflammatory pro-reparative subtype) expressed high levels of Cd163 and Mrc1; their proportion was restored 2 weeks after modeling, and they exhibited upregulated Eef1a1l1, Ccl6, and Alox5ap expression (adjusted P value < 0.05; Figure 3E), indicating the transition from the inflammatory phase to the repair phase. Peak macrophage recruitment occurred in the first week of modeling during the acute inflammatory phase; the macrophages recruited at this stage were mainly M1 proinflammatory cells. At this time, Cd74-expressing macrophages were activated, significantly upregulating the expression of inflammation-related genes and rapidly increasing the macrophage pool in response to injury. After 2 weeks of modeling, the demand for proinflammatory macrophages for hemorrhoidal tissue repair decreased, leading to a reduction in the recruitment of monocyte-derived macrophages, whereas the demand for M2 macrophages in the hemorrhoidal region increased to promote tissue repair.

Imbalance of endothelial cell subsets and vascular remodeling: Endothelial cells were clustered into 5 subsets. Modeling induced changes in the proportions of different subsets, with an enrichment of circulating endothelial cells and vascular endothelial cells (VECs), both of which expressed high levels of core endothelial markers such as Fabp4 and Plvap. In contrast, lymphatic endothelial cells and fibroblast-like endothelial cells were depleted, and a small number of alveolar-derived endothelial-like cells were present. Notably, VECs exhibited the specific upregulation of RT1-T24-4, Serpina3a, and Mt-atp6 expression 1 week after modeling (adjusted P value < 0.05; Figure 3F), suggesting that this subset plays a core role in early vascular remodeling and energy changes in individuals with hemorrhoidal disease. Combined with the functional annotations, the upregulation of Mt-atp6 reflects compensatory energy metabolism in VECs; Serpina3a may be involved in the inhibition of local inflammation and maintenance of vascular stability, and as a major histocompatibility complex class I molecule, the nonclassical function of RT1-T24-4 has been confirmed to mediate vascular endothelial cell migration and angiogenesis[23]. Its specific upregulation in VECs was highly consistent with the pathological characteristics of hemorrhoidal vascular remodeling, suggesting that RT1-T24-4 may be a key molecule that regulates the chemotaxis and recruitment of vascular-related cells to lesion sites. These results indicated that hemorrhoid modeling could induce an imbalance in vascular endothelial cell subsets, with the number of VECs peaking 1 week after modeling.

Phenotypic transition of vascular smooth muscle cells: Vascular smooth muscle cells were clustered into 5 subsets: (1) Fibroblast-smooth muscle cells (SMCs) (fibroblast-like transformed subtype) expressed Dcn, Col1a1, Col3a1, and Mmp2 at high levels, and presented sustained upregulation of Timp1, Prss23, and Col6a3 expression 1 week after modeling (adjusted P value < 0.05; Figure 3G); (2) Contractile-SMCs (contractile functional subtype) expressed Cnn1, Acta2, Tagln, and Myh11 at high levels, while the expression of contraction-related genes (Btg2, Fos, and Tcap) was consistently downregulated after modeling (adjusted P value < 0.05; Figure 3H), indicating impaired vascular contractile function and persistent abnormal dilation of hemorrhoidal blood vessels; (3) Synthetic-SMCs (synthetic proliferative subtype) expressed Mgp, Col1a1, and Col3a1 at high levels, presented sustained upregulation of synthesis- and proliferation-related genes (Mt-atp6, S100a6, and Kcnab1) after modeling (adjusted P value < 0.05; Figure 3I), and participated in vascular wall remodeling; (4) Ptn-SMCs (signal regulatory subtype) expressed Ptn and Fn1 at high levels and presumably mediated signal transduction; and (5) Cacna1c-SMCs (calcium channel regulatory subtype) expressed high levels of Cd44 and presumably regulated calcium channel function. During modeling, the proportion of contractile-SMC subsets with contractile function decreased, whereas the proportions of synthetic-SMCs with synthetic functions and fibroblast-SMCs that underwent a fibroblast-like transformation increased.

Fibroblasts constituted an important foundation of the hemorrhoidal tissue microenvironment, among which Igfbp5-fibroblasts and protomyofibroblasts were synergistically activated in the early stage and served as core fibroblast subsets that drive the fibrotic process. GO and KEGG enrichment analyses revealed that the core fibroblast subsets (Igfbp5-fibroblasts, protomyofibroblasts, and Fascia-fibroblasts) were significantly enriched in functions such as ECM remodeling, the phosphatidylinositol 3-kinase (PI3K)-protein kinase B (Akt) signaling pathway, actin cytoskeleton regulation, angiogenesis, and cell adhesion (all adjusted P value < 0.05) (Figure 4A-C). The KEGG enrichment analyses suggested that Cd74-expressing macrophages regulate actin polymerization and that M2 macrophages are significantly enriched in oxidative phosphorylation, antigen processing and presentation, and phagocytosis (all adjusted P value < 0.05) (Figure 4D). The KEGG enrichment analysis suggested that the contractile and fibroblast-like-transformed subtypes were significantly enriched in vascular smooth muscle contraction, ECM-receptor interactions, cell adhesion, actin regulation, and the PI3K-Akt signaling pathway (all adjusted P value < 0.05) (Figure 4E).

Figure 4 Plots of the results of the Gene Ontology and Kyoto Encyclopedia of Genes and Genomes enrichment analyses. The Y-axes indicate the names of the Gene Ontology (GO) terms or Kyoto Encyclopedia of Genes and Genomes (KEGG) pathways; the X-axes indicate the Rich factors (higher values indicate stronger enrichment). The dot size represents the number of genes; the dot color corresponds to the range of adjust P values. Adjust P < 0.05 was considered to indicate statistical significance. A: GO enrichment analysis of upregulated genes in Igfbp5-fibroblast [1 week of modeling (Model-1w) vs control group (Ctrl)]; B: GO enrichment analysis of upregulated genes in Protomyofibroblast (Model-1w vs Ctrl); C: KEGG pathway enrichment analysis of upregulated genes in three fibroblast subsets (Igfbp5-fibroblast, protomyofibroblast, and fascia-fibroblast) (Model-1w vs Ctrl); D: KEGG pathway enrichment analysis of upregulated genes in three macrophage subsets (M2 macrophage, Cd74-expressing macrophage, and monocyte-derived macrophage) (Model-1w vs Ctrl); E: KEGG pathway enrichment analysis of upregulated genes in three SMC subsets (fibroblast-smooth muscle cell (SMC), contractile-SMC, and synthetic-SMC) (Model-1w vs Ctrl). GO: Gene Ontology; KEGG: Kyoto Encyclopedia of Genes and Genomes.

The pseudotime analysis revealed the fibroblast differentiation trajectory: Proliferate-fibroblasts served as initial precursors, which underwent functional divergence at branch points 1/2 to differentiate into protomyofibroblasts and Igfbp5-fibroblasts, respectively, and ultimately codifferentiate into Fascia-fibroblasts (Figure 5A). Moreover, Ccl11-fibroblasts and fibrotic subsets originated from common precursors, forming parallel functional branches of inflammatory chemotaxis and fibrosis, and synergistically participated in the development of an inflammatory microenvironment and fibrotic progression (Figure 5B). The pseudotime analysis revealed the macrophage differentiation trajectory: Monocyte-derived macrophages served as initial precursors, transitioned through the Cd74-expressing intermediate macrophage subtype, and ultimately polarized into M2 macrophages, forming a dynamic precursor recruitment-proinflammatory activation-anti-inflammatory repair polarization process. These findings confirm the core role of macrophage polarization in the crosstalk between the immune microenvironment and fibrosis (Figure 5C). Vascular smooth muscle cells in the hemorrhoidal region underwent a pathological shift from a “contractile phenotype” to a “synthetic, fibrotic phenotype”, with differences in the temporal patterns of different phenotypic transitions (Figure 5D). The fibroblast-like transformation was more prominent in the early stage, whereas the synthetic proliferative phenotype was more obvious in the late stage. The reduction in the proportions of contractile functional subsets may exacerbate vascular dysfunction, which was likely an important cause of the loss of hemorrhoidal vascular tone and persistent abnormal dilation. Differences in the temporal responses of different subsets were also observed (Figure 5E). Endothelial cell subsets originated from common precursors; after modeling, the differentiation direction shifted from the preferential differentiation to lymphatic endothelial cells under physiological conditions toward the VEC and circulating endothelial cell branches, resulting in a differentiation-oriented vascular–lymphatic imbalance. These findings provide a cellular-level explanation for vascular dilation and lymphatic drainage disorders in hemorrhoidal tissues.

Figure 5 Pseudotime plots. Each dot represents a cell arranged along a pseudotime trajectory. Color gradients indicate cell state transitions; numbers in black circles indicate differentiation nodes. A: Pseudotime trajectory analysis of four fibroblast subsets (fascia-fibroblast, Igfbp5-fibroblast, protomyofibroblast, and proliferate-fibroblast); B: Pseudotime trajectory analysis of three fibroblast subsets (Ccl11-fibroblast, Igfbp5-fibroblast, and protomyofibroblast); C: Pseudotime trajectory analysis of three macrophage subsets (Cd74-expressing macrophage, M2 macrophage, and monocyte-derived macrophage); D: Pseudotime trajectory analysis of three smooth muscle cell (SMC) subsets: Fibroblast-SMC, contractile-SMC, and synthetic-SMC; E: Pseudotime trajectory analysis of three endothelial cell subsets: Circulating endothelial cell, lymphatic endothelial cell, and vascular endothelial cell. SMC: Smooth muscle cell; CEC: Circulating endothelial cell; LEC: Lymphatic endothelial cell; VEC: Vascular endothelial cell.

ST reveals the localization of key pathological events

A spatial transcriptomic analysis was performed to map the in situ distribution of key pathological events during internal hemorrhoid progression (Figure 6). HE staining first confirmed the progressive histopathological changes in hemorrhoidal lesions, including submucosal vascular dilation, inflammatory cell infiltration and tissue structure destruction, in the Model-1w and Model-2w groups (Figure 6A). Unsupervised clustering and dimensionality reduction analysis revealed differences in the global transcriptomic profiles between the control and model groups (Figure 6B and F). The spatial deconvolution analysis further revealed the spatial distribution of major cell types and functional subsets. Activated vascular endothelial cells, macrophages and fibroblasts formed a unique “pathological petal” structure in the mucosal and submucosal layers of Model-2w hemorrhoids, corresponding to angiogenesis, inflammatory cell infiltration and matrix remodeling regions, respectively, which may drive the disordered remodeling of hemorrhoidal tissues (Figure 6C). Under physiological conditions, the fibroblast subsets showed a regular layered distribution; in Model-1w, Igfbp5⁺ fibroblasts and protomyofibroblasts expanded spatially in the lesion core, whereas extensive colocalization of proliferative fibroblasts, Igfbp5⁺ fibroblasts and fascial fibroblasts was observed in Model-2w, confirming that spatiotemporal remodeling of fibroblast subsets was the core spatial basis of hemorrhoidal fibrosis (Figure 6D). The bubble plot further verified the dynamic changes in core functional gene modules across disease stages, with inflammatory genes enriched in Model-1w and ECM remodeling genes enriched in Model-2w (Figure 6E). In the Ctrl, type I/III collagen-encoding genes (Col1a1 and Col3a1) were restricted to the local submucosa, and smooth muscle contractile genes (Acta2, Myh11, Tagln, and Cnn1) were specifically enriched in the vascular wall. At 1 week after modeling, inflammatory and neovascularization signals were highly concentrated in the core vascular plexus region; collagen gene expression showed diffuse expansion along with vascular dilation, and the regions with high expression of the fibroblast activation marker S100a4 and the ECM-degrading enzyme Mmp2 were spatially colocalized with Cd74⁺ macrophages inflammatory regions, as determined by scRNA-seq, indicating that the inflammatory microenvironment initiated early fibroblast activation and ECM remodeling. At 2 weeks after modeling, Col1a1 and Col3a1 were widely expressed throughout the hemorrhoidal interstitium, accompanied by a further disruption of the distribution of smooth muscle-related genes; S100a4 expression was extensively overlapped with that of collagen genes, confirming that activated Igfbp5⁺ fibroblasts and protomyofibroblasts were the core effector cells responsible for excessive ECM deposition (Figure 6G and H).

Figure 6 Spatial transcriptome profiling of internal hemorrhoids. A: Images of hematoxylin and eosin staining; B: Dot plots of the spatial transcriptomic results. The clustering results from Seurat are overlaid on tissue sections; each dot indicates one spot, and different colors indicate each cluster; C: Maps of cell type annotations; D: Maps of the spatial deconvolution analysis of fibroblast subsets; E: Bubble plot of the top 3 marker genes per group; F: Maps showing the results of the dimensionality reduction clustering analysis; G: Maps showing the spatial localization of highly expressed genes; darker red dots indicate genes with higher expression in the spots; H: Violin plots. Ctrl: Control group; Model-1w: 1 week of modeling; Model-2w: 2 weeks of modeling.

Multiomics integration reveals intercellular communication networks

The CellPhoneDB analysis revealed the interaction strength of different cells in hemorrhoid tissue (Figure 7A). Further analysis revealed significant heterogeneity in the degree of cross-talk among the fibroblast subsets: The interaction between the fascia-fibroblasts and protomyofibroblasts was the most active (total interaction number = 113, P < 0.001), whereas the Ccl11⁺ fibroblasts exhibited minimal cross-talk with the other subsets (Figure 7B), indicating distinct signal transmission preferences among the different fibroblast subsets. The significant activation of multiple profibrotic signaling pathways (TGF-β, vascular endothelial growth factor, and platelet-derived growth factor receptor) in hemorrhoidal lesions after modeling (P = 0.0009). Core profibrotic ligand–receptor pairs, including Tgfb1-Tgfbr3, Pdgfd-Pdgfr and Col family-integrin complexes, were mainly enriched between endothelial cells and fibroblasts, as well as between macrophages and fibroblasts (all P < 0.001; Figure 7C). Mechanistically, modeling-induced venous congestion activated and injured endothelial cells, which released chemokines to recruit macrophages and drive vascular remodeling; moreover, activated endothelial cells and macrophages secreted TGF-β and other profibrotic factors to promote fibroblast-to-myofibroblast transdifferentiation.

Figure 7 Intercellular communication networks and protein expression. A and B: The X-axis and Y-axis indicate cell types; the color gradient from dark blue to dark red indicates an increasing number of interactions; C and D: The X-axis indicates cell type interactions, and the Y-axis indicates ligand-receptor pairs. Larger dots indicate smaller P values (greater significance); the color represents the average expression level (indicating interaction strength); E: Western blots showing protein expression. Ctrl: Control group; Model-1w: 1 week of modeling; Model-2w: 2 weeks of modeling.

The regression analysis indicated that core ligand-receptor axes driving fibroblast activation, namely, Fn1-integrin receptor complexes (Fn1-integrin aVb5 complex and Fn1-integrin aVb1_complex), were the most significantly and highly expressed across all subsets (P = 0.0009), with Col15a1-integrin a11b1 and Tgfβ3-TGFβ receptor complexes also serving as key regulatory nodes (all P < 0.001; Figure 7D). These results indicate that fibroblast subsets form a synergistic functional network via the ECM-integrin pathway, which probably mediates core pathological processes, including ECM remodeling and fibroblast-to-myofibroblast differentiation, in hemorrhoidal tissues.

Protein-level and histopathological validation

Key differentially expressed genes identified by transcriptomics were selected for the detection of protein expression levels via western blotting. After modeling, the protein levels of key differentially expressed genes (ACTA2, COL1A1, S100A4, CD74, and CD163) were significantly upregulated (P < 0.05; Figure 7E), indicating the activation of fibroblasts and the early inflammatory response and the initiation of fibrosis and compensatory smooth muscle cell proliferation. In the late stage, inflammation was alleviated, but fibrosis persisted, with aggravated smooth muscle injury and the specific downregulation of the CNN1 protein (P < 0.05).

HE staining was performed to observe the tissue morphology and cell types, and the degree of inflammatory cell infiltration and pathological characteristics across different regions were compared. Masson’s trichrome staining was performed to evaluate the collagen fiber distribution, quantify the degree of matrix remodeling, and verify the changes in the expression of matrix-related genes using spatiotemporal omics. The results (Figure 8) revealed that compared with rectal tissues from the Ctrl, the rectal tissues from the model groups exhibited progressive pathological changes. The Model-1w group showed an irregular cell arrangement, with vascular dilation and moderate inflammatory cell infiltration; the Model-2w group exhibited further vascular distortion, exacerbation of the abnormal tissue morphology, and a more irregular cellular arrangement. Masson’s trichrome staining confirmed the progressive aggravation of collagen fiber proliferation and irregular cellular arrangement, consistent with the trends in the expression of ECM-related genes from the spatial transcriptomic data and reflects the dynamic evolution of hemorrhoidal disease from the activation of inflammation to matrix fibrosis.

Figure 8 Images of histopathological staining (hematoxylin and eosin and Masson’s trichrome staining). Ctrl: Control group; Model-1w: 1 week of modeling; Model-2w: 2 weeks of modeling; HE: Hematoxylin and eosin.

DISCUSSION

In this study, we integrated scRNA-seq and ST to construct a spatiotemporal dynamic atlas of internal hemorrhoids in rats and systematically deduced the core regulatory mechanisms. The development of internal hemorrhoids is not an isolated process involving a single cell type but rather a dynamic pathological process involving the synergistic actions of multiple cell types and ordered progression in spatiotemporal dimensions. This feature can be gradually clarified through cellular and molecular changes at different modeling stages. After 1 week of modeling, the lesion is dominated by initiated inflammation and tissue repair: The increases in the numbers of macrophages and fibroblasts indicate the activation of acute inflammatory responses and the initiation of repair programs, whereas the decrease in the number of epithelial cells may involve mucosal repair or the epithelial-mesenchymal transition to promote vascular endothelial homeostasis. As the lesion progresses to the 2-week modeling stage, the intensity of the inflammatory response weakens, but the immune response persists, with increased numbers of T cells and plasma cells, indicating the activation of adaptive immunity and an attempt to regulate the inflammatory process. Another key change at this stage is the increased activity of vascular endothelial cells, with endothelial factors driving angiogenesis. However, the number and contractile function of vascular smooth muscle cells tend to decrease throughout the development of hemorrhoidal disease. Combined with the pathological feature of microvascular dilation, these changes are presumably important causes of early vascular dilation in hemorrhoidal tissues.

A comparative analysis of the control and model groups revealed that immune cell infiltration, tissue remodeling, and vascular abnormalities are core biological processes involved in the progression of hemorrhoids. Compared with traditional anal cushion prolapse and varicose vein theories, this study, by analyzing gene transcription data, revealed “macrophage polarization, endothelial cell angiogenesis and smooth muscle cell contractile dysfunction, and fibroblast activation” regulatory axis, which may be one of the core mechanisms of internal hemorrhoid development. Further analysis of the upstream link that initiates the regulatory axis revealed that Cd74⁺ macrophages are key switches for early inflammation and that their activation is the initiating event of the inflammatory cascade in hemorrhoidal disease. CD74 signaling is strongly activated during intestinal inflammation to protect the host by promoting epithelial regeneration and healing and maintaining mucosal barrier integrity[24]. This study confirmed that on the one hand, Cd74⁺ macrophages recruit immune cells and amplify inflammatory signals through the secretion of inflammatory factors such as tumor necrosis factor and interleukin 6; on the other hand, they directly initiate endothelial cell proliferation and migration through the secretion of vascular endothelial growth factor and fibrotic factors such as TGF-β to drive fibroblast activation, thereby forming an upstream “inflammation initiation-fibrosis activation” core regulatory chain. Further tracking of the functional dynamics of macrophages revealed their polarization from the M1 to the M2 phenotype: M1 macrophages are proinflammatory and are responsible mainly for the phagocytosis of necrotic tissue and the release of inflammatory factors to amplify the inflammatory response, which is common in the acute stage of hemorrhoids (e.g., hemorrhoid swelling and pain); M2 macrophages are anti-inflammatory and proreparative and inhibit excessive inflammatory responses and promote tissue repair, vascular regeneration, and fibrotic processes, corresponding to the remission stage of hemorrhoids (e.g., reduced swelling and repair). This dynamic polarization of macrophages is highly synchronized with the temporal progression of hemorrhoidal disease. In the process of fibrosis, Igfbp5⁺ fibroblasts and protomyofibroblasts are core effector subsets. Due to the regulatory effect of the growth factor Igfbp5 and the myofibroblast differentiation potential of protomyofibroblasts, these subsets mediate the synthesis and remodeling of the ECM. Their differentiation trajectory is highly consistent with the temporal progression of hemorrhoidal disease, suggesting that these two cell subsets are key targets for fibrosis intervention. With respect to vascular abnormalities, RT1-T24-4⁺ VECs mediate vascular endothelial migration and neovascularization and cooperate with the phenotypic transition of vascular smooth muscle cells to jointly induce hemorrhoidal vascular dilation and structural remodeling, providing molecular and cellular evidence for the sinusoidal blood vessel theory of hemorrhoids.

The results of single-cell sequencing and spatial analyses revealed that the core regulatory network of hemorrhoid development is as follows: Cd74⁺ macrophages, Igfbp5⁺ fibroblasts, protomyofibroblasts, and RT1-T24-4⁺ endothelial cells are core cell subsets that drive lesion development, and the PI3K-Akt signaling pathway is a key molecular hub connecting the major “inflammation-vascular abnormality-fibrosis” pathological processes. The spatial analysis confirmed that activated macrophages, fibroblasts, and endothelial cells are closely spatially colocalized in the hemorrhoid region. This spatial distribution provides a structural basis for functional synergy between cells, collectively promoting hemorrhoid progression. As a classic inflammation-fibrosis signaling axis, genes related to the PI3K-Akt pathway are abnormally expressed in hemorrhoidal tissues and participate in imbalanced angiogenesis and fibrosis[25]. Lipopolysaccharide promotes collagen synthesis in lung fibroblasts through aerobic glycolysis via the activation of the PI3K-Akt-mechanistic target of rapamycin/6-phosphofructo-2-kinase/fructose-2,6-bisphosphatase 3 pathway[26]. Inflammatory monocyte-derived amphiregulin drives intestinal fibrosis in experimental colitis and Crohn’s disease models by activating intestinal fibroblasts and promoting their proliferation via the PI3K-AKT pathway[27]. Unlike fibrosis in other solid organs, which is characterized by chronic inflammation as the initial factor and organ diffuse lesions, the fibrosis of internal hemorrhoids involves the abnormal dilatation and remodeling of anal cushion vessels as the initial core events, and the lesions are strictly confined to the submucosal vascular plexus of the anal cushion, forming a unique pathological petal-like spatial structure. Targeting early inflammation initiated by Cd74⁺ macrophages and inhibiting the PI3K-Akt signaling pathway can simultaneously target vascular abnormalities and fibrosis, disrupting key ligand-receptor pairs in “pathological petal” intercellular communication, which may provide new ideas for the development of more precise internal hemorrhoid treatment regimens.

Notably, this study was performed using a rat internal hemorrhoid model. Although the overlap between rat and human transcriptional information has been demonstrated, more clinical samples are needed to strengthen the scientific validity of this study. We also objectively recognize that certain species differences exists between rats and humans in terms of the fine anatomical structure of the anal cushion and the course of chronic disease. In our subsequent work, we will further optimize the model to better approximate the chronic pathological process of human internal hemorrhoids. Additionally, ST has not yet achieved single-cell resolution, and further in vivo functional experiments are warranted to precisely define the cellular localization. This work focused on early- to mid-stage internal hemorrhoid lesions, but the pathological mechanisms underlying the late chronic fibrotic phase of hemorrhoids remain to be explored. Based on the present results, our future research will focus on three key directions: Transcriptome sequencing of human hemorrhoid tissues to verify the cross-species conservation of our core observations; systematic in vitro and in vivo assays to confirm the causal relationship of the identified regulatory axis; and the construction of a chronic rat hemorrhoid model to elucidate the pathological mechanisms driving late-stage chronic fibrosis.

CONCLUSION

This study constructed a comprehensive atlas to explore the spatiotemporal transcriptomic dynamics underlying internal hemorrhoid development in rats, providing novel insights into the pathological mechanisms driving the progression of internal hemorrhoids through vascular inflammation to ECM remodeling.

ACKNOWLEDGEMENTS

The authors gratefully acknowledge the members of the experimental team, including the laboratory technicians, for their technical support, and the rats that served as experimental subjects for their contribution to this research.