Transcriptomic profiling of feline teeth highlights the role of matrix metalloproteinase 9 (MMP9) in tooth resorption

Tooth resorption (TR) in domestic cats is a common and painful disease characterised by the loss of mineralised tissues from the tooth. Due to its progressive nature and unclear aetiology the only treatment currently available is to extract affected teeth. To gain insight into TR pathogenesis, we characterised the transcriptomic changes involved in feline TR by sequencing RNA extracted from 14 teeth (7 with and 7 without signs of resorption) collected from 11 cats. A paired comparison of teeth from the same cat with and without signs of resorption identified 1,732 differentially expressed genes, many of which were characteristic of osteoclast activity and differentiation, in particular matrix metalloproteinase 9 (MMP9). MMP9 expression was confirmed by qPCR and immunocytochemistry of odontoclasts located in TR lesions. A hydroxamate-based MMP9 inhibitor reduced both osteoclast formation and resorption activity while siRNA targeting MMP9 also inhibited osteoclast differentiation although had little effect on resorption activity. Overall, these results suggest that increased MMP9 expression is involved in the progress of TR pathogenesis and that MMP9 may be a potential therapeutic target in feline TR.

Tooth resorption (TR) in domestic cats is a common and painful disease characterised by the loss of mineralised tissues from the tooth. Due to its progressive nature and unclear aetiology the only treatment currently available is to extract affected teeth. To gain insight into TR pathogenesis, we characterised the transcriptomic changes involved in feline TR by sequencing RNA extracted from 14 teeth (7 with and 7 without signs of resorption) collected from 11 cats. A paired comparison of teeth from the same cat with and without signs of resorption identified 1,732 differentially expressed genes, many of which were characteristic of osteoclast activity and differentiation, in particular matrix metalloproteinase 9 (MMP9). MMP9 expression was confirmed by qPCR and immunocytochemistry of odontoclasts located in TR lesions. A hydroxamate-based MMP9 inhibitor reduced both osteoclast formation and resorption activity while siRNA targeting MMP9 also inhibited osteoclast differentiation although had little effect on resorption activity. Overall, these results suggest that increased MMP9 expression is involved in the progress of TR pathogenesis and that MMP9 may be a potential therapeutic target in feline TR.
Tooth resorption (TR) of permanent teeth is the pathological loss of dental hard tissues and has been reported in many animals including humans, dogs, horses, domestic cats and wild cats [1][2][3][4][5] . Affecting multiple teeth, resorption is common in aged domestic cats but is rarer in other animals [6][7][8] . Destructive and progressive lesions usually initiate on the external surface of the coronal tooth root and extend into the cementum, enamel, dentine and in advanced stages, the pulp 7,9,10 . TR largely classifies into three Types (Type 1-3) by the Nomenclature Committee of the American Veterinary Dental College (AVDC). Type 1 is characterised by appearance of focal or multifocal lesions of the tooth with normal periodontal ligament and inflammatory changes while Type 2 presents resorption lesions with involvement of periodontal ligament and non-inflammatory replacement resorption 11 . Type 3 exhibits both features of Type 1 and Type 2. TR is caused by tooth-resorbing odontoclasts, cells functionally similar to bone-resorbing osteoclasts 7 . Both cell types are derived from circulating blood or bone marrow myeloid lineage progenitors which differentiate into mature multinucleated cells in their respective microenvironments 10,[12][13][14][15] . The main physiological role of odontoclasts is to resorb hard tissues around deciduous teeth allowing tooth eruption to commence during tooth development and to resorb the deciduous teeth for permanent tooth replacement 16 . In some cats these cells become dysregulated and attack the permanent teeth later in life 9,17 . Due to recognition of this disease and investigations into the prevalence of TR between the 1990s and 2000s, the aetiology of TR has been widely studied 9,18 . The central axis of differentiation of osteoclasts and odontoclasts is the pathway comprising colony stimulating factor 1 (CSF-1) and the receptor activator of nuclear factor κ B (RANK), RANK ligand (RANKL), and osteoprotegerin (OPG), a decoy receptor, capturing RANKL 12,15 . Since the microenvironments of physiological and pathologic TR are different, there has been a focus on identifying co-stimulating factors in pathological conditions. The accumulated evidence suggests that the initiation and progression of TR is complex and multifactorial and is also mediated by local rather than systemic mechanisms. Many possible molecules and risk factors that directly or indirectly impact on the RANK/ RANKL/OPG axis pathway have been postulated, including pro-inflammatory cytokines 19,20 , diet 21,22 , thyroid Tissue processing and RNA extraction. Samples including teeth, mandibles and maxilla were snap frozen in liquid nitrogen and stored at − 80 °C until RNA extraction. Information of samples are described in Table 1. Teeth were extracted from the alveolar sockets using dental equipment or using bone cutters while maintaining cold temperatures by working over dry ice. Tissue processing and RNA extraction protocols were optimised based on a guanidinium thiocyanate-phenol-chloroform extraction method and manufacturer's instructions for the Qiagen RNeasy Mini Kit. This protocol has been previously documented 28 . cDNA library and RNA sequencing. Tested RNA samples were sent to our academic facility (Edinburgh Genomics, The University of Edinburgh, UK) for cDNA library production and RNA sequencing. Briefly, one microgram of RNA from each of the thirteen samples and 0.8 µg of RNA from one low quantity sample (N27C23TR+) were used to generate fourteen cDNA libraries using the Illumina TruSeq stranded mRNA sample preparation kit according to the manufacturer's instructions (Illumina, San Diego, CA, USA). Paired-end sequencing was performed using the Hiseq 4000 system (Illumina). The low quality reads were filtered by Phred quality score (Q score 30) and 3′ adapter were trimmed with cutadapt (version 1.8.3). All the raw reads have Mapping to the feline reference genome and generation of read counts. Alignments to the Felis catus genome (version 6.2) were performed using STAR (version 2.5) 29 . The read counts from each sample were generated using HTSeq (version 0.6.0.1) 30 with mode 'union' . Duplicate reads were found using picard tools (version 1.141). Generation of MDS plots were generated using plotMDS function from edgeR package (version 3.12.1) to visualise the level of similarity of individual cases of a dataset. From the 14 RNA-seq samples, two samples N27C23TR+ and N08C8TR+ were identified as outliers in MDS plots of RNA-seq data therefore they were excluded for final analysis (Fig. S1) 28 .
Differential expression analysis. 'TR − ve cats' are defined as teeth collected from a cat without any TR lesion and 'TR + ve cats' are teeth collected from a cat with at least one TR lesion. 'TR tooth' compares TR status of individual teeth, hence TR − ve teeth were collected from both TR free cats and TR diagnosed cats as long as the tooth was free of TR when assessed on dental radiographs. 'Paired TR' compares TR − ve teeth and TR + ve teeth from the same individuals, so a 'within cat' comparison. From the remaining 12 libraries, three comparisons were performed: (1) TR − ve cats (n = 3) and TR + ve cats (n = 5), (2) TR − ve teeth (n = 7) and TR + ve teeth (n = 5), and (3) paired TR, TR − ve (n = 3) and TR + ve (n = 3). Differentially expressed (DE) genes were identified through comparisons between TR − ve and TR + ve groups using the R package edgeR (version 3.12.0) 31 . The trimmed mean of M-values normalization method was used. To explore sample dispersion between TR − ve and TR + ve groups, distance plot of biological coefficient of variation samples were generated and outliers were excluded for DE gene analysis. Significantly DE genes were identified using the threshold of false discovery rate (FDR)-adjusted p value < 0.05.

Enrichment analyses.
We performed two enrichment analyses on the up-and downregulated DE genes: of KEGG (Kyoto Encyclopaedia of Genes and Genomes) metabolic pathways 32 and GO (Gene Ontology) terms 33 . GO term enrichment was assessed using the R/Bioconductor package topGO v2.38.1 34 . KEGG pathway enrichment was performed using the R/Bioconductor package SPIA (Signalling Pathway Impact Analysis) v2.38.0 35 .
Quantitative PCR. Real-time PCR was performed using the Takyon Low ROX SYBR 2X MasterMix kit (Eurogentec) in the Stratagene MX3000P qPCR system (Agilent Technologies). Each reaction contained 300 nM of gene of interest or reference gene primer, 2 × MasterMix, 4 μl of diluted cDNA and nuclease free water, producing a final volume of 20 μl in a 96 well plate format. All samples were performed in duplicate and cycling threshold (Ct) values for each gene were generated. Calculation of relative gene expression analysis was based on delta-delta method 36,37 . The primers used are described in Table S1.
Generation of feline osteoclasts in vitro. The protocol was modified from a published rodent osteoclast culture protocol 38 . The long bones were dissected from cats at post mortem and cut off at the proximal and distal ends. Bone marrow was flushed out in α-MEM containing 10% FBS followed by red cell lysis in Lysis Buffer (150 mM NH 4 Cl, 10 mM NaHCO 3 , 0.1 mM EDTA) for 5 min, on ice. The bone marrow cells were resuspended in α-MEM containing 10% FBS, 2 mM L glutamine, 100 IU/ml benzyl penicillin, 100 mg/ml streptomycin, and CSF-1 (recombinant mouse CSF-1, 10 ng/ml, R&D Systems) and plated on a T75 flask for 24 h at 37 °C in 5% CO 2 /95% air. On the following day, non-adherent cells were collected and plated on a dentine disc in each well of a 24 well plate (dentine discs were a generous gift from Professor Timothy R. Arnett, University College London, London, UK) or Corning Osteo Assay Surface multiple well plates (1 × 10 5 cells/cm 2 ) with CSF-1 (10 ng/ml) at 37 °C in 5% CO 2 /95% air. On day 4, 90% of the medium was changed and, replaced with a new α MEM containing CSF-1 (10 ng/ml) and RANKL (recombinant mouse RANKL, 3 ng/ml, R&D Systems). The culture media was subsequently changed twice weekly for up to 10-14 days.
Characterisation of cultured osteoclast culture. Cultured osteoclasts on a cell culture plastic or dentine discs were washed twice with phosphate buffered saline (PBS) and fixed in 2.5% glutaraldehyde for 5 min. Tartrate resistant acid phosphatase (TRAP) activity was determined using Acid Phosphatase, Leukocyte (TRAP) Kit (Sigma-Aldrich) according to the manufacturer's protocol with some modifications including reduced incubation time (30 min) and without haematoxylin counter staining. TRAP positive cells were imaged via light microscopy (Leica DM LB2, Milton Keynes, UK). Cells were considered to be osteoclast-like when TRAP positive and with more than two nuclei when manually counted. Resorption pits in dentine discs were visualised and counted under light microscopy after removal of cells by sonication for 5 min in 0.25 M NH 4 OH followed by 2 min toluidine blue staining.
Inhibition of MMP9 in osteoclastogenesis. The MMP9 inhibitor containing a hydroxamate (-CON-HOH) group was purchased from Calbiochem (Merck Millipore, UK) and was reconstituted in DMSO (1 mg/ ml stock solution) and stored at − 20 °C. Primary feline osteoclasts derived from bone marrow were treated with a range of concentrations from 5 nM to 1 µM diluted in α-MEM until termination of osteoclast differentiation. Vehicle controls were treated with DMSO concentrations corresponding to the highest drug concentration used in the experiment. Untreated control cells were also included in each experiment. www.nature.com/scientificreports/ Targeting of siRNA of MMP9. The siRNA sequence against the feline MMP9 was designed using an online design tool, i-Score designer (i-score web service program). MMP9 feline siRNA targeting the MMP9 mRNA sequence (CCA GGA GAC TTG CGA ACT A) was constructed using the Silencer siRNA Construction Kit (Ambion) according to the manufacturer's protocol. Osteoclast precursors were incubated for two days in α-MEM containing 10% FBS, 2 mM L glutamine, 100 IU/ml benzyl penicillin, 100 mg/ml streptomycin, and CSF-1 (10 ng/ml) at 37 °C/5% CO 2 . Cells were detached and resuspended for electroporation in 100 μl Mouse Macrophage Nucleofector solution and siRNA. The cell suspension-DNA mixture was transferred to Amaxa electrode cuvettes and electroporated in an Amaxa Nucleofector II (Lonza, UK) using program Y-001. Following electroporation, cells were immediately suspended in medium containing CSF1 (10 ng/ml) and RANKL (3 ng/ ml) and seeded into 96-well plates at 1 × 10 5 cells/well or onto dentine discs. Cells were incubated at 37 °C/5% CO 2 until osteoclast formation was observed in control wells. Media were changed every other day.
Immunohistochemistry of feline teeth. Feline tooth samples were collected and fixed by immersion in 10% neutral buffered formalin for 48 h before being processed. All fixed teeth were demineralised in EDTA pH 7.0 for 4-6 weeks at room temperature (RT). Demineralised samples were dehydrated and embedded in paraffin wax following standard procedures. Haematoxylin and Eosin (H & E) staining was performed using an autostainer (Leica Autostainer XL). For MMP9 immunostaining, slides were dewaxed and rehydrated. Antigen retrieval was performed in 0.1 M citrate buffer pH 6.0 for 90 min at 70 °C. Endogenous peroxidase was blocked using 1% hydrogen peroxidase in methanol for 30 min. A further blocking was carried out using a normal goat  Table 2). All expressed genes and significantly differentially expressed genes are listed in Table S2. Additionally, twenty seven differentially expressed genes were common for all three comparisons, including 25 up-regulated genes and 2 down-regulated genes (Table 3).
Pathway analysis revealed up-regulated genes that are involved in osteoclast differentiation and calcium signalling. For further characterisation of the transcriptomic changes involved in TR, enrichment analyses were performed on the set of 1,732 genes DE between paired TR groups. KEGG downstream pathway analysis identified a total of 28 metabolic pathways within which these genes were involved (Table S3). Of these, the set of upregulated genes were enriched among the genes comprising 25 pathways, with the set of downregulated genes enriched only among 3 pathways. The most biologically relevant pathways were osteoclast differentiation ( Fig. 2 and Table S3, KEGG pathway 04380; p = 0.03) and calcium signalling ( Fig. 3 and Table S3, KEGG pathway 04020; p = 1.1 × 10 -4 ) which promote the differentiation of osteoclast precursors and induce osteoclastic tooth resorption in teeth. Relative expression of each gene in each pathway was visualised with colour labels in Figs. 2 and 3. Identified pathways and their active status were revealed based on the expression of genes which constitute that pathway (Table S3). For instance, key genes in the osteoclast differentiation pathway were highly expressed in TR + ve teeth (Table S3, total number of genes in osteoclast differentiation pathway = 93, number of up-regulated differentially expressed genes = 15, down-regulated differentially expressed genes = 3). RANK is significantly up-regulated (Table 4, log2 fold change = 1.157, FDR = 0.0013) in TR + ve teeth. RANKL was not significantly up-regulated (Table S2, log2 fold change = 0.810, FDR = 0.253) but OPG was down-regulated (Table 4, log2 fold change = −1.326, FDR = 7.0 × 10 -8 ), which led to an overall increase in the RANKL/ OPG ratio. Table 4 highlighted differentially expressed genes in osteoclast differentiation and calcium pathway with putative roles. Gene ontology enrichment analysis presented many gene sets involved in muscle physiology including actin binding, myofibril, actin cytoskeleton, muscle contraction, heart contraction, muscle fibre development, skeletal muscle cell differentiation and muscle organ morphogenesis (Table S4). www.nature.com/scientificreports/ RT-qPCR validated the differential expression of the selected genes from RNA-seq result. Twelve candidate genes were functionally associated with osteoclast differentiation and activity (Tables 3 and 4) and were selected for validation by qPCR (  (Table 5). RANKL expression was unchanged but the gene  and P2X4R (Fold change = 2.28). CTSK was highly upregulated in TR + ve teeth by qPCR (p = 0.003) but not by RNA-seq (FDR = 0.059). Overall the gene expression profile of TR identified by RNA-seq was largely consistent with that of the qPCR analysis. Differential expression of MMP9 and P2X4R were re-confirmed by qPCR on a further six paired TR − ve and TR + ve teeth samples although this was not the case for SPI1, RANKL, OPG and RANKL/OPG (Fig. S2).

MMP9 is highly expressed in odontoclasts present in feline tooth resorption pits. MMP9 was
successfully localised in feline dental tissues by immunohistochemistry (Fig. 4). Multinucleated odontoclast-like cells based on H&E staining (Fig. 4D,E) and MMP9 immunolabelled odontoclasts were identifiable in TR + ve sections (Fig. 4E,F) and there was no odontoclasts observed in TR − ve sections (Fig. 4A,B). Quantitative analysis was not carried out due to the limited number of sections and number of active odontoclasts. Since aggregated gingival epithelium also expressed some MMP9 (Fig. 4H), only multinucleated cells within an obvious resorption pit were considered to be odontoclasts. High expression of MMP9 was observed in actively resorbing odontoclasts in TR lesions (Fig. 4F,G). Fibroblasts of the periodontal ligament also expressed MMP9 in both TR − ve (Fig. 4B,C) and TR + ve teeth (Fig. 4F,G). It was also noted that the majority of TR + ve teeth had disrupted periodontal ligament fibres, narrowed or complete loss of periodontal ligament space and replacement with bone like tissues demonstrating a degree of ankylosis.
MMP9 expression was confirmed in cultured feline osteoclasts and increased upon osteoclast formation. MMP9 mRNA expression was significantly up-regulated in the early stages of osteoclast differentiation, with a 7.8 ± 0.2 (p < 0.001) fold change on day 3 compared to day 0. MMP9 mRNA expression remained high until day 6 where after it increased further 8.9 ± 0.2 (p = 0.0387) fold change on day 8 (Fig. 5).
Treatment with a semi-selective MMP9 inhibitor reduced osteoclast differentiation and resorption activity. In vitro osteoclast cultures were treated with a semi-selective MMP9 inhibitor (hydroxamatebased MMP inhibitor). Inhibitory effects on osteoclast formation and resorption activity were observed in a  www.nature.com/scientificreports/ dose dependent manner (Fig. 6). Osteoclast formation was not significantly reduced at low dose (5 nM, the half maximal inhibitory concentration = 5 nM) compared to vehicle control (p = 0.802), but the inhibitor significantly prevented osteoclast formation at higher doses, whereby 50 ± 6.3% (p < 0.0001) and 48.9 ± 5.8% (p < 0.0001) of the number of osteoclasts formed were observed compare to vehicle at 20 nM and 1000 nM treatments, respectively (Fig. 6A). There was no difference in osteoclast numbers when increasing the concentration from 20 to 1000 nM (p = 0.888). All three concentrations of MMP9 inhibitor resulted in reduction of resorption pits compared to vehicle control (Fig. 6B). This inhibitory effect was dose dependent with 5, 20 and 1000 nM treatments producing 62.2 ± 7.4% (p = 0.002), 55.4 ± 8.8% (p = 0.0003) and 12.2 ± 2.5% (p = 0.001) of the number of resorption pits observed in the vehicle control, respectively.

MMP9 specific siRNAs causes a moderate reduction in mRNA levels and reduced the number of osteoclasts formed.
Feline osteoclast precursors were transfected with feline MMP9 siRNA or a scrambled siRNA as control and seeded on dentine discs to investigate the role of MMP9 during osteoclast formation. Following electroporation of the osteoclast precursors with the MMP9 siRNA, a 44.00 ± 0.03% reduction (p = 0.0032) in MMP9 mRNA levels was observed in the transfected cells as assessed by qPCR at 48 h post transfection when compared to negative control transfected precursors (Fig. S3). Electroporation caused some cell death and therefore untransfected cells formed a higher number of osteoclasts and more resorption pits (Fig. 7B,D) but all groups of transfected cells were also able to differentiate into mature osteoclasts. Osteoclast formation was reduced (74.8 ± 7.2%; p < 0.0001) (Fig. 7A,B) in MMP9 siRNA transfected cells in comparison to scrambled siRNA control cells. However, no statistically significant differences in resorption activity was observed between scrambled siRNA and MMP9 siRNA treated wells (p < 0.2298) (Fig. 7).

Discussion
Current understanding of the aetiologies of feline TR is incomplete. Dysregulation of odontoclasts are responsible for TR but the mechanisms underlying their dysregulation are largely unknown. The central aims of this study were to investigate the transcriptomic profile of TR affected cats, identify candidate genes related to feline odontoclast dysregulation and assess their potential as therapeutic targets. First, we performed a RNA-seq to study the whole dental transcriptome of feline teeth and identify differentially expressed genes correlated with progression of TR. No previous RNA-seq analysis of feline TR has been reported, but there are several dental transcriptomic studies using RNA-seq on rodent developing teeth and microarray studies on human teeth that suggested that different cells and tissues expressed unique transcriptomic www.nature.com/scientificreports/ profiles (e.g. between odontoblasts and pre-secretory ameloblasts or deciduous and permanent periodontal ligament and pulp) [55][56][57] .
In this study, we hypothesised that the transcriptome of TR affected teeth is locally different from TR free teeth and compared transcriptomic changes between TR − ve and TR + ve samples. There were only a few changes in gene expression that were common between all three comparisons of TR −/+ ve teeth, TR −/+ ve cat and paired TR − ve and TR + ve teeth and this may be due to large individual variations in the samples available. There were only limited genes identified in common between TR − ve and TR + ve teeth. The majority of remaining genes were uncharacterised or were not reported to be involved in tooth or osteoclast biology. Of these, carbonic anhydrase 4 (CA4) has been reported to be expressed in the plasma membrane of osteoclasts 58 and CA6 is known to be abundantly present in salivary glands and saliva and involved in maintenance of pH homeostasis in oral cavity 59,60 . An oral microenvironment high in carbonic anhydrase may favour osteoclast differentiation 58 . Some homeobox genes (e.g. ARX, HOXA4, LBX1, TLX1) has been reported in tooth development 61 . One of the downregulated genes was SOSTDC1 which has been reported to be largely related to tooth development as an inhibitor of Wnt signalling to control tooth number and morphogenesis 62,63 . The role of these genes in TR + ve teeth are unclear but may suggest that permanent TR + ve teeth may undergo some form of remodelling process 61,64 .
The paired TR comparison gave the most interesting results; a large number of transcripts from the TR + ve teeth differed from TR − ve teeth (n = 1732) which implies that these transcriptomic changes occur locally within the same dentition. The fact that there was transcriptomic similarity between teeth of the same phenotype and that the different phenotypes could clearly be separated reveals that the TR affected teeth undergo similar cellular and molecular changes. The transcriptome corresponding to the phenotype may suggest that transcriptomic changes are at the core of TR pathogenesis, and this knowledge has the potential to increase our understanding of the aetiology of TR.
Previous feline TR studies suggested that the local expression of inflammatory cytokines (e.g. IL1B, IL6, TNF and IFNG) and RANKL may be involved in feline TR, however expression of both pro-inflammatory and inhibitory-inflammatory cytokines increased in TR teeth than control by qPCR which makes it difficult to conclude and the local level of actual cytokines (protein) has however never been measured in TR 19,20 . Therefore overall influence of inflammation cannot be conclusively determined. In human root resorption, P2X7R was considered a stimulator of osteoclast differentiation 65 .
In this study, both RNA-seq data and qPCR result could not identify any significant changes in the expression of inflammatory cytokine genes between TR -ve and TR + ve teeth (e.g. IL1B, IL6, TNF, IFNG), P2X7R or VDR, www.nature.com/scientificreports/  Table 5. RNA-seq result validation by qPCR. Total n = 6 from TR -ve = 3, from TR + ve = 3, *; 0.05 < p < 0.1, **; p < 0.05 by student t-test. www.nature.com/scientificreports/ Fig. S4). It is possible that inflammation only appears at a particular stage of TR development so that the overall oral environment possesses a similar extent of inflammation within this study population. It has been reported that the prevalence of TR increased with age in random mixed breed cats although TR occurred regardless of age in pure breed cats 66,67 . In our study population, the age range was broad (from 4 to 16 years), however the mean age (7.2 years) was similar to those in the previous studies (6.2 and 6.5 years) 66,68 . To avoid bias due to age, young cats were excluded (< 3 years) in our study and there was no statistically difference of the mean age between TR − ve (7.6 ± 3.2) and TR + ve cats (6.7 ± 1.9) (p = 0.68). Another possible explanation for the lack of any changes in the expression levels of inflammatory cytokines could be due to a bias of TR type in this study. In this study, type 2 was the dominant type of TR (Table 1). It has been proposed that this type might be involved in non-inflammatory or idiopathic TR although other studies have not confirmed this observation 11,69,70 . The finding that there is no change in RANKL mRNA levels was consistent with previous feline TR studies 19,20 . It has been reported that RANKL is highly expressed in the periodontal ligament of actively resorbing primary human teeth 71,72 , whereas the level of RANKL mRNA is very low or not expressed in human adult teeth 73 . In cats, RANKL mRNA levels in adult teeth and bones are similar 19 . Since odontoclast activity is tightly controlled by the RANKL/OPG ratio, the level of OPG and RANKL mRNA expression should be considered together. In healthy feline permanent teeth, the level of OPG mRNA is significantly higher than femoral or alveolar bone, providing evidence that the role of OPG may be a key factor, rather www.nature.com/scientificreports/ than RANKL, in modulating odontoclast function in feline TR 19 . This is supported by RNA-seq and qPCR data in the comparison of TR -/ + teeth from the same cat; the level of OPG mRNA level was low in the TR + ve teeth thus an increased RANKL/OPG ratio was observed in TR + ve teeth. However, the observation was not verified when repeated in a larger number of samples, hence the RANKL/OPG involvement in feline TR needs to be further investigated. Based on the large number of DE genes, there were several interesting metabolic pathways identified. Specifically, DE genes were found which are reported to be involved in osteoclast differentiation (SPI1, C-FMS, RANK, RANKL, TREM2, DAP12, Nox1 and Nox2) 39,41,[48][49][50][51] or highly expressed in functional and mature osteoclast (OCSTAMP, ACP5, MMP9, CALCR and CTSK) 42,43,45,46 . SPI1, the gene encoding the haematopoietic transcription factor PU.1, is essential for myeloid and B-lymphoid cell development 39 . This was reported to be the earliest molecule which binds to sites in both the c-fms promoter or intronic regulatory element (FIRE enhancer) for commitment to osteoclast lineage in mice 40,74 . The increased expression of SPI1 and C-FMS found in TR + ve teeth www.nature.com/scientificreports/ may suggest that TR + ve teeth possess more odontoclast precursors in their environment which would be likely to form mature odontoclasts and induce tooth resorption. ACP5 and CTSK are recognised osteoclast markers in bone and teeth 43,46 . Up-regulation of these genes also suggest increased odontoclast number and activity in TR teeth. MMP9 is a matrix metallopeptidase enzyme, gelatinase B and type IV collagenase. Its primary role is to degrade and remodel the extracellular matrix in many tissues including bone and teeth 75 . Increased MMP9 activity has been implicated in many diseases including osteoclastomas, Paget's disease, bone fracture repair and dental pulp inflammation 76,77 . It has also been reported that it is highly expressed in early stages of differentiating osteoclasts, mature osteoclasts as well as in odontoclasts in many species including human, mouse, rabbit and cow [78][79][80] . Recently accumulated evidence suggests that MMP9 have an important role in the progression of caries, www.nature.com/scientificreports/ specifically in the pulpal and periapical inflammation tissue destruction process. Significant elevation of MMP9 expression was observed in apical periodontitis. MMP9 is synthesized and secreted by various stimulators such as physical agents or cell cytokines (IL1B and TNF-α) 77 . Interestingly, Zhang et al. (2020) suggested that MMP9 may have anti-inflammatory properties 81 . They showed that MMP9 knockout mice experimentally induced with apical periodontitis suffered a more severe inflammation than wild type mice. It was also reported that MMP9 induced the reduction of IL1B, TNF-α, RANK, RANKL, TLR2, and TLR4 and increased OPG expression in LPSstimulated osteoblasts. These findings suggest that the role of MMP9 in bone destruction is complex and diverse and emphasises the role of MMP9 as an inflammatory regulator by osteoblasts 81 . In our data, MMP9 expression was higher in TR + ve teeth but there was no difference in expression of inflammatory cytokines (IL1B, TNF-α) or RANK. Feline TR starts normally at the external root and neck area rather than within the pulp, and often results in ankylosis between tooth and alveolar bone rather than alveolar bone destruction. Feline TR is a unique phenomenon and might indicate slightly different pathogenesis from typical inflammation. Our study primarily focused on the role of the odontoclast, but it might be necessary to investigate further roles of other cell types comprising of the tooth (odontoblast, periodontal ligament cells). However, expression of MMP9 and its role in bone and tooth resorption has not previously been described in cats. One study reported high expression of MMP9 in feline cancer 82 . Another pathway identified, the calcium signalling pathway, controls multiple cellular functions in many cells via release of calcium ions from internal stores and its entry from the extracellular fluid 53 . In particular, calcium signals in osteoclasts are responsible for diverse cellular functions including differentiation, bone resorption and gene transcription 52 . Genes in calcium signalling might induce Ca 2+ influx thus contributing to osteoclast differentiation (P2X2R, P2X4R, P2X6R, PLCB4) [51][52][53] . These genes are known to encode for nonselective cation channels permeable to Ca 2+ at the cellular level and thereby contribute to osteoclast differentiation. Originally, purinergic receptors were thought to play a role in ATP production which would in turn induce the formation and activation of osteoclasts 51,83 . Purinergic receptors e.g. P2X7 or P2Y have been implicated as potent local inhibitors of matrix mineralisation 83 . Expression of P2X2 receptor has been reported in both osteoblasts and osteoclasts but its precise role in bone is not clear 51 . Up-regulation of P2X4 receptor has been attributed to inflammation via macrophage invasion 84 . Although expression of the P2X6 receptor has not been reported in osteoclasts, its expression was reported to modulate the differentiation and migration of human mesenchymal stem cells 85 . The roles of these genes (and their encoded proteins) in the dental microenvironment is largely unknown, so further studies are required to investigate their possible involvement in tooth resorption.

RNA-seq qPCR Expression in TR+ teeth Expression in TR+
In addition to osteoclast biology, gene ontology enrichment analysis showed abundant genes involved in muscle physiology. Genes involved in structural myofibril components including ACTA1, ACTN2, ATP2A1, CSRP3, MYH11, MYOZ2, PDLIM3, SYNC, TCAP, TNNI1, TNNT3 and TPM2 were all up-regulated, however, the exact role of these genes in bone and tooth structure/function is unclear. These genes have also been identified to play a role in the development and differentiation of osteocytes, mechanosensing cells that coordinate the remodelling process mediated by osteoblasts and osteoclasts in developing bone and tooth 86 . This manifestation is unexpected as adult teeth normally possess a limited ability to remodel. However, as TR + ve teeth showed a predominantly disrupted periodontal ligament, new bone formation (ankylosis) and external TR lesions in the root rather than internally, the noted changes in gene expression may reflect repair and some limited remodelling in the roots of the TR + ve teeth, possibly resulting in the level of ankylosis noted in the histological sections. Of the above DE genes, we verified RNA-seq results with qPCR for the genes SPI1, OPG, MMP9, ACP5, CTSK, P2X2R, P2X4R, P2X6R and PCLB4 which all might be involved in osteoclast differentiation and activity, based on previous reports 43,46,51,74,75,87-89 . Finally, an in vitro feline osteoclast model was established to test the potential role of candidate genes. Since our data showed MMP9 was one of the most up-regulated genes in TR + ve teeth, it was chosen as a target for further investigation. Inhibition of MMP9 by a semi-selective synthetic inhibitor caused both a reduction of osteoclast formation and resorption activity. Synthetic MMP inhibitors have been developed for the treatment of several serious pathologies including periodontitis, although clinical trials gave disappointing results and only a few inhibitors (e.g. Periostat™) has been approved by the FDA for the treatment of periodontal disease 90 . The synthetic inhibitor possesses a hydroxamic acid analogue which binds to the enzymatic sites of MMPs with high-affinity. However, members of MMPs family share a basic structure which comprises an auto-inhibitory pro-domain rendering them enzymatic latency, the catalytic domain, and the C-terminal hemopexin-like domain, for the recognition of MMP substrates 91 . Although the inhibition of osteoclastic activity by the MMP9 inhibitor is concentration dependent, it is not clear if the effect of the MMP9 inhibitor is specific to MMP9 or not. Therefore, MMP9 inhibitor experiments were confirmed using the gene specific siRNA approach. When inhibition of MMP9 mRNA was successfully performed by feline MMP9 siRNA, osteoclast formation was inhibited, but not resorption activity. The reduction of resorption activity by synthetic inhibitors but not by the more specific siRNA approach may be a result of off-target inhibition of other members of the MMP superfamily. These results suggest that the role for MMP9 in TR pathogenesis requires further investigation before any firm conclusions can be reached. Strong MMP9 expression both in osteoclasts and odontoclasts has been reported in bovine, rabbit and mice skeletal tissue 78,92 , and was seen here in the histological tooth sections. Therefore MMP9 may be a potential therapeutic target in feline tooth resorption. The underlying mechanism of inhibition of osteoclast formation by MMP9 needs to be further investigated. Other candidate genes verified in this study could also be potential therapeutic targets in feline TR, and further studies, potentially looking at co-targeting of multiple targets should be performed.