Substrate Stiffness and Stretch Regulate Profibrotic Mechanosignaling in Pulmonary Arterial Adventitial Fibroblasts

Pulmonary arterial adventitial fibroblasts (PAAFs) are important regulators of fibrotic vascular remodeling during the progression of pulmonary arterial hypertension (PAH), a disease that currently has no effective anti-fibrotic treatments. We conducted in-vitro experiments in PAAFs cultured on hydrogels attached to custom-made equibiaxial stretchers at 10% stretch and substrate stiffnesses representing the mechanical conditions of mild and severe stages of PAH. The expression of collagens α(1)I and α(1)III and elastin messenger RNAs (Col1a1, Col3a1, Eln) were upregulated by increased stretch and substrate stiffness, while lysyl oxidase-like 1 and α-smooth muscle actin messenger RNAs (Loxl1, Acta2) were only significantly upregulated when the cells were grown on matrices with an elevated stiffness representative of mild PAH but not on a stiffness representative of severe PAH. Fibronectin messenger RNA (Fn1) levels were significantly induced by increased substrate stiffness and transiently upregulated by stretch at 4 h, but was not significantly altered by stretch at 24 h. We modified our published computational network model of the signaling pathways that regulate profibrotic gene expression in PAAFs to allow for differential regulation of mechanically-sensitive nodes by stretch and stiffness. When the model was modified so that stiffness activated integrin β3, the Macrophage Stimulating 1 or 2 (MST1\2) kinases, angiotensin II (Ang II), transforming growth factor-β (TGF-β), and syndecan-4, and stretch-regulated integrin β3, MST1\2, Ang II, and the transient receptor potential (TRP) channel, the model correctly predicted the upregulation of all six genes by increased stiffness and the observed responses to stretch in five out of six genes, although it could not replicate the non-monotonic effects of stiffness on Loxl1 and Acta2 expression. Blocking Ang II Receptor Type 1 (AT1R) with losartan in-vitro uncovered an interaction between the effects of stretch and stiffness and angiotensin-independent activation of Fn1 expression by stretch in PAAFs grown on 3-kPa matrices. This novel combination of in-vitro and in-silico models of PAAF profibrotic cell signaling in response to altered mechanical conditions may help identify regulators of vascular adventitial remodeling due to changes in stretch and matrix stiffness that occur during the progression of PAH in-vivo.


Introduction
Pulmonary arterial adventitial fibroblasts (PAAFs) are important for vascular extracellular matrix (ECM) homeostasis and remodeling [1,2], and there is evidence that PAAFs are regulated by matrix stiffness [3][4][5], stretch [6,7], or overstretch injury [8,9] and hypoxia [1]. During injury, PAAFs are activated and differentiate into myofibroblast subtypes that remodel vascular wall properties by directly altering the expression, degradation, or cross-linking of ECM proteins including collagen, fibronectin, and elastin. Given that the ECM also serves as a substrate for cell adhesion and sends physical and chemical cues that determine cell phenotype [10], it has been suggested that matrix stiffening may signal tissue remodeling and be causative drivers of pulmonary hypertension [11].
Pulmonary arterial hypertension (PAH) is a vasculopathy manifested by sustained elevation of pulmonary arterial pressures, vascular constriction, and irreversible vascular acrylamide (#A9099 and #146072, Sigma Aldrich, St. Louis, MO, USA): 3% acrylamide and 0.06% bis-acrylamide were used for the construction of 0.5 kPa gels; 4% acrylamide and 0.3% bis-acrylamide were used for the constructions of 3 kPa gels; and 10% acrylamide and 0.1% bis-acrylamide were used for the construction of 10 kPa gels [15].
Custom-made circular stretchers were designed using computer-aided design and constructed with polycarbonate. Polydimethylsiloxane (PDMS) membranes were built by mixing the Sylgard 186 elastomer kit (#4026144, Dow, Midland, MI, USA), extruding onto a wafer, degassing in a vacuum chamber then curing in the oven. The PDMS membranes were treated with 10% benzophenone (#A10739, Alfa Aesar Thermo Scientific, Haverhill, MA, USA) for polyacrylamide gel adherence, as previously described by Herum et al. [10]. The polyacrylamide gels were constructed to be 25 mm in diameter, cross-linked through exposure to ultraviolet light for 25 min, attached to PDMS membranes, and surrounded by silicone grease to prevent cell migration and media leakage. The gels were equilibrated in 1× Phosphate Buffered Saline (PBS, # 10010023, Gibco Thermo Fisher Scientific, Waltham, MA, USA) overnight, then Corning ® collagen I (100 g/mL, #354236, Sigma Aldrich, St. Louis, MO, USA) was attached with 1-Ethyl-3-(3-dimethylaminopropyl) carbodiimide (#00050, Chemplex, Mahwah, MJ, USA) and N-hydroxysuccinimide (#A10312, Alfa Aesar Thermo Scientific, Haverhill, MA, USA) to facilitate cell adherence. The stretcher was assembled so that two full turns were equivalent to 10% static stretch, as previously done [16,17]. PAAFs were trypsinized from the tissue culture plates using 0.25% Trypsin-EDTA (#25200056, Gibco Thermo Fisher Scientific, Waltham, MA, USA) and seeded onto the gels at a density of 140,000 cells per gel. Cells were cultured at 37 • C, 5% CO 2 , 100% humidity for three days. Cells were changed into serum-free media for 24 h before being stretched. The stretch condition was applied for 24 h based on the increase in gene expression shown by Herum et al. in left ventricular cardiac fibroblasts [10].

Inhibition Studies
For the inhibition experiments, PAAFs were seeded onto 0.5 kPa and 3 kPa gels at a density of 40,000 cells per gel and cultured for three days as described above. The media was changed to serum-free media, and each gel slated for inhibition was preincubated for 4 h with 1 M losartan (#3798, Tocris Bioscience, Minneapolis, MN, USA). The dose of losartan was delivered according to work by Kim et al. in adventitial fibroblasts from 6-week-old Sprague-Dawley rats [18]. The cells were then stretched for 24 h as previously described in the Stretcher Preparation subsection. RNA Isolation of these cells was conducted as described below.

RNA Isolation
For RNA extraction of the normotensive pulmonary artery, the adventitial layer was sectioned into 6 pieces and submerged in TRIzol reagent (Invitrogen #15596026, Thermo Fisher Scientific, Waltham, MA, USA). Tissue was homogenized using a BeadBug homogenizer with 1.5 mm zirconium beads. (Benchmark Scientific, Sayreville, NJ, USA). For PAAF experiments, RNA isolation was carried out using TRIzol and 5PRIME Phase Lock Gel Heavy tubes (#2302830, Quantabio, Beverly, MA, USA) and RNA was extracted using either a standard TRIzol total RNA isolation protocol or the RNeasy ® Mini kit (#74104, Qiagen ® , Hilden, Germany), which was then reverse transcribed into cDNA using the NEB cDNA ProtoScript First Strand Kit (#E6300L, New England Biolabs, Ipswich, MA, USA) and then RNA was extracted using the RNeasy Mini kit (#74104, Qiagen, Hilden, Germany). Quantitative real-time PCR was performed using the StepOnePlus TM Real-time PCR machine (Thermo Fisher Scientific, Waltham, MA, USA) and KAPA SYBR Fast Universal qPCR kit (#KK4601, Roche, Basel, Switzerland) using primers targeting genes of interest listed in Supplementary Table S1 (produced by Integrated DNA Technologies, San Diego, CA, USA). Relative gene expressions were compared against housekeeping gene 18S ribosomal RNA unless otherwise noted.
The PAAF signaling network model was implemented as a system of logic-based ordinary differential equations that were integrated numerically using an explicit Runge-Kutta method. Briefly, each state variable yi that represents each of the species in the network is normalized to a value between 0 and 1 and follows a Hill-type activation function. The system of equations for the 69 nodes in the model takes the form: where ωi represents the weight for each node, τ is the time constant, and f is the sigmoidal Hill function of yi, where f(EC50, n, B) is a function of half-maximal activation EC50, cooperativity n, and constant B. The form of biochemical interactions between two species is specified by 'AND' or 'OR' Boolean operators that represent the biochemical interactions between upstream nodes expressed as continuous functions as described in Wang et al. [12]. , and 8 phenotypic outputs (green diamonds), modified from our previous work [12]. Activation is shown with arrows and inhibition is shown with blunt head arrows. Blue arrows indicate non-PAAF-based experiments. Magenta nodes indicate the Phosphoinositide 3-kinase (PI3K) pathway, orange nodes indicate the Reactive Oxygen Species (ROS) pathway, purple nodes indicate the mitogen-activated protein kinase (MAPK) pathway, blue nodes indicate the calcineurin pathway, and red nodes indicate the Hippo signaling pathway.
The PAAF signaling network model was implemented as a system of logic-based ordinary differential equations that were integrated numerically using an explicit Runge-Kutta method. Briefly, each state variable y i that represents each of the species in the network is normalized to a value between 0 and 1 and follows a Hill-type activation function. The system of equations for the 69 nodes in the model takes the form: where ω i represents the weight for each node, τ is the time constant, and f is the sigmoidal Hill function of y i , where f(EC 50 , n, B) is a function of half-maximal activation EC 50 , cooperativity n, and constant B. The form of biochemical interactions between two species is specified by 'AND' or 'OR' Boolean operators that represent the biochemical interactions between upstream nodes expressed as continuous functions as described in Wang et al. [12]. Inputs of the model (PDGF, TNF, TGF-β, hypoxia, AngII, FGF, stiffness, stretch) were initialized at 0.25. Default reaction weights ω i , Hill coefficient n, and half-maximal activation EC50 were set to 1, 1.4, and 0.6, respectively. The time constants τ for each different Cells 2021, 10, 1000 6 of 18 reaction followed those used previously [12]: 0.1 h for signaling reactions; 1 h for transcription; and 10 h for translation. Similar to the network model analysis by Tan et al., (2017), we chose a change in normalized model output values of 0.1 as the threshold for considering the output to have changed significantly by mechanical stimulation or for a significant response to have been significantly inhibited [39]. While Tan et al., (2017) used a threshold of 0.05, we chose a more stringent threshold of 0.1, but our conclusions were not affected by this difference. Parameters in the model were not optimized or fitted. Rather we chose equal parameters for all reactions using values from Zeigler et al. [32]. While the parameters EC 50 , weight, and Hill coefficient were set to be the same value across all reactions, the time constant τ was chosen according to the type of reaction. Furthermore, parameter values were tested for consistency with the mathematical constraints described in Cao et al. for this class of model [36]. In a previous comprehensive analysis of parameter uncertainty, we verified that model accuracy was robust to the choices of these parameter values [12].
The input weights of stretch and stiffness were both set to 0.25 to represent the softest matrix, 0.5 kPa, and no applied stretch. We increased the stiffness input weight to 0.7 and 0.9 to represent the effects of 3 kPa and 10 kPa substrates, respectively, and evaluated the model at t = 72 h to mimic the in-vitro experimental time course.
To numerically simulate the effects of stretch on PAAFs after 24 h and the changes in substrate stiffness for 72 h when inhibiting nodes, the model was evaluated at those time points (i.e., y i (t = 24) and y i (t = 72)) after the corresponding input weights of stiffness and stretch were increased from 0.25 to 0.7. To simulate the effects of inhibition, y i max corresponding to blocked nodes were set to 0, while the other parameters remained the same. The change for each gene was calculated with respect to each condition's control group.
To simulate the different conditions under which losartan inhibited AT 1 R, we conducted eight sets of simulations. Four sets of these simulations were evaluated on 0.5 kPa substrate stiffness for 72 h with parameters at baseline, and input weights stiffness and stretch set to 0.25. For simulations involving stretch but no inhibition, the input weight of stretch was increased from 0.25 to 0.7 and the model was evaluated at t = 24 h. For the unstretched and stretched inhibited conditions, y i max corresponding to the AT 1 R node was set to 0 before applying changes to the stretch input weight and was evaluated at t = 24 h. The same combinations were used for the other four set of simulations on 3 kPa substrate conditions, but with a stiffness input weight of 0.7.
The implementation of the model in Python 3.9.0 along with simulation data are available on Github.

Statistics
Descriptive statistics were performed using JMP Pro Statistical software (version 14, SAS Institute Inc., NC, USA) for group comparisons of relative gene expression. For normally distributed data, one-way analysis of variance (ANOVA) was used to test for differences in means of three different stiffnesses and gene expression of the normotensive pulmonary artery adventitial layer for all six genes, followed by the Dunnett's post-hoc test. Otherwise, the non-parametric Wilcoxon-Kruskal-Wallis statistic was used followed by Dunnett's post-hoc test. Effects of stiffness and stretch were tested using two-way ANOVA with stiffness and stretch as fixed factors. For normally distributed data, the Dunnett's posthoc test was used. Otherwise, the non-parametric Wilcoxon-Kruskal-Wallis statistic was used followed by Dunnett's post-hoc test. For the inhibition studies, three-way ANOVA was used to compare the effects of stiffness, stretch, and inhibition, followed by Sidak's post-hoc test. Data are expressed as means ± standard error of the mean, unless otherwise specified. Statistical significance was determined at a level of 0.05. Data were graphed in GraphPad Prism software (v8.4.3.686, San Diego, CA, USA) and Illustrator (Adobe, San Diego, CA, USA, v24.2.3).

PAAFs Upregulate Profibrotic Genes in Response to Increased Substrate Stiffness and Stretch
When we immunostained cultures for markers of endothelial cells (vWF), smooth muscle cells (MYH11), and myofibroblasts (vimentin) (Figure 2A-H), only 3% were positive for vWF and 0.2% were positive for MYH11 suggesting high enrichment of PAAFs in our cell cultures. Intact pulmonary artery tissue sections were stained and imaged with the same antibodies and imaging settings as positive controls for these markers ( Figure 2I-R).

PAAFs Upregulate Profibrotic Genes in Response to Increased Substrate Stiffness and Stretch
When we immunostained cultures for markers of endothelial cells (vWF), smooth muscle cells (MYH11), and myofibroblasts (vimentin) (Figure 2A-H), only 3% were posi tive for vWF and 0.2% were positive for MYH11 suggesting high enrichment of PAAFs in our cell cultures. Intact pulmonary artery tissue sections were stained and imaged with th same antibodies and imaging settings as positive controls for these markers ( Figure 2I-R). PAAFs expanded on plastic reverted from a myofibroblast to a fibroblast phenotyp after three days of culture on 0.5 kPa stiffness 6-well plates, as assessed by their low ex pression of Acta2 (Appendix A, Figure A1) and their rounded appearance in culture (Fig  ure 3B) compared with the more stellate shapes and higher Acta2 expression in cells grown on stiffer substrates ( Figure 3C,D). Messenger RNA levels of Col1a1, Col3a1, Eln, Fn1 Loxl1, and Acta2 genes in fibroblasts cultured on 0.5 kPa substrates were not significantly different from those obtained by extracting RNA from the pulmonary artery adventitia o a normotensive rat (p > 0.05 by one-way ANOVA, Figure 3A). This finding suggests tha cells cultured on a 0.5 kPa substrate may mimic expression of PAAFs in vivo with respec to the six genes studied in this paper.
Compared with mRNA levels in PAAFs cultured on 0.5 kPa substrates, all six gene were significantly upregulated in response to increased matrix stiffness (p < 0.05 by one way ANOVA and Dunnett's post-hoc test). The expression of Acta2 and Loxl1 was signif icantly higher on cells grown on 3 kPa matrices but not significantly higher on cells grown on 10 kPa matrices, while Col1a1, Col3a1, Eln, and Fn1 were significantly upregulated on 10 kPa substrates (comparable to arterial stiffness in advanced PAH [14]), compared with PAAFs cultured on 0.5 kPa matrices ( Figure 3A). Interestingly, Acta2 and Loxl1 expression PAAFs expanded on plastic reverted from a myofibroblast to a fibroblast phenotype after three days of culture on 0.5 kPa stiffness 6-well plates, as assessed by their low expression of Acta2 (Appendix A, Figure A1) and their rounded appearance in culture ( Figure 3B) compared with the more stellate shapes and higher Acta2 expression in cells grown on stiffer substrates ( Figure 3C,D). Messenger RNA levels of Col1a1, Col3a1, Eln, Fn1, Loxl1, and Acta2 genes in fibroblasts cultured on 0.5 kPa substrates were not significantly different from those obtained by extracting RNA from the pulmonary artery adventitia of a normotensive rat (p > 0.05 by one-way ANOVA, Figure 3A). This finding suggests that cells cultured on a 0.5 kPa substrate may mimic expression of PAAFs in vivo with respect to the six genes studied in this paper.
Compared with mRNA levels in PAAFs cultured on 0.5 kPa substrates, all six genes were significantly upregulated in response to increased matrix stiffness (p < 0.05 by one-way ANOVA and Dunnett's post-hoc test). The expression of Acta2 and Loxl1 was significantly higher on cells grown on 3 kPa matrices but not significantly higher on cells grown on 10 kPa matrices, while Col1a1, Col3a1, Eln, and Fn1 were significantly upregulated on 10 kPa substrates (comparable to arterial stiffness in advanced PAH [14]), compared with PAAFs cultured on 0.5 kPa matrices ( Figure 3A). Interestingly, Acta2 and Loxl1 expression exhibited non-monotonic responses, with significant upregulation of gene expression on  Examining the transcriptional responses of the six genes to 10% equibiaxial stre for 24 h in PAAFs ( Figure 4A-F), Col1a1, Col3a1, Eln, Loxl1, and Acta2 were significan upregulated compared with unstretched cells independent of the substrate stiffness ( 0.05 based on group comparisons made using a two-way ANOVA). Although Fn1 expr sion did not significantly change after 24 h of stretch ( Figure 3D), it was significantly regulated after 4 h on all gel stiffnesses (from 1.62 ± 0.34 to 5.07 ± 1.22 on 0.5 kPa gels, 3 ± 0.26 to 6.93 ± 1.65 on 3 kPa gels and 5.08 ± 1.16 to 8.04 ± 2.32 on 10 kPa gels, p = 0.000 = 6). On the other hand, Col1a1 was only significantly upregulated after 24 h of stretch, not after 4 h. This suggests Fn1 is transiently induced by a short period of stretch, wh the Col1a1 response to stretch is much slower. This finding is consistent with reports id tifying Fn1 as an early response gene [40]. No significant interaction effects between s strate stiffness and stretch were found in the expression of any of the six genes (Table Examining the transcriptional responses of the six genes to 10% equibiaxial stretch for 24 h in PAAFs ( Figure 4A-F), Col1a1, Col3a1, Eln, Loxl1, and Acta2 were significantly upregulated compared with unstretched cells independent of the substrate stiffness (p < 0.05 based on group comparisons made using a two-way ANOVA). Although Fn1 expression did not significantly change after 24 h of stretch ( Figure 3D), it was significantly upregulated after 4 h on all gel stiffnesses (from 1.62 ± 0.34 to 5.07 ± 1.22 on 0.5 kPa gels, 3.36 ± 0.26 to 6.93 ± 1.65 on 3 kPa gels and 5.08 ± 1.16 to 8.04 ± 2.32 on 10 kPa gels, p = 0.0002, n = 6). On the other hand, Col1a1 was only significantly upregulated after 24 h of stretch, but not after 4 h. This suggests Fn1 is transiently induced by a short period of stretch, while the Col1a1 response to stretch is much slower. This finding is consistent with reports identifying Fn1 as an early response gene [40]. No significant interaction effects between substrate stiffness and stretch were found in the expression of any of the six genes (Table 1).   Increased ECM stiffness significantly upregulated protein expression of Collagen III and Smooth Muscle Actin (SMA) from a baseline of 0.5 kPa at both 3 kPa and 10 kPa based on a post-hoc test ( Figure 5). There was no significant difference between 3 kPa and 10 kPa protein expression for either Collagen III or SMA. This is consistent with the relative expression of RNA for Collagen III (Col3a1) and SMA (Acta2), shown in Figure 4B,F, which significantly increased from 0.5 kPa to 3 kPa, without further significant increase from 3 kPa and 10 kPa substrates.  Increased ECM stiffness significantly upregulated protein expression of Collagen III and Smooth Muscle Actin (SMA) from a baseline of 0.5 kPa at both 3 kPa and 10 kPa based on a post-hoc test ( Figure 5). There was no significant difference between 3 kPa and 10 kPa protein expression for either Collagen III or SMA. This is consistent with the relative expression of RNA for Collagen III (Col3a1) and SMA (Acta2), shown in Figure 4B,F, which significantly increased from 0.5 kPa to 3 kPa, without further significant increase from 3 kPa and 10 kPa substrates.

PAAF Network Model Simulates Gene Expression Activated by Stiffness and Stretch
A threshold change of 0.1 in the normalized variable representing each of the six genes was used to classify the change in each gene as significant. The model predicted significant upregulation of all six genes in response to an increase in substrate stiffness from 0.5 kPa to 3-or 10 kPa ( Figure 6). These model predictions matched our experimental observation that all six genes were significantly upregulated in cells grown on stiffer matrices. The model also predicted upregulation of the gene expression of Col1a1, Col3a1, Loxl1, and Acta2 and the return to baseline of Fn1 expression 24 h after induction by 10%

PAAF Network Model Simulates Gene Expression Activated by Stiffness and Stretch
A threshold change of 0.1 in the normalized variable representing each of the six genes was used to classify the change in each gene as significant. The model predicted significant upregulation of all six genes in response to an increase in substrate stiffness from 0.5 kPa to 3-or 10 kPa ( Figure 6). These model predictions matched our experimental observation that all six genes were significantly upregulated in cells grown on stiffer matrices.

PAAF Network Model Simulates Gene Expression Activated by Stiffness and Stretch
A threshold change of 0.1 in the normalized variable representing each of the six genes was used to classify the change in each gene as significant. The model predicted significant upregulation of all six genes in response to an increase in substrate stiffness from 0.5 kPa to 3-or 10 kPa ( Figure 6). These model predictions matched our experimental observation that all six genes were significantly upregulated in cells grown on stiffer matrices. The model also predicted upregulation of the gene expression of Col1a1, Col3a1, Loxl1, and Acta2 and the return to baseline of Fn1 expression 24 h after induction by 10% The model also predicted upregulation of the gene expression of Col1a1, Col3a1, Loxl1, and Acta2 and the return to baseline of Fn1 expression 24 h after induction by 10% equibiaxial stretch. However, the model predicted Eln expression to be downregulated with stretch while experimental results showed upregulation. We investigated whether the inhibitory effect of JNK1/2 on Eln may have outweighed the activating effect of PKC and found that decreasing the weight of the inhibition of JNK1/2 on Eln [41] by 50% allowed the model to predict the observed upregulation of Eln (Figure 4). While the model was in qualitative agreement with the data, it did not recapitulate the non-monotonic responses of Loxl1 and Acta2 ( Figure 3A), which were significantly upregulated by 3 kPa matrix stiffness (compared with 0.5 kPa) but not by 10 kPa substrates.

Angiotensin II Receptor Inhibition Unmasks an Interaction Between Stiffness and Stretch on Fibronectin Gene Expression
Based on a sensitivity analysis [12], we simulated the effects of stretch and increased stiffness in the presence of inhibitors of three mechanosensitive nodes in the model (AT 1 , and TGF-β, and MST1/2). Table 2 shows the effects of inhibiting AT 1 , TGF-β, and MST1/2 on changes in gene expression due to an increase in substrate stiffness from 0.5 and 3 kPa and due to stretch on 0.5 kPa stiffness matrices. Here, model-predicted differences in the normalized mRNA variable due to inhibitor treatments were considered significant if they exceeded a threshold of 0.1. Table 2. Changes in gene expression due to inhibition of Angiotensin II Receptor Type I (AT 1 ), Transforming Growth Factor-β (TGF-β), and Macrophage Stimulating 1 or 2 (MST1/2) receptors in response to stiffness and stretch. Numbers in bold indicate activity changes greater than a threshold of 0.1.

Effect of Stiffness on 3 kPa
Effect of Stretch on 0.5 kPa Genes 14 0 From the model simulations, the induction of Loxl1 expression by increased substrate stiffness is specifically regulated by MST1/2 signaling, whereas the responses of the other five genes to stiffness were all significantly inhibited by blocking TGF-β receptor. AT 1 receptor inhibition significantly attenuated the stiffness-dependent induction of Col1a1, Col3a1, Fn1, and Acta2, but had no significant effect on Eln or Loxl1, and the magnitude of inhibition was noticeably less than when TGF-β receptors were blocked. Blocking angiotensin signaling in the model with increased substrate stiffness downregulated the collagens by 20% and blocking TGF-β signaling downregulated the collagens by 30%, while blocking angiotensin downregulated Acta2 by 17% and blocking TGF-β downregulated it by 86% (Table 2).
Blocking TGF-β signaling in the model while applying stretch stimulation suppressed the upregulation of Acta2 by 28% and reduced the downregulation of Eln by 11%. Stretch induction of Col1a1 and Col3a1 was shown to be reduced by inhibition of MST1/2 (by 46%) and angiotensin II signaling (by 26%), while Loxl1 regulation by stretch was affected only by inhibiting MST1/2 (by 100%). Fn1 expression, which was not significantly altered by stretch, remained unchanged in the presence of all three inhibitors. This is in contrast to its response to substrate stiffness, where inhibiting AT 1 and TGF-β receptors had a significant effect ( Table 2).
Since the model suggested a significant role for angiotensin II signaling in regulating Fn1 expression in response to increased stiffness but not stretch (Table 2 and Figure 7B), we treated cultured PAAFs with 1 M losartan, an AT 1 receptor blocker. Losartan abrogated the induction of fibronectin mRNA expression by 3 kPa substrates compared with 0.5 kPa, and Fn1 expression by PAAFs grown on 0.5 kPa matrices remained unresponsive to stretch after 24 h ( Figure 7A). Stretch significantly upregulated the expression of Fn1 by PAAFs grown on stiffer 3 kPa substrates when angiotensin II signaling was blocked.
Model simulations of Fn1 expression in response to increased stiffness, stretch, and AT1 receptor inhibition corresponding to each experimental condition in Figure 7A are shown in Figure 7B. The model recapitulated the increase in fibronectin mRNA under control conditions, when substrate stiffness increased from 0.5 kPa to 3 kPa, and the qualitative response of Fn1 on 3 kPa matrices to induction by stretch. The model also correctly predicted no effect of losartan on Fn1 mRNA levels under stretched and unstretched conditions on 0.5 kPa substrates 24 h after stretch. However, while we observed in vitro that losartan inhibited Fn1 upregulation by increased stiffness in unstretched but not stretched PAAFs, the model could not reproduce this observation.

Discussion
In-vitro experiments in pulmonary arterial adventitial fibroblasts (PAAFs) were used to investigate the differential effects of equibiaxial stretch and increased substrate stiffness on six genes of a new mathematical model of PAAF cell signaling [12]. While both physical stimuli occur in PAH, these stimuli are thought to occur at different stages of the disease, in part because increased vascular fibrosis leads to ECM stiffening that in turn opposes the increase in arterial wall strain caused by increased wall stress. In this study, we used a novel combination of in-vitro and in-silico models to investigate how PAAFs respond to changes in ECM stiffness and strains representative of those associated with adventitial remodeling in PAH. While PAAFs are exposed to cyclic loading in vivo, we used static stretch as a model of the chronically elevated mean hemodynamic load (Herum KM et al., 2017) during PAH rather than acute phasic vascular loading, in part because cyclic stretch systems cannot recapitulate the high physiological cardiac frequencies (6 Hz) in rats (Layland J et al., 1995) [10,42]. Although there are no existing studies examining how PAAFs respond to cyclic stretch, Wu et al. (Wu J et al., 2013) reported that 10% cyclic stretch for 36 h led to 2-3 fold increases in Col1a1 and Col3a1 expression in mouse aortic fibroblasts, Losartan abrogated the induction of fibronectin mRNA expression by 3 kPa substrates compared with 0.5 kPa, and Fn1 expression by PAAFs grown on 0.5 kPa matrices remained unresponsive to stretch after 24 h ( Figure 7A). Stretch significantly upregulated the expression of Fn1 by PAAFs grown on stiffer 3 kPa substrates when angiotensin II signaling was blocked.
Model simulations of Fn1 expression in response to increased stiffness, stretch, and AT 1 receptor inhibition corresponding to each experimental condition in Figure 7A are shown in Figure 7B. The model recapitulated the increase in fibronectin mRNA under control conditions, when substrate stiffness increased from 0.5 kPa to 3 kPa, and the qualitative response of Fn1 on 3 kPa matrices to induction by stretch. The model also correctly predicted no effect of losartan on Fn1 mRNA levels under stretched and unstretched conditions on 0.5 kPa substrates 24 h after stretch. However, while we observed in vitro that losartan inhibited Fn1 upregulation by increased stiffness in unstretched but not stretched PAAFs, the model could not reproduce this observation.

Discussion
In-vitro experiments in pulmonary arterial adventitial fibroblasts (PAAFs) were used to investigate the differential effects of equibiaxial stretch and increased substrate stiffness on six genes of a new mathematical model of PAAF cell signaling [12]. While both physical stimuli occur in PAH, these stimuli are thought to occur at different stages of the disease, in part because increased vascular fibrosis leads to ECM stiffening that in turn opposes the increase in arterial wall strain caused by increased wall stress. In this study, we used a novel combination of in-vitro and in-silico models to investigate how PAAFs respond to changes in ECM stiffness and strains representative of those associated with adventitial remodeling in PAH. While PAAFs are exposed to cyclic loading in vivo, we used static stretch as a model of the chronically elevated mean hemodynamic load (Herum KM et al., 2017) during PAH rather than acute phasic vascular loading, in part because cyclic stretch systems cannot recapitulate the high physiological cardiac frequencies (6 Hz) in rats (Layland J et al., 1995) [10,42]. Although there are no existing studies examining how PAAFs respond to cyclic stretch, Wu et al. (Wu J et al., 2013) reported that 10% cyclic stretch for 36 h led to 2-3 fold increases in Col1a1 and Col3a1 expression in mouse aortic fibroblasts, which are comparable to the 3-fold increase in Col1a1 and 2to 4-fold increase in Col3a1 that we observed after 24 h of static stretch [43]. While the cells were maintained at a 10% static stretch for 24 h, measurements of cell area in cardiac fibroblasts using the same circular custom stretchers (Herum KM et al., 2017) showed that after cell area initially increased during static stretch, they returned to their original size within 1 h, well before the 24-h time point at which gene expression was measured [10].

Stiffness and Stretch Differentially Affect Expression of Six Profibrotic Genes
Stretch and increased substrate stiffness were both able to upregulate five out of the six profibrotic genes we investigated. However, while increasing stiffness from 0.5 kPa significantly induced all six genes, fibronectin expression was transiently upregulated by stretch at 4 h but was not significantly altered by stretch at 24 h. There was also a non-monotonic response to the two levels of increased substrate stiffness in the expression of Loxl1 and Acta2, which were both upregulated compared with 0.5 kPa substrates on 3 kPa matrices (similar to vessel walls during mild PAH), but the expression of both was not significantly altered compared with 0.5 kPa substrates on 10 kPa matrices (which are comparable in stiffness to vessel walls during severe PAH). Unlike observations in cardiac fibroblasts [10], we found no statistical interaction effects between the stretch and stiffness conditions in these six genes. These results suggest that the expression of Col1a1, Col3a1, and Eln could be expected to rise early in vivo as elevated pulmonary arterial pressure increases vascular wall strain and remain elevated as fibrosis increases adventitial ECM stiffness, even though this stiffening would also reduce arterial strain. In contrast, Loxl1 and Acta2 expression may initially rise but eventually return to baseline as wall stiffening becomes severe, and Fn1 mRNA may be induced only after the ECM has remodeled and stiffened.

Model Modifications to Investigate Differential Regulation by Stretch and Stiffness
By allowing stiffness and stretch to be separate inputs to the model, we investigated the pathways regulating the expression of six mechanosensitive genes (Col1a1, Col3a1, Eln, Fn1, Loxl1, Acta2) in response to each stimulus. While there is published evidence that TGF-β is activated by stretch in cardiac [10] and lung [44] fibroblasts, we only found experimental evidence of TGF-β activation by substrate stiffness in PAAFs [23,24]. Based on ample published data in other cell types, we refined the model of the MAPK signaling cascade in the original version of our model so that ERK1/2, p38, and JNK1/2 could be independently activated, and we updated the model to include the effects of stretch activated TRP channels observed by Yue and Suzuma et al. [37,38].
Comparing the predictions of this revised model against the same independent experimental data from rat and human PAAFs that we used to test our original implementation [12], we found no significant changes in model validation accuracy from what we reported previously [12]. Comparing predictions of the revised model with in-vitro PAAF experiments conducted here on the effects of stretch and stiffness on gene expression, the model correctly predicted the upregulation of all six ECM genes by increased stiffness though not the subsequent return to baseline levels on the stiffest matrices for Loxl1 and Acta2. The model also correctly predicted the observed upregulation of four ECM genes and the lack of response to stretch in Fn1 expression at 24 hrs. However, while we observed an increase in Eln mRNA after stretch, the model incorrectly predicted a decrease. Examining the regulation of elastin gene expression in the network, we found that halving the weight of JNK1/2 inhibition on Eln mRNA while leaving the activating weight of PKC on Eln the same reversed this result. Hence it is possible that the activating effect of PKC dominates the inhibiting effect of JNK1/2 in the regulation of elastin gene expression by stretch.
The equations in the model were formulated using studies from both in-vivo and in-vitro experiments. While we used the rat PAAFs to validate the gene expression in response to stimuli such as mechanical stretch or substrate stiffness, this approach allows us to predict how phenotypic outputs respond to mechanical load. It is reported that mechanical stretch may increase the stiffness of the substrate, which in turn decreases the stretch. However, the interactions between them have not yet been classified. Through this work, we can model the interactions by adjusting time parameters and the activated reactions to represent beneficial versus maladaptive remodeling in fibrosis. Furthermore, the model can simulate a high number of experimental designs and make corresponding predictions that would be difficult to reproduce experimentally. This feature also enables the model to predict effects of specific drugs through simulating the activation or inhibition of any target species in the network.

Crosstalk between TGF-β and Angiotensin II
Using this model to predict the effects of inhibiting key mechanoresponsive nodes in the network, we found that blocking AT 1 R in the model significantly decreased expression of Fn1 in response to stiffness but did not significantly decrease expression of Fn1 in response to stretch. Experimentally, we confirmed that blocking the AT 1 receptor with losartan inhibited the significant upregulation of Fn1 expression when substrate stiffness is increased and had no effect on the response to stretch on 0.5 kPa substrates. However, losartan unmasked a response to stretch in PAAFs grown on the stiffer matrices that was not seen in untreated cells or predicted by the model. These findings show that angiotensin II signaling is required for the Fn1 response to increased stiffness, and that Fn1 expression can be stretch-regulated by an angiotensin-independent pathway on stiffer matrices when the saturating effects of higher stiffness are blocked. The requirement for angiotensin receptors to be activated before fibronectin mRNA can be induced by elevated substrate stiffness may be related to angiotensin-mediated conversion of latent TGF-β to the active state [45,46]. The network already represents this feedback based on reports that cyclic strain activates AT 1 R to cause activation of TGF-β in rat cardiac fibroblasts and human fibroblasts [45,47]. However, because stiffness also directly activates TGF-β in the model, angiotensin signaling could be blocked in the model without preventing stiffness from inducing Fn1 expression. This suggests that increased stiffness only activates TGF-β signaling after angiotensin has activated latent TGF-β in our experiments. Angiotensin receptor inhibition also unmasked angiotensin-independent Fn1 expression in response to stretch at a higher substrate stiffness, but in the model, fibronectin mRNA can only be induced by stretch directly via the angiotensin receptor. Hence, our experimental results suggest that there is another stretch-activated pathway that regulates a transient response to stretch in fibronectin and is more active in cells grown on stiffer matrices independent of angiotensin II. This indicates a currently unknown angiotensin-independent stretch-activated pathway responsible for an initial rapid upregulation of fibronectin gene expression, possibly the STAT3 pathway, which is involved in stretch induction of fibronectin in renal epithelial cells though not in fibroblasts [48]. An angiotensin-dependent pathway that downregulates fibronectin expression after 4 h should also be added to the model. One way to further examine this potential crosstalk would be to treat PAAFs with a TGF-β blocker based on the high magnitude of inhibition predicted by the model simulations (Table 2). When blocking TGF-β in the model, expression of Col1a1, Col3a1, Eln, Fn1, and Acta2 was significantly inhibited in response to increased stiffness, with no significant changes observed on inhibiting stretch effects except for Acta2.

Limitations
We used rat PAAFs because of the detailed biomechanical measurements of ECM stiffness pulmonary arterial strain in the sugen-hypoxia rat model of PAH and normotensive control rats. However, human PAAF cell lines have been used to study fibrotic signaling in response to increased ECM stiffness [11], where they showed that ten out of twelve genes studied were differentially expressed when stiffness increased from 1 to 12 kPa. Their analysis identified an miR-130/301-PPARγ signaling network regulated by ECM stiffness and associated with ECM remodeling in human PAH. Studies of mechanosignaling in human PAAF cell lines would enable us to generate a similar model of profibrotic mechanosignaling in human cells that could include these networks. ECM remodeling depends on protein synthesis, post-translational modifications, and cell-mediated matrix assembly [2]. One limitation of this study is that we focused primarily on gene expression, but we did find that changes in collagen III and smooth muscle actin protein abundances in response to increased ECM stiffness were consistent with changes in their mRNA expression. Finally, while the model was able to predict the response of fibronectin gene expression on soft gels, it did not replicate all of the observed responses to AT 1 receptor inhibition on stiff gels or the transient response of Fn1 expression to stretch. These model limitations can nevertheless be used to identify candidate pathways and reactions that need to be added to the network.

Conclusions
In-vitro experiments using hydrogel substrates of various stiffnesses coating elastic membranes in cell-stretch devices showed that expression of profibrotic genes by PAAFs is differentially regulated by cell stretch and extracellular matrix stiffness. No interaction effects between stretch and stiffness were observed for the six genes studied here; however, the AT 1 receptor blockade uncovered an angiotensin-independent activation of Fn1 expression by stretch in PAAFs when grown on stiff but not soft substrates. A novel combination of in vitro and in silico models of PAAF profibrotic cell signaling in response to altered mechanical conditions may help identify regulators of the vascular adventitial remodeling that results from the changes in stretch and matrix stiffness occurring during the progression of PAH in vivo.

Supplementary Materials:
The following are available online at https://www.mdpi.com/article/ 10.3390/cells10051000/s1, Table S1: Primers used to run RT-PCR are listed, as are protocols for conversion of RNA to cDNA, the PCR setup, raw data from the StepOnePlus TM Real-time PCR machine, and analyzed data for all RT-PCR experiments in this paper, DVJ 2021 PAAF stiffness Gel Stretch RT-PCR.xlsx; Table S2: Imaging data for Collagen III and SMA, IntDenCTCFActa2Col3.xlsx. Funding: This research was funded by the USA National Heart, Lung and Blood Institute (NHLBI) grant number 1R01 HL155945-01, 1R25HL145817-01, 1T32HL105373, 1R01HL137100, and the National Science Foundation (NSF) 2046259.

Institutional Review Board Statement:
The study was conducted in accordance with the Institutional Animal Care and Use Committee of the University of California San Diego (protocol number S17237 approved on 17 February 2021).