Kinetics of Gene Expression Changes in Equine Fetal Interzone and Anlagen Cells Over 14 Days of Induced Chondrogenesis

Within developing synovial joints, interzone and anlagen cells progress through divergent chondrogenic pathways to generate stable articular cartilage and transient hypertrophic anlagen cartilage, respectively. Understanding the comparative cell biology between interzone and anlagen cells may provide novel insights into emergent cell-based therapies to support articular cartilage regeneration. The aim of this study was to assess the kinetics of gene expression profiles in these skeletal cell lines after inducing chondrogenesis in culture. Interzone and anlagen cells from seven equine fetuses were isolated and grown in a TGF-β1 chondrogenic inductive medium. Total RNA was isolated at ten time points (0, 1.5, 3, 6, 12, 24, 48, 96, 168, and 336 h), and gene expression for 93 targeted gene loci was measured in a microfluidic RT-qPCR system. Differential transcriptional responses were observed as early as 1.5 h after the initiation of chondrogenesis. Genes with functional annotations that include transcription regulation responded to the chondrogenic stimulation earlier (1.5–96 h) than genes involved in signal transduction (1.5–336 h) and the extracellular matrix biology (3–336 h). Between interzone and anlagen cell cultures, expression levels of 73 out of the 93 targeted genes were not initially different at 0 h, but 47 out of the 73 genes became differentially expressed under the chondrogenic stimulation. While interzone and anlagen cells are both chondrogenic, they display clear differences in response to the same TGF-β1 chondrogenic stimulation. This study provides new molecular insight into a timed sequence of the divergent developmental fates of interzone and anlagen cells in culture over 14 days.

Introduction

During early fetal skeletogenesis, articular cartilage develops from interzone located at the presumptive sites of synovial joints within the cartilaginous anlagen of bones. Therefore, two different types of cartilage differentiate in close proximity. While anlagen cartilage is transient, progressing through endochondral ossification, articular cartilage remains stable and functions throughout life to facilitate biomechanical load distribution and low friction movement between adjoining bones (1, 2). Despite the important functional properties of articular cartilage, its intrinsic ability to restore structural defects is limited in mature mammals (3) in stark contrast to the potential for bone tissue regeneration. Fractured bones repair quite well by recapitulating endochondral ossification; provided the fracture ends are brought together, stabilized, and not compromised by infection or loss of blood supply. Thus, research on the comparative cell biology between interzone and anlagen cells, as well as their developmental chondrogenic pathways may provide novel information relevant to improving mammalian articular cartilage regenerative treatments.

In an earlier study of these two cell types, chondrogenic potential was measured and compared after 21 days in three-dimensional pellet culture and continuous stimulation with a transforming growth factor beta 1 (TGF-β1) chondrogenic induction medium (4). The results demonstrated that while both interzone and anlagen cultures on 21 d showed increased expression of cartilaginous marker genes (aggrecan core protein, ACAN; type II collagen, COL2A1), the cell pellets displayed distinguishing histological characteristics, including proteoglycan amount and distribution, as well as cellular morphology and arrangement. Differences in gene expression and the kinetics of these changes between 0 and 21 d, which presumably would reflect molecular events responsible for differential histological properties of the two cell pellet cultures were not determined. In studies from other groups, the protocol for in vitro chondrogenic differentiation has also been reported after 7–28 days, with the expression of extracellular matrix (ECM) genes measured at the mRNA or protein level used as targeted functional outcomes (5, 6). However, gene expression changes induced by TGF-β, a well-established chondrogenic factor, start as early as 30 min−1 h after treatment (7, 8). Therefore, an important knowledge gap is what happens at earlier time points and the kinetics of gene expression differences in response to chondrogenic induction in these two cell types.

By comparing a timed sequence of the cellular response to TGF-β1 induced chondrogenesis between interzone and anlagen cell cultures, the present study investigated how gene expression patterns change over time and molecular details of chondrogenic divergence in these two cell types. The experiments were designed to test the hypothesis that divergent chondrogenic differentiation pathways in interzone and anlagen cultures will be evident within the first 24 h after in vitro chondrogenic induction. The aim of the study was to identify critical time points where interesting responses to the chondrogenic stimulation occur in the cell cultures over 14 days.

Materials and Methods

Cell Culture and Sample Collection

Primary equine fetal interzone cells, anlagen cells, and dermal fibroblasts (a negative control) were prepared from seven 45-days-old fetuses as previously described [(4); University of Kentucky, Institutional Animal Care and Use Committee approval number, 2014-1215], and stored frozen at passage 2 (P2). For the current experiments, cells were thawed and cultured in high glucose Dulbecco's modified Eagles medium (DMEM; 10569044, Gibco) supplemented with 10% (v/v) fetal bovine serum (heat inactivated; S11150H, Atlanta Biologicals) and 1% (v/v) penicillin/streptomycin (15070063, Gibco). When cell monolayers reached ~80% confluence, they were lifted using 0.25% Trypsin-EDTA solution (25200056, Gibco) and split into new flasks (seeding density of 10,000 cells/cm²). Cell viability at passage of >95% was confirmed by the trypan blue dye exclusion test.

When P4 monolayers reached ~80% confluence, one flask was used to collect total RNA. The cells were harvested in a guanidinium thiocyanate solution (1ml/T-75 flask; QIAzol Lysis Reagent; 79306, Qiagen), immediately snap-frozen, and stored at −80°C until total RNA isolation. The additional P4 monolayers were to establish cell pellets at P5 as previously described (4). Each cell pellet (500,000 cells/pellet) was maintained in 1 ml of chondrogenic inductive medium [high glucose DMEM supplemented with 1% penicillin/streptomycin; bovine serum albumin (12.5 mg/ml); ascorbic-2-phosphate (A8960-5G, Sigma, 50 μg/ml); TGF-β1 (human recombinant; T7039-50UG, Sigma-Aldrich, 10 ng/ml); 1% insulin-transferrin-selenium-sodium pyruvate; dexamethasone (100 nM); and 1% non-essential amino acid (11140-050, Gibco)] for the experimental period. The medium was changed every 3 days.

Pellet cultures were collected at ten time points: baseline (0 h), 1.5, 3, 6, 12, 24, 48, 96, 168, and 336 h after the initiation of chondrogenic induction. These time points were selected based on a review of literature relevant to the kinetics of gene expression after TGF-β treatment (7–14). At each time point, the collected pellets were washed with phosphate-buffered saline, snap-frozen, and stored at −80°C until total RNA isolation.

Total RNA Isolation

Thawed cell monolayer lysates and the pellets (3 pellets/1 ml of the guanidinium thiocyanate solution) were homogenized, then total RNA isolated and purified using spin-columns (RNeasy, Qiagen, 74106) and ethanol precipitation. RNA quantity was determined fluorometrically (Qubit, Q10211, Life Technologies) and purity spectrophotometrically (Nano Drop, Thermo Fisher Scientific). Finally, RNA structural integrity was determined with a Bioanalyzer 2100 (Agilent Technologies) using an Agilent RNA 6000 Pico kit (5067-1513, Agilent Technologies). Monolayer RNA samples resulted in 260/280 ratios of 2.0–2.1, 260/230 ratios of 2.3–2.7, and RINs of 8.7–10, except for one sample with a RIN of 6.8. A substantial majority of cell pellet samples had 260/280 ratios of 1.8–2.1, 260/230 ratios of 1.8–2.7, and RINs of 6.4–10. Out of 210 cell pellet RNA samples, ten samples had 260/230 ratios outside of this range, but gene expression patterns were consistent with experimental group averages, so the data were retained. Any genomic DNA contamination was removed with dsDNase during reverse-transcription using the Maxima First Strand cDNA Synthesis Kit (K1672, Thermo Fisher Scientific).

Quantitative Polymerase Chain Reaction

Targeted Gene Loci

Ninety-three genes of interest (Supplementary Table 1) were selected based on previously generated equine cartilaginous tissue RNA-seq data (15) and a literature review. These genes were either (1) differentially expressed between interzone and anlagen tissue samples at three developmental ages (day-45 gestation, day-60 gestation, and neonatal foals), (2) functionally annotated for fetal development, (3) functionally annotated for chondrogenic differentiation, and/or (4) established components or regulators of TGF-β signaling pathways.

Including three putative endogenous controls, B2M, GUSB, and RPLP0 (16), steady state mRNA levels for a total of 96 gene loci were assessed using commercially available (59 assays) and customized (37 assays) equine-specific TaqMan^® primer-probe sets (Thermo Fisher Scientific; Supplementary Table 1). Primer-probe sets designed to span exon junctions (88 assays) were prioritized, with eight assays designed within a single exon.

Positive Control RT-qPCR Assessments

Prior to the main microfluidic RT-qPCR analysis, a preliminary assessment was conducted (robotic ViiA™ 7 Real-Time PCR System, Thermo Fisher Scientific) in order to verify (1) amplification of two endogenous controls (GUSB and RPLP0) in all 231 cDNA samples (7 biological replicates × 3 cell lines × 11 time points, including P4 monolayer samples) and (2) amplification of all 96 targeted gene loci by the primer-probe sets in a positive control sample. The positive control was designed to maximize the number of expressed protein-coding genes by pooling equal parts of (1) a pooled total RNA sample composed of 43 different equine tissue/cell sources (17), and (2) a 35-day equine fetus. All 231 samples (10 ng of cDNA/reaction) expressed GUSB [cycle threshold (Ct) 21.44 ± 0.06] and RPLP0 (Ct 18.21 ± 0.06; data not shown). The positive control demonstrated amplification in 95 of the targeted gene loci (Ct 17.68–30.85) with the exception of COL10A1 (data not shown).

Microfluidic RT-qPCR

The cDNA samples (10 ng/μl) were arrayed across three 96 well-microfluidic chips (96.96 Fluidigm Dynamic Array). Seven 3-fold dilutions (125, 41.67, 13.89, 4.63, 1.54, 0.51, and 0.17 ng/μl) of the positive control sample were included on each chip to evaluate PCR efficiency and provide an inter-plate control. After 14-cycles of pre-amplification, steady state mRNA levels for the targeted 96 gene loci were measured in a microfluidic RT-qPCR system (Biomark HD high throughput amplification system, Fluidigm) using manufacturer-recommended protocols (18). The data were processed with Fluidigm Real-Time PCR Analysis software.

Data Analyses and Statistics

The two endogenous controls (GUSB and RPLP0) exhibiting the most stable expression across the sample set (Figure 1) were used for normalization within a sample (ΔCt = Ct of a gene of interest – average Ct of GUSB and RPLP0). Then, the ΔCt of each target gene was calibrated with ΔCt of the same target gene in the positive control sample to yield ΔΔCt. Finally, ΔΔCt values were converted to relative quantity [RQ = 2^−ΔΔCt; fold changes based on expression in the positive control; (19)]. To determine statistical differences between the datum points (targeted gene × cell type × time point), the fold change data were log-transformed and analyzed using SAS statistical software, version 9.3 (SAS Institute Inc.). One-way multivariate analysis of variance was conducted with Tukey's honest significance test for multiple comparison adjustments. The significance threshold was defined as P < 0.05.

FIGURE 1

Figure 1. Steady state mRNA levels (Ct, mean ± SEM; n = 7) of three prospective endogenous controls (B2M, GUSB, and RPLP0) across the 336-h experimental period measured by microfluidic RT-qPCR. GUSB and RPLP0 displayed more stability across the three cell types and time course.

Results

Of the 93 targeted genes of interest, six loci were not processed for further data analyses; five genes (ARHGEF15, CHODL, NPY, RET, and STAB1) had minimal expression in experimental samples (RQs <0.02), and one gene (IGF2) had low fluorescent intensity across 5 biological replicates. Microfluidic qPCR data (RQ = 2^−ΔΔCt values) for all gene loci investigated are available in Supplementary Table 2. The remaining 87 genes were categorized into three groups based on established functional annotation (Supplementary Table 3): (1) 15 genes regulating transcription, (2) 51 genes involved in signal transduction, and (3) 23 genes involved in ECM biology. Two genes (MASP1 and NEFL) were not categorized into any of the three annotation groups. Four genes were included in two of the three groups; ENTPD1, ENTPD2, and LEF1 in both the transcription and signaling groups, and THBS4 in both the signaling and ECM groups.

Kinetics of Gene Expression Changes Within Cell Type

Time point differences in steady state mRNA levels for individual gene loci were assessed by comparing values at the given time point to the 0 h baseline within a cell type. Significant upregulation and downregulation events were noted. At the first time point (1.5 h), the significant changes observed were all upregulation. Time points from 3 h on had instances of both upregulation and downregulation. Among the total 261 gene × cell type combinations (87 genes × 3 cell types), 110 combinations showed upregulation (Supplementary Figure 1A), 83 downregulation (Supplementary Figure 1B), and 22 combinations had mixed patterns of upregulation and downregulation (Supplementary Figure 1C). In the other 46 combinations, steady state mRNA levels were stable across all time points (Supplementary Figure 1D).

Monolayer vs. 0 h

To assess the effect of trypsin digestion and centrifugation required for establishing cell pellets from monolayers, steady state mRNA levels in P4 cell monolayers were compared to that of P5 0 h cell pellets. Among the 87 targeted loci, only four genes (APLNR, GDF5, S1PR3, and TLR2) had significantly lower expression levels in monolayers compared to 0 h cell pellets in one or more cell types (P < 0.05; Supplementary Figure 2). These genes were all in the signal transduction category.

Timing of Initial Differential Expression Relative to 0 h Within a Cell Type

Eighty-six of the 87 targeted genes displayed at least one significant change in response to chondrogenic induction. The lone exception was GALNT14, which did not change significantly at any time point for the three cell types (Supplementary Figure 1D). Within each functional annotation category, percentages of loci displaying their first onset of change (gene × cell type combinations) were measured (Figure 2; numeric values are available in Supplementary Table 4).

FIGURE 2

Figure 2. Histogram of percentages of the first response to the chondrogenic stimulation in new gene × cell type combinations within each of the three functional annotation groups.

Initial responses to chondrogenic induction were detected at the first collection time point (1.5 h) and involved gene loci with transcription regulation or signal transduction functional annotation. No ECM annotated genes were different at 1.5 h in any cell cultures. The majority of first responses were observed within the first 24 h in all three annotation groups, with peaks at 3 h in the transcription group, 6 h in the signaling group, and 12 h in the ECM group. Most (95.4%) of the initial changes in the transcription regulatory category occurred within the first day (0–24 h), while about 40% of the first responses in the signaling group and the ECM group occurred after 24 h.

In the transcription regulation group, the first half (58.1% of the total) of gene × cell type combinations showed their first chondrogenic responses within 6 h, with the mode observed at 3 h (34.9% of the total). The other 41.9% of gene × cell type combinations first responded to chondrogenic stimulation between 12 and 96 h. Only 4.6% of gene × cell type combinations responded for the first time between 48 and 96 h, and these were observed solely in the negative control, non-chondrogenic fibroblast cultures. In interzone and anlagen cell cultures, steady state mRNA levels for all of the genes with transcription regulation functional annotation had significantly changed within the first 24 h (Table 1).

TABLE 1

Table 1. Timing of the first response to the chondrogenic induction protocol within functional annotation group in each cell type (% of gene loci with significant changes).

In comparison to the transcription regulation group, genes involved in signal transduction showed slightly delayed responses. Roughly half (53.3%) of the gene × cell type combinations changed significantly within the first 12 h. Response to chondrogenic induction was even more delayed for genes involved in ECM biology. Indeed, no differences were observed in this annotation group at 1.5 h. The first half (50.9%) of gene × cell type combinations were observed between 3 and 12 h (Figure 2).

In the negative control dermal fibroblasts, responses to chondrogenic induction were delayed compared to interzone and anlagen cells. Transcriptional regulatory genes all showed their first responses within 24 h in the chondrogenic cells, compared to some initial changes delayed until 96 h in fibroblast cultures. In the signal transduction group, 8.3 and 2.1% of the first responses were recorded at the last time point (336 h) in interzone and anlagen cell cultures, respectively, while 15.4% of those was recorded at 336 h in fibroblast cultures. Also, the majority of first reactions in the ECM group were at 3 h (23.5%) and 12 h (23.5%) in interzone cultures, 12 h (28.6%) in anlagen cultures, and 96 h (23.5%) in fibroblast cultures (Table 1).

Relative Differences Between Cell Types at Each Time Point

Steady state mRNA levels at each time point were compared in pairwise comparisons between cell types, and the results categorized into four groups (Table 2; Supplementary Figure 3). While 14–18 genes were already differentially expressed at 0 h in the comparisons between cell types, most of the targeted loci (69–73 of the 87) did not show differences initially.

TABLE 2

Table 2. Four patterns of differential gene expression before and after inducing in vitro chondrogenesis.

Genes With Differential Expression Levels Prior to the Chondrogenic Induction

Across all three cell type comparisons, there were three common genes (DLX5, LEF1, and RUNX3) that had different baseline expression levels at 0 h, while also displaying cell type-specific responses to chondrogenic induction. Specific to comparisons between interzone and anlagen cultures, five genes (GDF6, MGP, OMD, PDLIM1, and RUNX2) started with different expression levels and also retained differential profiles following the chondrogenic induction (Table 2; Supplementary Figure 3A).

If initial mRNA levels of a gene were different in the comparisons to fibroblast cultures, the gene also showed different expression levels after inducing chondrogenesis. However, between interzone and anlagen cell cultures, three genes (COL2A1, COMP, and DIO2) with different initial mRNA levels lost differences across all post-chondrogenic induction time points (Table 2; Supplementary Figure 3B).

Genes With No Differential Expression Levels Prior to the Chondrogenic Stimulation

Among the 69–73 genes that initially did not have different expression levels in pairwise cell-type comparisons, 41–47 genes developed significant differences in steady state mRNA levels after inducing chondrogenesis. Between interzone and anlagen cultures, there were 47 genes in this category, but only nine (ALPK3, ASS1, AQP1, CTGF, FAM20A, FGF1, ITGA7, PLVAP, and TLR4) were specific to the comparison between these two chondrogenic cell lines (Table 2; Supplementary Figure 3C).

On the other hand, seven genes showed no difference in any cell type comparison: ADGRG2, BMP2, COL10A1, KCNJ8, LOC100630171, PTCH2, and SPARCL1. Focusing only on interzone and anlagen cultures, however, twelve additional genes responded to chondrogenic induction without cell-type differences: ADGRG1, ALPL, APLNR, DCN, ENTPD1, FAM132A, FZD1, GDF5, IGFBP5, IGFBP7, MMP2, and PLAT (Table 2; Supplementary Figure 3D).

Individual Gene Loci Demonstrate Clear Examples of Cell Type-Specific Differences

Cartilage biomarkers (COL2A1 and COMP) were upregulated in the chondrogenic cell lines over time (Figures 3A,B). While their steady state mRNA levels were also increased in fibroblast cultures toward the end of the experimental period, the levels were consistently lower compared to the chondrogenic cell cultures (Supplementary Figures 3A,C).

FIGURE 3

Figure 3. Steady state mRNA levels of example cell type-specific expression profiles across the experimental period. Relative quantities were calculated based on the positive control sample (pooled equine adult tissue and fetus RNA). Mean ± SEM (n = 7); Classic cartilaginous biomarkers: (A) COL2A1 and (B) COMP; Example genes that were differentially regulated in interzone cell cultures compared to the other cell types: (C) ABI3BP was downregulated and (D) PRKG2 was minimally expressed in interzone cell cultures; Example gene that was differentially regulated in anlagen cell cultures compared to the other cell types: (E) PANX3 was upregulated only in anlagen cell cultures toward later time points.

There were examples of genes that displayed unique expression patterns in interzone or anlagen cultures compared to the other two cell types. For interzone cells, changes in steady state mRNA levels of ABI3BP showed a clear example of downregulation (Figure 3C), while PRKG2 is an example of minimal expression (Figure 3D). In contrast, PANX3 was upregulated only in anlagen cell cultures toward the later time points, while its mRNA levels were consistently lower or not detected in the other two cell cultures (Figure 3E).

Discussion

Relationship Between Early and Delayed Gene Expression Changes

The results from the present study demonstrate a sequence whereby genes encoding proteins functionally annotated in transcription, signaling events, and ECM are expressed in roughly that order. The model consistent with this observation is that transcription regulating genes are primary responders, followed by secondary changes in downstream signaling events alter effector genes involved in ECM biology. To further elucidate these relationships and define “(early) primary” and “(delayed) secondary” genes from a perspective of cell biology, an experiment designed to include protein synthesis inhibitors will be necessary (20).

It is interesting to note that all changes in mRNA levels observed at 1.5 h were upregulation compared to the baseline. On the other hand, downregulation started to be detected by 3 h. This may confirm that degradation of transcripts requires more time than de novo synthesis (21). Regulatory genes may have relatively faster decay rates in order to rapidly alter signal transduction; 16 out of the 18 genes that showed their first responses as downregulation between 3 and 6 h were involved in either transcription or signaling events.

While the present study suggests some interesting time points where major regulatory responses to the chondrogenic stimulation likely occur, additional data on a transcriptome level will be needed to reveal important regulatory mechanisms responsible for the divergent chondrogenic pathways between interzone and anlagen cells. This knowledge would provide new information relevant to therapeutic efforts to improve articular cartilage regeneration.

Differential Gene Regulation Between Cell Types

Differences Between Chondrogenic and Non-chondrogenic Cell Cultures

Consistent with their established functional annotations, COL2A1 and COMP, classic cartilaginous biomarker genes (22, 23), were significantly upregulated in interzone and anlagen cell cultures confirming their chondrogenic potential. Overall, delayed responses in all functional annotation categories were observed in the fibroblasts compared to the chondrogenic cell lines. While this is likely due to the chondrogenic potential of interzone and anlagen cells, it is important to note that the targeted genes were selected based in part on established annotations related to chondrogenesis and expression in skeletal tissues. A full transcriptome analysis would help to determine if the relationship to chondrogenesis is a key variable.

Differences Between Interzone and Anlagen Cell Cultures

Interzone and anlagen tissues in equine fetuses at day-45 of gestation are readily distinguished morphologically with a dissecting microscope. Not surprisingly some gene expression differences were noted, with 14 out of the 87 targeted gene loci starting at different mRNA levels. Among the 14 genes, 11 genes still expressed different profiles after inducing chondrogenesis, and five of them showed these differential patterns only between interzone and anlagen cultures, but not in the comparisons to dermal fibroblasts. OMD and RUNX2 were upregulated in anlagen cell cultures after the chondrogenic induction, and their known annotation is related to ECM production in the bone (24, 25). Thus, the results are consistent with anlagen cells having developmental pathways leading to bone formation. Type X collagen (COL10A1) is a well-established biomarker of hypertrophic chondrocytes (26, 27). This was the single gene locus from the panel that was not detected in the positive control sample, a finding explained by the fact that the pool of RNAs used did not include hypertrophic cartilage. However, with the transcript pre-amplification steps used in the Fluidigm protocol, RQ values for COL10A1 could be calculated but did not exhibit significant cell type differences (Supplementary Figure 4). Subsequent RNA-seq analyses on five of the earlier time points (0, 1.5, 3, 12, and 96 h) confirmed no cell type differences (data not shown).

On the other hand, MGP and GDF6 were upregulated in interzone cultures, and their reported functional annotations are inhibiting ectopic tissue calcification (28) and diarthrodial joint formation (29), respectively. In addition, mRNA levels of PDLIM1 were greater in interzone cell cultures at the latest time points. From a previous RNA-seq dataset (Adam et al., in preparation), this gene was not differentially expressed between interzone and anlagen tissue lineages from 45-day-old equine fetuses and neonatal foals, but its expression was 3.41-fold greater (P < 0.0001) in interzone tissue compared to cartilaginous anlagen from 60-day-old fetuses. These observations may suggest that this transcription regulatory gene might become upregulated in interzone at a later stage during articular cartilage development. Taken together, expression patterns of these genes in the culture model are consistent with the articular and hypertrophic cartilaginous tissue outcomes of interzone and anlagen cells, respectively.

As noted above, steady state mRNA levels for 73 of the 87 targeted gene loci were not significantly different at 0 h, and the expression of these genes either diverged or changed in similar ways during the time course. Forty-seven genes differentially responded to the chondrogenic induction, but only nine of these were specific to the comparison between interzone and anlagen pellets. FGF1, a gene expressed in proliferative and hypertrophic chondrocytes (30), was upregulated early (6 h) in anlagen cultures and justifies further investigation. For the other 26 genes, the changes were concordant after chondrogenic induction. GDF5 is an established interzone biomarker (2, 31), but in the current in vitro model using primary cells derived from 45-day-old equine fetuses, its mRNA profiles were not different between interzone and anlagen cultures at both pre- and post-chondrogenic induction time points. The discrepancy may be due to the important variable of fetal age (developmental stage) between the two sample sets, but cell signaling and species-specific parameters cannot be excluded. Similarly, while ALPL is known to be involved in bone mineralization (32), its mRNA levels were not different between the two chondrogenic cell cultures. This may not be surprising given that the culture system was not intended to model osteogenesis.

Cell Type-Specific Expression Profiles

Several gene loci displayed clear cell type-specific expression profiles. Examples in interzone cell cultures were ABI3BP and PRKG2. ABI3BP is a novel ECM gene and is required to switch the cellular status from proliferation to differentiation in mesenchymal stem cells, including during chondrogenesis (33). The mRNA levels of this gene were initially greater in interzone cell cultures, started to decrease by 96 h, and became as low as the other cell lines at the last time point. In one study, this gene was found to be differentially upregulated in articular cartilage compared to hypertrophic growth plate both in vivo and in vitro (34). In contrast, mRNA levels of PRKG2 were minimal and remained low in interzone cultures, but were upregulated in the other two cell lines at later time points. PRKG2 is known to be involved in mammalian skeletal development, and its null mutation resulted in 23–30% decreased length of limb bones (35) due to impaired chondrocyte hypertrophy (36). Finally, PANX3 was specifically upregulated in anlagen cell cultures. This gene is expressed in cartilaginous anlagen of a developing limb (37), and when PANX3 was knocked down, hypertrophic differentiation was delayed and attenuated (38).

In conclusion, the data supported the hypothesis that divergent chondrogenic pathways in interzone and anlagen cell cultures will be evident within an early time frame—within the first 24 h after inducing in vitro chondrogenesis. The results confirmed the chondrogenic potential of interzone and anlagen cells, but also documented distinct functional responses to the same TGF-β1 chondrogenic induction signal. Several gene expression profiles were broadly consistent with well-established developmental fates of interzone and anlagen cells within limb buds. Overall, transcription regulatory responses preceded changes in genes with functional annotation related to signal transduction or ECM biology. Further investigation at a transcriptome level, together with functional experiments, should help identify key regulators and gene ontology relationships relevant to the mechanisms underlying divergent chondrogenesis between interzone and anlagen cells.

Data Availability Statement

The original contributions presented in the study are included in the article/Supplementary Material, further inquiries can be directed to the corresponding author/s.

Author Contributions

JM: secured funding. CM and JM: study design, data analysis, and drafting and revision. CM: data acquisition. All authors contributed to the article and approved the submitted version. All authors contributed to the article and approved the submitted version.

Funding

Financial support for this study was received from the Morris Animal Foundation (D16EQ-016), USDA Hatch (KY014050, Accession number 1006558), and the Lourie Foundation (Spy Coast Farm).

Conflict of Interest

The authors declare that the research was conducted in the absence of any commercial or financial relationships that could be construed as a potential conflict of interest.

Publisher's Note

All claims expressed in this article are solely those of the authors and do not necessarily represent those of their affiliated organizations, or those of the publisher, the editors and the reviewers. Any product that may be evaluated in this article, or claim that may be made by its manufacturer, is not guaranteed or endorsed by the publisher.

Acknowledgments

The authors graciously acknowledge efforts by Drs. Emma Adam, Rashmi Dubey, and Parvathy Thampi associated with independent studies in the laboratory that generated the primary cell lines used in the present experiments and/or data considered in selecting some of the targeted gene loci selection. The microfluidic RT-qPCR assay was conducted by Dr. Mark Band at the Roy J. Carver Biotechnology Center, a genomics core at the University of Illinois at Urbana-Champaign. Statistical consultation was provided by Applied Statistics Laboratory at the University of Kentucky.

Supplementary Material

The Supplementary Material for this article can be found online at: https://www.frontiersin.org/articles/10.3389/fvets.2021.722324/full#supplementary-material

References

1. Pitsillides AA, Ashhurst DE. A critical evaluation of specific aspects of joint development. Dev Dyn. (2008) 237, 2284–94. doi: 10.1002/dvdy.21654

PubMed Abstract | CrossRef Full Text | Google Scholar

2. Decker RS, Koyama E, Pacifici M. Genesis and morphogenesis of limb synovial joints and articular cartilage. Matrix Biol. (2014) 39, 5–10. doi: 10.1016/j.matbio.2014.08.006

PubMed Abstract | CrossRef Full Text | Google Scholar

3. Alford JW, Cole BJ. Cartilage restoration, part 1:Basic science, historical perspective, patient evaluation, and treatment options. Am J Sports Med. (2005) 33, 295–306. doi: 10.1177/0363546504273510

PubMed Abstract | CrossRef Full Text | Google Scholar

4. Adam EN, Janes J, Lowney R, Lambert J, Thampi P, Stromberg A, et al. Chondrogenic differentiation potential of adult and fetal equine cell types. Vet Surg. (2019) 48, 375–87. doi: 10.1111/vsu.13183

PubMed Abstract | CrossRef Full Text | Google Scholar

5. McCarthy HE, Bara JJ, Brakspear K, Singhrao SK, Archer CW. The comparison of equine articular cartilage progenitor cells and bone marrow-derived stromal cells as potential cell sources for cartilage repair in the horse. Vet J. (2012) 192, 345–51. doi: 10.1016/j.tvjl.2011.08.036

PubMed Abstract | CrossRef Full Text | Google Scholar

6. Rakic R, Bourdon B, Demoor M, Maddens S, Saulnier N, Galéra P. Differences in the intrinsic chondrogenic potential of equine umbilical cord matrix and cord blood mesenchymal stromal/stem cells for cartilage regeneration. Sci Rep. (2018) 8:13799. doi: 10.1038/s41598-018-28164-9

PubMed Abstract | CrossRef Full Text | Google Scholar

7. Franco DL, Mainez J, Vega S, Sancho P, Murillo MM, de Frutos CA, et al. Snail1 suppresses TGF-β-induced apoptosis and is sufficient to trigger EMT in hepatocytes. J Cell Sci. (2010) 123, 3467–77. doi: 10.1242/jcs.068692

PubMed Abstract | CrossRef Full Text | Google Scholar

8. Aomatsu K, Arao T, Sugioka K, Matsumoto K, Tamura D, Kudo K, et al. TGF-β induces sustained upregulation of SNAI1 and SNAI2 through Smad and non-smad pathways in a human corneal epithelial cell line. Invest Ophthalmol Vis Sci. (2011) 52, 2437–43. doi: 10.1167/iovs.10-5635

PubMed Abstract | CrossRef Full Text | Google Scholar

9. Ranganathan P, Agrawal A, Bhushan R, Chavalmane AK, Kalathur RKR, Takahashi T, et al. Expression profiling of genes regulated by TGF-beta: Differential regulation in normal and tumour cells. BMC Genomics. (2007) 8:98. doi: 10.1186/1471-2164-8-98

PubMed Abstract | CrossRef Full Text | Google Scholar

10. Rudini N, Felici A, Giampietro C, Lampugnani M, Corada M, Swirsding K, et al. VE-cadherin is a critical endothelial regulator of TGF-beta signalling. EMBO J. (2008) 27, 993–1004. doi: 10.1038/emboj.2008.46

PubMed Abstract | CrossRef Full Text | Google Scholar

11. Koyama N, Miura M, Nakao K, Kondo E, Fujii T, Taura D, et al. Human induced pluripotent stem cells differentiated into chondrogenic lineage via generation of mesenchymal progenitor cells. Stem Cells Dev. (2013) 22, 102–13. doi: 10.1089/scd.2012.0127

PubMed Abstract | CrossRef Full Text | Google Scholar

12. Yokota J, Chosa N, Sawada S, Okubo N, Takahashi N, Hasegawa T, et al. PDGF-induced PI3K-mediated signaling enhances the TGF-β-induced osteogenic differentiation of human mesenchymal stem cells in a TGF-β-activated MEK-dependent manner. Int J Mol Med. (2014) 33, 534–42. doi: 10.3892/ijmm.2013.1606

PubMed Abstract | CrossRef Full Text | Google Scholar

13. Nejadnik H, Diecke S, Lenkov OD, Chapelin F, Donig J, Tong X, et al. Improved approach for chondrogenic differentiation of human induced pluripotent stem cells. Stem Cell Rev Rep. (2015) 11, 242–53. doi: 10.1007/s12015-014-9581-5

PubMed Abstract | CrossRef Full Text | Google Scholar

14. Yamazaki S, Nakano N, Honjo A, Hara M, Maeda K, Nishiyama C, et al. The transcription factor Ehf is involved in TGF-β-induced suppression of FcεRI and c-Kit expression and FcεRI-mediated activation in mast cells. J Immunol. (2015) 195, 3427–35. doi: 10.4049/jimmunol.1402856

PubMed Abstract | CrossRef Full Text | Google Scholar

15. Adam EN, Zhang Y, Liu X, Coleman S, Flight R, Lambert J, et al. Gene expression patterns that distinguish equine articular cartilage from other cartilaginous tissues. Presentation Abstracts, 2017 ACVS Surgery Summit. (2017) 46:E1–75. doi: 10.1111/vsu.12710

CrossRef Full Text | Google Scholar

16. Mienaltowski MJ, Huang L, Stromberg AJ, MacLeod JN. Differential gene expression associated with postnatal equine articular cartilage maturation. BMC Musculoskelet Disord. (2008) 9:149. doi: 10.1186/1471-2474-9-149

PubMed Abstract | CrossRef Full Text | Google Scholar

17. Hestand MS, Kalbfleisch TS, Coleman SJ, Zeng Z, Liu J, Orlando L, et al. Annotation of the protein coding regions of the equine genome. PLoS ONE. (2015) 10:e0124375. doi: 10.1371/journal.pone.0124375

PubMed Abstract | CrossRef Full Text | Google Scholar

18. Fluidigm. Real-Time PCR Analysis PN, 68000088 N1 User Guide. South San Francisco CA: Fluidigm (2018).

Google Scholar

19. Livak KJ, Schmittgen TD. Analysis of relative gene expression data using real-time quantitative PCR and the 2–ΔΔCT method. Methods. (2001) 25, 402–8. doi: 10.1006/meth.2001.1262

CrossRef Full Text | Google Scholar

20. Yamamoto KR, Alberts BM. Steroid receptors: elements for modulation of eukaryotic transcription. Annu Rev Biochem. (1976) 45, 721–46. doi: 10.1146/annurev.bi.45.070176.003445

PubMed Abstract | CrossRef Full Text | Google Scholar

21. Eser P, Demel C, Maier KC, Schwalb B, Pirkl N, Martin DE, et al. Periodic mRNA synthesis and degradation co-operate during cell cycle gene expression. Mol Syst Biol. (2014) 10:717. doi: 10.1002/msb.134886

PubMed Abstract | CrossRef Full Text | Google Scholar

22. Bi W, Deng JM, Zhang Z, Behringer RR, De Crombrugghe B. Sox9 is required for cartilage formation. Nat Genet. (1999) 22, 85–9. doi: 10.1038/8792

CrossRef Full Text | Google Scholar

23. Akiyama H, Chaboissier MC, Martin JF, Schedl A, de Crombrugghe B. The transcription factor Sox9 has essential roles in successive steps of the chondrocyte differentiation pathway and is required for expression of Sox5 and Sox6. Genes Dev. (2002) 16, 2813–28. doi: 10.1101/gad.1017802

PubMed Abstract | CrossRef Full Text | Google Scholar

24. Ninomiya K, Miyamoto T, Imai, J.-,i, Fujita N, Suzuki T, et al. Osteoclastic activity induces osteomodulin expression in osteoblasts. Biochem Biophys Res Commun. (2007) 362, 460–6. doi: 10.1016/j.bbrc.2007.07.193

CrossRef Full Text | Google Scholar

25. Mevel R, Draper JE, Lie-a-Ling M, Kouskoff V, Lacaud G. RUNX transcription factors: orchestrators of development. Development. (2019) 146:dev148296. doi: 10.1242/dev.148296

PubMed Abstract | CrossRef Full Text | Google Scholar

26. Kwan AP, Cummings CE, Chapman JA, Grant ME. Macromolecular organization of chicken type X collagen in vitro. J Cell Biol. (1991) 114, 597–604. doi: 10.1083/jcb.114.3.597

PubMed Abstract | CrossRef Full Text | Google Scholar

27. Schmid TM, Cole AA, Chen Q, Bonen DK, Luchene L, Linsenmayer TF. Assembly of Type X collagen by hypertrophic chondrocytes. In: Yurchenco PD, Birk DE, Mecham RP, editors. Extracellular Matrix Assembly and Structure. San Diego: Academic Press (1994). p. 171–206.

Google Scholar

28. Luo G, Ducy P, McKee MD, Pinero GJ, Loyer E, Behringer RR, et al. Spontaneous calcification of arteries and cartilage in mice lacking matrix GLA protein. Nature. (1997) 386, 78–81. doi: 10.1038/386078a0

PubMed Abstract | CrossRef Full Text | Google Scholar

29. Settle SH, Rountree RB, Sinha A, Thacker A, Higgins K, Kingsley DM. Multiple joint and skeletal patterning defects caused by single and double mutations in the mouse Gdf6 and Gdf5 genes. Dev Biol. (2003) 254, 116–30. doi: 10.1016/S0012-1606(02)00022-2

PubMed Abstract | CrossRef Full Text | Google Scholar

30. Krejci P, Krakow D, Mekikian PB, Wilcox WR. Fibroblast growth factors 1, 2, 17, and 19 are the predominant FGF ligands expressed in human fetal growth plate cartilage. Pediatr Res. (2007) 61, 267–72. doi: 10.1203/pdr.0b013e318030d157

PubMed Abstract | CrossRef Full Text | Google Scholar

31. Hyde G, Dover S, Aszodi A, Wallis GA, Boot-Handford RP. Lineage tracing using matrilin-1 gene expression reveals that articular chondrocytes exist as the joint interzone forms. Dev Biol. (2007) 304, 825–33. doi: 10.1016/j.ydbio.2007.01.026

PubMed Abstract | CrossRef Full Text | Google Scholar

32. Golub EE, Boesze-Battaglia K. The role of alkaline phosphatase in mineralization. Curr Opin Orthop. (2007) 18, 444–8. doi: 10.1097/BCO.0b013e3282630851

CrossRef Full Text | Google Scholar

33. Hodgkinson CP, Naidoo V, Patti KG, Gomez JA, Schmeckpeper J, Zhang Z, et al. Abi3bp is a multifunctional autocrine/paracrine factor that regulates mesenchymal stem cell biology. Stem Cells. (2013) 31, 1669–82. doi: 10.1002/stem.1416

PubMed Abstract | CrossRef Full Text | Google Scholar

34. Hissnauer TN, Baranowsky A, Pestka JM, Streichert T, Wiegandt K, Goepfert C, et al. Identification of molecular markers for articular cartilage. Osteoarthr Cartil. (2010) 18, 1630–8. doi: 10.1016/j.joca.2010.10.002

PubMed Abstract | CrossRef Full Text | Google Scholar

35. Pfeifer A, Aszódi A, Seidler U, Ruth P, Hofmann F, Fässler R. Intestinal secretory defects and dwarfism in mice lacking cGMP-dependent protein kinase II. Science. (1996) 274:2082–6. doi: 10.1126/science.274.5295.2082

PubMed Abstract | CrossRef Full Text | Google Scholar

36. Kawasaki Y, Kugimiya F, Chikuda H, Kamekura S, Ikeda T, Kawamura N, et al. Phosphorylation of GSK-3β by cGMP-dependent protein kinase II promotes hypertrophic differentiation of murine chondrocytes. J Clin Invest. (2008) 118:2506–15. doi: 10.1172/JCI35243E1

PubMed Abstract | CrossRef Full Text | Google Scholar

37. Bond SR, Lau A, Penuela S, Sampaio AV, Underhill TM, Laird DW, et al. Pannexin 3 is a novel target for Runx2, expressed by osteoblasts and mature growth plate chondrocytes. J Bone Miner Res. (2011) 26:2911–22. doi: 10.1002/jbmr.509

PubMed Abstract | CrossRef Full Text | Google Scholar

38. Oh SK, Shin JO, Baek JI, Lee J, Bae JW, Ankamerddy H, et al. Pannexin 3 is required for normal progression of skeletal development in vertebrates. FASEB J. (2015) 29:4473–84. doi: 10.1096/fj.15-273722

PubMed Abstract | CrossRef Full Text | Google Scholar