Transcriptomic and proteomic studies of condylar ossification of the temporomandibular joint in porcine embryos

Abstract Background The ossification mechanism of the temporomandibular joint (TMJ) condyle remains unclear in human embryo. The size and structure of TMJ, shape of articular disc and the characteristics of omnivorous chewing in the pig are similar to those of humans. The pig is an ideal animal for studying the mechanism of ossification of the TMJ condyle during the embryonic period. Method In a previous study by our group, it was found that there was no condylar ossification on embryonic day(E) 45, but the ossification of condyle occurred between E75 and E90. In this study, a total of 12 miniature pig embryos on E45 and E85 were used. Six embryos were used for tissue sections (3 in each group). The remaining six embryos were used for transcriptomic and proteomic studies to find differential genes and proteins. The differentially expressed genes in transcriptome and proteomic analysis were verified by QPCR. Results In total, 1592 differential genes comprising 1086 up‐regulated genes and 506 down‐regulated genes were screened for fold changes of ≥2 to ≤0.5 between E45 and E85. In the total of 4613 proteins detected by proteomic analysis, there were 419 differential proteins including 313 up‐regulated proteins and 106 down‐regulated proteins screened for fold changes of ≥2 to ≤0.5 between E45 and E85. A total of 36 differential genes differing in both transcriptome and proteome analysis were found. QPCR analysis showed that 14 of 15 selected genes were consistent with transcriptome analysis. Conclusion Condylar transcriptome and proteomic analysis during the development of TMJ in miniature pigs revealed the regulatory genes/proteins of condylar ossification.


| INTRODUC TI ON
The temporomandibular joint (TMJ) is a synovial, double, and bicondylar joint, which is composed of mandibular condyle, articular disc, articular fossa, superior and inferior articular cavities and associated muscles and tendons. 1,2 The mandibular condyle which is mainly formed by endochondral ossification is the growth and development center of the mandible. The condyle plays an important role in the development of the cranio-maxillofacial region, with abnormal development leading to underdevelopment or overdevelopment of the mandible. Its growth and development consists of the proliferation, maturation, ossification and new bone deposition of the condylar cartilage. 3 The TMJ condylar cartilage is different from the cartilage of long bone joint. It is secondary cartilage affected by both local factors and growth factors. The mechanism of secondary cartilage growth and cartilage formation depends largely on the external environment. The cell source of condylar cartilage is also different from that of long bone articular cartilage. TMJ condylar cartilage is fibrocartilage, developed from embryonic ectoderm and characterized by type I collagen, while long bone joint cartilage comes from the mesoderm and is hyaline cartilage, mainly type II collagen and mucopolysaccharide. 4,5 Studies of embryonic morphological development of TMJ began as early as 1952, but studies of the ossification mechanism of TMJ condylar cartilage lagged behind. There are few studies on the molecular regulation mechanism of TMJ development in the embryonic stage. These studies were based on rodents. Some classical signal pathways and transcription factors were discovered such as Wnt/ ß-catenin, IHH, TGF-ß, SOX9 and RUNX2. [6][7][8][9][10][11][12][13][14][15] However, there are great differences in TMJ structure between rodents and humans.
For example, compared with humans, TMJ articular fossa in rodents is shallow and flat, there are no articular nodules, the volume and function of the lateral pterygoid muscle are smaller, in the process of development the superior articular cavity is formed first, followed by the inferior articular cavity, and the TMJ articular disc rarely becomes fibrocartilage with increasing age. These differences have undoubtedly delayed progress in TMJ-related research. 16 In recent years, the whole genome of the pig has been sequenced. It has been found that the genomic composition of the pig is highly homologous to that of humans and the physiological and anatomical characteristics of pig are very similar. For example, the size, structure of TMJ, shape of the articular disc and characteristics of omnivorous chewing are very similar to those of humans. Pig epiblast stem cells are excellent models for discovering the properties of human pluripotent stem cells. Chondrogenesis of the TMJ originates from the ectoderm and pig ectodermal stem cells can mimic human ectodermal cell development. 17 As a model animal for TMJ studies, the pig has obvious advantages over rodents (rat and mice), carnivores (dog and cat) and herbivores (cattle, sheep and rabbits). 18 However, the signal molecules and transcription factors related to chondroossification in porcine TMJ embryo development have not been reported.
In previous experiments of this study, involving histological analysis of pig TMJ development at embryonic days (E) 35, E45, E55, E75, E90 and postnatal day (P) 1, it was found that day E45 and the days between E75 and E90 were the key periods of condylar cartilage development and ossification. 19 In this study, transcriptome and proteomic analysis were used to find the key regulatory factors of condylar ossification and the signal pathways related to condylar cartilage differentiation.

| ME THODS
The experiments were conducted according to the guidelines of the Ethical Review of Laboratory Animal Welfare (GB/T358922018).

| Histological staining
The embryos were dissected to fully expose the TMJ. The whole TMJ (including condyle, articular disc, articular fossa, superior and inferior articular cavities) was obtained from the zygomatic arch to the mandible. Twelve TMJs obtained from E45 and E85 embryos (n = 3 for each group) were fixed in 10% formalin for 2 days. E85 embryos (n = 3) were decalcified in 0.5 mol L −1 ethylenediaminetetraacetic acid solution (Solarbio Biological Company) for 14 days and rinsed with tap water for 10 min. Temporomandibular joint tissues from E45 and E85 embryos were embedded in paraffin for analysis.
Semi-serial sagittal sections (3 μm) of the TMJ were prepared using a microtome (HM 325; Thermo Fisher Scientific, Waltham, MA, USA). The samples were stained with Hematoxylin-Eosin (HE) and Safranin O-Fast green for analysis. The HE staining was performed using a Leica AutoStainer XL (Leica Biosystems, Wetzlar, Germany).
Safranin O-Fast green staining was performed following the instructions for the kit. The iview6.3.6 software was used to analyze the tissue sections.

| Transcriptome and proteomic analysis
Twelve TMJs were obtained from E45 and E85 embryos (n = 3 for each group). The condyle was excised between the lowest point of the sigmoid notch and the highest point of the condyle. The condyle sample from each embryo was divided into two parts.
The first portion of tissue was used to extract RNA (RNeasy Mini Kit, Cat#74106, Qiagen), construct RNA library and Illumina sequencing (Shanghai Huaying Biopharmaceutical Technology Co., Ltd). Genomic mapping analysis of the sequencing result was performed using HISAT2 software. The reference genome version is Sus_scrofa.Sscrofa10. The second portion of tissue was used for protein extraction. The protein samples obtained were redissolved, subjected to reductive alkylation, enzymatic hydrolysis, desalting and tandem mass spectrometry. The proteome Discoverer 2.4 software was used to search the database to identify the proteome. Analysis of protein expression used the relative quantitative method of Label free (Shanghai Huaying Biopharmaceutical Technology Co., Ltd.). Student's t test was used to calculate the p value of the two groups. The intergroup ratio of fold changes ≥2 or ≤0.5 was selected as indicating for differential gene/protein. The genes with differences in transcriptome and proteology levels were screened as candidate genes.
Differential genes/proteins between groups were analyzed by GO analysis and KEGG analysis.

| Verification of differential genes
Fifteen different genes were selected for QPCR verification. The existing literature shows that these 15 differential genes were associated with bone development, regeneration and pathogenesis of osteoarthritis. They may also be associated with ossification of the temporomandibular joint condyle during the embryonic stage.
Total RNA was extracted (RNeasy Mini Kit, Qiagen) from condylar tissue. The RNA was reverse transcribed into cDNA using the PrimeScriptTM RT regent Kit with gDNA Eraser (Takara, RR047A).
The mRNA sequence information of the 15 candidate genes was downloaded from NCBI, and the PCR primers were designed using the BLAST module in the NCBI database (

| Histological staining
The crown-rump lengths (CRL) of E45 and E85 embryos were respectively 6.5 cm and 14.5 cm. The complete TMJ structures including condyle, articular disc, articular fossa, superior and inferior articular cavities were observed by histological anatomy (Figure 1).
Histological staining showed that the condyle is separated from the temporal bone (Figure 2A,B,a,b) at E45. The articular disc is not fully developed (Figure 2A,B,f) at this time. The inferior articular cavity is formed (Figure 2A,B,c), but the superior articular cavity is not fully formed -only cracks are visible at both ends of the articular disc ( Figure 2A

| Transcriptome analysis
The sequencing results of the E85 and E45 condylar tissue transcript groups were to analyze porcine genomic mapping using HISAT2 software. A total of 25 322 genes were detected, and 19 584 differential genes between E85 and E45 were found. Of these, 1592 genes were consistent with intergroup ratio fold changes of ≥2 and ≤0.5. These genes were screened and a volcanic map was drawn ( Figure 3A). Among the 1592 differential genes, 1086 genes were upregulated at E85 compared with E45. A total of 506 genes were downregulated (the top 15 genes with upregulation and downregulation differences are shown in Table 2). The genes with similar expression patterns are clustered together to draw the cluster analysis map of different genes ( Figure 3B).
The differential genes were classified by GO analysis as cell com- Similar to the GO classification statistics, the number of differentially expressed genes in each pathway category of KEGG was counted ( Figure 3E). The differential genes were mainly concentrated in the signal transduction, metabolism and immune system categories. The differentially expressed genes were analyzed by pathway using the Fisher precise test. p values were obtained using the Fisher accurate test calculation. Pathway with p < 0.05 values indicated a significantly enriched pathway ( Figure 3F). The KEGG enrichment analysis showed that differential genes were mainly enriched in cell phagosome (ssc04145), cell adhesion factor (ssc04514), rheumatoid arthritis pathway (ssc05323) and osteoclast differentiation (ssc04380) categories.

| Proteomic analysis
A total of 4613 proteins were obtained from E45 and E85 condylar tissue by protein spectrum analysis. The normalized signal mean of all samples in the two groups was calculated. The inter-group ratio fold change was calculated, and Student's t test was used to calculate the p value of the two groups. Proteins with ratios of fold changes of ≥2 and ≤0.5 were regarded as differential proteins. There were 419 differential proteins including 313 proteins that were upregulated between E85 and E45 and 106 proteins that were downregulated.
The top 15 differential proteins expressed are shown in Table 3.
A volcanic diagram and a cluster analysis map of the 419 differential proteins were drawn ( Figure 4A,B), and the 10 items with the smallest p values were selected to draw Barplot map (less than 10 items not shown, Figure 4C). GO enrichment analysis of differential proteins before and after condylar ossification showed that differential proteins were mainly enriched in the process of small molecule metabolism in BP, cytoplasm in CC, and catalytic activity and ion binding in MF. KEEG pathway enrichment analysis ( Figure 4D) showed that the differential proteins were mainly concentrated in the metabolic pathway (ssc01100) and oxidative phosphorylation pathway (ssc00190) categories.

| Cross-analysis in transcriptome and proteomics
Thirty-seven differential genes including 29 upregulated genes/ proteins and 8 downregulated genes/proteins differing in both the transcriptome and proteomic analysis were found. It is interesting that the expression of IHH and TNMD genes was upregulated in the proteomics, but downregulated in transcriptome. The results of the cross-analysis of transcriptome and proteomics results are shown in Table 4.

| Quantitative real-time PCR
Fifteen differential genes were selected and verified by QPCR. The

| DISCUSS ION
The TMJ condylar cartilage and its subchondral bone formation is an important and understudied topic in dental research. The transition from chondrocytes to osteocytes in articular cartilage of long bone has been reported. 20 There are still contradictory views on the biomolecular regulation mechanism that controls condylar cartilage replaced by bone tissue in the process of natural growth. At present, there are two main theories about the cause of cartilage osteogenesis. One explanation is that bone tissue is formed by programmed cell death (apoptosis) of hypertrophic chondrocytes, followed the MMp-13 and MMP-9 absorb the laterally mineralized matrix, finally potential bone marrow and blood vessel formed.
The invading capillaries bring osteogenic progenitor cells and bone marrow stem cells that differentiate into osteoblasts. [20][21][22] Another theory is that hypertrophic chondrocytes are transformed into osteoblasts without apoptosis. 23 In this study, we used transcription and proteomics analysis to detect the molecular mechanism of condylar ossification in order to find the key regulatory molecules of condylar ossification.
In this study, histological observations showed that the condyle at E45 was occupied by a large number of chondrocytes, and scattered hypertrophic chondrocytes could be seen at the bottom of condyle, but no mineralized bone tissue was found. In contrast, at E85 a layer of mineralized bone tissue was formed at the base of the condyle. This difference in bone morphology is accomplished under a specific molecular regulation mechanism. In order to find these important ossification regulatory molecules, we used transcriptome analysis to identify 1592 differential genes (including 1086 upregulated and 506 downregulated genes) during the development of TMJ between E45 and E85. Among these differential genes, the expression of the top 15 genes was more than 32-fold different between E45 and E85. ( Table 2). The differential genes between E85 and E45 may play important roles in bone formation, but there are many biological processes such as splicing and modification in gene expression that prevent these genes from being translated into proteins to perform biological functions.
Therefore, proteomic analysis was additionally performed in this study to more accurately identify target genes/proteins. In total, 419 differential proteins were found by proteomic analysis, and the expression of the top 15 differential proteins was more than 10 times ratio at E85 than at E45. (Table 3). Cross-analysis of all the differential genes and proteins showed that there were a total of 37 differential genes ( Table 4).
The molecular regulation mechanism of embryonic condyle ossification has not previously been reported. This preliminary study of the regulation mechanism of TMJ condyle ossification identified 37 genes which may play an important role in the formation of condylar

TA B L E 3
The top 15 differentially expressed proteins at E85 and E45.
promotes osteoblast differentiation, but also inhibits osteoblast maturation to maintain normal bone homeostasis. Gli1 + cells gradually produce osteoblasts in all bone sites, and can mediate theTGF-ß signaling pathway and play an important role in endochondral bone formation after fracture. 28 Some studies found that Gli1 + cells possess mechanical stress, which may regulate changes in intracellular F I G U R E 4 Proteomic analysis of condylar tissue at embryonic day (E) 85 and E45. (A) Differential protein volcano map. The ordinate is -log10 (p value) value, and the abscissa is log2 (fold change) value. Red dots, up-regulated differential proteins; green dots, down-regulated differential proteins. Other proteins are represented by red and black dots. (B) Differential protein cluster analysis map. The abscissa represents the sample name between the groups, the ordinate represents the differential genes, red represents a high expression value of the differential genes in the grouping sample, and blue represents a low expression value of the differential protein in the grouping sample. (C) GO enrichment analysis of differential proteins. The vertical axis is numbers of genes, and the horizontal axis is the entry. (D) KEGG Pathway analysis. The horizontal axis is the GeneRatio. The vertical axis is the entry. The color of the points indidcates the p value. The size of the point indicates the number of genes.
Ca 2+ concentration through the inositol 1,4,5-trisphosphate receptor (IP3R) to respond to mechanical force. 29,30 Compared with the above four genes, CLEC3B an CRABP2 have been less well studied in osteogenesis. CLEC3B is a major structural protein of the ocular lens and is abundant in the kidney, central nervous system, skeletal muscle and heart. 31 In a study of osteoarthritis (OA), it was found TA B L E 4 Genes differing in transcriptome and proteome analysis.  34 In this study, it was found that the expression of CRYAB, CLEC3B, SPP1 and ENPP1 at E85 were significantly higher than that at E45 and the development of condyle, and the amount of mineralized bone tissue in E85 was also significantly higher than that in E45. Gli1 was significantly downregulated at E85, which was inconsistent with the results of postnatal Gli1 in osteoblast proliferation and differentiation. It may be that Gli1 may be regulated by different signal molecules in the embryonic stage than in the postnatal stage, and its regulatory role in the development of TMJ at the embryonic stage remains to be further researched.
In addition to genes related to osteogenesis, genes related to chondrocytes and osteoclasts were also found, such as Matn1,  37,38 In our study the expression of the gene/protein at E85 was significantly higher than that at E45, and the protein expression was more than 40 times higher, suggesting that G6P and Pi exchange was active and metabolism was enhanced in E85 condylar cells. ACP5 is involved in the production of reactive oxygen species, normal bone development, osteoblast regulation and macrophage function. 39 Medical studies have found that mutations in ACP5 can cause metaphyseal dysplasia, autoimmune regulation disorders and nerve damage. [40][41][42] However, the role of ACP5 in TMJ osteogenesis has not been reported. The role of ATP6V0d2 in osteoclast formation has been studied. It is not only regulates cell fusion in osteoclast differentiation, but also plays an important role in osteoclast specific proton pump mediated by bone resorption and extracellular acidification. 43 The upregulation of Atp6v0d2 is beneficial to the formation of osteoclasts. In line with these findings, mice with ATP6V0d2 gene deletion showed increased bone mass, but knockout of ATP6V0d2 gene in ovariectomized mice did not inhibit the occurrence of osteoporosis. 44 In our study it was found that the expression of SLC37A2, Acp5 and ATP6V0d2 at E85 were upregulated compared with E45, especially in the proteomic analysis, which showed that both SLC37A2 and Acp5 were upregulated more than 40-fold, indicating that a large number of osteoclasts were produced compared with E45. The osteoblasts and osteoclasts were associated with each other to maintain the dynamic balance of bone metabolism.
Studies on the molecular regulation mechanism of TMJ condyle ossification in large mammal embryos have not been reported. Here we report a preliminary study of the key genes regulating condyle ossification of TMJ in pig embryo. In this study, 15 differential genes associated with bone development, regeneration and pathogenesis of osteoarthritis through literature survey were selected for QPCR verification. The other 22 genes were not verified by QPCR, but the future study based on the validation results for the 15 genes will verify other genes. In follow-up studies, all the 37 differential genes will be verified by QPCR, immunohistochemistry and western blot in TMJ of minipigs of different gestational ages. Insufficient sample size of embryos at E45 and E85 is another limitation of this study. In future research, our research group will expand the sample size to exclude differences between different litters. We believe that with further in-depth studies, the discovery of molecules regulating the ossification of the TMJ condyle during the embryonic period will provide a basis for the treatment of temporomandibular joint diseases and tissue engineering research.

CO N FLI C T O F I NTER E S T S TATEM ENT
The authors have no potential conflicts of interest to declare with respect to the research or the authorship and/or publication of this article.

EH I C S S TATEM ENT
This study does not contain any studies with human subjects. It was approved and supervised by the Institutional Animal Care and Use Committee of Chinese PLA General Hospital (Beijing, China; approval document no. 2020-X16-109).