Comparative proteomics of paired vocal fold and oral mucosa fibroblasts.

Injuries of the vocal folds frequently heal with scar formation, which can have lifelong detrimental impact on voice quality. Current treatments to prevent or resolve scars of the vocal fold mucosa are highly unsatisfactory. In contrast, the adjacent oral mucosa is mostly resistant to scarring. These differences in healing tendency might relate to distinct properties of the fibroblasts populating oral and vocal fold mucosae. We thus established the in vitro cultivation of paired, near-primary vocal fold fibroblasts (VFF) and oral mucosa fibroblasts (OMF) to perform a basic cellular characterization and comparative cellular proteomics. VFF were significantly larger than OMF, proliferated more slowly, and exhibited a sustained TGF-β1-induced elevation of pro-fibrotic interleukin 6. Cluster analysis of the proteomic data revealed distinct protein repertoires specific for VFF and OMF. Further, VFF displayed a broader protein spectrum, particularly a more sophisticated array of factors constituting and modifying the extracellular matrix. Conversely, subsets of OMF-enriched proteins were linked to cellular proliferation, nuclear events, and protection against oxidative stress. Altogether, this study supports the notion that fibroblasts sensitively adapt to the functional peculiarities of their respective anatomical location and presents several molecular targets for further investigation in the context of vocal fold wound healing. BIOLOGICAL SIGNIFICANCE: Mammalian vocal folds are a unique but delicate tissue. A considerable fraction of people is affected by voice problems, yet many of the underlying vocal fold pathologies are sparsely understood at the molecular level. One such pathology is vocal fold scarring - the tendency of vocal fold injuries to heal with scar formation -, which represents a clinical problem with highly suboptimal treatment modalities. This study employed proteomics to obtain comprehensive insight into the protein repertoire of vocal fold fibroblasts, which are the cells that predominantly synthesize the extracellular matrix in both physiological and pathophysiological conditions. Protein profiles were compared to paired fibroblasts from the oral mucosa, a neighboring tissue that is remarkably resistant to scarring. Bioinformatic analyses of the data revealed a number of pathways as well as single proteins (e.g. ECM-remodeling factors, transcription factors, enzymes) that were significantly different between the two fibroblast types. Thereby, this study has revealed novel interesting molecular targets which can be analyzed in the future for their impact on vocal fold wound healing.

Introduction

Modern society is strongly dependent on successful communication skills, making a healthy voice more important than ever before. Voice production is inextricably linked with physiological vibration of the vocal folds, which in turn is enabled by the finely tuned actions of laryngeal muscles, as well as by the sophisticated microarchitecture of the vocal fold mucosa (VFM) [1]. Detrimental changes of this delicate structure lead to pathological vibratory behavior and thus are often linked to voice problems, such as hoarseness and an unstable voice [2]. Indeed, the majority of patients seeking treatment for voice disorders present with organic vocal fold (VF) changes [3]. Scarring of the VF (due to chemical or physical traumata, e.g. surgery of the VF) is one prominent VF pathology that is characterized by a suboptimal reconstitution of the VF lamina propria by its main resident cells, the VF fibroblasts (VFF) [4]. Presence of a VF scar can be accompanied by a lifelong impairment of voice quality, yet current treatment modalities are largely unsatisfactory. Next-generation medical strategies could comprise transplantation of ex vivo generated VFM [5,6] or pharmacological approaches, i.e. the application of anti-fibrotic drugs [7,8] in order to resolve VF scars or even prevent their development. However, there is still a marked paucity of knowledge about VFM cells, which hampers the development of novel clinical applications.

In contrast to the VFM, the oral mucosa (OM) is endowed with a markedly scar-resistant phenotype. Studies that compared healing of oral wounds and scar-prone adult dermis have suggested several factors that contribute to superior healing of the OM [9]. These findings are largely in agreement with comparisons between adult dermis and fetal dermis, the latter being one of the most intensively studied tissues capable of scarless healing [10,11]. In general, wounds of scar-resistant tissues were found to contain reduced pro-inflammatory and pro-fibrotic mediators, lower blood vessel density, as well as fast-proliferating epithelial cells and fibroblasts [9,12]. Despite the evident risk of VF injuries to heal with scar formation, few of these and other healing-related features have been investigated in the VFM in a comparative manner. Experiments in rats have identified significant differences in the expression of the cytokines transforming growth factor beta 1 (TGF-β1) and transforming growth factor beta 3 (TGF-β3) in VFM versus OM and skin, which likely are related to the distinct healing tendencies of these tissues [13]. Further, a recent study has for the first time shed light on global transcriptional changes in the VFM after injury [14]. While the classical division of wound healing into the three phases inflammation, proliferation, and remodeling might in principle also apply to the VFs, the study's results suggest a considerable overlap of the cell division/proliferation and inflammatory activity phases. Further, previously undescribed mRNAs with dynamic expression during VF wound healing were identified [14], some of which could eventually become therapeutic targets.

An alternative way to better define the weak spots of VF wound healing is the direct comparison of VF cells to cells arising from tissues with superior healing properties, which is the focus of the present study. Although VF cells have been benchmarked against other cells with respect to small sets of genes/proteins of interest (e.g. related to the extracellular matrix (ECM) [15] or surface antigens [16,17]), it should be noted that, to date, no comprehensive comparison between a distinct VF cell type and its counterpart in other tissues has been performed. To close this gap, we first established cell cultures of VFM and OM to obtain paired, near-primary populations of the fibroblasts homing to the respective tissue. Second, we performed a basic cellular characterization of the fibroblast populations, covering the aspects of marker gene expression, size, morphology, proliferation, migration, 3D culture characteristics and response to a pro-fibrotic stimulus. Third, fibroblast pairs were subjected to label-free quantitative proteomics to provide a comprehensive comparison between cells of the VFM and the OM, including significantly overrepresented gene sets and corresponding pathways.

Materials and methods

Standard 2D cell culture

Paired tissue samples from sheep VFM and OM (Vestibulum oris) were obtained from the local slaughterhouse (approximately 1 h post mortem) and transported to the laboratory in Dulbecco's Modified Eagle's Medium (DMEM; Life Technologies, Carlsbad, CA) supplemented with 1× l-Glutamine–Penicillin–Streptomycin solution (GPS; Sigma-Aldrich, St. Louis, MO). Subsequent steps were performed in a laminar flow hood. For decontamination, tissues were incubated in DMEM supplemented with 10× GPS at room temperature for 10 min, followed by dipping in 70% Ethanol for 10 s and extensive rinsing in Phosphate buffered saline (PBS, Life Technologies). VF and OM tissues were then pinned on a sterile cork plate with the epithelium downwards. This enabled the preparation of small tissue pieces (approximately 2 mm × 2 mm × 2 mm) consisting mainly of connective tissue and only minimal residual epithelium. Tissue pieces were placed in 12 well plates (Greiner Bio-One, Kremsmünster, Austria) and allowed to adhere by 5–10 min of incubation (without lid). Finally, tissue pieces were carefully overlaid with serum-reduced growth medium (GM), consisting of DMEM/F12 Nutrient Mix (Life Technologies) (3:1), 5% Fetal Bovine Serum (FBS; Biochrom, Berlin, Germany), 1× Insulin-Transferrin-Selenium (Life Technologies), 200 nM 3,3′,5-Triiodo-l-thyronine, 100 ng/mL human Insulin-like Growth Factor-I, 12.5 ng/mL human Fibroblast growth factor-2 (all Sigma), and 100 μg/mL Normocin (Invivogen, San Diego, CA). Tissue pieces were cultivated in a humidified atmosphere at 37 °C and 5% CO2 and GM was partially (50% v/v) replaced by fresh GM every 2–4 days. After robust outgrowth of fibroblastoid cells from adherent tissue pieces (but before the appearance of a confluent cell layer) the first passaging was performed using 0.5% trypsin (Sigma). Cells were subsequently passaged every 3–4 days at passaging ratios between 1:2 and 1:8. For downstream analyses, seeding of cells into culture plates at defined cell densities was based on manual counting of trypsinized cell suspensions using a hemocytometer and subsequent calculation of cell concentrations.

Real time cell growth assay

Cell proliferation was measured using the xCELLigence (OLS OMNI Life Science Bremen, Germany) Real-Time Cellular Analysis system. Briefly, the background impedance was measured following the addition of 100 μL of GM to the 16-well E-plates (OLS). Cell suspensions containing 5000 cells in 100 μL GM were seeded into the wells and attachment and proliferation were monitored using the xCELLigence Real-Time Cellular Analysis system. Impedance changes expressed as cell index (CI) were automatically calculated and correlated with cell growth over a period of 100 h.

3D spheroid cell culture

20,000 cells resuspended in GM were seeded in hanging drops (20 μL) under the lids (Lactan, Graz, Austria) of cell culture plates. The dish was carefully inverted to keep the drops intact. 10 mL PBS were added to the bottom of the dish to prevent dehydration. Drops were cultured at 37 °C in 5% CO2 without medium change for up to three days. To investigate attachment and cell spreading, spheroids were subsequently pipetted into 12 well plates with 1 mL GM per well.

2D migration assay

Pairs of VFF and oral mucosa fibroblasts (OMF) with similar passage numbers were seeded into 2-chamber silicone inserts which were placed in a 24 well plate (ibidi, Martinsried, Germany). Depending on preliminary experiments, 20,000–30,000 cells resuspended in 70 μL GM were pipetted into each chamber and incubated overnight, which resulted in confluent cell layers. Subsequently, silicone inserts were removed with sterile forceps, cells were washed 3 times with PBS and 1 mL of standard medium (SM), consisting of DMEM/F12 Nutrient mix (3:1) supplemented with 5% FBS and 100 μg/mL Normocin, was added per well. Closure of the 500 μm gap between the 2 regions of cells was monitored using a Cell-IQ system (Chipman Technologies, Tampere, Finland) with automated phase contrast imaging (4–6 regions per well) performed every 30 min over a period of 72 h. Semi-automatic segmentation of the image series and quantitative analysis was performed using the Cell-IQ software. For each cell type, at least 3 wells (assayed in parallel) were averaged as technical replicates.

Reverse transcription quantitative real-time PCR (RT-qPCR)

Cells were seeded in GM into 12 well plates (Greiner) at 5000–10,000 cells/cm2. After 4 days, the medium was changed to SM. The next day, medium was again changed to SM ± 5 ng/mL TGF-β1 (R&D Systems, Minneapolis, MN). Cells were lysed in 700 μL QIAzol after 24 h and 72 h. RNA isolation was performed using the miRNeasy Mini kit (Qiagen, Hilden, Germany) according to the manufacturer's protocol. Concentrations of RNA samples were determined on a NanoDrop 2000c UV–Vis spectrophotometer (Thermo Scientific, Waltham, MA) and 1 μg total RNA was used as input for cDNA synthesis using the QuantiTect Reverse Transcription kit (Qiagen) according to the manufacturer's instructions. For RT-qPCR, cDNA samples were diluted with nuclease-free water to 2.5 ng RNA/μL, and 4 μL of diluted solutions were combined with 5 μL GoTaq® qPCR Master Mix (Promega, Mannheim, Germany) and 1 μL of the respective primer pair mix (each primer: 2 μM) in 384 well plates. Sequences of primers are listed in Supplemental Table 1 in the supporting information (SI). Each combination of cDNA and primer pair was assayed in technical triplicates. RT-qPCR was performed on a LightCycler® 480 system (Roche, Vienna, Austria) using the following program: 2 min/95 °C (denaturation), 45 cycles of 10 s/95 °C and 1 min/60 °C (amplification), ramping at 2.5 °C/min from 55 °C to 95 °C (melting curve analysis). CT values were calculated using the AbsQuant/2nd Derivative Max method of the LightCycler® 480 software. Relative quantification of expression levels of the genes of interest was further performed according to the 2−ΔΔCT method [18] using UXT as internal reference RNA.

LC-MS/MS analysis

Pairs of VFF and OMF with similar passage numbers were seeded into T75 cell culture flasks (Greiner). At ~50% optical confluence, old medium was removed and cells were washed twice with PBS before addition of fresh DMEM/F12 Nutrient Mix (3:1). After incubation for 6 h, cells were trypsinized and the resulting cell suspension was centrifuged at 170 ×g for 7 min. Supernatants were discarded and cell pellets were resuspended in PBS, followed by determination of cell concentration and centrifugation as before. Supernatants were thoroughly removed, cell pellets were snap frozen in liquid nitrogen and stored at − 80 °C. One million cells were lysed in 100 μL lysis buffer (100 mM Tris pH 8.0, 10 mM TCEP, 40 mM 2-chloroacetamide and 1% SDS) by sonication for 2 min on ice. The lysates were heat denatured for 10 min at 90 °C. The supernatants were collected after removal of cell debris by centrifugation at 1000 ×g for 10 min at room temperature. Protein was determined by BCA-RAC assay (Thermo Scientific). Protein patterns a were visualized by SDS-PAGE on 4–12% NuPAGE Bis-Tris gels (Thermo Fisher Scientific), Krypton (Thermo Fisher Scientific) staining and fluorescence laser scanning on a Bio-Rad FX Pro Plus Imager. For LC-MS/MS analysis 20 μg protein were subjected to filter aided sample preparation as described previously [19] with minor modifications using 3 kD cut off ultracentrifugation filters (Millipore). Buffer was exchanged three times with 8 M Urea, 100 mM TrisHCl, pH 8.0. The samples were predigested with rLysC (Promega, enzyme/protein 1:100) for 4 h at 37 °C, diluted to 2 M Urea with 100 mM ammonium bicarbonate and digested with trypsin (Promega, enzyme/protein 1:50) overnight. Peptides were collected by centrifugation at 14000 ×g for 45 min at room temperature. Samples were acidified with formic acid (final concentration of 0.1%). 500 ng protein digest of each sample yielding similar total ion chromatograms (cf. Supplemental Fig. S1 in the SI) were injected and analyzed by nano-HPLC (Dionex Ultimate 3000) equipped with a C18, 5 μm, 100 Å, 5 × 0.3 mm, enrichment column and an Acclaim PepMap RSLC nanocolumn (C18, 2 μm, 100 Å, 500 × 0.075 mm) (all Thermo Fisher Scientific). Samples were concentrated on the enrichment column for 2 min at a flow rate of 5 μL/min with 0.1% heptafluorobutyric acid as isocratic solvent. Separation was carried out on the nanocolumn at a flow rate of 300 nL/min at 60 °C using the following gradient, where solvent A is 0.1% formic acid in water and solvent B is acetonitrile containing 0.1% formic acid: 0–6 min: 4% B; 6–150 min: 4–25% B; 150–155 min: 25–95% B, 155–165 min: 95% B; 165–165.1 min: 4% B; 165.1–180 min: 4% B. The maXis II ETD mass spectrometer (Bruker, Vienna, Austria) was operated with the captive source in positive mode with following settings: peptide fragmentation: CID, mass range: 200–2000 m/z, 2 Hz, capillary 1300 V, dry gas flow 3 L/min with 150 °C, nanoBooster 0.2 bar, precursor acquisition control top17.

Processing of LC-MS/MS data

The LC-MS/MS data were analyzed by the Data analysis software (Bruker), using the Sum Peak algorithm, and by MaxQuant by searching the public Uniprot database with taxonomy Ovis aries (downloaded on 13.06.2016, 13,148,561 residues, 27,448 sequences) and common contaminants. Carbamidomethylation on Cys was entered as fixed modification, oxidation on methionine as variable modification. Detailed search criteria were used as follows: trypsin, max. missed cleavage sites: 2; search mode: MS/MS ion search with decoy database search included; precursor mass tolerance ± 0.006 Da; product mass tolerance ± 40 ppm; acceptance parameters for identification: 1% PSM false discovery rate (FDR); 1% protein FDR. In addition a label free quantitation was performed using MaxQuant [20] requiring a minimum of 2 ratio counts of quantified razor and unique peptides. Total LFQ intensities of individual samples were as follows: VFF.BR1: 8.615.297.702; OMF.BR1: 9.078.542.438; VFF.BR2: 8.810.529.124; OMF.BR2: 8.819.990.495; VFF.BR3: 8.697.640.277; OMF.BR3: 9.160.585.501; VFF.BR4: 8.713.845.133; OMF.BR4: 8.890.439.610. The mass spectrometry proteomics data have been deposited to the ProteomeXchange Consortium via the PRIDE [21] partner repository with the dataset identifier PXD005106. Further processing, data visualization, and statistical testing were performed with R (www.r-project.org). Briefly, identified common contaminant peptides were removed, and UniProtKB AC/IDs of the remaining hits were mapped to their respective gene symbols using the “Retrieve/ID mapping” function of www.uniprot.org. Spectra that were associated to different gene symbols were reduced to a single gene symbol based on the following criteria: (i) highest number of peptide counts (razor + unique); (ii) highest sequence coverage. Hierarchical Clustering was performed using the “complete linkage” method (hclust function of the R stats package). Gene Ontology (GO) term over-representation analysis was performed using the “Gene List Analysis” tool [22] from www.pantherdb.org with default settings. For proteins identified in both cell types, log2-transformed fold changes of label-free quantification (LFQ) intensity values were calculated for each VFF-OMF pair (pFC). Proteins were filtered for differential abundance (|mean(pFC)| > 0.58) and a paired t-test was conducted for hits where 3–4 pFC could be calculated, followed by applying the Benjamini-Hochberg correction for multiple testing. Gene Set Enrichment Analysis (GSEA) [23,24] was conducted using the javaGSEA desktop application with preset values of gene set sizes. Gene sets with a FDR-corrected p-value (FDR-q) < 0.1 were considered as significantly enriched.

Results and discussion

Establishment of paired VFF and OMF populations from sheep

As every surgical procedure on the VF bears the risk of scar development and thus lifelong voice disorder (dysphonia), obtaining non-pathological VFM samples from humans is rarely feasible. To circumvent this problem, we used sheep as a source of VFM as well as OM. This animal model has previously been demonstrated to be comparable to the human VF with respect to anatomical dimensions, phonation indices, and tissue histochemistry [25]. Outgrowth of fibroblastoid cells from VFM and OM tissue pieces started two to seven days after start of incubation in GM (Supplemental Fig. S2A in the SI). Occasionally, these primary cultures also contained cells with epithelial or endothelial morphology (Supplemental Fig. S2BC in the SI); however, these cells disappeared after the first passages. Further, expression of marker genes for skeletal muscle (ACTN3), epithelium (CDH1), and endothelial cells (VWF) was either close to background, or undetectable at passages 4–6 (Supplemental Fig. S2D in the SI). Thus, sufficiently pure populations of paired VFF and OMF could be successfully established from the corresponding sheep tissues. These near-primary cells were used for subsequent analyses.

Basic cellular characterization of paired VFF and OMF

Under standard cell culture conditions, VFF had a more elongated shape and were significantly larger than OMF (Fig. 1AB). Proliferation assays revealed that OMF exhibited faster growth than VFF (Fig. 1C), which was in line (i) with previous studies comparing human adult dermal fibroblasts and OMF [26]; and (ii) with the fact that passaging of OMF had to be performed at higher split ratios to prevent the formation of a confluent cell layer. Both cell types could be cultivated in a three-dimensional setting, yet with morphological differences: while VFF-derived spheroids appeared rather encapsulated, spheroid boundaries were less defined for OMF-derived constructs. Correspondingly, OMF started to spread from the spheroid earlier than VFF (Fig. 1D; Supplemental Fig. S3 in the SI). Analysis of cell migration in a 2D migration assay revealed that two OMF populations had faster migratory characteristics than the respective VFF; however no cell-specific differences in migration could be measured for the other biological replicates (Supplemental Fig. S4 in the SI). TGF-β1 is a well-known pro-inflammatory cytokine that has been identified in various wounds (including VF wounds) [27] and that is also implicated in fibrosis in many tissues, e.g. via promoting myofibroblast differentiation as well as the expression of other pro-inflammatory mediators [28,29]. Twenty-four hours after exposure to TGF-β1 we identified a significant increase in the levels of the myofibroblast marker alpha smooth muscle actin (α-SMA, ACTA2) in OMF and a similar trend in VFF (Fig. 2A, left panels). For both cell types, significantly higher α-SMA mRNA levels were measured 72 h after TGF-β1 treatment (Fig. 2A, right panels). Thus, although we could confirm the classical induction of α-SMA by TGF-β1 for our in vitro models, the results are distinct from previous reports describing a decreased responsiveness for OMF [30]. Interleukin 6 (IL6) has been identified as a cytokine favoring scar development [31]. While TGF-β1 promoted a long-lasting significant increase of IL6 in VFF, this effect was not observed for OMF, where IL6 expression had returned to baseline levels after 72 h (Fig. 2B). It can therefore be hypothesized that successful VF wound healing might be hampered at least partly due to the prolonged IL6 production of VFF, as opposed to a more short-lived response in the OM. Hyaluronic acid (HA) is an important component of the physiological VF ECM [32,33]. We therefore analyzed the expression of hyaluronan synthases, the key enzymes in HA biosynthesis. We found HAS1-3 to be expressed in all cell populations, yet only HAS1 mRNA to be significantly induced in VFF 24 h after TGF-β1 exposure (Supplemental Fig. S5 in the SI). Collectively, the basic cellular characterization of our newly established fibroblast populations revealed several phenotypical differences prompting a more comprehensive comparison of their repertoire of expressed proteins.

Fig. 1

Comparative cellular characteristics of paired vocal fold fibroblasts (VFF) and oral mucosa fibroblasts (OMF). (A) Representative phase contrast microscopy images from two biological replicates (BR1, BR2); scale bar: 100 μm. (B) Mean cell diameter ± s.e.m. (n = 4) was obtained from trypsinized (i.e. spherical) cells. Data was analyzed by paired t-test; *… P ≤ 0.05. (C) Growth curves of fibroblast pairs. (D) Phase contrast microscopy images of 3D spheroid cultures (BR2) were acquired 6 h and 24 h post re-plating; scale bar: 100 μm.

Fig. 2

Cellular response to pro-fibrotic stimulus. Paired vocal fold fibroblasts (VFF) and oral mucosa fibroblasts (OMF; n = 4) were exposed to TGF-β1 (5 ng/mL) or vehicle (“C”) and expression of (A) α-smooth muscle actin (ACTA2) and (B) interleukin 6 (IL6) was analyzed 24 h and 72 h later. TGF-β1 effects were analyzed by paired t-test.

Comparison of VFF and OMF proteomes

As a first qualitative analysis, we separated VFF- and OMF-derived protein lysates by SDS-PAGE. Subsequent protein staining revealed highly similar band patterns within the four biological replicates of each cell type (Supplemental Fig. S6 in the SI). We therefore proceeded with trypsin digestion of the samples and subjected the resulting peptides to LC-MS/MS analysis. Altogether, 1609 unique proteins with association to a known gene were confidently identified (Supplemental Table S2 in the SI; cf. Materials and Methods). Of these, 53% (858 proteins) were detected in both VFF and OMF, while 35% (564 proteins) were found only in biological replicates of VFF and 12% (187 proteins) were identified only in biological replicates of OMF (Fig. 3A). In line with these overall values, the number of detected proteins was higher in every VFF biological replicate compared to the respective OMF pair (Fig. 3B). To further estimate congruencies in the protein repertoire, we performed pairwise comparisons of samples. The resultant average “overlap index” (i.e. the fraction of proteins detected in both samples, normalized to the total number of proteins detected in the respective combination) was 0.571 for comparisons across the two distinct cell types; while it was slightly higher (0.586 for VFF; 0.623 for OMF) for comparisons within the same cell type (Supplemental Fig. S7 in the SI). In order to obtain subsets of high-confidence hits, we filtered the total data for proteins that were repeatedly detected in independent samples. Interestingly, unsupervised clustering of the resulting data clearly distinguished the two fibroblast types (Fig. 3C). Collectively, these results revealed that protein profiles exist which are specific for VFF and OMF, respectively. This is in line with a previous study that found site-specific gene expression profiles of human adult dermal fibroblasts [34]. It also supports the notion that, although basic traits of distinct fibroblast populations might be similar, adaptations of the “fibroblast phenotype” meeting the functional peculiarities of the respective anatomical location exist.

Fig. 3

Proteomic profiling of paired vocal fold fibroblasts (VFF) and oral mucosa fibroblasts (OMF). (A) Venn diagram depicting the number of proteins detected solely in VFF, solely in OMF, and in both cell types. (B) Number of robustly detected proteins (≥2 peptide counts (razor + unique)) for each biological replicate (BR1–BR4). (C) Heat map showing log10-transformed normalized LFQ intensity values for three groups: group 1, 520 proteins detected in both cell types; group 2, 45 proteins detected only in VFF, group 3, 9 proteins detected only in OMF. Within each group, proteins are sorted in a decreasing manner according to their mean abundance. Clustering tree resulting from unsupervised clustering of samples is shown on top.

Proteins with differential abundance

We next focused on the subset of 45 proteins that were reliably (i.e., present in 3–4 biological replicates) detected in VFF, yet never in OMF (“VFF-specific”; Fig. 4A). In line with cells being the starting material of our analysis, the general GO terms “cell”, “cell part”, “intracellular”, “intracellular part”, and “cytoplasmic part” were significantly overrepresented. Interestingly, also GO terms relating to the extracellular space and vesicles exhibited significant overrepresentation; the subset of corresponding proteins included the alpha 2 chain of fibrillar collagen type I (COL1A2), the alpha 1 chain of the fibril-associated collagen type 12 (COL12A1), enzymes responsible for collagen maturation (prolyl 4-hydroxylase subunit alpha 1 (P4HA1), procollagen-lysine,2-oxoglutarate 5-dioxygenase 2 (PLOD2)), as well as proteins involved in ER-Golgi and vesicular trafficking (SEC23A, SEC23IP, FLOT2, ATP6V1E1). Further, also the GO terms “anchoring junction” and “adherens junction” were significantly overrepresented due to the presence of several actin-binding and -remodeling proteins (ITGA6, CNN1, DOCK7, SWAP70). Thymosin beta 4 (THYMB4X) is a critical regulator of F-actin assembly [35] and was the most abundant VFF-specific protein in our screen. THYMB4X is also present in blood and wound fluid [36], and has been described to block expression of pro-fibrotic genes [37] and accelerate dermal wound healing [38]. With respect to the scarless healing properties of OM, it was therefore surprising that we did not detect synthesis of THYMB4X by OMF. As for VF wounds, THYMB4X mRNA showed a non-significant increase in the early phase of wound healing in rat (analysis of publicly available microarray data (Gene Expression Omnibus; GSE62204) from [14]). However, studies investigating a potential anti-fibrotic action of THYMB4X in VF wound healing have not been performed to date. Unexpectedly, we identified fibroblast growth factor 2 (FGF2), well-known for its mitogenic and angiogenic activities [39,40], as another VFF-specific protein. FGF2 accelerates wound healing [41] and has recently been shown to ameliorate VF scarring [42], possibly via reshaping the synthesis of ECM components by VFF [43]. Our result thereby implies that (i) the faster growth of OMF might be due to other factors than FGF2, and that (ii) the amount of FGF2 produced by VFF in vivo might be too limited to enable scarless wound healing of the VF.

Fig. 4

Cell type-selective protein sets and their significantly enriched Gene Ontology terms. Groups of proteins detected solely in (A) vocal fold fibroblasts (VFF) or (B) oral mucosa fibroblasts (OMF) are depicted as heat map with proteins specified by gene symbol. Right part of panels depicts results of statistical overrepresentation tests performed for the three main Gene Ontology (GO) aspects molecular function, cellular component, and biological process. Grey fields indicate significant (P < 0.05 after Bonferroni correction for multiple testing) overrepresentation of the respective GO term (specified at the bottom).

The converse subset of proteins that were reliably detected in OMF, but never in VFF (“OMF-specific”) was considerably smaller (9 proteins) and thus, the only significantly enriched GO term was “protein kinase B binding” (Fig. 4B). The corresponding proteins were Girdin (CCDC88A), which interacts with actin, is activated by Akt/Protein kinase B (PKB) and promotes DNA synthesis and cell motility [44,45], as well as adaptor protein containing pleckstrin homology domain, phosphotyrosine binding domain and leucine zipper motif (APPL1), a multidomain protein that mediates insulin as well as adiponectin signaling [46]. Notably, the beneficial effects of several growth factors and small molecule pharmacological agents in dermal wound healing (or specifically on dermal fibroblasts) have been linked to activation of Akt/PKB signaling [47-50]. In addition, protective effects of adiponectin have recently been demonstrated in lung fibrosis and might be mediated at least partly by direct action on the resident fibroblasts [51]. In conjunction with our results, these findings imply the future evaluation of pharmacological approaches involving Akt/PKB activation and/or adiponectin in VF wounds.

Among the large subset of proteins that were detected in both cell types (Fig. 3C), we aimed to select those that exhibited a significantly differential abundance (cf. Experimental Procedures). Twenty such proteins were identified (Fig. 5), corresponding to ~4% of the total of proteins that were reliably detected in both cell types; this proportion was similar as in a recent study that compared gingival and dermal fibroblasts [52]. Candidates classified as VFF-enriched comprised again actin-binding proteins such as nexilin F-actin binding protein (NEXN), actinin alpha 1(ACTN1), and cysteine and glycine rich protein 1 (CSRP1), which has also been described as a critical transcriptional coregulator in smooth muscle cell differentiation [53]. Conversely, candidates classified as OMF-enriched comprised several proteins with an established protective function against oxidative stress (glucose-6-phosphate dehydrogenase (G6PD), peroxiredoxin (PRDX1), transaldolase 1 (TALDO1)) and DNA damage (replication protein 1 (RPA1)). It can therefore be speculated that the superior wound healing of the OM might be partially due to a superior ensemble of stress-protective mechanisms in OMF. Indeed, the premature senescence of dermal fibroblasts from chronic wounds has been proposed to be due to a reduced ability to cope with oxidative stress [54]. Conversely, favorable effects of antioxidative agents have already been demonstrated for lung fibroblasts [55]. As reactive oxygen species are also known to be heavily induced in the early phase of VF wound healing [56], assessing the potential of antioxidative agents in VF injury is a promising future prospect. In analogy to the analysis of VFF- and OMF-specific protein sets, we also subjected the significantly VFF- or OMF-enriched protein sets to overrepresentation analysis. While no GO term was significantly overrepresented within the VFF-enriched set, the GO term “cytoskeleton” was significantly overrepresented within the OMF-enriched set. The subset of corresponding proteins included the two non-motor actin binding proteins Plectin (PLEC) and Coronin-1C (CORO1C). As these proteins also interact with intermediate filaments and microtubules [57,58], this result might reflect the distinct mechanical stresses implied on fibroblasts in OM versus VFM.

Fig. 5

Vulcano plot. Candidate proteins with significantly differential protein abundance (P < 0.05 after Benjamini-Hochberg correction for multiple testing) between paired vocal fold fibroblasts (VFF) and oral mucosa fibroblasts (OMF) are shown in red and specified by gene symbol.

Gene set enrichment analysis

We eventually subjected the entire dataset of reliably detected proteins to GSEA in order to identify biological categories that are most distinct between the two investigated fibroblast types. In total, eight gene sets were found to be significantly VFF-enriched (Table 1). Interestingly, these included gene sets of the “matrisome”, which was recently defined (based on a combination of in silico and experimental approaches) as the ensemble of architectural ECM proteins and associated proteins [59]. Further, Kyoto Encyclopedia of Genes and Genomes (KEGG) pathways related to cell-cell (tight junction) and cell-ECM (focal adhesion) interactions were also identified, as was a hallmark gene set defining epithelial to mesenchymal transition. Collectively, these results likely reflect the necessity of VFF to adapt to an environment of exceptional mechanical forces, i.e. the periodic vibration due to phonation events. It may be further postulated that the larger repertoire of ECM proteins produced by VFF is a prerequisite for the sophisticated structure of the VF lamina propria, as opposed to the less complex connective tissue layer of the OM. Such an inherent disproportion in the fibroblasts' ECM-synthesizing abilities could partially explain the suboptimal results that have been achieved in the past when using autologous oral mucosa grafts as a surrogate structure for resected VFM [2,60]. Procollagen C-endopeptidase enhancer (PCOLCE) was found as one of the VFF-enriched proteins contributing to the “matrisome” signature. As PCOLCE promotes the maturation of fibrillar procollagens into mature collagen fibers [61,62], future studies in which the function of this enzyme is antagonized could be of particular interest in the context of VF wound healing. Likewise, antagonism of PLOD2, a “matrisome” component that was identified as VFF-specific (Fig. 4A), might be a promising therapeutic approach, as the collagen cross-linking action of this enzyme is a key feature of fibrosis [63]. GSEA further revealed that potential direct targets of the miR-30 family of microRNAs were significantly enriched in the VFF proteome. We therefore queried miRTarBase, a compendium of experimentally validated miRNA-mRNA interactions [64], and found that indeed 9 of the 15 proteins constituting the respective core enrichment (Table 1) have been validated as direct miR-30 targets. It would be interesting to analyze in the future whether miR-30 family members are indeed more abundant in OMF versus VFF, and if so, whether increasing miR-30 microRNAs in VFF might convey anti-fibrotic properties, as recently described for hepatic stellate cells [65] as well as for an in vitro model of kidney fibrosis [66].

Table 1

Gene sets significantly enriched in vocal fold fibroblasts versus oral mucosa fibroblasts.

Gene sets database	Gene sets significantly enriched in VFF versus OMF	FDR q-value	Gene symbols of proteins within core enrichment
Canonical pathways	NABA_MATRISOME	0.0633	P4HA1, PLOD2, FGF2, COL1A2, P4HA2, S100A2, COL6A1, SERPINH1, P3H1, ANXA3, COL12A1, PCOLCE, COL15A1, COL6A2
	NABA_MATRISOME_ASSOCIATED	0.0689	P4HA1, PLOD2, FGF2, P4HA2, S100A2, SERPINH1, P3H1, ANXA
KEGG	KEGG_REGULATION_OF_ACTIN_CYTOSKELETON	0.0440	MYH10, FGF2, ACTN1, ITGA6, MYH9, PPP1R12A, VCL, ACTN4, ARPC3, RDX, ROCK1, CRK, IQGAP1, ROCK2, ACTB, ARPC4
	KEGG_FOCAL_ADHESION	0.0507	COL1A2, ACTN1, ITGA6, COL6A1, PPP1R12A, VCL, FLNB, ACTN4, PARVA, COL6A2, ROCK1, FLNA, CRK, COL6A3, TLN1, ROCK2, ACTB
	KEGG_TIGHT_JUNCTION	0.0657	MYH10, ACTN1, TJP1, EPB41L2, PPP2CB, MYH9, ACTN4
Reactome	–	–	–
Hallmark	HALLMARK_APICAL_JUNCTION	0.0076	MYH10, NEXN, ACTN1, TJP1, EPB41L2, MYH9, ACTA1, CNN2, VCL, ACTN4, PARVA
	HALLMARK_EPITHELIAL_MESENCHYMAL_TRANSITION	0.0084	PLOD2, FGF2, COL1A2, TPM1, FERMT2, SERPINH1, TAGLN, DPYSL3, CALU, COL12A1, PCOLCE, COL6A2, FLNA, GJA1
Transcription factor targets	–	–	–
miRNA binding	TGTTTAC,MIR-30A-5P,MIR-30C,MIR-30D,MIR-30B,MIR-30E-5P	0.0806	PTGFRN, MYH10, VAT1, P4HA2, DOCK7, ITGA6, SEC23A, PPP1R12A, HSPA5, CALU, LPP, GFPT2, LIN7C, VAT1L, GJA1
GO biological process	–	–	–
GO cellular component	–	–	–
GO molecular function	–	–	–

Distinct gene set databases (column 1) were queried to identify gene sets (column 2) which are significantly (false discovery rate (FDR) q-value <0.1, column 3) VFF-enriched when comparing the protein expression profiles of VFF with OMF. Proteins (specified by gene symbol) detected in the VFF/OMF samples and having the strongest contribution to the enrichment result (“core enrichment”, as defined by the Gene Set Enrichment Analysis (GSEA) method) are listed in column 4.

GSEA also revealed that 20 gene sets were significantly enriched in OMF versus VFF (Table 2), a conspicuously larger number than the eight gene sets identified as VFF-enriched. This discrepancy is surprising given the clearly smaller repertoire of proteins expressed by OMF (Fig. 3). A possible explanation could be that molecular characterization of cells from the VF has been performed in less depth so far. Notably, a considerable fraction of OMF-enriched gene sets was related to the nucleus and biological processes within this compartment, most importantly processing of RNA. This finding is in line with the significant upregulation of two proteins constituting the nuclear lamina, lamin A/C (LMNA) and lamin B1 (LMNB1; Fig. 5). However, it is also partly surprising as e.g. upregulation of mRNA splicing would intuitively suggest a greater divergence in the OMF protein fingerprint, which is not the case. Alternatively, it can be hypothesized that OMF express a more sophisticated array of non-coding RNAs, which are likewise often subjected to ribonucleoprotein-dependent processing/maturation steps. Future comparative studies employing state-of-the-art RNA sequencing techniques should shed light on this issue. Another interesting finding was the significant enrichment of the GO term “cell proliferation”, which is in line with the faster growth characteristics of OMF (Fig. 1C). This result was corroborated by the significant enrichment of two hallmark gene sets categorized as targets of the transcription factors MYC and E2F, both implicated in cell cycle progression [67]. In the case of E2F, further evidence is provided as several OMF-enriched proteins harbor consensus binding sites for the E2F transcription factor family in their promoter-associated regions. It is thus tempting to speculate that the OMF's fast-proliferating nature is at least partly founded in transcriptional activities of distinct E2F isoforms, as previously described for other fibroblast types [68,69].

Table 2

Gene sets significantly enriched in oral mucosa fibroblasts versus vocal fold fibroblasts.

Gene sets database	Gene sets significantly (FDR q-value <0.1) enriched in OMF versus VFF	FDR q-value	Gene symbols of proteins within core enrichment
Canonical pathways	–	–	–
KEGG	KEGG_SPLICEOSOME	0.02613	DDX39B, PUF60, DDX5, PRPF19, SRSF6, HNRNPU, HNRNPM, SRSF1, SNRNP70, HNRNPK, HNRNPC, SF3B1, PCBP1, SF3B2
Reactome	REACTOME_PROCESSING_OF_CAPPED_INTRON_CONTAINING_PRE_MRNA	0.01592	DHX9, SRSF6, HNRNPU, HNRNPM, HNRNPL, SRSF1, SNRNP70, HNRNPK, HNRNPC, SF3B1, PCBP1, TPR, SF3B2, PTBP1
	REACTOME_MRNA_SPLICING	0.03449	DHX9, SRSF6, HNRNPU, HNRNPM, HNRNPL, SRSF1, SNRNP70, HNRNPK, HNRNPC, SF3B1, PCBP1, SF3B2, PTBP1
	REACTOME_MRNA_PROCESSING	0.01884	DHX9, SRSF6, HNRNPU, HNRNPM, HNRNPL, SRSF1, SNRNP70, HNRNPK, HNRNPC, SF3B1, PCBP1, TPR, SF3B2, PTBP1
	REACTOME_METABOLISM_OF_LIPIDS_AND_LIPOPROTEINS	0.08798	ACLY, HMGCS1, ME1, DHCR7, LBR
Hallmark	HALLMARK_E2F_TARGETS	<10-E4	CSE1L, PA2G4, SYNCRIP, RANBP1, SRSF1, PCNA, STMN1, PAICS, RBBP7, SMC3, MCM3, LMNB1, MCM6, H2AFX, RPA1, LBR
	HALLMARK_G2M_CHECKPOINT	<10-E4	SFPQ, HNRNPU, NCL, SYNCRIP, SRSF1, BUB3, STMN1, NUMA1, G3BP1, SMC2, MCM3, LMNB1, MCM6, H2AFX, LBR
	HALLMARK_MYC_TARGETS_V1	0.05351	HSP90AB1, CCT4, PSMC6, YWHAE, RRM1, CCT7, PSMA4, CCT3, RAN, RPS5, KPNB1, PABPC4, FAM120A, RUVBL2, IMPDH2, DUT, PSMD3, HNRNPU, CCT5, TARDBP, PA2G4, PSMA1, SYNCRIP, RANBP1, SRSF1, SERBP1, PCNA, HNRNPC, BUB3, CBX3, CNBP, PCBP1, G3BP1, MCM6
Transcription factor targets	V$E2F1_Q4_01	0.00835	KPNB1, RPS19, HIST1H1D, IMPDH2, NCL, RANBP1, SRSF1, SERBP1, PCNA, STMN1, CNBP, SMC3, SMC2, MCM3, DDX17, PTMA, MCM6, UCHL1
	V$E2F_Q3_01	0.03749	KPNB1, RPS19, HIST1H1D, NCL, RANBP1, SRSF1, SERBP1, PCNA, STMN1, RALY, SMC3, NUMA1, SMC2, MCM3, DDX17, PTMA, MCM6, UCHL1
	V$E2F_Q4_01	0.03923	KPNB1, RPS19, HIST1H1D, IMPDH2, NCL, RANBP1, SRSF1, SERBP1, PCNA, STMN1, CNBP, SMC3, SMC2, MCM3, DDX17, PTMA, MCM6, UCHL1
	V$USF_01	0.07438	SAE1, RANBP1, STMN1, NUDC, EIF4B, SMC3, PTMA, RPA1
miRNA binding	–	–	–
GO biological process	CELL_PROLIFERATION_GO_0008283	0.02559	FSCN1, CSE1L, PA2G4, PCNA, BUB3, NUDC, RBBP7, KHDRBS1, CCDC88A, APPL1, PRDX1
	RNA_SPLICING	0.03236	SFPQ, SRSF6, KHSRP, SYNCRIP, SRSF1, SNRNP70, HNRNPC, NONO, IVNS1ABP, SF3B2, PTBP1
GO cellular component	NUCLEAR_PART	<10-E4	HNRNPU, PSMA1, HNRNPM, HNRNPL, SNRNP70, HNRNPK, HNRNPC, RALY, IVNS1ABP, SF3B2, PTBP1
	RIBONUCLEOPROTEIN_COMPLEX	0.00929	HNRNPU, PSMA1, HNRNPM, HNRNPL, SNRNP70, HNRNPK, HNRNPC, RALY, IVNS1ABP, SF3B2, PTBP1
	NUCLEUS	0.01046	DUT, DDX19B, HNRNPU, RPS3A, TARDBP, PARP1, THRAP3, IPO5, NCL, CSE1L, TNPO1, PA2G4, PSMA1, SYNCRIP, RANBP1, HNRNPM, HNRNPL, SNRNP70, HNRNPK, PTRF, HNRNPC, MATR3, NONO, CBX3, EHD4, TPR, RALY, SMC3, NUMA1, IVNS1ABP, G3BP1, SMC2, SF3B2, KHDRBS1, MCM3, DDX17, LMNB1, PTBP1, PTMA, H2AFX, RPA1, DHCR7, ACTL6A, APPL1, NUP160, LBR
	NUCLEAR_ENVELOPE	0.01540	DDX19B, PARP1, IPO5, MATR3, TPR, LMNB1, DHCR7, NUP160, LBR
	NUCLEAR_LUMEN	0.01985	SRP68, LYAR, RPS19, RUVBL2, PARP1, THRAP3, NCL, HNRNPL, HNRNPK, SMC3, IVNS1ABP, PTBP1, RPA1, ACTL6A, APPL1
GO molecular function	DNA_BINDING	0.02572	HNRNPU, TARDBP, PARP1, PA2G4, CNBP, PCBP1, KHDRBS1, MCM3, RPA1, LBR

Distinct gene set databases (column 1) were queried to identify gene sets (column 2) which are significantly (false discovery rate (FDR) q-value <0.1, column 3) OMF-enriched when comparing the protein expression profiles of VFF with OMF. Proteins (specified by gene symbol) detected in the VFF/OMF samples and having the strongest contribution to the enrichment result (“core enrichment”, as defined by the Gene Set Enrichment Analysis (GSEA) method) are listed in column 4.

Conclusion

We have introduced near-primary fibroblasts from sheep VFM and OM as two novel in vitro models for the basic study of (patho)physiological aspects related to the respective tissue. This is especially important in the realm of laryngology, as in most instances ethical considerations prohibit the acquisition of VF biopsies from healthy humans. We further provide the first comprehensive proteomic comparison of VFF with another cell type, i.e. OMF, realized in a paired experimental design that precludes any identified difference being the result of genetic variation. Our proteome study has revealed that VFF and OMF are indeed distinct cell types, with VFF having a broader spectrum of expressed proteins, particularly a more sophisticated array of factors constituting and modifying the ECM. Last, comparative analyses of the VFF and OMF protein fingerprints have brought to light several interesting molecular pathways that could be utilized to tackle the problem of VF scarring from unprecedented angles.

Supplementary Material

Appendix A. Supplementary data

Supplementary data to this article can be found online at http://dx.doi.org/10.1016/j.jprot.2017.01.010.