Cellular Contact Guidance Emerges from Gap Avoidance

Summary In the presence of anisotropic biochemical or topographical patterns, cells tend to align in the direction of these cues—a widely reported phenomenon known as “contact guidance.” To investigate the origins of contact guidance, here, we created substrates micropatterned with parallel lines of fibronectin with dimensions spanning multiple orders of magnitude. Quantitative morphometric analysis of our experimental data reveals two regimes of contact guidance governed by the length scale of the cues that cannot be explained by enforced alignment of focal adhesions. Adopting computational simulations of cell remodeling on inhomogeneous substrates based on a statistical mechanics framework for living cells, we show that contact guidance emerges from anisotropic cell shape fluctuation and “gap avoidance,” i.e., the energetic penalty of cell adhesions on non-adhesive gaps. Our findings therefore point to general biophysical mechanisms underlying cellular contact guidance, without the necessity of invoking specific molecular pathways.


INTRODUCTION
Cell organization plays a crucial role in the micro-architecture of tissues, dictating their biological and mechanical functioning. 1 In vivo, cells are embedded in the extracellular matrix (ECM), which generally comprises a network of organized micrometer-scale fibers that provide cells with anisotropic geometrical cues. [2][3][4] Cells typically align in the direction of the cues-a phenomenon known as cellular contact guidance. Contact guidance has been widely observed for several decades 5,6 and has been shown to affect various downstream cell behaviors, including survival, motility, and differentiation. 7-9 Uncovering the underlying mechanisms is critical for a better understanding of tissue and organ morphogenesis and regeneration. 10 Numerous studies have demonstrated and examined contact guidance effects, typically by culturing cells on microfabricated anisotropic substrates consisting of microgrooves, adhesive lines, or fibers (see Tamiello et al. 11 for an extensive review). Collectively, these studies have shown that contact guidance intriguingly occurs on substrate patterns of a wide range of sizes. For example, it was shown that contact guidance at length scales smaller than focal adhesions (FAs) arises as a result of constrained FA alignment and maturation in the direction of anisotropy. 9,12 More recently, we showed that cells constrained on single lines of fibronectin align in the line direction not only when the line width is smaller than the cell body, thereby In the present study, we investigated the emergence of cellular contact guidance resulting from anisotropic cues at these intermediate length scales, from the typical size of FAs to the typical cell length. We found that myofibroblasts exhibited two regimes of cellular alignment: one at small length scales of the anisotropic cue, where cell alignment is induced by multiple patterns, and one at large-length scales, as a result of spatial confinement of the whole cells. Interestingly, these alignments occurred in the absence of spatially constrained FA alignment, suggesting that contact guidance can result from an alternative mechanism. To understand this, we extended a recently established statistical mechanics framework for living cells, 14 enabling the simulation of cells on substrates with a heterogeneous distribution of ligands. Quantitative comparison between our experimental and computational data indicates that contact guidance arises from the minimization of cellular adhesions on non-adhesive regions and that non-adhesive gaps play a decisive role for cell alignment.

Cellular Alignment and Intracellular Organization Are Controlled by Substrate Anisotropy at Micron-to Cell-Size Scale
To mimic the fibrillar nature of the ECM and examine the effects of substrate anisotropy on cellular alignment, we created micropatterns consisting of parallel lines of fibronectin with defined widths and inter-line spacings using microcontact printing and seeded human myofibroblasts on these microcontact printed substrates. This experimental approach allowed us to systematically probe cellular contact-guidance response in a purely planar setup, without the need to consider the influence of the third dimension as in studies employing microfabricated grooves and ridges. The line width w and inter-line spacing s ranged from 2 to 200 mm and were initially chosen to be equal (i.e., w = s) to maintain a constant cell-substrate contact area. The lower limit of w = 2 mm was motivated by the size of typical mature FAs, 15,16 whereas the upper limit of w = 200 mm is large enough to allow cell spreading on one single line (cf. individual myofibroblast length on fibronectin-coated substrate; Figure 1C).
Images obtained 24 h after seeding showed that cell morphology and orientation are strongly influenced by the width of the lines ( Figure 1A). On the thinnest lines (w = 2 mm), there was only a weak cell-alignment effect due to the micropatterned lines. When w was increased up to 20 mm, cells increasingly elongated and aligned parallel to the lines. The trend inverted when w was further increased up to 200 mm. We quantified the changes in cell shape and orientation for more than 600 cells on the substrates using an automated morphometric analysis of the immunofluorescence images. 17 Briefly, we fitted an ellipse to the cell outline and defined the orientation angle q as the angle between the major axis of the best-fitted ellipse and direction of the lines ( Figure 1B). The analysis revealed that, with increasing w, the distribution of q around 0 strongly narrowed up to w = 20 mm, indicating enhanced cell alignment, and then broadened back toward an isotropic distribution, resembling that on homogeneous substrates (control; Figure 1D). Furthermore, we noted that w = 20 mm also demarcates the transition from the situation where cells adhered on more than one line to the situation where cells fitted within single lines.

OPEN ACCESS
These data demonstrate that the cellular orientation response can be divided into two regimes: regime I for w < 20 mm, where cell alignment was induced by multiple lines, and regime II for w R 20 mm, where cell alignment was influenced by the spatial confinement within single lines. In regime II, as we previously observed, 13 cell alignment decreased with increasing w. In regime I, cell alignment decreased with decreasing w. The two regimes were also captured by the order parameter Q (see Equation 1 in the Experimental Procedures and Figure S1), clearly showing that the order transitions at w = 20 mm.
Cell alignment on the micropatterned lines was accompanied by changes in cell morphology, as quantified through aspect ratio and adhesion area ( Figure S2A). The largest degree of alignment was observed when cells were most elongated. Because cell morphology is governed by the dynamics and properties of intracellular structural components, such as the nucleus and individual actin fibers, 18 we hypothesized that the two observed regimes should also reflect in the organization of these components. Indeed, the actin fibers and nuclei showed similar trends as for the cell orientations, including the two regimes of alignment (Figures 1D and (D) The cell, actin fiber, nucleus, and FA orientation, where 0 represents the direction of the lines. The boxes of the boxplots represent the quartiles of the distributions, with the whiskers indicating the outliers in the experiments and the 5 th and 95 th percentiles of the distributions. Note that, with this data representation, the median is at 0 and the box ranges from -45 to 45 when the distribution of cell orientation is perfectly isotropic. The data reported are results from three independent samples; at least 60 cells were considered per condition.

OPEN ACCESS
Cell Reports Physical Science 1, 100055, May 20, 2020 3 Article S2B). This confirms that micrometer-scale variations in the extracellular patterns can tune cell orientation as well as intracellular organization.
Contact Guidance Does Not Require Constrained Alignment of FAs On patterns with dimensions smaller than FAs, it has been previously reported that spatially constrained alignment of FAs play a central role in contact guidance. 9 Generally, when cells adhere to substrates homogeneously coated with fibronectin, they first form small, nascent FAs (0-2 mm long), which can either disappear or develop into 2-to 6-mm-long, mature FAs. 19 Because a minimum length of 2 mm has been shown to be required for local contacts to establish adhesions 20 and the FAs of our cells can grow to lengths much larger than 2 mm ( Figure S3), we reasoned that the lines with w = 2 mm would provide an area for FA maturation and orientation only in the direction of the lines and that FAs could have a wide distribution of orientations on 5-and 10-mm lines, where FA maturation is not spatially constrained. By examining the FA and cell orientation on these line widths, we therefore can test the requirement of FA alignment for cellular contact guidance. Indeed, Ray et al. 9 have recently suggested that relaxing spatial constraints on FA maturation and consequently FA organization diminishes cellular alignment. Surprisingly, we obtained the opposite result: on lines of w = 2 mm, a wide distribution of orientation angles of FAs was observed, whereas the FAs became slightly more aligned at larger line width w ( Figure 1D). This trend is similar to, but weaker than, the orientation response of the cells. Therefore, our data show that, at length scales larger than FA size, increasing the adhesive area for FAs leads to the counterintuitive increase of FA and cell alignment in the direction of the lines. This suggests that contact guidance at these length scales does not arise from spatially constrained alignment of FAs, which is an underlying mechanism of contact guidance at smaller scales. 9,12 To further confirm this observation, we investigated in more detail the morphology and organization of FAs in regime I. The analysis showed that lines of w = 2 mm were able to constrain FA maturation. In particular, vinculin staining showed that FAs were primarily formed on the adhesive fibronectin lines and that the distributions of the FAs were markedly influenced by the lines (Figure 2A). Specifically, for w = 5 and 10 mm, we observed rows of FAs lining up along the edge of the lines (Figure 2A, insets). For cells on the 2-mm-wide lines, we found small, non-aligned adhesions, in contrast to the large, non-aligned adhesions on homogeneously coated substrates of fibronectin. Moreover, FAs were significantly smaller for 2-mm lines relative to the homogeneous substrate ( Figure S2C), consistent with the notion that lines of 2-mm width provided a geometrical constrain for FA maturation. Increasing the line width w resulted in more aligned, elongated FAs in the direction of the lines (Figure 2A). We further characterized the size and shape of individual FAs by determining their length and aspect ratio, respectively. The length and aspect ratio of FAs on 2-mm lines were significantly smaller compared to those on homogeneous substrates, whereas the long axes of FAs on 5-and 10-mm lines were equal to those on the homogeneous substrates ( Figures 2C and 2D). The aspect ratio of FAs on 5-mm lines was smaller compared to FAs on 10-mm lines, meaning that FAs on 5-mm lines are long and wide, although FAs on 10-mm lines are thinner.
Given these differences in the FA morphology, we classified the FAs into three categories: ''round'' FAs are those with aspect ratio smaller than 1.6; ''non-aligned'' FAs are non-round FAs whose orientation angle q is larger than 20 ; and ''aligned'' FAs are non-round FAs with q smaller than 20 . 17 We found that the fraction of aligned FAs increases and the fraction of non-aligned FAs decreases for increasing line widths ( Figure 2E). Interestingly, for 2-mm lines, only a small amount of FAs aligned in direction ll OPEN ACCESS of the lines (11%), comparable to that on a homogeneous substrate (13%). These results strongly suggest that constrained FA maturation on 2-mm lines was not enough to induce a cell alignment distribution different than for cells on homogeneous substrates.
To evaluate further whether the orientation of FAs has a role in contact guidance for linear patterns larger than FAs, we constructed 2D density plots showing the distribution of the orientations of cells and FAs. The data show that strong cell alignment (i.e., contact guidance) does not correspond to a similarly strong alignment of mature FAs. Specifically, the density plots demonstrate that the distribution of orientation angles of FAs was scattered, showing all possible orientations (Figure 2B). This was again comparable to the orientations of FAs on homogeneous substrates (control), despite the increasing cell alignment for increasing w.
Together, these results indicate that, at the investigated length scales, a mechanism independent of FA alignment determines contact guidance of myofibroblasts. In an earlier study, we demonstrated that contact guidance for cells constrained on single lines (i.e., in regime II) emerges from an entropy-mediated mechanism, driven by the non-thermal shape fluctuations of cells. 13 We therefore asked whether such a mechanism may in fact be responsible for contact guidance over both regimes that we observed in the present study. To test this, we extended the statistical homeostatic mechanics framework developed by Shishvan et al. 14 to model cell behavior. This framework was adopted because it captures an essential statistical feature of cell behavior: cells exhibit large fluctuations in terms of single-cell geometrical parameters (e.g., area, orientation, and aspect ratio), with clear trends that arise only when the statistics of these observables are compared for different experimental conditions. This is evident in our experimental observations (e.g., Figure 1D). The main hypothesis is that cells continually change their shape, orientation, and intracellular components by (1) minimizing their free energy over the short-term period (seconds), while (2) maintaining an average homeostatic free energy over the long-term period (minutes). Recent reports strongly suggest that these hypotheses hold for different cell types. For example, Suresh et al. 21 adopted an analogous framework to predict stem cell differentiation.
Shishvan et al. 14 computed the free energy associated with each cell configuration on a substrate, f, by accounting for the energy contribution of the stress fiber cytoskeleton f cyto and the passive elasticity of other cellular components F elas . By simulating a large number of possible cell configurations, each with its free energies, previous studies have shown that this modeling framework can give critical insights into the ensemble behavior of adherent cells cultured on homogeneous substrates of varying stiffnesses 14 or on single adhesive lines of different width. 13 To model cell behavior on substrates with heterogeneous distribution of ligands as in the micropatterned substrates used in the present study, we extended this framework by considering the free-energy contribution of cell adhesions f adh , which are approximated as linear springs with spring constants k n and k a on non-adhesive and adhesive areas, respectively (see Experimental Procedures for a detailed description). As a consequence of these energy contributions, cells elongate up to a certain extent to minimize the contribution of both the stress fibers and passive elastic components to the cell free energy, while at the same time avoid forming adhesions on non-adhesive areas of the substrates because they are associated with high values of free energy.
A comparison between the computational and experimental results for myofibroblasts on substrates with lines of width w is shown in Figures 3 and S6. Remarkably, the simulations were able to model the two distinct regimes of cellular alignment observed experimentally ( Figure 3A). Similar to the experimental observations, in regime I, cells increasingly aligned and elongated in the direction of the lines with increasing w up to a maximum at w = 20 mm ( Figures 3A and 3B). It is worth emphasizing that these two regimes of cellular alignment emerged solely from the inclusion of f adh , without adjusting any parameters for the other free-energy contributions. Furthermore, dynamic cellular protrusions were not simulated, leading to the expected underestimation of cell area compared to experimental data ( Figure 3C). The simulations predicted that cells bridge a decreasing number of adhesive lines with increasing w ( Figure 3D). This was in qualitative agreement with the experiments and indicated that cells thin along the direction perpendicular to the lines, thereby resulting in a reduced number of gaps they are in contact with. In regime II, cells spread on single adhesive lines ( Figure 3D). At the onset of regime II (i.e., ll OPEN ACCESS on 20-mm-wide lines), cells fit on single lines by thinning along the perpendicular direction, while strongly elongating and orienting in the direction of the lines ( Figures  3A and 3B). With increasing w, the perpendicular thinning and longitudinal elongation necessary for cells to fit on single lines are smaller and therefore cells have increasing freedom to vary their main orientation, consistent with our earlier study. 13 These results indicate that the hypotheses of free-energy minimization over the short timescale (seconds) and homeostasis over the long timescale (minutes) can explain the variation of cell morphologies observed for cells on substrates with alternating (non-)adhesive lines.
Given the agreement between the model and experimental observations, we proceeded to analyze the model features that are responsible for cell alignment in the simulations. The probability density functions of the total Gibbs free energy of cells in regime I (i.e., w < 20 mm) exhibited smaller interquartile ranges compared to cells in regime II ( Figure 3E), suggesting that adhesion to multiple lines forces cells to explore less free-energy states and morphologies. In addition, we found increasing absolute values of the elastic Gibbs free energy for increasing line widths ( Figure 3F), corresponding to an increase in cell deformation, which was compensated by a corresponding decrease in the cytoskeletal Gibbs free energy ( Figure 3G). The Gibbs free energy of cell adhesions was relatively large for cells on multiple lines compared to that for the other line widths ( Figure 3H). This suggests that cells reduce the adhesion energy by minimizing the contact with the non-adhesive lines by (1) reducing their total spread area ( Figure 3C) or by (2) thinning in the direction perpendicular to the lines, resulting in overall cell alignment ( Figures 3A and 3B). Interestingly, cells on 50 mm and 200 mm show almost identical probability density functions for the cellular free energies, although a much higher alignment was observed for the 50-mm-wide lines compared to the 200-mm ones ( Figures 3A and 3E). Thus, it appears that, in regime II, a relatively high degree of alignment can be achieved without substantially perturbing the cell free energy.
Notably, cells on 20-mm-wide lines have energy profiles that are remarkably different from all the others; the total and cytoskeletal Gibbs free energies are much lower, with wider interquartile ranges compared to the other line widths ( Figures 3E and  3G). The Gibbs free energy of cells was affected by forcing cells to spread along single lines that are small compared to the size of the cell, inducing a very high degree of alignment. Taken together, the statistical thermodynamics framework suggests that contact guidance results from the tendency of cells to spread to minimize their cytoskeleton free energy while avoiding the formation of cell adhesions on non-adhesive lines. In other words, contact guidance may emerge from the maximization of actomyosin polymerization into stress fibers and minimization of FA formation on non-adhesive areas.

The Non-adhesive Gap Size Determines the Degree of Cell Alignment
The insights from the model imply that cellular contact guidance may in fact be induced by the cells' reluctance to bridge across non-adhesive gaps. The gaps between adhesive patterns have been previously suggested to act as a barrier to the formation of stable adhesions and protrusions. 22 Live-cell imaging of myofibroblasts stained for the actin cytoskeleton (silicon rhodamine [SiR] actin) and FAs (GFP-Talin) for 24 h indicated that, while spreading, cells formed numerous protrusions and adhesions at both ends of the cell body, and smaller adhesions were also formed on the side of the cells, even after complete cell spreading and alignment. Moreover, cells were observed to form exploratory protrusions extending perpendicular to the lines, bridging the non-adhesive area, but these were often short lived and eventually retracted without initiating further cell spreading in that direction (see an example in Figure S3 and Video S1). As a result, over time, the cells elongated and oriented in the direction of the lines. This corroborates the notion that non-adhesive gaps between the adhesive lines hinder stable adhesions and thereby promote anisotropic cell spreading.
An implication of these experimental and computational findings is that contact guidance at length scales from the FA to the cell size is induced not by the width  Article in the cell aspect ratio, which was much higher when cells stayed within single lines ( Figure 4D). When s > 20 mm, the myofibroblasts failed to reach the neighboring adhesive line and became highly elongated, resulting in a pronounced alignment in direction of the lines ( Figure 4E). This was accompanied by intracellular remodeling, such as increasing organization of the actin fibers and elongation of the nucleus. Remarkably, the transition at s = 20 mm coincides exactly with the transition from regime I to regime II that we saw earlier ( Figure 1D), indicating that the latter is caused by the myofibroblasts' (in)ability to bridge multiple lines separated by large gaps.
We also checked the effect of varying line widths w, but not the inter-line spacing s, by comparing patterns of 2 3 5 mm with 5 3 5 mm and 2 3 10 mm with 10 3 10 mm lines. Consistent with our hypothesis, varying the line width w alone had little effect on the orientation of cells and FAs ( Figures 5B, 5D, and 5E). By contrast, inter-line spacing s has a clear effect on cell morphology, as quantified by their area and aspect ratio ( Figure 5C). For increasing s, cells spread more; however, the FAs were not affected. By comparing the distribution median of the absolute

OPEN ACCESS
value of cell angles (using the independent-samples median test), we observed that the difference between cells on substrates with different values of w and equal values of s was not statistically significant ( Figure 5E). In contrast, changing values of s generally led to statistically different median values (p < 0.05 for all comparisons except for 2 3 5 mm and 2 3 10 mm, p = 0.54).
Taken together, our experiments combined with computational modeling support the idea that cell alignment induced by anisotropic cues larger than FAs and smaller than the typical cell length can result from the minimization of cellular adhesions on non-adhesive gaps-a phenomenon we define as ''gap avoidance.''

DISCUSSION
The ability of cells to align in the direction of topographical or biochemical linear patterns, i.e., contact guidance, is critical in various physiological contexts. [23][24][25] Previous studies on contact guidance have shown that anisotropic features created by microcontact printing can direct cell alignment. 16,[26][27][28][29][30] Although numerous studies have already extensively characterized the ability of a variety of cells to respond to anisotropic features present in their environment, a universal mechanism underlying this phenomenon remains elusive. Particularly, although spatial constriction of FAs 9,12,31 or actin fiber structures (e.g., stress fibers, filopodia, lamellipodia) 6,30,32,33 have been proposed to play a role in this process ( Figure 5D), the contribution of each of these intracellular components to the emergence of contact guidance at intermediate length scales is still unclear.
Our combined experimental and computational work demonstrates the existence of two distinct regimes of cellular alignment, governed by the length scale of the anisotropic cue ( Figure 1D). These distinct regimes have not been observed previously, as previous studies have focused on limited ranges of line widths, in some cases using microfabricated substrates that present confounding additional cues in the third dimension as well, such as grooves and ridges (see Tamiello et al. 11 for a review). For w <20 mm, myofibroblasts aligned parallel to the lines, bridging multiple lines. Decreasing w results in reduced cellular alignment, and at w = s = 2 mm, contact guidance effect is greatly decreased. A recent study by Ramirez-San Juan et al. 30 has observed a much stronger alignment of NIH 3T3 fibroblasts on parallel fibronectin lines with w = s = 2 mm than that of the myofibroblasts in our study. The discrepancy between this previous study and our observations might result from the differences in cell size. Indeed, Ray et al. 9 showed that the degree to which different cell types respond to anisotropic substrates can vary significantly, and 3T3 fibroblasts are typically much smaller than human myofibroblasts. 34 The model presented in the present study suggests that the energy penalty of a gap size has in proportion a stronger effect if the cell size is smaller. Thus, the stronger alignment of 3T3 cells compared to that of myofibroblasts in response to 2 3 2 mm lines might derive from the different sizes of these two cell types. To test whether this model prediction is worth considering, we analyzed the alignment of human cardiomyocyte progenitor cells (CMPCs) (average spread length 100 mm) on our fibronectin lines. We found that CMPCs showed a very similar contact-guidance response to that of myofibroblasts (average spread length 180 mm), as expected ( Figure S4). Moreover, the transition from cell alignment on single lines to alignment on multiple lines for CMPCs occurs at s = 10 mm, smaller than for myofibroblasts (s = 20 mm), in qualitative agreement with our model prediction. Article systematically varying spacing between adhesive lines could further elucidate the role of cell type and size on contact guidance.
We observed that the two regimes of cellular orientation were also reflected in the orientation of the individual actin fibers and nucleus, although the two regimes were less clear for FAs. Previous studies demonstrated that FA dynamics are affected by the geometric dimensions of anisotropic substrates. 9,12,31,35 In addition, higher levels of aligned mature adhesions were observed along lines or microgrooves. However, a closer examination of the FAs of cells in regime I revealed that many FAs were not fully aligned in direction of the lines, although the cells showed a clear increase in alignment for increasing line widths ( Figures 1A and 2B). Given that the typical FA size is larger than 2 mm ( Figure 2D), we expected that, on 2mm-wide lines, FAs would orient and elongate in direction of the lines. However, we found that FAs were non-aligned, smaller, and less elongated for w = 2 mm compared to on homogeneously coated substrates ( Figures 2C and 2D). These results indicate that, at intermediate-length scales, contact guidance arises from a mechanism other than constrainment of FA alignment.
In a seminal study on contact guidance, Dunn and Heath 6 postulated that this phenomenon arises when non-adhesive gaps are wide enough to cause a mechanical restriction for the formation of actin protrusions, causing cells to align in direction of the lines. Although the alignment of filopodia and lamellipodia 22,30,33 has been proposed to play a role in this process, the exact contribution of these protrusions is not fully clear. To determine the effects of non-adhesive gaps on cell alignment, we varied s and showed that cells can still form protrusions to bridge the non-adhesive gaps of length scales much larger than the characteristic length of filopodia ($5 mm; Figures 4A and 4B). 36 This confirms that the formation of actin protrusions was not restricted and mechanical restriction of the formation of actin protrusions is not responsible for cellular contact guidance. Additionally, previous studies have also found that blocking filopodia formation did not affect contact guidance. 37 Figure 5B shows that cellular alignment is not significantly affected by increasing w but is enhanced by increasing s, indicating that the interline spacing s is a more important parameter in guiding cellular alignment than line width w.
To understand the mechanisms determining contact guidance, several computational models have been developed that account for the (re)organization of intracellular structural components, i.e., FAs, actin stress fibers, and actin protrusions. The computational model of Loosli et al. 38 successfully predicted experimental observations reported in the literature, 27 but the hypotheses at the basis of the model were phenomenological. Similarly, Vigliotti et al. 39 could predict contact guidance for cells on nanogrooved substrates 40 with a cooperative feedback model for the interplay between FA formation and actin stress fibers development. However, this approach can only be used for groove sizes that are very small compared to the cell size. Recently, an increasing number of computational models were established for cell spreading on substrates, based on the hypothesis that free energy minimization drives this biological process. [41][42][43][44] This hypothesis successfully explained the spreading of cells on smooth substrates 42 or patterned islands. 43 However, these models were unable to account for the full range of experimental observations, except by including a fitting (phenomenological) parameter for the cell fluctuations. Shishvan et al. 14 recently postulated that the fluctuations in cell morphology are constrained by the fact that cells remodel their shape and intracellular structures to maintain an overall homeostatic state over the long-term period.

OPEN ACCESS
Here, we extended this framework by taking the adhesive properties of the substrate into account. With this inclusion, the experimentally observed regimes of cellular alignment and aspect ratio were successfully captured.
The proposed computational framework does not include the directionality of FA formation or the dynamics of actin-rich protrusions. Still, the model was able to fully predict the two regimes of cellular alignment, clearly suggesting that alignment of FA or protrusions are not necessary ingredients for the alignment of cells in response to anisotropic adhesive lines, in agreement with our experimental results. The model further suggests that the tendency of cells to align in the direction of the lines is a direct consequence of minimizing the number of non-adhesive gaps to bridge in order to reduce the contribution of the cellular adhesion energy to the free energy. This is consistent with the observation of Romsey et al., 29 who demonstrated that the number of lines a cell contacts affects the fidelity of contact guidance. Indeed, the average number of adhesive gaps to bridge or number of adhesive lines that each cell contacted decreased with increasing s ( Figure 3D). Altogether, our study showed that cellular alignment is determined by the adhesion properties between the cells and the non-adhesive gaps. Future studies should be directed to dissecting in greater detail the role of dynamic processes, such as actin protrusions, in this gap-avoidance behavior. Moreover, it will be instructive to explore the role of these cellular processes in the transition between adhesion-mediated and gap-mediated contact guidance, for example, by systematically monitoring cellular responses above and below 2 mm, which is technically impossible using our current experimental approach.
To summarize, we have shown that alignment of myofibroblasts can be induced by anisotropic geometrical cues ranging from micron to hundreds-of-micron scale. Our experimental and computational results suggest that cell alignment at these length scales emerges from the tendency of cells to elongate and maximize actin polymerization, while avoiding the formation of cellular adhesions on non-adhesive gaps (''gap avoidance''). This is a generic biophysical mechanism underlying the morphological fluctuations of cells that does not require specific biochemical regulation or molecular pathway. Thus, this understanding not only offers an attractive, alternative explanation for how substrate anisotropy regulates the cellular orientation response but can also be relevant for devising new experimental strategies for directing tissue morphogenesis and regeneration.

EXPERIMENTAL PROCEDURES Design and Fabrication of Stamps
Microcontact printing was used to pattern adhesive lines ranging from 2 mm to 200 mm. The micropatterned substrates were fabricated via standard photolithography techniques, according to previous protocols. 45 Briefly, the desired features were generated onto a silicon master by deep reactive-ion etching (Philips Innovation Services, Eindhoven, the Netherlands) from a chromium photomask (Toppan Photomask, Corbeil Essonnes, France). The silicon surface was passivated with a fluorosilane, and microstamps were obtained by molding the silanized silicon master with polydimethylsiloxane (PDMS) (Sylgard 184; Dow Corning) and curing agent (10:1), which was cured at 65 C overnight. The cured PDMS stamps containing the desired features were then peeled off from the master.

Microcontact Printing
For microcontact printing, the stamps were cleaned by sonicating in 70% ethanol for 30 min and dried using compressed air. Article stamps were incubated for 1 h at room temperature with a 50 mg/mL rhodamine fibronectin (FN) solution (Cytoskeleton). The substrates for printing were flat PDMS-coated glass coverslips that were oxidized in a UV/Ozone cleaner (PDS UV-ozone cleaner; Novascan, Ames, IA) for 8 min just before use. The FN-coated stamps were dried under sterile air flow and gently deposited on the substrates for 15 min at room temperature. Uncoated regions were blocked by immersing the micropatterned coverslips for 5 min in a 1% solution of Pluronic F-127 (Sigma-Aldrich). Finally, the coverslips were three times washed with phosphate-buffered saline (PBS) and stored in PBS at 4 C before use. As a control substrate, a homogeneous fibronectin coating was obtained using a flat PDMS stamp.

Cell Culture
Human vena saphena cells (HVSCs) were harvested from the vena saphena magna obtained from patients according to Dutch guidelines of secondary used material and have previously been characterized as myofibroblasts. 46 The myofibroblasts were cultured in advanced Dulbecco's modified Eagle's medium (Invitrogen, Breda, the Netherlands) supplemented with 10% fetal bovine serum (Greiner Bioone), 1% penicillin/streptomycin (Lonza, Basel, Switzerland), and 1% GlutaMax (Invitrogen). Only cells with a passage lower than 7 were used in this study. To avoid effects on cell alignment deriving from cell-cell contacts, the microcontact printed substrates were seeded with a cell density of 2,000 cells/cm 2 . The myofibroblasts were cultured for 24 h at 37 C and 5% CO 2 on these substrates.

Immunofluorescence Labeling
After culture on the micropatterned substrates, the myofibroblasts were washed with PBS and fixed with 3.7% formaldehyde in PBS (Sigma-Aldrich) for 15 min at room temperature. The cells were permeabilized with 0.5% Triton X-100 (Merck) in PBS for 10 min at room temperature.  were transduced with Talin-GFP CellLight BacMam 2.0 (Life Technologies) according to manufacturer recommendations to stain for the FAs. Culture flasks were incubated with approximately 30 BacMam particles/cell for at least 24 h. Then, the cells were trypsinized and resuspended in 100 nM SiR-actin (Cytoskeleton) for 45 min at 37 C. This cell suspension was centrifuged for 7 min at 1,000g and resuspended in medium. We added 100 nM of SiR-actin to the medium and seeded the cells at 2,000 cells/cm 2 on the printed substrate. After 30 min incubation to allow cell adhesion, the sample was placed in the microscopy incubator under controlled temperature (37 C) and CO 2 (5%) conditions. Data were collected approximately 90 min after cell seeding by imaging at 512 3 512 pixels every 15-20 min over 24 h with a confocal microscope (Leica TCS SP 5 confocal microscope equipped with a Leica MC170 HD camera) using a 203/0.7 (0.32 mm/pixel) objective.

Image Analysis
The cellular, nuclear, and FA morphological parameters were assessed from Phalloidin-Atto, DAPI, and vinculin-stained images, respectively, that were analyzed using a custom-built script publicly available either in Mathematica (Wolfram Research, Mathematica, Version 11.1, Champaign, IL, USA, 2017) 17 or in MATLAB (The Mathworks, Natick, MA, USA; https://gitlab.tue.nl:443/stem/sfalab/tree/v0. 01). Briefly, the images of cells, nuclei, and FAs were binarized and fitted with an ellipse using a least-square algorithm. This was used to quantify the orientation, which was defined as the angle q between the major axis of the ellipse and the direction of the patterned lines, where q = 0 represents the direction of the lines. The cell orientation order can be characterized by the order parameter Q: 49 where <: > denotes ensemble averaging over all measurements. The values of Q range from 0, representing random alignment, to 1, representing perfect alignment. To represent the shape of the cells, nuclei, and FAs, we computed the aspect ratio (i.e., the ratio of the major axis to the minor axis) and area. For the quantification of the FA orientation, only elongated FAs (aspect ratio larger than 1.6) were taken into account. Actin fiber orientation was quantified using a fiber orientation algorithm, based on the work by Frangi et al. 50

Statistical Analysis
The data were obtained and pooled from at least three independent experiments (substrates), each containing multiple (4-8) micropatterns. For each condition, we analyzed at least 60 cells unless indicated otherwise. To assess differences between the different line widths on the morphological features of the FAs, oneway ANOVA with a Bonferroni post hoc test was used.

Modeling
To get a better understanding of the mechanisms determining the response of myofibroblasts to multiple adhesive (fibronectin) lines alternated with non-adhesive lines, we extended a recently developed statistical thermodynamics framework for cells. 14  Article A system comprising a single cell in contact with a substrate and surrounded by a nutrient bath was analyzed. This is an open system because the cell exchanges molecular species with the surrounding nutrient bath. Over the short-term period of seconds, this exchange is negligible and the system can be assumed as closed. Therefore, the state of the intracellular molecules over this short-term period can be predicted by Gibbs free-energy minimization as in classical statistical mechanics. Over the long-term period of minutes to hours, the exchange of molecules through the cellular membrane is not negligible. Over this time span, the cell is then out of equilibrium in terms of classical statistical mechanics and exhibits fluctuations in terms of its observables. However, it is known that cells, via intracellular biochemical processes, actively strive to maintain an average homeostatic number of molecular species within their body. 5 This, in terms of Gibbs free energy, translates to enforcing the constraint that the average Gibbs free energy of the system is constant over the long-term period. 14 By enforcing that the molecular species within the cell body minimize their Gibbs free energy over the short-term period (seconds), while the average Gibbs free energy of the system is constant over the long-term period (minutes to hours), the probability distributions of the cell orientation, aspect ratio, etc. can be obtained. Due to the cell movement and exchange of nutrients, these observables will exhibit fluctuations much larger than thermal fluctuations but proportional to a homeostatic temperature, a concept that was introduced in Shishvan et al. 14 Similar to previous studies, 13,14 the Gibbs free energy of a cell configuration ðjÞ was defined as where V cell denotes the cell volume and f is the specific Helmholtz free energy of the cell. This last term was computed as where the terms f cyto and F elas , respectively associated with the stress fiber and cell passive components, were defined as in Shishvan et al. 14 Compared to previous studies, here, the term f adh was added to consider the Helmholtz free energy of cell adhesions. Specifically, cell adhesions on the substrate were approximated as linear springs, and it was assumed that springs associated with non-adhesive parts of the substrate have a very low stiffness k n compared to the stiffness k a of springs on adhesive areas (k n << k a ). From this, it follows that the Helmholtz free energy of a cell adhesion at a cell material point x can be computed as where FðxÞ is the magnitude of the force (in material point x) that the linear spring has to exert to equilibrate the residual forces resulting from the interaction between stress fibers and cell passive components. We observe that, from k n << k a , it follows ðFðxÞÞ 2 2kn [ ðFðxÞÞ 2 2ka . Consequently, cell configurations without cell adhesions on non-adhesive points have much lower values of Helmholtz (and Gibbs) free energy and are considerably more likely to occur. Therefore, the approach efficiently models the preference of cells to form and mature adhesions only on adhesive lines.

Material Parameters for Myofibroblasts
A single set of material parameters was adopted and unchanged for all simulations of this study. Most of the model parameters describe inherent properties of stress fiber proteins that are not expected to differ among cell types. Differences among cell types are recapitulated in terms of the stress fiber protein volume fraction F 0 , the passive elastic properties described by the shear modulus of the cell (m), the in-plane bulk modulus of the cell k, and m (a material constant that describes the nonlinearity of the deviatoric elastic cell response). Overall, increasing the term F 0 causes a widening of the Gibbs free-energy distribution and an increase of the homeostatic area of a cell on a homogeneous substrate, although an increase of the cellular elastic properties has the opposite effect. In addition to these parameters, the homeostatic area can also change because of the inherent size of a cell type. In the model, this parameter is represented by the scaling parameter l w . For more details on the model parameters, we refer the reader to the Supplemental Information and to previous studies. 14,51 The parameters adopted in the present study are based on validations performed in previous works. 14,52 In particular, all parameters except for k a and k n were taken in agreement with Shishvan et al. 14 (k a and k n were not present in this previous study). Shishvan et al. demonstrated that, with these parameters, the model can predict morphological features (e.g., area and aspect ratio) of smooth muscle cells spreading on substrates with different stiffness. Myofibroblasts have very similar mechanical properties as smooth muscle cells, and therefore, the same parameters were adopted in the present study. The parameter k a was then chosen such that the focal adhesion energy computed for cells on homogeneous and stiff substrates was similar to the energy reported in McEvoy et al., 52 who analyzed the focal adhesion free energy for cells on substrates with different rigidities. Finally, k a /k n = 100 was chosen such that adhesions formed on non-adhesive areas were hardly observed in the simulations. A list of the material parameters is given in Table S1.

Computational Simulation
The cell in the initial, undeformed (circular) configuration was discretized with three-noded triangular elements. Different configurations ðjÞ were then obtained by varying the displacements u ðjÞ , which uniquely identify the strains in the cell. From these values, by considering the chemical equilibrium in the cell body, f cyto and F elas are computed, as well as the passive stress within the cell and the stress fiber stress. Given these passive and active stress values, the residual traction forces Tðx ðjÞ Þ can be identified at each node, together with the consequential values of Fðx ðjÞ Þ = À Tðx ðjÞ Þ. From these forces, the Helmholtz free energy f adh and the total Gibbs free energy can be computed. For obtaining the statistics of cell shapes on a substrate, the homeostatic Gibbs free energy and the associated homeostatic temperature were first computed iteratively using Metropolis algorithm, analyzing 2,000,000 cell-shape perturbations. The same procedure was then repeated to simulate cells on microcontact printed substrates. This was achieved without changing the material parameters but by changing the homeostatic temperature of the system to enforce the condition that the average Gibbs free energy of a cell in every substrate is equal to the homeostatic energy found for the free-standing cell. See Supplemental Experimental Procedures for a complete description of the simulations.
The morphologies of representative cells from the simulations and experiments are compared in Figure S6. In general, the cellular shapes obtained with the Article computational simulations are comparable to those experimentally observed. For each cell shape resulting from the simulations, the cellular orientation and aspect ratio were obtained by fitting the outline of the cell with the best-fitting ellipse (in a least-square sense), as in the experiments. The statistics of these simulations are reported in Figure S5.

Data and Code Availability
Codes used to analyze the experimental images are publicly available at https:// gitlab.tue.nl:443/stem/sfalab/tree/v0.01 and https://doi.org/10.1371/journal. pone.0195201. Codes used to perform the computational modeling of cells are publicly available at https://zenodo.org/record/3710132. All data associated with the study are included in the paper and the Supplemental Information or from the lead author upon reasonable request.

Modeling framework
To predict the shape and orientation of cells on multiple adhesive lines of fibronectin, we extended the statistical framework of Shishvan et al. [1], which is here briefly described. In what follows, we refer to: "the cell", to indicate a generic cell to model; "the substrate", as a generic substrate with alternating (non-)adhesive lines; "cellular configuration", to indicate the mapping of the cell material points on the substrate; and "the system", as the system composed of the cell and the substrate.

Characterization of a configuration
The cell was approximated as a two-dimensional body in the x 1 − x 2 plane, with through-thickness stress Σ 33 = 0. X indicates the coordinates of each material point of the cell in the undeformed configuration, omitting the index of all different material points for brevity. The cell in the undeformed configuration was assumed cylindrical, with radius R 0 and thickness b 0 = R 0 /5. The other configurations were characterized by the displacements of each material point u (j) (X) ≡ u (j) , such that x (j) = X + u (j) are the coordinates of the material points of the cell in the configuration with index (j).

The probability of a configuration
By assuming a rigid and purely elastic substrate not deformed by cells, the Gibbs free-energy of a configuration (j) was defined as: where V cell denotes the volume of the cell and f is the specific Helmholtz free-energy of the cell. As described in Shishvan et al. [1], by assuming Gibbs free-energy minimization over the short-term period (seconds) and homeostasis over the long-term period (minutes), it can be derived that a configuration (j) has probability where Z = j exp −ζG (j) is the partition function and ζ is a constant such that Here,Ḡ is the homeostatic Gibbs free-energy, corresponding to the free-energy of a free-standing cell (i.e. isolated, not interacting with any substrate). The average magnitude of the fluctuations of the energy of a cell on a specific substrate is proportional to 1/ζ, which is referred to as homeostatic temperature [1]. This quantity is a characteristic of the interaction between the cell and the substrate, and was computed for each different substrate.

Definition of the specific Helmholtz free-energy f
In Shishvan et al. [1], the Helmholtz free-energy f was computed as f = f cyto + Φ elas , where f cyto and Φ elas are the specific Helmholtz free-energies associated with stress fibers and other cellular passive components (e.g. cell membrane), respectively. Here, we added a term f adh corresponding to the specific Helmholtz free-energy of cell adhesions to account for the heterogeneity of the substrate in exam, where cells can easily form adhesions on the adhesive lines, while the formation of cell adhesions on the non-adhesive lines is obstructed: The terms f cyto and Φ elas were computed as in [1], by using the model of Vigliotti et al. [2] to compute f cyto , and a two-dimensional adaptation of the Ogden hyperelastic model [3] to calculate Φ elas . In what follows, the computation of these terms is summarized, followed by the definition of the term f adh introduced in the present study.

The specific Helmholtz free-energy of stress fibers f cyto
Vigliotti et al. [2] envisioned stress fibers as formed by a number of functional units with reference length l 0 . The quantities describing stress fibers are defined in a representative volume element (RVE) that, in its undeformed configuration, is defined as a cylinder with the axis along the thickness of the two-dimensional cell, radius n R l 0 /2, and thickness b 0 (such as the cell). Stress fibers initiate from the center of this RVE, oriented at angles φ (−π/2 ≤ φ < π/2) with respect to the axis x 1 . In the undeformed state, n R functional units are present in the stress fibers in each direction while, along the thickness, n s layers of identically distributed stress fibers are present. Given an infinitesimal area of the cylinder of the RVE oriented at an angle φ, the number of stress fibers crossing this area is dΠ = ηdφ, such that η(φ) is defined as the angular stress fiber concentration. When a nominal strain ε nom (φ) is applied along the normal of the surface dA, the stress fibers in that direction deform and, in response, remodel. In particular, due to the deformation of stress fibers and their functional units, the number of these functional units increases (in case of tensile strain) or decreases (for compressive strain) till their deformed length reaches an optimal value l ss . It follows that, in addition to the nominal strain of the RVE, we can define a nominal strainε nom (φ) for the stress fiber functional units, relative to l 0 . A relationship between ε nom (φ) andε nom (φ) can be found by observing that n 0 ≡ n R [1 + ε nom ] functional units of length l 0 are necessary to cover the deformed stress fiber length n R l 0 [1 + ε nom ]. In contrast, if n functional units are present along such length, each of them has the length n R l 0 [1 + ε nom ] /n, and thus At steady-state each unit has length l ss , thus and are the nominal strain and number of stress fiber functional units along a direction φ at steady-state.
In the time-scale analyzed in this study, the total number of functional units within the cell normalized over the cell volume was assumed to be constant and it was denoted with N 0 ≡ V cell N T dV /V 0 , where V cell and V 0 represent the volume of the cell in the deformed and undeformed configurations, respectively. Given the conservation of this value, normalized quantities were defined When the cell is at steady-state,ε nom =ε nom = 0, and from the assumption that stress fibers exert tension with mechanisms similar to striated muscle, it follows that in these conditions they exert the maximum stress σ max . For these conditions, the chemical potential associated with each functional unit that is part of a stress fiber is where the first term on the right-hand side accounts for the enthalpy µ b of n R functional units, while the second term considers the entropy of mixing between packets of unbound and bound stress fiber proteins with sites for bound proteins, and the entropy of mixing between unbound proteins and lattice sites. For a complete derivation of this chemical potential, we refer the reader to [1]. In Eq. (9),η max represents the maximum value ofη such that all available lattice sites N L are occupied, whileN L ≡ N L /N 0 . The enthalpy µ b of bound proteins can be quantified as where µ b0 is the internal energy of n R bound stress fiber functional units and Ω ≡ A 0 n R l 0 is the volume of n R functional units of length l 0 . The chemical potential associated with every unbound functional unit is with µ u the internal energy of unbound functional units potentially forming n R stress fiber units. Now, given ρ 0 ≡ N 0 /V cell , the specific Helmholtz free-energy accounting for the stress fiber cytoskeleton is This quantity can be computed, given the previous definitions, once calculatedη andN u . Note that f cyto has to be computed for equilibrium morphological microstates, which means dG (j) = 0. It can be demonstrated [1] that this requirements translates to the constraints of chemical and mechanical equilibria X b = X u and Σ ij = 0, where Σ ij is the total Cauchy stress, characterized as a summation between the stress fiber stress, the stress from the cellular adhesions, and the passive stress of the other passive cellular components. From the chemical equilibrium X b = X u and Eqs. (9) and (11), an expression forη in terms ofN u follows: The normalized valueN u can be computed by recalling that N T = N u + N b is the total number of functional units in a material point x i , while N 0 V 0 /V R is the total number of functional units that can be formed within the cell. These two quantities are related as The integral V cell N u dV can be easily computed because, from the chemical equilibrium, it follows that X u is constant across the cell and therefore, due to Eq. (11), also N u is constant. Now, from the definitions (14), we have from whichN u can be computed, since it is the only unknown (given Eq. (13)).

The specific Helmholtz free-energy of cellular passive components Φ elas
The specific Helmholtz free-energy of passive components of the cell was modeled with a modified version of the Ogden hyperelastic strain energy density function, with the modification necessary to specialize the formulation to a two-dimensional setting. In particular, where µ is the shear modulus, m a material constant, λ I and λ II are the principal stretches, κ is the in-plane bulk modulus, whilek and J c are material parameters introduced to account for the penalty associated with cellular configurations with an areal stretch λ I λ II below a threshold J c . For a complete derivation and motivation of Eq. (16), we refer the reader to [1].

The specific Helmholtz free-energy of cell adhesions f adh
Cell adhesions on the substrate were approximated as linear springs, and it was assumed that springs associated with non-adhesive parts of the substrate have a very low stiffness k n compared to the stiffness k a of springs on adhesive areas (k n k a ). From this it follows that , if x (j) is adhesive, where F (x (j) ) is the magnitude of the force (in material point x (j) ) that the linear spring has to exert to equilibrate the residual forces resulting from the interaction between stress fibers and cell passive components. From k n k a it follows (F (x (j) )) 2 2kn (F (x (j) )) 2 2ka ; therefore, the approach efficiently models the preference of cells to form and mature adhesions only on adhesive lines (Eqs. (1)-(17)).

The magnitude of the force of cellular adhesions
In the present study, we assumed that cellular adhesions exert forces of a magnitude equal to the forces necessary to mechanically equilibrate the active and passive forces of the cell. This translates to assume that, given the total Cauchy stress Σ ij as a summation between the stress fiber stress, the stress from cellular adhesions, and the passive stress of the other passive cellular components, we have with σ sf ij stress fiber stress, σ adh ij stress due to cellular adhesions, and σ p ij stress from other passive cellular components. It follows that, to quantify the magnitude of the force of cellular adhesions, computing σ sf ij and σ p ij is sufficient. In agreement with [1], the stress fiber stress was computed as with F 0 ≡ n s (A 0 l 0 )ρ 0 , while φ * is the angle with respect to x 1 , in the deformed configuration, that stress fibers originally oriented at an angle φ have in the deformed configuration. Note that Eq. (19) was determined with the assumption the cell is incompressible. This assumption, leads also to the determination of the Cauchy stress of the passive components of the cell: where σ p k ≡ λ k ∂Φ elas ∂λ k is the principal passive stress, with associated p

Numerical method for computing the Gibbs free-energy of a configuration
The cell in the undeformed configuration was discretized with three-noded triangular elements of size R 0 /10. Different configurations (j) were then obtained by varying the displacements u (j) , which uniquely identify the strains in the cell. From these values, by considering the chemical equilibrium in the cell body, f cyto and Φ elas can be computed, as well as the passive stress within the cell and the stress fiber stress. Given these passive and active stress values, the residual traction forces T (x (j) ) can be identified at each node, together with the consequential values of F (x (j) ) = −T (x (j) ). From these forces, the Helmholtz free-energy f adh (Eq. (17)) and the total Gibbs free-energy (Eqs. (1) and (4)) can be computed.
For the calculation of the statistics of cell shapes on a substrate (Eq. (17)), the values ofḠ and ζ needs to be determined. First, the homeostatic energyḠ of a free-standing cell was computed. Only one equilibrium configuration is possible for a free-standing cell: the cell needs to be traction free and in chemical equilibrium without the contribution of the adhesions with the substrate. The total Gibbs free-energy can be computed for this unique configuration: eq G (j) = P 1 eq G 1 = G 1 , where the index 1 refers to the unique configuration of the free-standing cell (which has P 1 eq = 1). An iterative method was pursued to calculate ζ. Given a hypothesized value ζ 0 , the corresponding average Gibbs free-energy G 0 = 1 Z j G (j) exp −ζ 0 G (j) was approximated with the Metropolis algorithm [4]. If G 0 =Ḡ (with a maximum absolute error equal to 2%Ḡ), ζ 0 was accepted; otherwise, the procedure was repeated with a new value ζ 1 . For the random walks among the possible configurations required for the Metropolis algorithm, M points were chosen among the material points of the cell. Hypothetical and random deformations U

Material parameters for myofibroblasts
The material parameters were chosen similar to Buskermolen et al. [6], who simulated myofibroblasts with the same statistical framework, omitting the specific Helmholtz free-energy associated with cell adhesions and considering the contribution of the cell nucleus. The parameter differences between the two studies derive from these dissimilarities. The parameters for the cell adhesions were chosen by considering the assumption that k n k a . A list of the material parameters is given in Table S1.

Orientation, aspect ratio, and area of a configuration
For each configuration (j), as in experiments, the cellular orientation and aspect ratio were obtained by fitting the outline of the cell with the best-fitting ellipse (in a least square sense), calculated using a least-square algorithm. The area of the cell was calculated by summing the areas of the deformed elements composing the cell in the deformed configuration (j).

Scaling of the computational results
We first performed computational simulations of myofibroblasts cultured on substrates with alternating adhesive and non-adhesive lines equal in widths, where the size of these widths is normalized to the diameter of the undeformed cell. In what follows, we describe the procedure to find the scaled version of these widths. We first simulated the experiments for a large range of different line widths, spanning from lines that were small enough to enable the cell to spread on multiple lines, to lines wide enough to accommodate the cell on a single line. The computational results were presented with boxplots, with each boxplot corresponding to a different width of adhesive and non-adhesive lines. The value of the line widths is reported on the horizontal axis normalized to the diameter of the undeformed cell. Hereafter, we indicate this normalized value with l w .
The simulations predicted that cells generally tend to co-align with adhesive lines, with the variability of alignment dependent on l w . In particular, the alignment increased with l w when l w ≤ 0.3 and decreased with l w for l w ≥ 0.7. Interestingly, l w values between 0.4 and 0.6 corresponded to a change of trend, with alignment increasing from 0.3 to 0.5 and again decreasing for l w from 0.5 to 0.7.
To compare the computational results with the experiments, the normalized line widths l w have to be scaled with experimental measures. In our experiments, 20 µm seems a characteristic size for myofibroblasts, corresponding to the smallest line width where cells can spread on single lines (Figure 1 of the main text). Moreover, cells on these lines exhibited the maximum degree of alignment. Therefore, the normalized line widths were scaled by identifying the value of l w that corresponds to 20 µm. For brevity, we indicate this value with l ref w . To determine l ref w we observe that, in the computational simulations, cells start spreading on single lines when l w ≥ 0.6 ( Figure S4C). In addition, cells exhibit very high degrees of alignment for values of l w in between 0.6 and 0.7. Consequently, given the similarities with the experimental results for 20 µm lines, we can assume that 0.6 < l ref w ≤ 0.7. To make the scaling more specific, another requirement for l ref w can be deduced by observing the experimental results of cells on large lines (such as 100 µm and 150 µm); cells on such lines still exhibit a relatively high parallel alignment. However, in the computational simulations, large values of l w correspond to decreased alignment. Therefore, large experimental measures of line widths should correspond to relatively small l w , which translates to the requirement that l ref w should be relatively small. Considering these two requirements, as a first approximation, we assumed that 20 µm corresponds to l ref w = 0.625 (such that l w = 1 becomes equivalent to 32 µm). The computational simulations were then repeated for values of l w corresponding with the other line widths analyzed in the experiments. The results of these computational simulations are reported in Figure 5 of the main text. In Figure ??, morphologies of representative cells from the simulations and experiments are compared. In general, the cellular shapes obtained with the computational simulations are similar to the ones that were experimentally observed.      Material constant modulating the nonlinearity of the deviatoric elastic response of the cellular passive components k 1 GPa Parameter modulating the elastic penalty as a result of a large reduction of cell area J c 0.6 Penalty parameter associated with a large reduction of cell area σ max 240 kPa Maximum stress-fibre contractile stress ρ 0 3 · 10 6 µm −3 Density of stress-fibre proteins F 0 0.032 Stress-fibre protein volume fractioñ ε ss nom 0. 35 Steady state stress-fibre functional unit strain µ b0 − µ u 1kT Difference between bound and unbound stress-fibre potential Ω 10 −7.1 µm 3 Volume of the reference number of functional units that are within a stress-fibre in the undeformed representative volume element η max 1 Concentration of bound stress-fibre proteins b 0 /R 0 0.2 Ratio between the thickness and radius of an undeformed cell k n /(µR 0 ) 1.061 · 10 2 Nondimensional representation of the stiffness of linear springs approximating focal adhesions on non-adhesive lines k a /k n 10 2 Ratio between the stiffness of linear springs approximating focal adhesions on adhesive areas and the value for non-adhesive areas