Abstract Unlike highly regenerative animals, such as axolotls, humans are believed to be unable to counteract cumulative damage, such as repetitive joint use and injury that lead to the breakdown of cartilage and the development of osteoarthritis. Turnover of insoluble collagen has been suggested to be very limited in human adult cartilage. The goal of this study was to explore protein turnover in articular cartilage from human lower limb joints. Analyzing molecular clocks in the form of nonenzymatically deamidated proteins, we unmasked a position-dependent gradient (distal high, proximal low) of protein turnover, indicative of a gradient of tissue anabolism reflecting innate tissue repair capacity in human lower limb cartilages that is associated with expression of limb-regenerative microRNAs. This association shows a potential link to a capacity, albeit limited, for regeneration that might be exploited to enhance joint repair and establish a basis for human limb regeneration.

INTRODUCTION Some animals, such as zebrafish, bichir, and axolotl, with a high regenerative capacity, regulate limb regeneration by a circuit of microRNA (miRNA) conserved across species (1). Unlike these animals, humans cannot regenerate whole limbs. Humans are also believed to be unable to counteract the cumulative damage of repetitive joint use or one substantial, usually sports- or trauma-related, injury that leads to the breakdown of cartilage and the development of osteoarthritis (OA). Although humans have a limited natural regenerative capacity, evidenced by the ability to regrow distal portions of amputated digits during childhood (2), by the carbon-14 (14C) bomb pulse method, turnover of insoluble collagen (residual after hyaluronidase and trypsin digestion) was suggested to be very limited in human adult cartilage and Achilles tendon (3, 4). We have focused on endogenous molecular clocks based on spontaneous posttranslational protein modification by deamidation—the nonenzymatic hydrolysis of the amide group on the side chains of asparagine (Asn, N) and glutamine (Gln, Q) to yield aspartate (Asp, D) and glutamate (Glu, E) (5), respectively. Until a series of studies demonstrated the deamidation of cytochrome c in vivo (6–8), the deamidation of proteins was considered an in vitro artifact introduced by the protein purification process. Because both nondeamidated (“young”) and deamidated (“old”) epitopes can be simultaneously quantified by mass spectrometric analysis due to a mass shift increase of 0.984 Da with deamidation (9–11), posttranslational modification due to deamidation is a particularly attractive means for monitoring protein turnover, defined as the balance between protein synthesis and protein degradation. Deamidation may be considered irreversible since there are no known extracellular repair mechanisms for these spontaneous nonenzymatic modifications (12), and this reversal has not been observed (5). Because newly synthesized proteins eventually replace posttranslationally modified proteins, monitoring the accumulation of deamidated proteins by mass spectrometry (MS) can be a way of estimating protein turnover and tissue repair. The turnover of cartilage proteins is particularly suited to study using protein molecular clocks because many cartilage proteins are long lived and therefore susceptible to accumulation of many types of age-related posttranslational nonenzymatic protein amino acid modifications (13–15). Deamidation rates of Asn and Gln vary from hours to decades. Each site of nonenzymatic amino acid deamidation serves as its own intrinsic molecular clock due to the deamidation rate varying according to the two-dimensional (2D; pentapeptide sequences) (5, 16, 17) and 3D (contextual) features (18) of a protein epitope. This method of determining biological half-life is theoretically applicable to any protein with Asn and/or Gln residues. However, to date, reliable prediction of deamidation rates based on 2D and 3D protein features is only available for Asn. Moreover, median deamidation rates of Asn are much faster than those of Gln (19), allowing estimations of half-life for even short-lived proteins. For all these reasons, we focused specifically on Asn deamidation for our analyses. For those animals with high limb-regenerative capacity, blastema formation at the site of limb injury is critical (20, 21). Genetic defects in blastema formation can block limb/appendage regeneration (22, 23). Most recently, an important blastema miRNA regulatory circuit shared by three highly regenerative animal systems has been shown to control limb regeneration (1). We hypothesized that miRNA, critical for blastema formation during limb regeneration across species (20–23), would play a role in articular cartilage protein turnover. To assess whether human cartilage shares any of the evolutionarily conserved miRNA regulatory circuitry of the blastema, essential for limb regeneration in highly regenerative animals (1), we extracted total RNA from human ankle, knee, and hip cartilages for quantification of blastema-relevant miRNA (miR-21, miR-31, and miR-181c) expression in cartilage. To evaluate protein turnover by joint site, disease state, and cartilage depth, we used proteomic tools to quantify the amounts of native and deamidated forms of cartilage proteins. Last, we assessed the association of these cartilage proteins with regenerative miRNA. To our knowledge, this is the first study to comprehensively profile the regional patterns of protein turnover and miRNA expressions in human adult limb joint cartilage in health and disease. The results of this study can provide insight into the potential regenerative capacity within human postnatal cartilages and lay the foundation for identifying the underlying mechanisms determining and regulating this regenerative capacity.

RESULTS A determination of protein turnover based on spontaneous, nonenzymatic Asn deamidation of peptides requires crystal structure for the protein domain containing the peptide of interest (16). Actual crystal structure was available for cartilage oligomeric matrix protein (COMP) and fibronectin (FN1); in silico–predicted protein structure could be generated using the Iterative Threading Assembly Refinement server (I-TASSER) for COMP and FN1, and an additional five cartilage proteins, including clusterin (CLU), prolargin (PRELP), collagen type III (COL3A1), collagen type II (COL2A1), and aggrecan core protein (G1 and G3) domains (24). The abundance of the deamidated (Asp containing) relative to the nondeamidated (Asn containing) forms of these proteins revealed distinct spatial and disease state–related differences in the turnover of these proteins, with the highest protein turnover rates (lowest ratio of deamidated to nondeamidated protein) in ankle joints (compared with knee and hip; Fig. 1A), at the cartilage surface (compared with deep; Fig. 1B), and in OA cartilage (compared with healthy; fig. S1B). To evaluate the independent effects of joint site, disease state, age, and gender and account for the repeated measures from cartilage depth, we used a multivariable mixed model. A significant gradient of protein turnover, indicative of tissue anabolism reflecting an innate tissue repair capacity by joint site (highest ankle, intermediate knee, and lowest hip), was evident for the majority of cartilage proteins, including the aggrecan G1 domain (aggrecan N terminus), COMP, PRELP, CLU, and COL2A1 (Fig. 1A). We also observed a uniformly higher turnover rate at the cartilage surface that was generalizable to multiple proteins, including FN1, CLU, COL3A1, and aggrecan G1 and G3 domains (Fig. 1B). We developed a high-sensitivity, sandwich enzyme-linked immunosorbent assay (ELISA) specific to the deamidated form of the N terminus of COMP (14). The abundance of deamidated COMP (DCOMP) in cartilage extracts by sandwich ELISA associated strongly with the abundance quantified by MS (Pearson r = 0.77, P = 0.01) (fig. S1A). These results suggested that artifactual deamidation during tissue processing for MS was minimal. Fig. 1 Spatial and cartilage depth-related distributions of deamidated proteins in human adult cartilage. Protein turnover in human lower limb cartilage was based on Asn deamidation quantified by MS. Counting only proteins identified in at least two biological replicates of one joint site, a total of 469 proteins were identified in cartilage, of which 285 proteins (61%) had at least one deamidated residue. We observed significant independent differences in the deamidation protein ratio in articular cartilage (n = 54 samples in total) stratified by (A) joint site and (B) cartilage depth. A significant lower mean deamidated protein ratio, indicative of higher tissue anabolism reflecting a higher innate tissue repair capacity, was consistently observed (A) in ankle cartilage compared with knee and hip cartilage, and (B) at the surface zone of cartilage compared with the deep zone. Abundance of deamidated protein forms is represented as the mean standardized value of the deamidation ratio (z score derived from the ratio of deamidated Asn/total Asn); the error bars represent the SD. P values were derived from multivariable mixed models. In contrast to these proteins, the highest turnover of the aggrecan G3 (C-terminal) domain was in the hip, as evidenced by the lowest abundance of deamidated aggrecan G3 in the hip compared with the knee and ankle joints. This may indicate the retention of the G3 domain in the ankle after cleavage by metalloproteinase and ADAMTS (A Disintegrin and Metalloproteinase with Thrombospondin motifs) enzymes (25–28). In fact, on the basis of our MS analyses in the same samples, the abundance of aggrecan G3 was linked to the abundance of tenascin-C (Pearson r = 0.45, P = 0.01) (fig. S1B), a high-affinity ligand for aggrecan G3 (29, 30). Tenascin-C was significantly enriched in the ankle compared with the knee and hip cartilage. This suggests that the cleaved G3 domain is preferentially retained in ankle cartilage due to its interaction with extracellular matrix (ECM) proteins enriched in ankle cartilage. A higher turnover with OA, reflected by less abundant mean deamidation ratio compared with healthy cartilage, was observed in all proteins, but only statistically significant for COMP (P = 0.012) after controlling for all other factors, namely, joint site and cartilage depth. Prior efforts to estimate the half-lives of cartilage proteins based on Asp racemization (31–33) required specific protein purification and were restricted to very long-lived proteins, such as cartilage collagen, due to the very slow rate of Asp racemization (7.7 × 10−4 to 8.3 × 10−4 per year). In contrast, mass spectrometric analysis of deamidated peptides does not require specific protein purification. Because of the fact that MS can distinguish protein identity based on peptide sequence and the rates of Asn deamidation are commensurate with human biology, ranging from hours to years (34), the estimation of protein half-life can potentially be applied to any Asn-containing protein epitope in any tissue. We confirmed that predictive modeling is feasible and valid for half-life determinations when crystal structures are not available based on the strong association of half-lives derived from protein crystal structure and in silico I-TASSER–predicted models for two proteins (Pearson r = 0.93, P < 0.0001 for COMP; Pearson r = 0.76, P < 0.0001 for FN1). As in silico–predicted structures were available for all proteins, we used these to compare half-lives among proteins. The estimated protein half-lives based on Asn deamidation analyses ranged from days to years (Fig. 2A). Within healthy cartilage, the shortest mean half-life, 3.4 days [95% confidence interval (95% CI), 1.8 to 4.9], was observed for the C-terminal G3 domain of aggrecan. The N-terminal G1 domain had a much longer mean half-life, 87.9 days (95% CI, 63 to 112.7), than the G3 domain (a mean 84.5-fold greater for matched samples, P < 0.0001) (Fig. 2B), consistent with the known retention of the N terminus of aggrecan in the cartilage matrix after metalloproteinase and ADAMTS proteolysis (25–28). In contrast to insoluble COL2A1—reported to have no significant turnover in adults, even with the occurrence of disease (3), soluble COL2A1 had a mean half-life of only 2420 days. These data demonstrate that this is a powerful approach for obtaining data simultaneously on multiple individual proteins and even epitopes within a protein, as demonstrated here for the aggrecan G1 and G3 domains, without the need to purify the specific protein of interest. The major determinant of half-life differences was protein identity; the next major determinant of protein half-life was joint site (cartilage protein half-lives of ankle < knee < hip; Fig. 2C). Shorter protein half-lives were observed in ankle cartilages for aggrecan G1, COMP, CLU, and PRELP. Fig. 2 Estimated protein half-lives based on protein deamidation. (A) The proof of concept for estimating protein half-lives by deamidation was demonstrated by the analysis of abundant cartilage ECM proteins. The half-lives varied from days to years (range, 2.8 to 2420 days) depending on the protein. The symbols indicate the 5th to 95th percentile. (B) Within healthy cartilage, the N-terminal G1 domain had a much longer mean half-life, 87.9 days (95% CI, 63 to 112.7) compared with the G3 domain (3.4 days), representing a mean 84.5-fold greater half-life for G1 than G3 for matched samples (paired t test P < 0.0001). (C) Shorter protein half-lives were consistently observed in aggrecan G1 domain, COMP, PRELP, and CLU in ankle compared with knee and hip cartilage. Although the half-life of the COL2A1 protein was shorter in ankle cartilage, it was not significant after accounting for all independent factors. The symbols indicate the mean half-lives, and the error bars indicate the SEM. P values were derived from multivariable mixed models. To assess whether miRNA, critical for blastema formation during limb regeneration across species (20–23), would play a role in the differential cartilage protein turnover we observed in a gradient pattern of the human lower limb and in the up-regulation of cartilage protein turnover with OA, we quantified the blastema-relevant miRNA (miR-21, miR-31, and miR-181c) expression in cartilage. miR-21, the most highly expressed miRNA of the three, was differentially expressed in the human limb cartilages in patterns matching the turnover of proteins based on deamidation. Evaluated by a multivariable mixed model, including the independent effects of joint site, disease state, cartilage depth, age, and gender, we found a significant gradient of miRNA expression by joint site (highest in ankle) for miR-21 and miR-181c (P = 0.0158 and 0.018, respectively; Fig. 3A). A significant gradient of miRNA expression (highest in the superficial layer of cartilage) was also observed for miR-21 and miR-31 (P < 0.0001; Fig. 3B). We also observed that miR-181c was enriched in OA compared with healthy cartilage (P = 0.044; Fig. 3C). Abundance of miR-21 was inversely associated with proportions of deamidated cartilage ECM proteins, including aggrecan, COMP, FN1, PRELP, CLU, COL3A1, and COL2A1 (Fig. 3D). These results are consistent with a positive association of miR-21 with cartilage tissue anabolism. Higher miR-31 abundance was also linked to lower deamidation ratios of COMP, PRELP, and CLU (fig. S2). Stratifying cartilage by disease state, miR-21 was significantly associated with proportions of all seven deamidated proteins in OA cartilage but only linked to PRELP and CLU in healthy cartilage, suggesting a large effect of this miRNA on ECM production induced under the stress of OA (Fig. 3E). Fig. 3 miRNAs that regulate whole-limb regeneration in animals are enriched in the cartilage of ankle joints, OA-affected cartilage, and cartilage superficial layer. (A to C) miRNA distribution in human postnatal articular cartilage stratified by (A) joint site, (B) disease state, and (C) cartilage depth (n = 48 in total). Blastema miRNAs (miR-21, miR-31, and miR-181c) are presented as fold change from the mean expression of internal control (ctrl) miRNAs (miR-423, miR-191, and miR-U6), all quantified using miScript SYBR Green PCR Kit. P values were derived from multivariable mixed models. (D) Association of miR-21 with total protein and deamidated protein abundance. Higher miR-21 expression was consistently associated with lower deamidated protein ratio but had a less consistent association with total protein abundance. (E) The association of miR-21 with deamidated protein ratio was observed consistently in OA cartilage but less consistently in healthy cartilage. P values were derived from Pearson correlations. r = correlation. All P and r values are reported, but regression (green) lines are only shown for correlations with P < 0.05. Sup, superficial cartilage; Mid, middle depth cartilage. We then investigated the association of these miRNAs with the most abundant 193 cartilage proteins that were identified by MS in at least 50% of samples. On the basis of cluster analysis, all three miRNAs yielded similar overall patterns of association with the cartilage proteins in OA cartilage (Fig. 4A; distance equal to 0.16 to 0.17). In contrast, the patterns of association with cartilage proteins were different for the three miRNAs in healthy cartilage, especially miR-181c (Fig. 4A). The pattern for miR-181c in healthy non-OA cartilage was inversely linked to the pattern in OA cartilage (P = 0.0005). These findings suggest entrainment or synchronization of cartilage protein expression by all three miRNA in the context of the stress of OA. Only 19% of cartilage proteins that were significantly correlated with miR-21 were identified to be ECM proteins (fig. S3). This indicates that miR-21 effects are not limited to the matrix but rather that miR-21 and the other regenerative miRNA orchestrate a generalized anabolic response in human OA cartilage. Nearly 50% of the cartilage proteins that associated significantly with miR-21 were identified exclusively in OA cartilage (Fig. 4B). This subset of unique OA cartilage proteins was highly enriched in collagen proteins (Fig. 4C). Associations of these collagen proteins and the three miRNAs were highly variable in healthy cartilage but markedly similar in OA cartilage (Fig. 4D). We believe that collagens are indirect rather than direct targets of miRNAs. Using a representative example of miR-21 and collagen proteins, we now provide evidence in support of the ability of miRNA to indirectly activate collagen proteins. In this study, miR-21 gene expression was inversely associated with gene expression of its known target, transforming growth factor-β receptor 2 (TGFBR2) (r = −0.66, P = 0.0002). In turn, TGFBR2 gene expression was inversely associated with the protein expression of numerous collagen proteins (COL1A1, COL1A2, COL2A1, COL5A2, COL9A1, COL9A3, and COL16A1; range of P values, 0.042 to 0.0001; range of r values, −0.41 to −0.73) quantified by MS. The net effect is an up-regulation of expression of collagen proteins by miR-21. Together, these findings suggest that major cartilage proteins, such as collagens, are coordinately regulated by evolutionarily conserved blastema miRNAs in response to the stress of OA. Fig. 4 miRNAs and cartilage proteins respond collectively to the stress of OA. (A) Heatmap representing the associations of miR-21, miR-31, and miR-181c with the most abundant 193 cartilage proteins identified in at least 50% of samples in healthy (H) and OA cartilage. The distance (d) listed between each paired dataset was derived from the cluster analysis; distance equal to 1 means two paired sets were unassociated. The associations of miRNAs and cartilage proteins were ordered according to the pattern produced by miR-21 in the healthy cartilage as the reference. The most similar association patterns were evident for the miRNA in OA cartilage (distance equal to 0.16 to 0.17), suggesting entrainment or synchronization of cartilage protein expression by all three miRNA in the context of the stress of OA. The association of miR-181c with cartilage proteins changed the most from non-OA to OA cartilage (distance equal to 1 versus 0.16); in fact, the overall patterns for miR-181c in non-OA and OA were inversely associated (table below graphic). (B) The expression of miR-21 was significantly linked to 106 cartilage proteins: 50 of these proteins were identified exclusively in OA cartilage (among them 12 collagen components); 29 of these proteins were identified in both OA and healthy cartilage (included 3 collagen components); and 27 of these proteins were identified exclusively in healthy cartilage (no unique collagen components identified). (C) String analysis identified protein-protein interactions in OA cartilage that were significantly associated with miR-21 expression. Of note, a tightly integrated network of protein members of the collagen family was highly associated with miR-21 and exclusively identified in OA cartilage (B and C). (D) Heatmap illustrating the associations of miRNA and abundance of collagen proteins ordered according to the pattern produced by miR-21 in the healthy cartilage. Similar to results for the total measured proteome, the most similar association patterns were evident for the miRNA in OA cartilage (distance equal to 0.14 to 0.21), suggesting entrainment or synchronization of collagen protein expression by all three miRNA in the context of the stress of OA.

DISCUSSION It is generally believed that human articular cartilage has limited regenerative capacity. Still, studies have discovered progenitor cells residing in mature cartilage (35) and stem cells in Ranvier’s groove of cartilage and other intra-articular tissues (36). Intrinsic regenerative capacity of articular cartilage has also been suggested by clinical studies observing an apparent increase in radiographic joint space width and cartilage thickness as a result of joint distraction (37). Together, these studies indicate the potential of cartilage repair capability in postnatal (mature) articular cartilage. Our study demonstrated a position-dependent gradient (distal high, proximal low) of protein turnover in human lower limb cartilage. We also demonstrated that the expression of blastema miRNA is notably associated with the protein turnover gradient. These findings reveal a dynamic anabolic effect in human limbs, which reflect a potential innate, albeit limited, regenerative capacity in human cartilage. These results are consistent with studies showing increased expression of matrix proteins in OA knee cartilage (38–40). The lack of repair capacity in proximal joints may explain, at least in part, the higher prevalence of hip and knee OA compared with ankle OA, particularly the rarity of severe ankle OA (41–43), and the lower prevalence of surface fibrillation of ankle cartilage with age compared with knee cartilage (44). Together, these data suggest a role for regenerative miRNA in cartilage homeostasis, turnover, and intrinsic repair capacity. Successful regeneration in limb-regenerating animals is a position-dependent response (45). There are now several lines of evidence emerging to suggest that this anabolic repair or limited regenerative capacity in humans is primarily dictated by location rather than joint shape or loading. First, as mentioned above, in vivo evidence demonstrates a higher prevalence of hip and knee OA compared with ankle OA. Second, on the basis of principal component analysis and independent of age and OA disease status, a recent study showed a major and statistically significant difference in DNA methylation of knee compared with hip joint cartilages (46). Several other recent studies investigating genome-wide methylation profiles in cartilage have also identified joint site differences (comparing knee and hip cartilages) in methylation patterns that are independent of disease status (47–51). Third, a study comparing synovial fibroblasts from various joints, including upper and lower human limb joints, demonstrates a topographical difference in transcriptomic patterns (52). There are several limitations of this study. ECM proteins were extracted with guanidine-HCl, rather than by enzyme digestion or strong acid solution, thereby reflecting a subset of cartilage proteins and, in particular, only a subset of collagens and collagen-associated proteins (53). The COL2A1 extracted here was likely only loosely associated with the collagen framework. This extractable COL2A1 likely represents a portion of the total collagen framework that has a much higher turnover rate than the insoluble fraction containing cross-linked collagens. There are also limitations based on the use of in silico–predicted models. Both COL2A1 and COL3A1 have triple helices, but the I-TASSER protein structure prediction program could not accommodate this structure. Therefore, for these proteins, the number of interstrand hydrogen bonds was likely underestimated, which would have the effect of underestimating the deamidation rate and collagen half-life estimates. Nevertheless, the half-life of COL2A1 estimated here was the longest among the proteins we monitored. Through novel use of in vivo protein deamidation molecular clocks, we found a distal-proximal gradient of protein turnover, indicative of innate tissue anabolism in human postnatal lower limb cartilages that is notably associated with the expression of evolutionarily conserved limb-regenerative (blastema) miRNAs. We also identified a coordinate regulation of cartilage protein expression in association with the expression of regenerative miRNAs in the context of the stress of OA. These findings indicate a hitherto unappreciated innate, albeit limited, regenerative capacity of human cartilage that is linked to miRNAs and up-regulated in OA. These data also suggest that anabolic treatments may be required in addition to anticatabolic treatments, especially for hip joints, to prevent or slow the progression of OA. Future functional studies on the effects of these miRNAs would be required to further elucidate their regulatory role in cartilage repair. Mechanisms associated with the mammalian regenerative response involve the regulation of pathways that are intricately linked to a large number of human disease states (45) and as shown here, OA. Further insights into the parallels of animal limb regeneration and human limb tissue repair and regeneration could inform future therapeutic approaches to arthritis. Moreover, injection of key regenerative miRNA in a joint, singly or in combination, could potentially enhance endogenous repair and resist degeneration of joint tissues in arthritis of all types including after traumatic injury or insult. Findings such as these, related to innate repair, could also be applied to tissue-engineered constructs to enhance exogenous tissue generation prior to transplantation. Identification of key “missing” (expressed in limb-regenerating organisms but not in mammals) factors, including miRNA, could lead to future use of a “molecular cocktail” for attempted recapitulation of blastema-mediated limb regeneration in humans.

MATERIALS AND METHODS Clinical cartilage specimen sources Under the Duke Institutional Review Board approval, all articular cartilages were collected from Duke University Hospital as waste surgical specimens. Full-thickness cartilage was collected from perilesional regions of the load-bearing area of the hip, knee, and ankle joints from patients with end-stage OA who had total arthroplasty surgery. Full-thickness healthy non-OA cartilage was collected at the time of surgery for acute trauma; the absence of OA was determined by the surgeon and confirmed by macroscopic inspection upon acquisition of the sample in the laboratory. A total of 18 individual specimens, including three types of joints (hip, knee, and ankle), two types of disease state (healthy and OA), and three biological replicates of each type with matched age range, were collected. The mean age of the healthy non-OA patients was 58.8 years (range, 30 to 82); the mean age of the patients with OA was 59.8 years (range, 42 to 87). Within 2 hours of surgical acquisition, all specimens were snap frozen in tubes on dry ice and stored in a −80°C freezer until extracted. Tissue dissection and intracellular proteome isolation Frozen cartilage specimens without subchondral bone were embedded in Tissue-Tek O.C.T. (Sakura, Alphen aan den Rijn, The Netherlands) for cryosectioning. Serial transverse frozen sections (12-μm thickness) were collected starting from the cartilage surface; 20 adjacent sections were collected at different distances from the surface to represent the superficial (0 to 240 μm), middle (480 to 720 μm), and deep (960 to 1200 μm) layers. Twenty sections between each layer were skipped to avoid cross-contamination. In all, a total of 54 samples were acquired, representing each joint site, disease state, and cartilage depth. Sample preparation for MS analysis Extraction of the cartilage with guanidine-HCl, representing the ECM-enriched proteome, was performed as previously described (54). Briefly, for mass spectrometric analyses, frozen cartilage sections were extracted using 4 M guanidine-HCl extraction buffer with 0.2% RapiGest (Waters Corporation, Milford, MA) for 24 hours on an orbital shaker at 4°C. Protein extracts were reduced with 4 mM dithiothreitol, alkylated with 16 mM iodoacetamide, and then precipitated with ethanol (9:1) overnight at 4°C. Samples were digested with 2 μg of trypsin, filtered through a 30-kDa filter (Pall Life Sciences, Port Washington, NY) and then a reversed-phase C18 column (The Nest group, Southborough, MA) to remove residual polysaccharides and salts. MS analysis Discovery experiments (nontargeted MS) were performed on a quadrupole Orbitrap benchtop mass spectrometer (Q Exactive) (Thermo Fisher Scientific, Waltham, WA) equipped with an EASY-nLC 1000 system (Thermo Scientific, Waltham, MA) as previously described (55). An MS scan [400 to 1200 mass/charge ratio (m/z)] was recorded in the Orbitrap mass analyzer set at a resolution of 70,000 at 200 m/z, 1 × 106 automatic gain control target, and 100-ms maximum ion injection time. Database searching Protein identification was performed using the Homo sapiens taxonomy (20,200 sequences) setting of the Swiss-Prot database (SwissProt_2015_06) with Proteome Discoverer 2.1 (version 2.1.1.21, Thermo Scientific) as previously described (53). Quantification of MS1 precursor ions was analyzed using the Proteome Discoverer 2.1 (version 2.1.1.21, Thermo Scientific). A relative quantification approach was used for deamidation measurement as follows: The ratio of a deamidated peptide/total peptide was used to represent the abundance of the deamidated peptide. The ratio of a specific deamidated peptide across the samples was standardized (z score); the mean standardized values of the peptides derived from one protein were used to represent the level of deamidated protein abundance. miRNA quantification by real-time polymerase chain reaction Matched frozen cartilage tissue from patients used for proteomic analysis was dissected transversely from the cartilage surface. Total RNA was extracted using the miRNeasy Mini Kit (Qiagen) following the manufacturer’s instruction with modifications. Reverse transcription was conducted with miScript II RT Kit (Qiagen) following kit instructions. Quantification of blastema miRNAs (miR-21, miR-31, and miR-181c) and internal control miRNAs (miR-423, miR-191, and miR-U6) was assessed by quantitative real-time polymerase chain reaction (PCR) using the miScript SYBR Green PCR Kit (Qiagen). The ΔCt values for blastema miRNAs, representing the miRNA normalized to the mean of the internal control miRNAs, were used to evaluate the correlation of miRNA and cartilage ECM protein expression; these ΔCt values were computed as the Ct of the blastema miRNA minus the mean Ct of the internal control miRNAs. Deamidation rate prediction and protein half-life estimation Half-life determinations require an estimate of deamidation rate constants based on the local 2D and 3D protein structure. Protein crystal structures were available through the Protein Data Bank (www.rcsb.org/pdb) for portions of some proteins (FN1 and COMP) or generated using the I-TASSER (http://zhanglab.ccmb.med.umich.edu/I-TASSER) (24, 56). Swiss-PdbViewer (v4.10) (57) was used to visualize the structures and manually acquire the steric covariate values for estimating deamidation rate (16, 34). Deamidation, like racemization, is a first-order process (17). Therefore, accumulation of deamidated residues can be described by two linear processes: the rate constants of deamidation (K id ) and the protein turnover (K T ). The algorithm can be described as d C A s n ( D / D + N ) ) / d t = K i d C A s n N / ( D + N ) - K T C A s n N / D / ( D + N ) where C Asn is the molar concentration of Asn in the protein, and D/(D + N) and N/(D + N) represent the fractions of the deamidated and nondeamidated native proteins. The measured amount of deamidated residues in the sample is the net result of accumulated deamidated residue accumulation and their loss due to protein turnover. The original estimation of deamidation rates was reported in days (17). Therefore, the deamidation rates derived from our algorithm also roughly correspond to units of days. Quantification of DCOMP by immunoassay Proteins extracted, as described above from cryosections (but without reduction, alkylation, or trypsin digestion), were precipitated with ice-cold EtOH (95% ethanol with 50 mM NaAc). DCOMP was quantified using a previously described sandwich ELISA assay (14). Data analysis Since our experimental design contains repeated measures from three types of cartilage depth per individual specimen, multivariable mixed models for continuous outcome variables (e.g., deamidated protein ratio and miRNA abundance) include multiple independent factors, joint site (hip/knee/ankle), disease state (healthy/OA), depth (superficial, middle, and deep), age, and gender of the cartilage samples. The first-order autoregressive [AR(1)] was used as a variance-covariance matrix to account for within- and between-specimen correlations. The primary predictors of interest were joint site and cartilage depth. Similarity measurement of correlations of miRNA and cartilage proteins was performed using Cluster 3.0 (58). Protein-protein interactions were visualized using STRING (59). Statistical significance of each factor and the overall model was reported at the nominal significance level (P < 0.05). The multivariable mixed-model analyses were performed using JMP Pro 13 (SAS, Cary, NC). All the graphs were prepared in Microsoft Excel, PowerPoint 2013, or GraphPad Prism 8.

SUPPLEMENTARY MATERIALS Supplementary material for this article is available at http://advances.sciencemag.org/cgi/content/full/5/10/eaax3203/DC1 Supplementary Materials and Methods Fig. S1. Disease-related distributions of deamidated proteins in human adult cartilage. Fig. S2. The correlation of miR-31 with total protein and deamidated cartilage protein abundance. Fig. S3. Classification of cartilage proteins correlated with miR-21.

This is an open-access article distributed under the terms of the Creative Commons Attribution-NonCommercial license, which permits use, distribution, and reproduction in any medium, so long as the resultant use is not for commercial advantage and provided the original work is properly cited.

Acknowledgments: We thank B. Donald and S. Endo-Streeter for helping with the initial protein structure modeling. We also thank Y.-J. Li, the Associate Director of the Biostatistics Core of the Duke Translational Medicine Institute, for advice related to the statistical modeling for this study. Funding: This study was supported by an OARSI Collaborative Scholarship to M.-F. H., a Collaborative Exchange Award from the Orthopaedic Research Society to V.B.K., and NIH/NIA P30-AG-028716. Mass spectrometers were funded by the Crafoord Foundation and the Inga-Britt and Arne Lundberg Foundation, and other financial support was obtained from the Swedish Research Council (2014-3303) to P.Ö. Author contributions: M.-F.H., P.Ö., and V.B.K. were involved in the conception and design of the study as well as the interpretation of the data. M.E.E. and M.P.B. provided the cartilage specimens. M.-F.H. and P.Ö. conducted the MS experimental work. M.-F.H. conducted the analysis and drafted the manuscript. All authors critically revised the article and approved the final version for submission. Competing interests: The authors declare that they have no competing interests. Data and materials availability: All data needed to evaluate the conclusions in the paper are present in the paper and/or the Supplementary Materials. Additional data related to this paper may be requested from the authors.