Next Article in Journal
RNA-Seq Reveals the Roles of Long Non-Coding RNAs (lncRNAs) in Cashmere Fiber Production Performance of Cashmere Goats in China
Previous Article in Journal
Morphological and Molecular Bases of Male Infertility: A Closer Look at Sperm Flagellum
 
 
Search for Articles:
Title / Keyword
Author / Affiliation / Email
Journal
Article Type
 
 
Advanced Search
 
Section
Special Issue
Volume
Issue
Number
Page
 
Journals
Genes
Volume 14
Issue 2
10.3390/genes14020385
genes-logo
Submit to this Journal Review for this Journal Propose a Special Issue
► Article Menu

Article Menu

share Share announcement Help format_quote Cite question_answer Discuss in SciProfiles thumb_up
...
Endorse textsms
...
Comment
first_page
settings
Order Article Reprints
Font Type:
Arial Georgia Verdana
Font Size:
Aa Aa Aa
Line Spacing:
Column Width:
Background:
Article

Multiple Genetic Loci Associated with Pug Dog Thoracolumbar Myelopathy

by
Gustaf Brander
1,2,
Cecilia Rohdin
3,
Matteo Bianchi
1,
Kerstin Bergvall
3,
Göran Andersson
4,
Ingrid Ljungvall
3,
Karin Hultin Jäderlund
5,
Jens Häggström
3,
Åke Hedhammar
3,
Kerstin Lindblad-Toh
1,2,*,† and
Katarina Tengvall
1,*,†
1
Science for Life Laboratory, Department of Medical Biochemistry and Microbiology, Uppsala University, 751 23 Uppsala, Sweden
2
Broad Institute of MIT and Harvard, Cambridge, MA 02142, USA
3
Department of Clinical Sciences, Swedish University of Agricultural Sciences, 750 07 Uppsala, Sweden
4
Department of Animal Breeding and Genetics, Swedish University of Agricultural Sciences, 750 07 Uppsala, Sweden
5
Department of Companion Animal Clinical Sciences, Norwegian University of Life Sciences, N-1432 Ås, Norway
*
Authors to whom correspondence should be addressed.
†
These authors contributed equally to this work.
Genes 2023, 14(2), 385; https://doi.org/10.3390/genes14020385
Submission received: 21 December 2022 / Revised: 24 January 2023 / Accepted: 27 January 2023 / Published: 1 February 2023
(This article belongs to the Section Animal Genetics and Genomics)
Browse Figures
Review Reports Versions Notes

Abstract

:
Pug dogs with thoracolumbar myelopathy (PDM) present with a specific clinical phenotype that includes progressive pelvic limb ataxia and paresis, commonly accompanied by incontinence. Vertebral column malformations and lesions, excessive scar tissue of the meninges, and central nervous system inflammation have been described. PDM has a late onset and affects more male than female dogs. The breed-specific presentation of the disorder suggests that genetic risk factors are involved in the disease development. To perform a genome-wide search for PDM-associated loci, we applied a Bayesian model adapted for mapping complex traits (BayesR) and a cross-population extended haplotype homozygosity test (XP-EHH) in 51 affected and 38 control pugs. Nineteen associated loci (harboring 67 genes in total, including 34 potential candidate genes) and three candidate regions under selection (with four genes within or next to the signal) were identified. The multiple candidate genes identified have implicated functions in bone homeostasis, fibrotic scar tissue, inflammatory responses, or the formation, regulation, and differentiation of cartilage, suggesting the potential relevance of these processes to the pathogenesis of PDM.

1. Introduction

The pug has been bred to promote characteristics such as small stature, flat face, large head, bulging eyes, wrinkled forehead, and soft, thin ears. The selection for specific traits and a small population size has contributed to a strong predisposition to several disorders. The 2020 Breed Health and Conservation Plan for Pugs included spinal problems in their list of priority health issues [1], with a two-fold increase in the risk of spinal cord disorder in pugs compared to other breeds [2]. The orthopedic foundation of America has reported that pelvic limb abnormalities are the second most common health concern in the pug breed. In addition, the insurance data from the biggest animal insurance company in Sweden, Agria Pet Insurance, has shown a seven-fold increase in the mortality rate in pugs with ataxia and paresis, which are the clinical consequences of a spinal cord problem. Pug dogs with thoracolumbar myelopathy (PDM) commonly present with a chronic, progressive clinical course of ataxia and paraparesis, often accompanied by incontinence [3,4,5,6,7]. Although no exact prevalence of PDM has been established, gait abnormalities (including, but not exclusively, PDM) are a common problem in pugs [4].
PDM has a late age of onset and usually exhibits neurological deficits at the age of around seven [5,7,8]. In Sweden, the affected pugs are usually euthanized within a year of the onset of pelvic limb ataxia [8]. The histopathologic characteristics of PDM include meningeal fibrosis with associated spinal cord destruction, often accompanied by neighboring vertebral column pathology. In addition, central nervous system (CNS) inflammation has been reported in a considerable number of pugs with PDM [8]. Differences in meningeal and vertebral column pathologies have led to various proposed etiologies, however shared pathological findings suggest they belong to a spectrum of the same disorder. Due to the breed-specific appearance of PDM, there are likely genetic risk factors contributing to the condition, although no genes associated with the disorder have hitherto been reported.
In this study, we analyzed 89 pugs (51 cases and 38 controls) with the aim of identifying the genes associated with PDM. First, we applied a Bayesian model adapted to complex trait mapping, BayesR, a well-suited method when defining multiple genetic risk factors with smaller effects, which is common for complex traits. BayesR assesses the effect size variances of all genetic variants simultaneously, resulting in a lower risk of false negatives. Moreover, by treating the effect sizes as random effects, the corresponding estimates are unbiased [9]. Second, we performed a cross-population extended haplotype homozygosity analysis (XP-EHH), which aims to detect candidate regions under recent selection by comparing the integrated extended haplotype homozygosity (EHH) in two subpopulations [10,11]. This methodology was employed to explore the possibility that the increased risk for PDM derives from pleiotropic or hitchhiking effects resulting from artificial selection. The selection of the specific pug characteristics could have created large genetic regions with high linkage disequilibrium (LD) and, if it is associated with PDM risk, this could be detected as an extended haplotype homozygosity differing between cases and controls. In total, we identified 19 associated loci with BayesR and three candidate regions under selection. These loci and regions harbor candidate genes with implicated functions in bone homeostasis, cartilage, fibrotic scar tissue, and inflammatory responses, and are therefore potentially relevant to the pathogenesis of PDM.

2. Materials and Methods

2.1. Ethical Statement

The study was approved by the Local Ethical Committees in Uppsala and Stockholm, Sweden. Signed informed consent was obtained from all of the dog owners.

2.2. Study Population

All of the dogs included in this study were privately owned and sampled in collaboration with veterinary clinics in Sweden, between 2010–2020. In addition, samples were included from phenotyped dogs that had donated blood to the biobank at the Swedish University of Agricultural Sciences, Uppsala, Sweden.
The study was comprised of affected pugs with signs of a myelopathy localized to the thoracolumbar spinal cord and control pugs without these signs. The affected pugs presented with a history of more than one month of a bilateral pelvic limb gait abnormality, described by the owner as incoordination and weakness. In addition, the author, C. Rohdin, performed a neurological examination of potential cases to confirm pelvic limb ataxia and paraparesis, suggestive of a myelopathy localized to the thoracolumbar spinal cord. Magnetic resonance imaging was performed, as previously described in Rohdin et al. [8], to confirm the focal spinal cord pathology in 37 affected pugs. We have previously shown that pugs presenting with a history of more than one month of bilateral pelvic limb gait abnormality alongside a neurological examination confirming pelvic limb ataxia and paraparesis, all had spinal cord pathology involving the meninges and spinal cord [8]. A subset of the recruited dogs was also included in two previous studies investigating the presence of vertebral malformations and pathology, respectively, of pugs with PDM [4,8].
The control pugs were at least eight years old with no history or signs of a neurological disorder. Specifically, to be included as a control, no signs of incoordination, weakness or incontinence should be present; nor should they exhibit worn-down nails or skin on the dorsum of the paws in the pelvic limbs.
Before quality control, the affected pugs consisted of 27 females (40.3%) with an average age at onset of 81.5 months (SD 23.7, median 84, interquartile range 68–96.5), and the control pugs consisted of 28 females (60.9%) with an average age at phenotyping of 132 months (SD 18.9, median 129, interquartile range 116–145).
Blood samples were collected in EDTA tubes and stored in −70 °C until analysis. One subset of the sample set, which included 31 pugs (10 cases/21 controls), was genotyped using Illumina 230k CanineHD BeadChip (Illumina, San Diego, CA, USA), and the other subset, 86 pugs (60 cases/26 controls), was whole-genome sequenced at low coverage using Gencove’s low-pass sequencing followed by imputation against a reference panel of 676 dogs (Pipeline Dog low pass v2; Gencove Inc, New York, NY, USA). The chromosomal coordinates were based on the dog CanFam3.1 genome assembly [12].

2.3. Imputation

The Illumina genotyped pugs underwent genetic imputation using the 86 Gencove low-pass sequenced pugs as the reference panel. To remove low confidence genotypes, variants that were not present in the reference panel were excluded prior to imputation, as were variants with minor allele frequency (MAF) < 0.05 and call rate < 95%. One case dog from the Illumina dataset with a <95% call rate was excluded. The remaining 30 dogs were phased using Shapeit v4.1.3 [13] and the imputation was performed using Impute2 v2.3.2 [14]. To validate the imputation, a selection of 5% random variants (n = 4965) were removed from the Illumina genotype dataset before imputing a second time. For this set of variants, the concordance rate between the imputed genotypes and the Illumina genotypes was 96.7%.

2.4. Quality Control

Quality control (QC) was performed using PLINK v1.9. [15]. After imputation, all individuals and variants had a call rate of >95%. The Hardy-Weinberg Equilibrium threshold, set to 1e−10 for both the cases and controls, removed 475 variants, and 24,745 variants with MAF below 0.01 were removed. The multidimensional scaling revealed no clusters or outliers. To correct for any batch effect deriving from the different genotyping technologies used in the two sample sets, we contrasted these two groups and found 6229 variants with a statistical difference at p < 1e−03, which we subsequently discarded. After QC, 2,140,239 variants remained. All of the dogs showed concordance between the reported and genotype predicted sex, and no dogs had an outlying (3 SD) average heterozygosity rate.
Using KING v2.2.7 to test for relatedness [16], five duplicates were identified within or between the genotyping technique groups. The duplicates were confirmed from the official records of individual registration numbers. An additional 22 pairs of pugs with a relatedness score above 0.177 (corresponding to a kinship between first- and second-degree relatives) were identified. From the duplicate pairs and the pairs above the relatedness threshold of 0.177, the individual with the lower call rate was excluded, resulting in 89 pugs (51 cases and 38 controls) in the final dataset.

2.5. Genetic Analyses

We used R v4.0.2 [17] and the following R-packages: GWASTools v1.40.0 [18] and GENESIS v2.24.0 [19] for handling the genotype data, and SNPRelate v1.28.0 [20] for defining the PCs and kinship. All of the figures were generated with ggplot2 v.3.3.5 [21].
The BayesR software (v. 1, updated 1 April 2021) implements a Bayesian mixture model for the analysis of complex traits, assuming that the variant effects belong to four different variance distributions: zero effect, 0.01σ2 g, 0.001σ2 g, and 0.0001σ2 g. Markov Chain Monte Carlo sampling is applied to arrive at the posterior inference about the variant effects based on the four distributions. Random effects, rather than fixed effects, provide less biased estimates than the traditional GWAS models and, as the effect of all variants are assessed simultaneously, the analysis is not subjected to multiple testing [9]. Before running the BayesR model, we LD-pruned the imputed QCed dataset using Plink v1.9 (window size 25 kb, step size 5 kb, correlation threshold [r2] 0.999). The Illumina genotyped variants that had been deleted in the LD-pruning were then merged back in. The final dataset consisted of 266,687 variants. The model, adjusted for genotyping technique (Illumina vs. Gencove) and sex, was run using 300,000 iterations and 100,000 burn-in samples, with every tenth sample saved for posterior inference. This was repeated five times to assess the convergence of the model, and the absolute values of the average variant effects were reported as the final results. The top 50 effect size variants were defined as effect variants and the associated loci were defined by variants in LD at r2 > 0.8. The sum of the risk alleles at the top effect variant per associated locus (risk index) was compared between cases and controls using the two-tailed Welsh Two Sample T test, and the phenotypic variance explained by the associated loci and sex, was calculated using ANOVA (R package stats v. 4.1.2).
To investigate whether the artificial selection for desired pug characteristics has resulted in the accumulation of risk or protective factors for PDM, we performed a selection signature analysis, comparing cases and controls. First, we used fastPHASE v.1.4.8 [22] with random starts set to 10, to identify the haplotypes for cases and controls and then the R package rehh v3.2.1 [23,24] to identify EHH in cases and controls separately. EHH detects the transmission of an extended haplotype without recombination [10] and XP-EHH compares the integrated EHH between two populations at each variant position [11]. XP-EHH detects regions with an extended haplotype around the alleles at a position in one population, but with a more rapid decay of LD in the other. We compared cases to the controls; thus, the positive XP-EHH scores indicate the selection in cases, whereas the negative scores indicate the selection in controls. The candidate regions under selection were defined by the variants with −log10(p) XP-EHH > 5.

2.6. Data Availability

The genotype datasets (plink files: bed, bim, fam, pheno), after quality control and relatedness filtering, have been uploaded onto the SciLifeLab data repository (DOI 10.17044/scilifelab.21948521).

3. Results

The final sample set consisted of 89 pugs (51 PDM cases and 38 controls). Thirty-two of the affected pugs (62.7%), and fourteen of the control pugs (36.8%) were males, as expected given the higher prevalence of PDM in males. The mean onset of clinical signs for the affected pugs was 78 months (standard deviation [SD] = 23), and the mean age at phenotypic confirmation of the control pugs was 132 months (SD = 19). An LD-pruned subset of the imputed variants was used (n = 266,687) for the Bayesian analysis, while the complete set of imputed variants (n = 2,140,239) was used for the XP-EHH analysis.

3.1. Bayesian Analysis Defined Nineteen Genetic Loci

From the 50 (0.02%) top variant effects in the BayesR model, we identified 19 associated loci, defined by variants in LD r2 > 0.8 with the top variant, harboring, in total, 67 genes (Table 1 and Figure 1). The associated loci with the highest variant effects and with more than one effect variant per region were chr4:11 Mb (intragenic in BicC Family RNA Binding Protein 1 [BICC1]), chr12:65 Mb (intragenic in SEC63 Homolog, Protein Translocation Regulator [SEC63]), chr17:53 Mb (top effect variant intragenic in Calsequestrin 2 [CASQ2]), and chr33:31 Mb (intragenic in Large 60S Subunit Nuclear Export GTPase 1 [LSG1]).
When comparing the risk index, i.e., the sum of the risk alleles from the 19 associated loci, between cases and controls, we found a statistically significant difference (p = 3.7e−12; Figure 2), with the risk index explaining 45.6% of the PDM variance (p = 4.3e−16;) and sex explaining 15.3% (p = 1.1e−07).

3.2. Selection Signature Analysis Defined Three Regions

The XP-EHH analysis identified three candidate regions under selection. The most significant region (p = 2.3e−7) was located on chromosome 15, with the top two selection variants mapping to positions chr15:28,234,106 and chr15:28,234,130. This signal included 1635 selection variants (i.e., variants with −log10(p) XP-EHH > 5) encompassing ~1.2 Mb (chr15:27,882,677–29,065,882) and overlapping the gene MGAT4 family member C (MGAT4C). The XP-EHH value was negative (−5.18) and the integrated haplotype homozygosity (iHH) for allele A was higher in controls (iHH = 7.3 Mb) than in cases (iHH = 2.0 Mb), thus indicating selection in the controls. The second candidate region under selection was detected on chromosome 17. The top selection variant mapped to position chr17:57,069,701, and the signal included 19 selection variants, encompassing ~7 Kb (chr17:57,064,222–57,071,120). No genes have been annotated within this region, however, the closest genes were SEC22 Homolog B, Vesicle Trafficking Protein (SEC22B) and NOTCH receptor 2 (NOTCH2), located 23 Kb and 44 Kb from the top selection variant, respectively. The third candidate region under selection was located on chromosome 30, with the top variant mapping to position chr30:8,951,408. This region harbored 79 selection variants, encompassing ~20 Kb (chr30:8,933,970–8,954,088), which overlapped the gene EH domain containing 4 (EHD4). Both the regions on chromosome 17 and chromosome 30 indicated signals of selection in the cases (Figure 3).

4. Discussion

This is the first study to identify genetic loci associated with PDM, a disorder so far only reported in pugs. The BayesR analysis identified 19 associated loci harboring 67 genes, out of which more than half (n = 34) were putative candidate genes for PDM (Table 1). Of these 34 genes, 14 had implicated functions in more than one of the below defined categories, with possible connections to PDM pathogenesis. The XP-EHH analysis identified three candidate regions under selection, defining four candidate genes in total, within or next to the selection signals, with potential functions related to PDM. The candidate genes, potentially associated with the PDM phenotype, are implicated in bone homeostasis, cartilage, fibrosis, and inflammation.
The stability and structure of the spine is maintained by the vertebrae (bone), the vertebral endplates, the facet joints, and the intervertebral discs (cartilage). Bone is a multifunctional, dynamic, mineralized connective tissue undergoing considerable change over time. It exhibits different types of cells, including immature osteoblasts, that differentiate into long-lived osteocytes, and bone-marrow derived osteoclasts [80,81]. Osteocytes are the most common cellular components of bone and are essential for bone mass regulation. The potential loss of the structural vertebral column integrity through spinal instability has been implicated as the cause for the clinical signs of PDM. Vertebral anomalies are suggested to interfere with spinal biomechanics and to result in instability and repeated spinal cord injury [3,6,82,83,84]. Pugs have been identified to possess the progressive loss of vertebral column integrity and an increase in spinal kyphosis (the abnormal curvature of the spine) [85,86]. From the BayesR analyses, we found 24 candidate genes (50.7% of the total number of genes) with potential functions in bone homeostasis. Nine of these have previously been reported to be involved in osteoblast metabolism. For instance, Tripartite Motif Containing 38 (TRIM38) and TGFB Induced Factor Homeobox 1 (TGIF1) play critical roles in bone remodeling and are involved in osteoblast and osteoclast differentiation [41,76]. Furthermore, BICC1, Nescient Helix-Loop-Helix 2 (NHLH2), and Protein Kinase AMP-Activated Non-Catalytic Subunit β 1 (PRKAB1) have been shown to affect bone mineral density [31,57,63]. In addition to being a tumor suppressor gene and being involved in inflammation, Programmed Cell Death 4 (PDCD4) is evolutionarily highly conserved and has a role in regulating osteogenic differentiation and bone defect repair [68].
Skeletal development, maintenance, and remodeling is regulated by osteoclasts [87]. Of the genes identified in our study, ten (14.9%) were specifically associated with the functions in osteoclasts, including TRIM38 and TGIF1, discussed above. Ras Homolog Family Member U (RhoU) is an encoding member of the Rho family, a family of small GTPases proteins active in the organization of the actin cytoskeleton [88] and involved in migration, epithelial cell morphogenesis, and osteoclastogenesis [39,89,90,91,92]. RhoU has been identified in the differentiation of macrophages into osteoclasts (osteoclastogenesis), and in decreasing bone resorption by its role in the adhesion structures of osteoclasts [89]. In addition, the disruption of Mitochondrial transcription factor A (TFAM) has been implicated in increased bone resorption through its presence in osteoclasts [32]. Of the candidate genes from the XP-EHH analysis, NOTCH2 is of particular interest, as it is implicated in osteoclastogenesis [93]. NOTCH2 expression in osteoblasts has been observed to regulate the cancellous bone volume and microarchitecture [94].
Cartilage is a tissue comprised of extracellular matrix components, including compact collagen, which is a fibrous protein that adds to tissue strength and structure. It is present in the vertebral column in the intervertebral discs and vertebral endplates, and in the articular surfaces of bone [95]. Soft cartilage is a desirable characteristic in pugs; for instance, kennel clubs (e.g., the American Kennel Club and United Kennel Club; [96,97]) define the ear characteristics of the standard pug breed as thin, small, and soft, similar to black velvet. However, the selection for this trait may have affected the supportive structure in more parts of the body, including the cartilage in the spine, thereby increasing the risk for PDM. The candidate region under selection on chromosome 30, with suggested selection pressure acting in the cases, includes the gene EHD4, which has been implicated in the cartilage functions [98]. Fourteen (20.9%) candidate genes from the associated loci have been involved in cartilage-related processes. Six of them are directly associated with cartilage, and ten through their association with osteoarthritis, which can be described as the degradation of cartilage in the joints at the end of bones. TGIF1, in addition to its role in bone remodeling described above, has been directly implicated in controlling cartilage as it encodes a transcription regulator described to inhibit differentiation into cartilage [43]. TGIF1 has also been described as playing a role in directing the differentiation of mesodermal cells (during TGFβ signaling) toward fibrogenesis instead of following chondrogenic differentiation [99]. From the XP-EHH analysis, NOTCH2 has been observed to mediate chondrogenesis differentiation in cartilage progenitor/stem cells [100] and to play a role in chondrocyte maturation [93], and MGAT4C has presented with upregulated mRNA expression in osteoarthritis [36,101]. Interestingly, the XP-EHH values of the selection variants overlapping MGAT4C were negative, i.e., indicating selection in the controls, which may suggest that breeders are aware of the disorder and are actively breeding against it.
Seven (10.4%) of our identified candidate genes have previously been associated with the fibroblast functions, which can be linked to the excessive scar tissue of the meninges, observed in PDM. One of these genes is TFAM, which shows that a reduction in protein expression, and the associated mitochondrial damage, translates into an enhanced sensitivity of fibroblasts to profibrotic stimuli [38]. Fibrosis, or scarring, is a form of tissue repair in which connective tissue replaces the original parenchyma. When this process is disturbed, as in chronic diseases, it often leads to increased fibrosis [102]. The meninges are the fibrous coverings of the CNS and consist of a variety of cell types, including CNS fibroblasts and specialized immune cells. Fibroblasts from the meninges are major drivers of fibrotic scar formation following CNS injury [103].
Nine (13.4%) candidate genes from the associated loci have previously been associated with inflammatory responses, the activation of a highly coordinated immunological response specific for the initial stimulus [104]. This represents a clear link to the CNS inflammation observed in PDM. PDCD4, already described for its role in osteogenic differentiation above, is widely expressed in the immune cells of humans and other animals. It plays an important role in the macrophage function and exhibits both inflammatory and anti-inflammatory functions [105]. PDCD4-deficient mice, immunized with myelin oligodendrocyte glycoprotein to induce experimental autoimmune encephalomyelitis, have shown resistance to autoimmune encephalomyelitis and developed a reduced degree of spinal cord inflammation [106]. Interestingly, PDCD4 promotes chronic inflammation and could therefore be relevant for the gradual and protracted onset of PDM. Furthermore, SEC22B, identified in the XP-EHH analysis, has been implicated in phagosome formation, crucial for the defense against pathogens [107]. In addition to being implicated in fibrosis and bone resorption, TFAM has also been shown to induce pro-inflammatory and cytotoxic responses of microglia [40].
We identified 19 associated loci in BayesR, defined by the 50 top effect variants; a cutoff previously used in dogs [108]. While the ratio between false negatives and false positives has been shown to be favorable for BayesR, with lower effect sizes, the risk of false positives increases [9]. A stricter definition, using the top 10 effect variants, would result in seven associated loci and 26 candidate genes, out of which 13 have potential implications in PDM and covering all of the above-mentioned biological functions. Being aware of the higher risk of false positives among the 19 associated loci, we still defined the definition of the top 50 effect variants to avoid losing the true positive associations. Given the large proportion of PDM variance explained (45.6%) by the risk index, the associated loci are indicated to be of high relevance to the development of PDM. These results will hopefully help us understand the etiology of the disease in dogs better. While there is no known human disorder corresponding to PDM, the rare disorder adhesive arachnoiditis in human [109,110] has shown similarities with PDM in terms of meningeal fibrosis and subsequent spinal cord destruction. However, no comparative studies have yet been performed. Signals on chromosome 17 were identified both in the BayesR analysis and the XP-EHH analysis. They were, however, approximately 4 Mb apart and neither the LD region nor the candidate region of selection overlapped with the other, pointing to them being different signals. Future studies should include larger sample sets to confirm our current findings and validate the lowest effect loci in particular. Analyses of the gene and protein expression differences in the relevant tissues from pugs would be of high interest to explore the functions of the candidate genes in PDM.

5. Conclusions

Taken together, this study suggests that the genes implicated in bone homeostasis, fibrotic scar tissue, inflammatory responses, and formation, regulation, and differentiation of cartilage are likely to be implicated in the underlying etiology of PDM. Even though the details on how the candidate genes may be involved in these biological processes remain to be fully elucidated, we believe that by utilizing the advanced genome-wide mapping methods, we have increased the knowledge about this devastating disorder.

Author Contributions

Conceptualization, K.L.-T., G.A. and Å.H.; methodology, K.L.-T., Å.H., M.B., G.A. and K.T.; software, G.B. and K.T.; validation, G.B., K.T. and M.B.; formal analysis, G.B.; investigation, G.B. and C.R.; resources, K.L.-T., C.R. and K.B.; data curation, G.B. and K.T.; writing—original draft preparation, G.B., K.T. and C.R.; writing—review and editing, G.B., C.R., M.B., K.B., G.A., I.L., K.H.J., J.H., Å.H., K.L.-T. and K.T.; visualization, G.B.; supervision, Å.H., I.L., K.H.J., J.H., K.T. and K.L.-T.; project administration, K.L.-T.; funding acquisition, K.L.-T. and C.R. All authors have read and agreed to the published version of the manuscript.

Funding

This research was funded by a Distinguished Professorship from the Swedish Research Council, with additional funding from the Thure F and Karin Forsberg’s Research Foundation.

Institutional Review Board Statement

The study was approved by the Local Ethical Committees in Uppsala and Stockholm, Sweden.

Informed Consent Statement

Informed consent was obtained from all subjects involved in the study.

Data Availability Statement

The genotype datasets (plink files: bed, bim, fam, pheno), after quality control and relatedness filtering, have been uploaded onto the SciLifeLab data repository (DOI 10.17044/scilifelab.21948521).

Acknowledgments

Computations and data handling were enabled by resources in projects, SNIC 2017/7-384, SNIC 2017/7-385, and SNIC 2021/5-579, provided by the Swedish National Infrastructure for Computing (SNIC) at UPPMAX, partially funded by the Swedish Research Council through grant agreement no. 2018-05973.

Conflicts of Interest

The authors declare no conflict of interest.

References

  1. Brachycephalic Working Group. The 2020 Pug Breed Health and Conservation Plans (BHCPs). 2020. Available online: http://www.ukbwg.org.uk/wp-content/uploads/2020/05/Pug_BHCP_2020.pdf (accessed on 9 August 2022).
  2. O’Neill, D.G.; Sahota, J.; Brodbelt, D.C.; Church, D.B.; Packer, R.M.A.; Pegram, C. Health of Pug dogs in the UK: Disorder predispositions and protections. Canine Med. Genet. 2022, 9, 4. [Google Scholar] [CrossRef]
  3. Fisher, S.C.; Shores, A.; Simpson, S.T. Constrictive myelopathy secondary to hypoplasia or aplasia of the thoracolumbar caudal articular processes in pugs: 11 cases (1993–2009). J. Am. Vet. Med. Assoc. 2013, 242, 223–229. [Google Scholar] [CrossRef]
  4. Rohdin, C.; Häggström, J.; Ljungvall, I.; Nyman Lee, H.; De Decker, S.; Bertram, S.; Lindblad-Toh, K.; Hultin Jäderlund, K. Presence of thoracic and lumbar vertebral malformations in pugs with and without chronic neurological deficits. Vet. J. 2018, 241, 24–30. [Google Scholar] [CrossRef] [PubMed]
  5. Alisauskaite, N.; Cizinauskas, S.; Jeserevics, J.; Rakauskas, M.; Cherubini, G.B.; Anttila, M.; Steffen, F. Short- and long-term outcome and magnetic resonance imaging findings after surgical treatment of thoracolumbar spinal arachnoid diverticula in 25 Pugs. J. Vet. Intern. Med. 2019, 33, 1376–1383. [Google Scholar] [CrossRef]
  6. Driver, C.J.; Rose, J.; Tauro, A.; Fernandes, R.; Rusbridge, C. Magnetic resonance image findings in pug dogs with thoracolumbar myelopathy and concurrent caudal articular process dysplasia. BMC Vet. Res. 2019, 15, 182. [Google Scholar] [CrossRef] [PubMed]
  7. Lourinho, F.; Holdsworth, A.; McConnell, J.F.; Gonçalves, R.; Gutierrez-Quintana, R.; Morales, C.; Lowrie, M.; Trevail, R.; Carrera, I. Clinical features and MRI characteristics of presumptive constrictive myelopathy in 27 pugs. Vet. Radiol. Ultrasound 2020, 61, 545–554. [Google Scholar] [CrossRef] [PubMed]
  8. Rohdin, C.; Ljungvall, I.; Häggström, J.; Leijon, A.; Lindblad-Toh, K.; Matiasek, K.; Rosati, M.; Wohlsein, P.; Hultin Jäderlund, K. Thoracolumbar meningeal fibrosis in pugs. J. Vet. Intern. Med. 2020, 34, 797–807. [Google Scholar] [CrossRef]
  9. Moser, G.; Lee, S.H.; Hayes, B.J.; Goddard, M.E.; Wray, N.R.; Visscher, P.M. Simultaneous Discovery, Estimation and Prediction Analysis of Complex Traits Using a Bayesian Mixture Model. PLoS Genet. 2015, 11, e1004969. [Google Scholar] [CrossRef]
  10. Sabeti, P.C.; Reich, D.E.; Higgins, J.M.; Levine, H.Z.P.; Richter, D.J.; Schaffner, S.F.; Gabriel, S.B.; Platko, J.V.; Patterson, N.J.; McDonald, G.J.; et al. Detecting recent positive selection in the human genome from haplotype structure. Nature 2002, 419, 832–837. [Google Scholar] [CrossRef]
  11. Sabeti, P.C.; Varilly, P.; Fry, B.; Lohmueller, J.; Hostetter, E.; Cotsapas, C.; Xie, X.; Byrne, E.H.; McCarroll, S.A.; Gaudet, R. Genome-wide detection and characterization of positive selection in human populations. Nature 2007, 449, 913–918. [Google Scholar] [CrossRef] [Green Version]
  12. Hoeppner, M.P.; Lundquist, A.; Pirun, M.; Meadows, J.R.S.; Zamani, N.; Johnson, J.; Sundström, G.; Cook, A.; FitzGerald, M.G.; Swofford, R. An improved canine genome and a comprehensive catalogue of coding genes and non-coding transcripts. PLoS ONE 2014, 9, e91172. [Google Scholar] [CrossRef]
  13. Delaneau, O.; Marchini, J.; Zagury, J.F. A linear complexity phasing method for thousands of genomes. Nat. Methods 2012, 9, 179–181. [Google Scholar] [CrossRef]
  14. Howie, B.N.; Donnelly, P.; Marchini, J. A flexible and accurate genotype imputation method for the next generation of genome-wide association studies. PLoS Genet. 2009, 5, e1000529. [Google Scholar] [CrossRef]
  15. Chang, C.C.; Chow, C.C.; Tellier, L.C.A.M.; Vattikuti, S.; Purcell, S.M.; Lee, J.J. Second-generation PLINK: Rising to the challenge of larger and richer datasets. Gigascience 2015, 4, 7. [Google Scholar] [CrossRef]
  16. Manichaikul, A.; Mychaleckyj, J.C.; Rich, S.S.; Daly, K.; Sale, M.; Chen, W.M. Robust relationship inference in genome-wide association studies. Bioinformatics 2010, 26, 2867–2873. [Google Scholar] [CrossRef]
  17. Ihaka, R.; Gentleman, R. R: A Language for Data Analysis and Graphics. J. Comput. Graph. Stat. 1996, 5, 299–314. [Google Scholar] [CrossRef]
  18. Gogarten, S.M.; Bhangale, T.; Conomos, M.P.; Laurie, C.A.; McHugh, C.P.; Painter, I.; Zheng, X.; Crosslin, D.R.; Levine, D.; Lumley, T.; et al. GWASTools: An R/Bioconductor package for quality control and analysis of genome-wide association studies. Bioinformatics 2012, 28, 3329–3331. [Google Scholar] [CrossRef] [PubMed]
  19. Gogarten, S.M.; Sofer, T.; Chen, H.; Yu, C.; Brody, J.A.; Thornton, T.A.; Rice, K.M.; Conomos, M.P. Genetic association testing using the GENESIS R/Bioconductor package. Bioinformatics 2019, 35, 5346–5348. [Google Scholar] [CrossRef] [PubMed]
  20. Zheng, X.; Levine, D.; Shen, J.; Gogarten, S.M.; Laurie, C.; Weir, B.S. A high-performance computing toolset for relatedness and principal component analysis of SNP data. Bioinformatics 2012, 28, 3326–3328. [Google Scholar] [CrossRef] [PubMed]
  21. Wickham, H. Ggplot2: Elegant Graphics for Data Analysis; Springer-Verlag: New York, NY, USA, 2016; ISBN 978-3-319-24277-4. [Google Scholar]
  22. Scheet, P.; Stephens, M. A fast and flexible statistical model for large-scale population genotype data: Applications to inferring missing genotypes and haplotypic phase. Am. J. Hum. Genet. Genet. 2006, 78, 629–644. [Google Scholar] [CrossRef] [Green Version]
  23. Gautier, M.; Vitalis, R. rehh: An R package to detect footprints of selection in genome-wide SNP data from haplotype structure. Bioinformatics 2012, 28, 1176–1177. [Google Scholar] [CrossRef] [PubMed]
  24. Gautier, M.; Klassmann, A.; Vitalis, R. rehh 2.0: A reimplementation of the R package rehh to detect positive selection from haplotype structure. Mol. Ecol. Resour. 2017, 17, 78–90. [Google Scholar] [CrossRef]
  25. Coudert, A.E.; Del Fattore, A.; Baulard, C.; Olaso, R.; Schiltz, C.; Collet, C.; Teti, A.; De Vernejoul, M.C. Differentially expressed genes in autosomal dominant osteopetrosis type II osteoclasts reveal known and novel pathways for osteoclast biology. Lab. Investig. 2014, 94, 275–285. [Google Scholar] [CrossRef]
  26. Mosig, R.A.; Dowling, O.; DiFeo, A.; Ramirez, M.C.M.; Parker, I.C.; Abe, E.; Diouri, J.; Al Aqeel, A.; Wylie, J.D.; Oblander, S.A. Loss of MMP-2 disrupts skeletal and craniofacial development and results in decreased bone mineralization, joint erosion and defects in osteoblast and osteoclast growth. Hum. Mol. Genet. 2007, 16, 1113–1123. [Google Scholar] [CrossRef] [PubMed]
  27. Mosig, R.A.; Martignetti, J.A. Loss of MMP-2 in murine osteoblasts upregulates osteopontin and bone sialoprotein expression in a circuit regulating bone homeostasis. DMM Dis. Model. Mech. 2013, 6, 397–403. [Google Scholar] [CrossRef]
  28. Paic, F.; Igwe, J.C.; Nori, R.; Kronenberg, M.S.; Franceschetti, T.; Harrington, P.; Kuo, L.; Shin, D.G.; Rowe, D.W.; Harris, S.E.; et al. Identification of differentially expressed genes between osteoblasts and osteocytes. Bone 2009, 45, 682–692. [Google Scholar] [CrossRef]
  29. Jia, T.; Lao, J. Bioinformatics Analyses of Regulatory Network of Biomarkers in Chondrocytes from Patients with Osteoarthritis. Braz. Arch. Biol. Technol. 2022, 65. [Google Scholar] [CrossRef]
  30. Boer, C.; de Kruijf, M.; Broer, L.; Hofman, A.; Uitterlinden, A.; van Meurs, J. Novel susceptability loci for osteoarthritis of the hand: Variants in coding EN GENE regulatory regions. Osteoarthr. Cartil. 2015, 23, A196. [Google Scholar] [CrossRef]
  31. Mesner, L.D.; Ray, B.; Hsu, Y.-H.; Manichaikul, A.; Lum, E.; Bryda, E.C.; Rich, S.S.; Rosen, C.J.; Criqui, M.H.; Allison, M.; et al. Bicc1 is a genetic determinant of osteoblastogenesis and bone mineral density. J. Clin. Investig. 2014, 124, 2736–2749. [Google Scholar] [CrossRef] [PubMed]
  32. Miyazaki, T.; Iwasawa, M.; Nakashima, T.; Mori, S.; Shigemoto, K.; Nakamura, H.; Katagiri, H.; Takayanagi, H.; Tanaka, S. Intracellular and extracellular ATP coordinately regulate the inverse correlation between osteoclast survival and bone resorption. J. Biol. Chem. 2012, 287, 37808–37823. [Google Scholar] [CrossRef] [Green Version]
  33. Itzstein, C.; Coxon, F.P.; Rogers, M.J. The regulation of osteoclast function and bone resorption by small GTPases. Small GTPases 2011, 2, 117–130. [Google Scholar] [CrossRef]
  34. Ory, S.; Brazier, H.; Pawlak, G.; Blangy, A. Rho GTPases in osteoclasts: Orchestrators of podosome arrangement. Eur. J. Cell Biol. 2008, 87, 469–477. [Google Scholar] [CrossRef] [PubMed]
  35. Wang, Y.; Zhao, X.; Lotz, M.; Terkeltaub, R.; Liu-Bryan, R. Mitochondrial biogenesis is impaired in osteoarthritis chondrocytes but reversible via peroxisome proliferator-activated receptor γ coactivator 1α. Arthritis Rheumatol. 2015, 67, 2141–2153. [Google Scholar] [CrossRef] [PubMed]
  36. Zhang, X.; Bu, Y.; Zhu, B.; Zhao, Q.; Lv, Z.; Li, B.; Liu, J. Global transcriptome analysis to identify critical genes involved in the pathology of osteoarthritis. Bone Jt. Res. 2018, 7, 298–307. [Google Scholar] [CrossRef] [PubMed]
  37. Zhang, Y.; Zhu, T.; He, F.; Chen, A.C.; Yang, H.; Zhu, X. Identification of Key Genes and Pathways Associated with Differences of Subchondral Bone in Osteoarthritis. 2020; preprint. [Google Scholar] [CrossRef]
  38. Zhou, X.; Trinh-Minh, T.; Tran-Manh, C.; Gießl, A.; Bergmann, C.; Györfi, A.-H.; Schett, G.; Distler, J.H.W. Impaired Mitochondrial Transcription Factor A Expression Promotes Mitochondrial Damage to Drive Fibroblast Activation and Fibrosis in Systemic Sclerosis. Arthritis Rheumatol. 2022, 74, 871–881. [Google Scholar] [CrossRef]
  39. Ory, S.; Brazier, H.; Blangy, A. Identification of a bipartite focal adhesion localization signal in RhoU/Wrch-1, a Rho family GTPase that regulates cell adhesion and migration. Biol. Cell 2007, 99, 701–716. [Google Scholar] [CrossRef]
  40. Schindler, S.M.; Frank, M.G.; Annis, J.L.; Maier, S.F.; Klegeris, A. Pattern recognition receptors mediate pro-inflammatory effects of extracellular mitochondrial transcription factor A (TFAM). Mol. Cell. Neurosci. 2018, 89, 71–79. [Google Scholar] [CrossRef]
  41. Saito, H.; Gasser, A.; Bolamperti, S.; Maeda, M.; Matthies, L.; Jähn, K.; Long, C.L.; Schlüter, H.; Kwiatkowski, M.; Saini, V. TG-interacting factor 1 (Tgif1)-deficiency attenuates bone remodeling and blunts the anabolic response to parathyroid hormone. Nat. Commun. 2019, 10, 1354. [Google Scholar] [CrossRef]
  42. Eichholz, K.F.; Woods, I.; Riffault, M.; Johnson, G.P.; Corrigan, M.; Lowry, M.C.; Shen, N.; Labour, M.N.; Wynne, K.; O’Driscoll, L.; et al. Human bone marrow stem/stromal cell osteogenesis is regulated via mechanically activated osteocyte-derived extracellular vesicles. Stem Cells Transl. Med. 2020, 9, 1431–1447. [Google Scholar] [CrossRef]
  43. Chen, L.; Jiang, C.; Tiwari, S.R.; Shrestha, A.; Xu, P.; Liang, W.; Sun, Y.; He, S.; Cheng, B. TGIF1 gene silencing in tendon-derived stem cells improves the tendon-to-bone insertion site regeneration. Cell. Physiol. Biochem. 2015, 37, 2101–2114. [Google Scholar] [CrossRef]
  44. Kimura, M.Y.; Koyama-Nasu, R.; Yagi, R.; Nakayama, T. A new therapeutic target: The CD69-Myl9 system in immune responses. Semin. Immunopathol. 2019, 41, 349–358. [Google Scholar] [CrossRef] [PubMed]
  45. Vacher, J.; Bruccoleri, M.; Pata, M. Ostm1 from mouse to human: Insights into osteoclast maturation. International Journal of Mol. Sci. 2020, 21, 5600. [Google Scholar] [CrossRef] [PubMed]
  46. Chan, C.K.F.; Gulati, G.S.; Sinha, R.; Tompkins, J.V.; Lopez, M.; Carter, A.C.; Ransom, R.C.; Reinisch, A.; Wearda, T.; Murphy, M.; et al. Identification of the Human Skeletal Stem Cell. Cell 2018, 175, 43–56. [Google Scholar] [CrossRef] [PubMed]
  47. Jasenc, L.; Stražar, K.; Mihelič, A.; Mihalič, R.; Trebše, R.; Haring, G.; Jeras, M.; Zupan, J. In Vitro Characterization of the Human Skeletal Stem Cell-like Properties of Primary Bone-Derived Mesenchymal Stem/Stromal Cells in Patients with Late and Early Hip Osteoarthritis. Life 2022, 12, 899. [Google Scholar] [CrossRef]
  48. Tajuddin, S.M.; Schick, U.M.; Eicher, J.D.; Chami, N.; Giri, A.; Brody, J.A.; Hill, W.D.; Kacprowski, T.; Li, J.; Lyytikäinen, L.P.; et al. Large-Scale Exome-wide Association Analysis Identifies Loci for White Blood Cell Traits and Pleiotropy with Immune-Mediated Diseases. Am. J. Hum. Genet. 2016, 99, 22–39. [Google Scholar] [CrossRef]
  49. Kowanetz, M.; Valcourt, U.; Bergström, R.; Heldin, C.-H.; Moustakas, A. Id2 and Id3 Define the Potency of Cell Proliferation and Differentiation Responses to Transforming Growth Factor β and Bone Morphogenetic Protein. Mol. Cell. Biol. 2004, 24, 4241–4254. [Google Scholar] [CrossRef]
  50. Huynh, N.P.T.; Gloss, C.C.; Lorentz, J.; Tang, R.; Brunger, J.M.; McAlinden, A.; Zhang, B.; Guilak, F. Long non-coding rna graslnd enhances chondrogenesis via suppression of interferon type II signaling pathway. Elife 2020, 9, e49558. [Google Scholar] [CrossRef]
  51. He, P.; Zhang, Z.; Liao, W.; Xu, D.; Fu, M.; Kang, Y. Screening of gene signatures for rheumatoid arthritis and osteoarthritis based on bioinformatics analysis. Mol. Med. Rep. 2016, 14, 1587–1593. [Google Scholar] [CrossRef]
  52. Severa, M.; Coccia, E.M.; Fitzgerald, K.A. Toll-like receptor-dependent and -independent Viperin gene expression and counter-regulation by PRDI-binding factor-1/BLIMP1. J. Biol. Chem. 2006, 281, 26188–26195. [Google Scholar] [CrossRef] [Green Version]
  53. Cho, H.; Proll, S.C.; Szretter, K.J.; Katze, M.G.; Gale, M.; Diamond, M.S. Differential innate immune response programs in neuronal subtypes determine susceptibility to infection in the brain by positive-stranded RNA viruses. Nat. Med. 2013, 19, 458–464. [Google Scholar] [CrossRef]
  54. You, D.; Yang, C.; Huang, J.; Gong, H.; Yan, M.; Ni, J. Long non-coding RNA MEG3 inhibits chondrogenic differentiation of synovium-derived mesenchymal stem cells by epigenetically inhibiting TRIB2 via methyltransferase EZH2. Cell. Signal. 2019, 63, 109379. [Google Scholar] [CrossRef] [PubMed]
  55. Azuma, K.; Shiba, S.; Hasegawa, T.; Ikeda, K.; Urano, T.; Horie-Inoue, K.; Ouchi, Y.; Amizuka, N.; Inoue, S. Osteoblast-Specific γ-Glutamyl Carboxylase-Deficient Mice Display Enhanced Bone Formation with Aberrant Mineralization. J. Bone Miner. Res. 2015, 30, 1245–1254. [Google Scholar] [CrossRef] [PubMed]
  56. Nie, Z.; Chen, S.; Deng, S.; Long, L.; Peng, P.; Gao, M.; Cheng, S.; Cao, J.; Peng, H. Gene expression profiling of osteoblasts subjected to dexamethasone-induced apoptosis with/without GSK3β-shRNA. Biochem. Biophys. Res. Commun. 2018, 506, 41–47. [Google Scholar] [CrossRef]
  57. Swan, A.L.; Schütt, C.; Rozman, J.; del Mar Muñiz Moreno, M.; Brandmaier, S.; Simon, M.; Leuchtenberger, S.; Griffiths, M.; Brommage, R.; Keskivali-Bond, P.; et al. Mouse mutant phenotyping at scale reveals novel genes controlling bone mineral density. PLoS Genet. 2020, 16, e1009190. [Google Scholar] [CrossRef] [PubMed]
  58. Macsai, C.E.; Georgiou, K.R.; Foster, B.K.; Zannettino, A.C.W.; Xian, C.J. Microarray expression analysis of genes and pathways involved in growth plate cartilage injury responses and bony repair. Bone 2012, 50, 1081–1091. [Google Scholar] [CrossRef] [PubMed]
  59. Meng, J.; Ma, X.; Ma, D.; Xu, C. Microarray analysis of differential gene expression in temporomandibular joint condylar cartilage after experimentally induced osteoarthritis. Osteoarthr. Cartil. 2005, 13, 1115–1125. [Google Scholar] [CrossRef]
  60. Calender, A.; Rollat Farnier, P.A.; Buisson, A.; Pinson, S.; Bentaher, A.; Lebecque, S.; Corvol, H.; Abou Taam, R.; Houdouin, V.; Bardel, C.; et al. Whole exome sequencing in three families segregating a pediatric case of sarcoidosis. BMC Med. Genom. 2018, 11, 23. [Google Scholar] [CrossRef]
  61. Hershey, C.L.; Fisher, D.E. Mitf and Tfe3: Members of a b-HLH-ZIP transcription factor family essential for osteoclast development and function. Bone 2004, 34, 689–696. [Google Scholar] [CrossRef]
  62. Yao, J.; Wang, Y.; Cao, C.; Song, R.; Bi, D.; Zhang, H.; Li, Y.; Qin, G.; Hou, N.; Zhang, N.; et al. CRISPR/Cas9-mediated correction of MITF homozygous point mutation in a Waardenburg syndrome 2A pig model. Mol. Nucleic Acids 2021, 24, 986–999. [Google Scholar] [CrossRef]
  63. Zhang, J.G.; Tan, L.J.; Xu, C.; He, H.; Tian, Q.; Zhou, Y.; Qiu, C.; Chen, X.D.; Deng, H.W. Integrative analysis of transcriptomic and epigenomic data to reveal regulation patterns for BMD variation. PLoS ONE 2015, 10, e0138524. [Google Scholar] [CrossRef]
  64. Thaler, R.; Zwerina, J.; Rumpler, M.; Spitzer, S.; Gamsjaeger, S.; Paschalis, E.P.; Klaushofer, K.; Varga, F. Homocysteine induces serum amyloid A3 in osteoblasts via unlocking RGD-motifs in collagen. FASEB J. 2013, 27, 446–463. [Google Scholar] [CrossRef] [PubMed]
  65. Lind, T.; Melo, F.R.; Gustafson, A.-M.; Sundqvist, A.; Zhao, X.O.; Moustakas, A.; Melhus, H.; Pejler, G. Mast cell chymase has a negative impact on human osteoblasts. Matrix Biol. 2022, 112, 1–19. [Google Scholar] [CrossRef]
  66. Zhengquan, D.; Lei, W. Gene expression profiles in osteoarthritis: A bioinformatic analysis. Chin. J. Tissue Eng. Res. 2019, 23, 335. [Google Scholar] [CrossRef]
  67. Chen, T.C.; Lai, C.H.; Chang, J.L.; Chang, S.W. Mitomycin C retardation of corneal fibroblast migration via sustained dephosphorylation of paxillin at Tyrosine 118. Investig. Opthalmol. Vis. Sci. 2012, 53, 1539–1547. [Google Scholar] [CrossRef]
  68. Jiang, Y.; Li, S.; Zhou, Q.; Liu, S.; Liu, X.; Xiao, J.; Jiang, W.; Xu, Y.; Kong, D.; Wang, F.; et al. PDCD4 Negatively Regulated Osteogenic Differentiation and Bone Defect Repair of Mesenchymal Stem Cells through GSK-3β/β-Catenin Pathway. Stem Cells Dev. 2021, 30, 806–815. [Google Scholar] [CrossRef]
  69. Liu, L.; Feng, Y.; Hu, S.; Li, H.; Li, Y.; Ke, J.; Long, X. PDCD4 suppresses autophagy and promotes apoptosis via Akt in chondrocytes of temporomandibular joint osteoarthritis. Oral Dis. 2021, 27, 547–558. [Google Scholar] [CrossRef]
  70. Jo, S.H.; Kim, D.E.; Clocchiatti, A.; Dotto, G.P. PDCD4 is a CSL associated protein with a transcription repressive function in cancer associated fibroblast activation. Oncotarget 2016, 7, 58717–58727. [Google Scholar] [CrossRef]
  71. Selfors, L.M.; Schutzman, J.L.; Borland, C.Z.; Stern, M.J. soc-2 encodes a leucine-rich repeat protein implicated in fibroblast growth factor receptor signaling. Proc. Natl. Acad. Sci. USA 1998, 95, 6903–6908. [Google Scholar] [CrossRef]
  72. Zhu, Y.; Liu, L.; Hu, L.; Dong, W.; Zhang, M.; Liu, Y.; Li, P. Effect of Celastrus orbiculatus in inhibiting Helicobacter pylori induced inflammatory response by regulating epithelial mesenchymal transition and targeting miR-21/PDCD4 signaling pathway in gastric epithelial cells. BMC Complement. Altern. Med. 2019, 19, 91. [Google Scholar] [CrossRef] [Green Version]
  73. Alam, I.; Sun, O.; Koller, D.L.; Liu, L.; Liu, Y.; Edenberg, H.J.; Li, J.; Foroud, T.; Turner, C.H. Differentially expressed genes strongly correlated with femur strength in rats. Genomics 2009, 94, 257–262. [Google Scholar] [CrossRef]
  74. Gorrell, L.; Omari, S.; Makareeva, E.; Leikin, S. Noncanonical ER–Golgi trafficking and autophagy of endogenous procollagen in osteoblasts. Cell. Mol. Life Sci. 2021, 78, 8283–8300. [Google Scholar] [CrossRef]
  75. Liu, H.; Wang, W.; Shen, W.; Wang, L.; Zuo, Y. ARHGAP24 ameliorates inflammatory response through inactivating Rac1/Akt/NF-κB pathway in acute pneumonia model of rat. Ann. Transl. Med. 2020, 8, 1289. [Google Scholar] [CrossRef]
  76. Kim, K.; Kim, J.H.; Kim, I.; Seong, S.; Kim, N. TRIM38 regulates NF-κB activation through TAB2 degradation in osteoclast and osteoblast differentiation. Bone 2018, 113, 17–28. [Google Scholar] [CrossRef]
  77. Ribet, A.B.P.; Ng, P.Y.; Pavlos, N.J. Membrane Transport Proteins in Osteoclasts: The Ins and Outs. Front. Cell Dev. Biol. 2021, 9, 644986. [Google Scholar] [CrossRef]
  78. Wang, Q.; Notay, K.; Downey, G.P.; McCulloch, C.A. The Leucine-Rich Repeat Region of CARMIL1 Regulates IL-1-Mediated ERK Activation, MMP Expression, and Collagen Degradation. Cell Rep. 2020, 31, 107781. [Google Scholar] [CrossRef]
  79. Hu, M.-M.; Xie, X.-Q.; Yang, Q.; Liao, C.-Y.; Ye, W.; Lin, H.; Shu, H.-B. TRIM38 Negatively Regulates TLR3/4-Mediated Innate Immune and Inflammatory Responses by Two Sequential and Distinct Mechanisms. J. Immunol. 2015, 195, 4415–4425. [Google Scholar] [CrossRef]
  80. Downey, P.A.; I Siegel, M. Bone biology and the clinical implications for osteoporosis. Phys. Ther. 2006, 86, 77–91. [Google Scholar] [CrossRef]
  81. Florencio-Silva, R.; Sasso, G.R.D.S.; Sasso-Cerri, E.; Simões, M.J.; Cerri, P.S. Biology of Bone Tissue: Structure, Function, and Factors That Influence Bone Cells. BioMed Res. Int. 2015, 2015, 421746. [Google Scholar] [CrossRef]
  82. Parker, A.J.; Smith, C.W. Meningeal cyst in a dog. J. Am. Anim. Hosp. Assoc. 1974, 10, 595–597. [Google Scholar]
  83. Flegel, T.; Müller, M.K.; Truar, K.; Löffler, C.; Oechtering, G. Thoracolumbar spinal arachnoid diverticula in 5 pug dogs. Can. Vet. J. 2013, 54, 969. [Google Scholar]
  84. Chen, A.V.; Bagley, R.S.; West, C.L.; Gavin, P.R.; Tucker, R.L. Fecal incontinence and spinal cord abnormalities in seven dogs. J. Am. Vet. Med. Assoc. 2005, 227, 1945–1951. [Google Scholar] [CrossRef]
  85. Wyatt, S.; Gonçalves, R.; Gutierrez-Quintana, R.; de Decker, S. Outcomes of nonsurgical treatment for congenital thoracic vertebral body malformations in dogs: 13 cases (2009–2016). J. Am. Vet. Med. Assoc. 2018, 253, 768–773. [Google Scholar] [CrossRef]
  86. De Rycke, L.M.; Crijns, C.; Chiers, K.; van Bree, H.J.J.; Gielen, I. Late-onset wedge-shaped thoracic vertebrae in a six-month-old pug. Vet. Rec. Case Rep. 2016, 4, e000317. [Google Scholar] [CrossRef]
  87. Bar-Shavit, Z. The osteoclast: A multinucleated, hematopoietic-origin, bone-resorbing osteoimmune cell. J. Cell. Biochem. 2007, 102, 1130–1139. [Google Scholar] [CrossRef]
  88. Kjøller, L.; Hall, A. Signaling to Rho GTPases. Exp. Cell Res. 1999, 253, 166–179. [Google Scholar] [CrossRef]
  89. Brazier, H.; Pawlak, G.; Vives, V.; Blangy, A. The Rho GTPase Wrch1 regulates osteoclast precursor adhesion and migration. Int. J. Biochem. Cell Biol. 2009, 41, 1391–1401. [Google Scholar] [CrossRef]
  90. Chuang, Y.Y.; Valster, A.; Coniglio, S.J.; Backer, J.M.; Symons, M. The atypical Rho family GTPase Wrch-1 regulates focal adhesion formation and cell migration. J. Cell Sci. 2007, 120, 1927–1934. [Google Scholar] [CrossRef]
  91. Brady, D.C.; Alan, J.K.; Madigan, J.P.; Fanning, A.S.; Cox, A.D. The Transforming Rho Family GTPase Wrch-1 Disrupts Epithelial Cell Tight Junctions and Epithelial Morphogenesis. Mol. Cell. Biol. 2009, 29, 1035–1049. [Google Scholar] [CrossRef]
  92. Alan, J.K.; Berzat, A.C.; Dewar, B.J.; Graves, L.M.; Cox, A.D. Regulation of the Rho Family Small GTPase Wrch-1/RhoU by C-Terminal Tyrosine Phosphorylation Requires Src. Mol. Cell. Biol. 2010, 30, 4324–4338. [Google Scholar] [CrossRef]
  93. Zieba, J.T.; Chen, Y.-T.; Lee, B.H.; Bae, Y. Notch signaling in skeletal development, homeostasis and pathogenesis. Biomolecules 2020, 10, 332. [Google Scholar] [CrossRef]
  94. Zanotti, S.; Canalis, E. Notch1 and Notch2 expression in osteoblast precursors regulates femoral microarchitecture. Bone 2014, 62, 22–28. [Google Scholar] [CrossRef] [Green Version]
  95. Hunziker, E.B.; Lippuner, K.; Keel, M.J.B.; Shintani, N. An educational review of cartilage repair: Precepts & practice—Myths & misconceptions—Progress & prospects. Osteoarthr. Cartil. 2015, 23, 334–350. [Google Scholar] [CrossRef]
  96. American Kennel Club. Official Standard of the Pug. 2008. Available online: https://images.akc.org/pdf/breeds/standards/Pug.pdf (accessed on 9 August 2022).
  97. United Kennel Club. Official UKC Breed Standard: Pug. 2012. Available online: https://www.ukcdogs.com/docs/breeds/pug.pdf (accessed on 9 August 2022).
  98. Kuo, H.J.; Tran, N.T.; Clary, S.A.; Morris, N.P.; Glanville, R.W. Characterization of EHD4, an EH Domain-containing Protein Expressed in the Extracellular Matrix. J. Biol. Chem. 2001, 276, 43103–43110. [Google Scholar] [CrossRef]
  99. Lorda-Diez, C.I.; Montero, J.A.; Martinez-Cue, C.; Garcia-Porrero, J.A.; Hurle, J.M. Transforming growth factors β coordinate cartilage and tendon differentiation in the developing limb mesenchyme. J. Biol. Chem. 2009, 284, 29988–29996. [Google Scholar] [CrossRef]
  100. Moqbel, S.A.A.; Zeng, R.; Ma, D.; Xu, L.; Lin, C.; He, Y.; Ma, C.; Xu, K.; Ran, J.; Jiang, L.; et al. The effect of mitochondrial fusion on chondrogenic differentiation of cartilage progenitor/stem cells via Notch2 signal pathway. Stem Cell Res. Ther. 2022, 13, 127. [Google Scholar] [CrossRef]
  101. Zhang, Y.; Zhu, T.; He, F.; Chen, A.C.; Yang, H.; Zhu, X. Identification of Key Genes and Pathways in Osteoarthritis via Bioinformatic Tools: An Updated Analysis. Cartilage 2021, 13, 1457S–1464S. [Google Scholar] [CrossRef]
  102. Karsdal, M.A.; Kraus, V.B.; Shevell, D.; Bay-Jensen, A.C.; Schattenberg, J.; Rambabu Surabattula, R.; Schuppan, D. Profiling and targeting connective tissue remodeling in autoimmunity—A novel paradigm for diagnosing and treating chronic diseases. Autoimmun. Rev. 2021, 20, 102706. [Google Scholar] [CrossRef]
  103. Derk, J.; Jones, H.E.; Como, C.; Pawlikowski, B.; Siegenthaler, J.A. Living on the Edge of the CNS: Meninges Cell Diversity in Health and Disease. Front. Cell. Neurosci. 2021, 15, 703944. [Google Scholar] [CrossRef]
  104. Natoli, G.; Ghisletti, S.; Barozzi, I. The genomic landscapes of inflammation. Genes Dev. 2011, 25, 101–106. [Google Scholar] [CrossRef]
  105. Jiang, Y.; Jia, Y.; Zhang, L. Role of programmed cell death 4 in diseases: A double-edged sword. Cell. Mol. Immunol. 2017, 14, 884–886. [Google Scholar] [CrossRef]
  106. Hilliard, A.; Hilliard, B.; Zheng, S.-J.; Sun, H.; Miwa, T.; Song, W.; Göke, R.; Chen, Y.H. Translational Regulation of Autoimmune Inflammation and Lymphoma Genesis by Programmed Cell Death 4. J. Immunol. 2006, 177, 8095–8102. [Google Scholar] [CrossRef] [Green Version]
  107. Sun, W.; Tian, B.X.; Wang, S.H.; Liu, P.J.; Wang, Y.C. The function of SEC22B and its role in human diseases. Cytoskeleton 2020, 77, 303–312. [Google Scholar] [CrossRef]
  108. Baker, L.A.; Momen, M.; McNally, R.; Berres, M.E.; Binversie, E.E.; Sample, S.J.; Muir, P. Biologically Enhanced Genome-Wide Association Study Provides Further Evidence for Candidate Loci and Discovers Novel Loci That Influence Risk of Anterior Cruciate Ligament Rupture in a Dog Model. Front. Genet. 2021, 12, 593515. [Google Scholar] [CrossRef]
  109. Jurga, S.; Szymańska-Adamcewicz, O.; Wierzchołowski, W.; Pilchowska-Ujma, E.; Urbaniak, Ł. Spinal adhesive arachnoiditis: Three case reports and review of literature. Acta Neurol. Belg. 2021, 121, 47–53. [Google Scholar] [CrossRef]
  110. Anderson, T.L.; Morris, J.M.; Wald, J.T.; Kotsenas, A.L. Imaging appearance of advanced chronic adhesive arachnoiditis: A retrospective review. Am. J. Roentgenol. 2017, 209, 648–655. [Google Scholar] [CrossRef]
Genes 14 00385 g001 550
Figure 1. BayesR analysis resulted in four top associated loci and 19 loci in total. The top 50 effect variants are above the dotted red line defining multiple associated loci and variants with the highest effect sizes marked out with labels (A). Variants in LD (r2 > 0.8; red) with the top effect variant and annotated protein coding genes are presented for each of the four top associated loci (B).
Figure 1. BayesR analysis resulted in four top associated loci and 19 loci in total. The top 50 effect variants are above the dotted red line defining multiple associated loci and variants with the highest effect sizes marked out with labels (A). Variants in LD (r2 > 0.8; red) with the top effect variant and annotated protein coding genes are presented for each of the four top associated loci (B).
Genes 14 00385 g001
Genes 14 00385 g002 550
Figure 2. The risk index was higher in PDM cases compared to controls. The histogram shows the distribution of summed risk alleles at the associated loci (risk index) in cases and controls (A), and the boxplots indicate the median, first and third quartiles, and range (B).
Figure 2. The risk index was higher in PDM cases compared to controls. The histogram shows the distribution of summed risk alleles at the associated loci (risk index) in cases and controls (A), and the boxplots indicate the median, first and third quartiles, and range (B).
Genes 14 00385 g002
Genes 14 00385 g003 550
Figure 3. XP-EHH analysis defined three candidate regions under selection. XP-EHH p-values and XP-EHH values are presented across the whole genome (A). The candidate regions under selection are presented as EHH plots in cases and controls for chromosomes 15, 17, and 30 with haplotype lengths indicated by integrated extended haplotype homozygosity (iHH) for the respective alleles in cases and controls (B). XP-EHH p-values corresponding to the iHH regions (C) and zoom-ins around the top selection variants (shaded in grey) are presented with the annotated protein coded genes in (D).
Figure 3. XP-EHH analysis defined three candidate regions under selection. XP-EHH p-values and XP-EHH values are presented across the whole genome (A). The candidate regions under selection are presented as EHH plots in cases and controls for chromosomes 15, 17, and 30 with haplotype lengths indicated by integrated extended haplotype homozygosity (iHH) for the respective alleles in cases and controls (B). XP-EHH p-values corresponding to the iHH regions (C) and zoom-ins around the top selection variants (shaded in grey) are presented with the annotated protein coded genes in (D).
Genes 14 00385 g003
Table
Table 1. PDM associated loci identified with BayesR.
Table 1. PDM associated loci identified with BayesR.
Top VariantA1/A2Effect Size 1EAF
A/U		VAR EXP (%)Variants in Top 50 (N)LD Region 2 (Size in Kb)Genes Closest to Top VariantGenes in LD RegionBone HomeostasisCartilageOsteoarthritisFibrotic Scar TissueInflammatory
Response
chr2:60108854G/A5.36E−050.42/0.594.771chr2:60095511–60653394 (558)CES5A (intergenic 55 Kb); GNAO1 (intergenic 128 Kb)CAPNS2, CES1, CES5A, IRX6, LPCAT2, MMP2, SLC6A2CES1 [25], MMP2 [26] [27], IRX6 [28] GNAO1 [29]
chr3:71553918A/G−3.64E−050.58/0.3413.571chr3:71521315–71768914 (248)APBB2APBB2 APBB2 [30]
chr4:10805369C/T4.88E−050.14/0.3710.8317chr4:9844757–11012379 (1168)BICC1BICC1, CCSAP, RAB4A, RHOU, TFAM, UBE2D1BICC1 [31], TFAM [32], RHOU [33] [34] TFAM [35] [36], BICC1 [37]TFAM [38], RHOU [39]TFAM [40]
chr6:65661296T/C−4.09E−050.16/0.044.671chr6:65386176–65714252 (328)ADGRL2ADGRL2
chr7:70300935A/G−3.90E−050.22/0.091.661chr7:69701081–71152660 (1452)DLGAP1DLGAP1, MYL12B, MYOM1, TGIF1TGIF1 [41], MYL12B [42]TGIF1 [43]MYOM1 [37] MYL12B [44]
chr12:65015175C/T5.76E−050.21/0.413.802chr12:64968699–65110915 (142)SEC63SEC63, OSTM1OSTM1 [45]
chr12:66289122G/A4.18E−050.16/0.380.271chr12:66266562–66343853 (77)CD164 (intergenic 6 Kb); CCDC162P (intergenic 5 Kb)CCDC162P/C12H6orf183, CD164, PPIL6CD164 [46]CD164 [46]CD164 [47] CCDC162P [48]
chr17:5061849A/G3.90E−050.29/0.463.671chr17:4062154–5436291 (1374)RNF144A (intergenic 424 Kb); ID2 (intergenic 895 Kb)CMPK2, RNF144A, RSAD2ID2 [49]RNF144A (as GRASLND, i.e., RNF144-AS1; human transcript) [50] RSAD2 [51,52 and 53]
chr17:9192869G/C−3.62E−050.33/0.127.781chr17:8871659–9520837 (649)TRIB2 (intergenic 33 Kb)TRIB2 TRIB2 [54]
chr17:39587560A/G−4.43E−050.19/0.071.641chr17:39502199–39589843 (88)SH2D6 (intergenic 15K Kb), MAT2A(intergenic 55 Kb)GGCX, MAT2AGGCX [55]
chr17:53278611T/A8.22E−050.28/0.531.4411chr17:53239890–53619454 (380)CASQ2, SLC22A15CASQ2, MAB21L3, NHLH2, SLC22A15CASQ2 [56], NHLH2 [57]CASQ2 [58]CASQ2 [59]
chr17:54018787C/T−4.27E−050.54/0.300.061chr17:54018787–54024163 (5)IGSF3IGSF3 IGSF3 [60]
chr20:22144318A/G−5.05E−050.42/0.270.311chr20:21579384–22373451 (794)MITF (intergenic 271 Kb), FRMD48 (intergenic 258 Kb)FRMD4B, MITFMITF [61] MITF [62]
chr26:15777160C/T3.86E−050.01/0.145.002chr26:15350328–16204269 (854)CITBICDL1, CCDC60, CIT, GCN1, PRKAB1, PXN, RAB35, RPLP0, TMEM233PKAB1 [63], PXN [64], TMEM233 [65] PXN [66]PXN [67]
chr28:21986613C/T3.96E−050.44/0.570.031chr28:21986613–22282079 (295)RBM20BBIP1, PDCD4, RBM20, SHOC2PDCD4 [68] PCDC4 [69]PDCD4 [70], SHOC2 [71]PDCD4 [72]
chr28:30266337G/T−6.34E−050.33/0.240.471chr28:30041108–30695769 (655)ARF1 (intergenic 26 Kb), PLPP4 (intergenic 244 Kb)SEC23IP, PLPP4SEC23IP [73], ARF1 [74]
chr32:9378372C/T−4.51E−050.28/0.220.411chr32:9192329–9505692 (313)ARHGAP24ARHGAP24 ARHGAP24 [75]
chr33:31097826C/T−5.78E−050.10/0.040.635chr33:31088620–31160245 (72)LSG1LSG1
chr35:23618123T/C−3.65E−050.44/0.410.021chr35:23525268–23956076 (431)SCGNCARMIL1, HIST1H2AA, HIST1H2BA, SCGN, SLC17A1, SLC17A2, SLC17A3, SLC17A4, TRIM38TRIM38 [76], SLC17A2 (i.e., SLC34A1) and SLC17A1 [77]CARMIL1 [78] CARMIL1 [78]CARMIL1 [78], TRIM38 [79]
1 Effect size off A1; 2 LD r2 > 0.8; Abbreviations: EAF effect allele frequency; VAR EXP variance explained.
Disclaimer/Publisher’s Note: The statements, opinions and data contained in all publications are solely those of the individual author(s) and contributor(s) and not of MDPI and/or the editor(s). MDPI and/or the editor(s) disclaim responsibility for any injury to people or property resulting from any ideas, methods, instructions or products referred to in the content.

Share and Cite

MDPI and ACS Style

Brander, G.; Rohdin, C.; Bianchi, M.; Bergvall, K.; Andersson, G.; Ljungvall, I.; Hultin Jäderlund, K.; Häggström, J.; Hedhammar, Å.; Lindblad-Toh, K.; et al. Multiple Genetic Loci Associated with Pug Dog Thoracolumbar Myelopathy. Genes 2023, 14, 385. https://doi.org/10.3390/genes14020385

AMA Style

Brander G, Rohdin C, Bianchi M, Bergvall K, Andersson G, Ljungvall I, Hultin Jäderlund K, Häggström J, Hedhammar Å, Lindblad-Toh K, et al. Multiple Genetic Loci Associated with Pug Dog Thoracolumbar Myelopathy. Genes. 2023; 14(2):385. https://doi.org/10.3390/genes14020385

Chicago/Turabian Style

Brander, Gustaf, Cecilia Rohdin, Matteo Bianchi, Kerstin Bergvall, Göran Andersson, Ingrid Ljungvall, Karin Hultin Jäderlund, Jens Häggström, Åke Hedhammar, Kerstin Lindblad-Toh, and et al. 2023. "Multiple Genetic Loci Associated with Pug Dog Thoracolumbar Myelopathy" Genes 14, no. 2: 385. https://doi.org/10.3390/genes14020385

APA Style

Brander, G., Rohdin, C., Bianchi, M., Bergvall, K., Andersson, G., Ljungvall, I., Hultin Jäderlund, K., Häggström, J., Hedhammar, Å., Lindblad-Toh, K., & Tengvall, K. (2023). Multiple Genetic Loci Associated with Pug Dog Thoracolumbar Myelopathy. Genes, 14(2), 385. https://doi.org/10.3390/genes14020385

Note that from the first issue of 2016, this journal uses article numbers instead of page numbers. See further details here.

Article Metrics

Back to TopTop