Matrix Stiffness-Induced Transcriptome Alterations and Regulatory Mechanisms Revealed by RNA-Seq in Endothelial Cells ()
1. Introduction
Over the past two decades, it is widely accepted that mechanical signals elicit specific responses in cells, as biological components do. Numerous studies have shown that ECM stiffness, or elasticity, affects fundamental cellular processes, including spreading, growth, proliferation, migration, stem cell differentiation, organoid formation [1] and regeneration [2]. Furthermore, ECM stiffness is emerging as a prominent mechanical cue that precedes disease and drives its progression and is becoming a new potential therapeutic target [3]. ECM stiffness changes in many pathological conditions, such as atherosclerosis. It has been shown that the stiffness of normal arterial blood vessels is about 40 kPa, but in the early stage of atherosclerosis, the stiffness of blood vessels will be greatly reduced due to the accumulation of lipids in the intima [4] [5], thus may lead to abnormal function of the epithelium and even other cells beneath it.
Under physiological conditions, endothelial cells (ECs) play an important role in maintaining the integrity of the vessel wall, and are involved in the modulations of metabolic homeostasis (trophic functions), vascular hemodynamics (tonic functions), vascular permeability, coagulation, and cell extravasation (trafficking) [6]. Endothelial cells can respond to both physical and chemical signals. When pathological signals are delivered to endothelial cells, arterial remodeling may occur, which contributes to common diseases, such as atherosclerosis [7]. ECs cultured on matrices with different stiffness showed significant differences in cell morphology and cell mechanics [8]. On the one hand, under the activation of endogenous molecules, ECs are capable of acting like immune cells, secreting cytokines, chemokines, growth factors, etc. [6]. In addition, ECs exposed in vitro to TGF-β (tumor growth factor-β), TNF-α (tumor necrosis factor-α), and IFN-γ (interferon-γ) undergo Endothelial-to-Mesenchymal Transition [9]. On the other hand, ECs can sense physical cues such as fluid mechanical forces [10] [11] and matrix stiffness [12]. Studies have shown that matrix stiffness affects not only endothelial cell morphology but also its biological processes. Furthermore, matrix stiffness has a significant impact on the interactions between endothelial cells and monocytes. When ECs were cultured on matrices with different stiffness (8 kPa, 20 kPa, 40 kPa), more monocytes adhered to ECs cultured on 8 kPa and 40 kPa compared to those cultured on 20 kPa, which was correlated with the expression levels of vascular cell adhesion molecule 1 (VCAM-1) and intercellular adhesion molecule 1 (ICAM-1) [13]. In addition, stiffness also regulates autophagy of ECs [14], neuroblastoma angiogenesis [15], and endothelial monolayer integrity [16]. However, the underlying mechanisms remain to be elucidated, and the effect of stiffness change on the gene expression of ECs is not well understood.
ECs are localized at the interface between blood and tissue, which plays a critical role in maintaining multi-organ health and homeostasis [17]. Endothelial dysfunction is regarded as a hallmark of many cardiovascular diseases such as atherosclerosis. And it manifests as impaired endothelium-dependent vasodilation, elevated oxidative stress, chronic inflammation, leukocyte adhesion and hyperpermeability, and endothelial cell senescence [18]. As discussed above, physical signals such as matrix stiffness play an integral role in disease initiation and progression. Therefore, it is particularly important to further study how matrix stiffness influences endothelial cells. In this study, in order to discover the potential link between stiffness and ECs, we examined the ECs transcriptome expression profiles, which cultured on matrices with different stiffness (1 kPa and 40 kPa, respectively), and performed a series of bioinformatics analyses based on differentially expressed genes. We also aimed to find the potential markers for assessing the occurrence of early disease.
2. Materials and Methods
2.1. Preparation of Matrices with Different Stiffness
Matrices with various stiffness were prepared as previously described [5] [13]. The following experimental procedures were performed in a sterile environment. Briefly, 40% (w/v) acrylamide (Sangon Biotech, China) and 2% (w/v) bis-acrylamide (Sangon Biotech, China) were mixed in different ratios, and crosslinked with tetramethylethylenediamine (TEMED, Klamar, China) and 10% ammonium persulfate solution (Table 1). Before polymerization, 100 μL mixture was added onto a glass plate pre-treated with gel slick (Lonza, Switzerland), and a circular silanizedcoverslip with a diameter of 20 mm was put onto the gel solution. After 20 min of polymerization, the coverslip with the polymerized matrix was pried up from the glass plate, and the matrix on the coverslip was rinsed with 2 mL 50 mM HEPES (Thermo Fisher, USA) for 3 times. Then the matrix was treated with Sulfo-SANPAH (Thermo Fisher, USA) and coated with 1 mg/mL Collagen I (Corning, USA). Finally, the matrix was sterilized with 70% (v/v) ethanol solution for 15 min, and washed with Dulbecco’s Phosphate-Buffered Saline (DPBS). The prepared matrices were stored in a 4˚C refrigerator for up to about 2 weeks and UV sterilized for 30 min before cell culture.
![]()
Table 1. Composition of the PA gels with different stiffness.
2.2. Cell Culture
Human umbilical vein endothelial cells (HUVECs) were purchased from ScienCell Research Laboratories (ScienCell, USA) and cultured to passage 6 for experiments. HUVECs were cultured with EC medium supplemented with 5% fetal bovine serum, 1% EC growth factor and 1% penicillin/streptomycin (Sciencell, USA). In all experiments, HUVECs were cultured on matrices with different stiffness for more than 24 hours.
2.3. Reverse-Transcription and Real-Time Quantitative Polymerase Chain Reaction (RT-qPCR)
Total RNA was extracted using the TRIzol plus RNA Purification kit (Thermo Fisher Scientific, USA), and reverse transcribed using the FastQuant RT kit (with gDNase) (TIANGEN, China). The cDNA samples were analyzed by RT-qPCR (Thermo Fisher Scientific, USA) using SuperRealPreMix Plus (SYBR Green) (TIANGEN, China) and the specific primer pairs (Table 2). GAPDH was used as a normalization control.
2.4. RNA Isolation and Library Preparation
Total RNA for transcriptome sequencing was extracted using the TRIzol reagent according to the manufacturer’s protocol. RNA purity and quantification were
![]()
Table 2. Human mRNA primer sequences.
evaluated using the NanoDrop 2000 spectrophotometer (Thermo Scientific, USA). RNA integrity was assessed using the Agilent 2100 Bioanalyzer (Agilent Technologies, Santa Clara, CA, USA). Then the libraries were constructed using TruSeq Stranded mRNA LT Sample Prep Kit (Illumina, San Diego, CA, USA) according to the manufacturer’s instructions. The transcriptome sequencing and analysis were conducted by OE Biotech Co. Ltd. (Shanghai, China).
2.5. RNA Sequencing, Identification of Differentially Expressed Genes and Gene Enrichment Analysis
The libraries were sequenced on an Illumina HiSeq X Ten platform and 150 bp paired-end reads were generated. Raw data (raw reads) were processed using Trimmomatic [19]. The reads containing ploy-N and the low-quality reads were removed. The clean reads were mapped to the human genome (GRCh38) using HISAT2 [20]. FPKM [21] of each gene was calculated using Cufflinks [22], and the read counts of each gene were obtained by HTSeq-count [23]. Differential expression analysis was performed using the DESeq (2012) R package [24]. P value < 0.05 and fold change > 2 or fold change < 0.5 was set as the threshold for significantly differential expression. Hierarchical cluster analysis of differentially expressed genes (DEGs) was performed to demonstrate the expression pattern of genes in different groups and samples. The identified DEGs were used subsequently for Gene Ontology (GO, http://geneontology.org/) and Kyoto Encyclopedia for Genes and Genomes (KEGG, http://www.genome.jp/kegg/) pathway analysis, which performed respectively using R based on the hypergeometric distribution.
2.6. Protein-Protein Interaction Network Analysis of DEGs
In order to explore the interactions between DEGs, a protein-protein interaction (PPI) network was built using the online analysis database STRING (http://string-db.org/). Cytoscape (http://www.cytoscape.org/) software was used to visualize the network.
2.7. Analysis of Transcription Factors (TFs) and Their Target Genes
To explore the role of the TFs in mechanotransduction, we investigated the distribution of TF families. The relationship between TFs and their target genes was studied based on the hTFtarget database (http://bioinfo.life.hust.edu.cn/hTFtarget#!/tf).
2.8. Statistical Analysis
Graph Pad Prism 6.0 software was used for statistical analyses and graphing. For all bar graphs, data were expressed as means ± SD. Differential analyses of mRNA expression between different matrix stiffness groups were carried out using Student’s t-test. Statistically significance was indicated by p < 0.05.
3. Results
3.1. Transcriptome Sequencing and Identification of DEGs
Samples were categorized by matrix stiffness into soft matrix group (1 kPa) and stiff matrix group (40 kPa), and were sequenced on an Illumina HiSeq X Ten platform. After removing the adapter sequences and discarding low-quality reads, the amount of data of each sample was distributed in the range of 6.22 - 7.09 G, the Q30 base was distributed in the range of 93.27% - 94.02%, and the average GC content was 50.57%. By aligning reads to the reference genome, the mapping rates of all samples were greater than 95%. In summary, the data acquired in this study were eligible to be used for later analysis (Table 3).
Principal component analysis (PCA) was performed to show the linear combination of the original high-dimensional data after dimensionality reduction [25]. In this study, PCA results showed that the principal component 1 (PC1) and principal component 2 (PC2) had 69.5% and 16.5% variations, respectively, which indicated that there were obvious differences between different matrix stiffnesses and that the same matrix stiffness had good repeatability. Based on this, the samples were distinguished into 2 groups of different dimensions (Figure 1(a)). The RNA-seq results showed that gene expression was significantly different between the two groups. Therefore, the heatmap of the gene expression profile based on different gene expression values between the two groups was acquired (Figure 1(b)). Totally, there were 1775 DEGs identified, among which 450 DEGs were up-regulated and 1325 DEGs were down-regulated in HUVECs cultured on soft matrix (1 kPa) compared to those cultured on stiff matrix (40 kPa) (Figure 1(c)).
3.2. GO Enrichment Analysis of DEGs
Based on DEGs, GO analysis was conducted to evaluate the biological significance of matrix stiffness on ECs. In this study, we counted the number of differential genes in each GO term, and used the hypergeometric distribution algorithm to calculate their significance. We found that enriched GO terms were within three main categories: biological process, molecular function, and cellular component.
In the biological process category, most enriched genes were related to positive regulation of transcription by RNA polymerase, suggesting that changes in
![]()
Figure 1. Differential expression analysis of ECs samples between 1 kPa and 40 kPa. (a) Principal component analysis. Information about the separation of wildness training giant pandas and captive giant pandas still retained after dimensionality reduction, PC1 explains 69.5% of the variance, and PC2 explains 16.5% of the variance. (b) Heat map of DEGs. The red indicates relatively high expression, while blue indicates relatively low expression. (c) The MA plot shows the gene expression distribution of Log2EC and q-value for each gene. Red dots indicate genes of up-regulation, green dots indicate genes of down-regulation, and gray dots indicate genes with no significant change in expression.
matrix stiffness would exert extensive impact on many aspects of the ECs. In the molecular function category, most DEGs enriched among protein binding, ATP binding and metal ion binding. In the cellular component category, most DEGs were associated with nucleoplasm, cytosol, nucleus, centriole and membrane (Figure 2(a)).
Moreover, we also compared the gene distribution of GO level 2 containing 64
![]()
Figure 2. GO enrichment results of DEGs.
terms (Figure 2(b)). In each of the above category, the enriched terms were associated with cell-cell interaction and cell-ECM interaction, including biological adhesion and cell junction. In the molecular function category, binding showed distinct expression differences. These results implied that changes in matrix stiffness can affect cell junction, which in turn may disrupt vessel integrity and endothelial barrier function.
3.3. KEGG Pathway Enrichment Analysis of DEGs
KEGG pathway enrichment analysis was performed to further explore the signaling pathways involved in mediating the effects of matrix stiffness on cells. In this study, we show the top 20 signaling pathways with the highest enrichment scores after screening (Figure 3(a)). And, the gene expression involved in each pathway is shown in Figure 3(b). Among these pathways, focal adhesion and regulation of actin cytoskeleton, which play an important role in cell-ECM interactions and cell-cell communications, may indicate the potential targets of mechanotransduction [26] [27]. In addition, DEGs were also enriched in intercellular junctions (including adherens junction, tight junction and gap junction) and endocytosis, which play an important role in endothelial cell barrier function [28].
![]()
Figure 3. KEGG enrichment results of DEGs.
3.4. Expression Pattern of TFs and Its Target Genes
TFs play an indispensable role in regulating gene expression, are ubiquitous in all organisms and differences in their expression may provide important clues for understanding mechanotransduction. Therefore, the distribution of differential TF families was counted followed by their target genes (Figure 4). By comparing the distribution of TF families (up- and down-regulated), we found clear differences between the two groups (Figure 4(a)). Among these TF families,
![]()
Figure 4. TF families and its target genes distribution. (a) The x-axis and the y-axis represent different TF families and the number of genes, respectively. The ordinate is the number of genes; red represents all genes in the TF family; orange represents the differential TFs in the TF family; yellow and green represent the differences that are up-regulated and down-regulated TFs, respectively. (b) The x-axis and the y-axis represent different TFs and the number of its target genes, respectively. Red and orange represent the up-regulated and down-regulated genes, respectively.
zf-C2H2, chromosome remodeling factor, bHLH and ARID, are involved in mechanotransduction [29] [30] [31]. In order to further study the regulation of TFs on gene expression, we extracted the differential target genes corresponding to differential TFs according to the relationship list between TFs and target genes, and screened out the top 25 differential TFs (Figure 4(b)). In addition, we also found that some transcription factors like BRD4, BRD2, ARID1A and ARID2, which have been shown to be involved in mechanotransduction, were associated with a large number of DEGs [30] [32] [33].
3.5. Protein-Protein Interaction (PPI) Analysis
PPI analysis was performed based on the top 35 DEGs. In order to explore the interactions between DEG, we established a PPI network using an online tool STRING database (https://stringdb.org/), and then used Cytoscape (http://www.cy-toscape.org/) software to visualize the network (Figure 5). We
![]()
Figure 5. Protein-protein interaction analysis of DEGs. Node size represents the fold change, blue represents down-regulated genes and red represents up-regulated genes. Line width represents interaction score (>0.7).
used node size to represent the fold of differential expression, color to distinguish up- or down-regulated genes, and line width to indicate interaction score (>0.7). In this network, some genes encoding ECM proteins, which play a considerable role in maintaining the integrity of ECM, such as fibronectin 1 (FN1), collagen α2(IV) chain (COL4A2), matrix metalloproteinase 2 (MMP2) and matrix metalloproteinase 14 (MMP14), have obvious expression differences and interactions. Moreover, cadherin 5 (CDH5), which plays a critical role in endothelial adherens junction assembly and maintenance, presents the same variation. In addition, ITGA5 and ITGB5, which encode integrin subunits that mediate cell-ECM interactions and function as mechanosensor, have the same results.
3.6. Differential Gene Expression Validation by RT-qPCR
Based on the above bioinformatics analysis results, we found that the genes encoding some ECM proteins such as FN1 and COL4A, some genes encoding junction proteins such as CDH5 and F11R, the genes mediating cell-ECM communication such as ITGA5 had obvious differences in different matrix stiffness. ECM maintains the structural integrity of tissues and regulates cell and tissue functions. Dysregulation of the ECM composition and structure can lead to several human diseases [34]. The integrity of cell junctions is the basis for cells to function as a barrier. Therefore, six genes (including CDH5, F11R, FN1, COL4A2, MMP14 and ITGA5) were subject to PCR validation. As shown in Figure 6, the expression of these genes was reduced on soft matrix (1 kPa), which is consistent with the RNA-seq results. These results suggested that matrix stiffness may affect the endothelial cell barrier function as well as the structural integrity of ECM.
4. Discussion
The stiffness of the blood vessel is an important determinant of vascular homeostasis. The stiffness of the blood vessel is affected by aging and various pathological conditions, including hypertension, diabetes, hypercholesterolemia, chronic kidney disease, and unhealthy lifestyles [35]. For example, in the state of hypertension, a high shear stress environment can damage ECs, making them unable to regulate the entry and exit of substances such as lipids, which eventually leads to the changing of the stiffness of blood vessel and the formation of plaques [36]. A central function of vascular ECs, located at the interface between blood and tissue, is to serve as a barrier that regulates the entry and exit of solutes and cells, thereby maintaining homeostasis. ECs can directly sense the stimulation of mechanical signals and convert them into chemical signals [16] [37]. However, the understanding of how matrix stiffness affects EC function is still incomplete. So, in this study, we studied ECs cultured on matrices with different stiffness and compared their transcriptomes to determine any differences in DEGs using RNA-seq. A total of 1775 DEGs were identified, including 1325 up-regulated genes and 450 down-regulated genes on 1 kPa matrix compared
![]()
Figure 6. Comparison of expression level of selected genes using RT-qPCR. Gene expression is normalized to the GAPDH. The results are means of three biological replicates (±SE).
with 40 kPa matrix. PCA analysis and cluster analysis presented in Figure 1 indicated two distinct groups, suggesting that soften matrix stiffness (1 kPa) have great impact on EC function.
GO enrichment analysis further demonstrated that most of the DEGs were enriched to cellular components. In the biological process category, most enriched genes were related to positive regulation of transcription by RNA polymerase, suggesting that changes in matrix stiffness would affect gene transcription, thus exert extensive impact on many aspects of the ECs. In addition, the enrichment of proteins in nucleoplasm and nucleus (Figure 2(a)) suggested that matrix stiffness change may largely result in transcription level changes due to epigenetic modification and regulation of gene expression by binding of the TFs. We further compared the changes of the DEGs of each term at GO level 2. Cell junctions and ECM and proteins are among those with significant differences. Matrix stiffness is one of the important features of the ECM. These results suggest that stiffness change could affect cell-cell interactions as well as cell-ECM interactions [38] [39], which in turn may disrupt vessel integrity and endothelial barrier function.
The regulatory role of the TFs is crucial in all aspects of cellular functions, including mechanotransduction. For example, YAP and TAZ, best-known as mechanosensitive transcriptional regulators, whose nuclear translocation is influenced by ECM stiffness, function as important effectors of mechanotransduction to regulate cell behaviors [40]. Among all the analyzed the distribution of TF families, zf-C2H2 transcription family had the most differential genes. Then we screened out KDM2B, which belongs to the ARID family by counting the DEG each TF targeted. A recent study showed that KDM2B directly recruits RNA polymerase II to further initiate and promote Il6 transcription, which in turn facilitates inflammatory responses [41]. In addition, BRD4, which is similar to YAP, can also respond to changes in matrix stiffness and regulate pro-fibrotic gene expression [30]. Furthermore, at low mechanical stress, ARID1A binds to YAP/TAZ, thereby preventing the formation of TEAD-YAP/TAZ complex [42].
In the KEGG pathway enrichment analysis, DEGs were concentrated in several important signaling pathways such as Notch signaling pathway, cell-cell junctions, and focal adhesion (Figure 2). Notch signaling, which is an evolutionarily highly conserved signaling pathway, plays an important role in maintaining tissue homeostasis, development and disease occurrence [43]. He et al. found that Piezo1, which is known as a mechanosensitive ion channel, requires the Notch signaling cascade to mechanistically regulate the development of inflammation and renal fibrosis [44]. In addition, there are also significant differences in the proteins involved in adherens junction, gap junction, and tight junction, which play an important role in maintaining vascular barrier function and mediating cell-cell interaction [45]. We observed the down-regulation of claudin, nectin, VE-cadherin, vinculin and catenin, on 1 kPa matrix compared to 40 kPa matrix. Furthermore, there are also distinct differences in focal adhesion and ECM-receptor interaction that serve an important role in cell-ECM interaction and mechanotransduction. All these results implied that soft matrix caused by pathological conditions may directly attenuate vascular barrier function and homeostasis.
To further understand the biological significance of the DEGs in regulating EC function, PPI network analysis of DEGs was performed using STRING database, followed by visualization using Cytoscape (Figure 5). As shown in the results, some ECM proteins such as FN1, COL4A2 and MMP14, were all down-regulated on 1 kPa matrix compared to 40 kPa matrix. More importantly, the expression of CHD5 involved in adherens junction and F11R involved in tight junction, respectively, was lower on 1 kPa matrix. These finding were confirmed by RT-qPCR results (Figure 6). Both of them play key roles in maintaining ECM integrity, cell-ECM interactions, and transmitting mechanical signals [46]. Besides, ITGA5, which encodes integrin α5, also showed significant differences. Given the critical role of integrins in translating mechanical signals into biochemical signals, a large number of therapeutic approaches based on targeting integrins have been developed [47]. These results suggest that matrix stiffness may affect the integrity of ECM and endothelial cell barrier function through ITGA5.
5. Conclusion
In summary, we characterized the transcriptome profile of ECs cultured on matrices with different stiffness. The dynamic changes of DEGs of ECs contributed to understanding the molecular mechanisms of mechanotransduction. Our results suggested that matrix stiffness may affect ECM integrity, cell-ECM interactions, and cell-cell junctions, which in turn affect endothelial cell barrier function. The above findings provided unique insights for further studies on the mechanisms underlying how matrix stiffness affects cell behavior and could serve as the basis for the prevention and treatment of diseases.
Abbreviations
Endothelial cells (ECs);
Differentially expressed genes (DEGs);
Transcription factor (TF);
Protein-protein interaction (PPI);
Extracellular matrix (ECM);
Human umbilical vein endothelial cells (HUVECs);
Reverse-transcription and real-time quantitative polymerase chain reaction (RT-qPCR);
Fibronectin 1 (FN1);
Cadherin 5 (CDH5);
F11 receptor (F11R);
Integrin, alpha 5 (ITGA5);
Collagen α2(IV) chain (COL4A2);
Matrix metallopeptidase 14 (MMP14).