Integrating transcriptome-wide association study and mRNA expression profile identified candidate genes related to hand osteoarthritis
Arthritis Research & Therapy volume 23, Article number: 81 (2021)
Osteoarthritis (OA) is a common skeletal system disease that has been partially attributed to genetic factors. The hand is frequently affected, which seriously affects the patient’s quality of life. However, the pathogenetic mechanism of hand osteoarthritis (hand OA) is still elusive.
A genome-wide association study (GWAS) summary of hand OA was obtained from the UK Biobank dataset, which contains data from a total of 452,264 White British individuals, including 37,782 OA patients. The transcriptome-wide association study (TWAS) of hand OA was performed using FUnctional Summary-based ImputatiON (FUSION) with the skeletal muscle and blood as gene expression references. The significant genes identified by TWAS were further subjected to gene set enrichment analysis (GSEA) with the Database for Annotation, Visualization and Integrated Discovery (DAVID) tool. Furthermore, we compared the genes and gene sets identified by our TWAS with that of a knee OA mRNA expression profile to detect the genes and gene sets shared by TWAS and mRNA expression profiles in OA. The mRNA expression profiles of 18 normal knee cartilages and 20 OA knee cartilages were acquired from the Gene Expression Omnibus database (accession number: GSE114007).
TWAS identified 177 genes with P < 0.05 for the skeletal muscle, including ANKRD44 (P = 0.0001), RIC3 (P = 0.0003), and AC005154.6 (P = 0.0004). TWAS identified 423 genes with P < 0.05 for the blood, including CRIM1 (P = 0.0002), ZNF880 (P = 0.0002), and NCKIPSD (P = 0.0003). After comparing the results of the TWAS to those of the mRNA expression profiling of OA, we identified 5 common genes, including DHRS3 (log2fold = − 1.85, P = 3.31 × 10− 9) and SKP2 (log2fold = 1.36, P = 1.62 × 10− 8). GSEA of TWAS identified 51 gene ontology (GO) terms for hand OA, for example, protein binding (P = 0.0003) and cytosol (P = 0.0020). We also detected 6 common GO terms shared by TWAS and mRNA expression profiling, including protein binding (PTWAS = 2.54 × 10− 4, PmRNA = 3.42 × 10− 8), extracellular exosome (PTWAS = 0.02, PmRNA = 1.18 × 10− 4), and cytoplasm (PTWAS = 0.0183, PmRNA = 0.0048).
In this study, we identified 5 candidate genes and 6 GO terms related to hand OA, which may help to uncover the pathogenesis of hand OA. It should be noted that the possible difference in the gene expression profiles between hand OA and knee OA may affect our study results, which should be interpreted with caution.
Osteoarthritis (OA) is a kind of musculoskeletal system disease that occurs in many elderly individuals. The hand is one of the locations where OA occurs most often, especially in women. The main symptoms are joint stiffness, pain, and even swelling. OA affects up to 10% of men and 13% of women worldwide . Among women aged 60 to 70 in the USA, the proportion of women suffering from hand OA is as high as 75% . The hand flexion and extension activity of hand OA patients is seriously affected, and hand function and quality of life are reduced as a result. With the aging of society, an increasing number of hand OA patients may be seen in the future.
With the development of bioinformatics, an increasing number of studies have focused on the genetic underpinnings of OA. At the genetic level, there seem to be several mechanisms that may differ between hand OA and knee or hip OA . Through genome-wide association studies (GWASs), researchers have found several loci significantly associated with OA, such as chr12 near the matrix Gla protein (MGP) gene and chr12 near the CCDC91 gene . Several significant functional polymorphisms were also found between hand OA patients and normal subjects, such as ASPN D15 and CILP rs2073711 TT . Unfortunately, the genetic mechanism underlying hand OA remains unclear.
GWASs have been successfully applied for gene mapping of complex human diseases and traits. However, this approach is limited in evaluating disease risk, because most GWAS-identified SNPs are located in non-coding regions of the genome. Variant SNPs may play a role in gene expression levels . Expression quantitative trait loci (eQTL) analysis aims to identify genes related to variation in gene expression . Therefore, integrating GWAS and eQTL analysis helps researchers to identify genes with causal associations to disease more effectively. The unique advantages of transcriptome-wide association studies (TWASs) have gradually been highlighted in recent years. TWAS combines the precomputed gene expression weights with summary statistical data from GWAS to recognize novel causal genes of elusive diseases . In recent years, TWASs have been utilized by an increasing number of researchers to identify genetic loci associated with the disease. For example, a TWAS of approximately 3000 subjects found 69 novel genes significantly related to BMI, lipid profile, and height that have a bearing on obesity . A TWAS investigating chronic low-grade inflammation recognized 448 genes related to inflammatory biologic age .
In this study, we used a large GWAS dataset for hand OA in conjunction with cartilage mRNA expression profiling. A TWAS was performed first to find genetic loci that may be associated with hand OA. Then, gene ontology (GO) and pathway enrichment analysis were carried out for the notable genes screened by TWAS. To identify common genes and biological pathways, the significant genes identified by the TWAS were compared with the OA profile analysis results, as well as GO and pathway analyses (Fig. 1). It should be noted that the difference in the gene expression profiles between hand OA and knee OA may affect our study results. More biological studies should be conducted to confirm our study results.
Material and method
GWAS summary datasets for hand OA
The GWAS summary data for hand OA was derived from the UK Biobank samples (UK Biobank fields: 20002) . Briefly, the UK Biobank cohort contains 452,264 White British individuals, including 37,782 osteoarthritis patients . Information about various phenotypes has been collected from every participant, and blood samples were taken when the subjects visited a UK Biobank assessment center. DNA extraction and genotyping were performed in the Affymetrix Research Services Laboratory. There were 62,394 genotyped variants that passed quality control through the Applied Biosystems UK Biobank Axiom Array. This dataset contains 9,113,133 imputed variants after filtering. The IMPUTE4 program was used to perform the imputation (http://jmarchini.org/software/). Detailed information on the subjects, genotyping, imputation, and quality control can be found in the published study .
Gene expression profiles of cartilage
The mRNA expression profiles of knee OA were acquired from the Gene Expression Omnibus database (accession number: GSE114007). Briefly, the normal human knee cartilage tissues were collected from 5 females and 13 males (age 18–61, mean 38), who had never suffered from joint disease or trauma and processed within 24–48 h postmortem. The knee cartilage specimens affected by OA were gathered from 12 females and 8 males undergoing knee replacement surgery (age 52–82, mean 66). Body mass indices between the normal controls (BMI = 32.4 ± 8.0) and OA patients (BMI = 30.7 ± 8.1) were not significantly (p = 0.506) different . The original image data were transformed into sequence data by the Illumina Genome Analyzer Pipeline Software (Casava v1.8.2). Genes designated as expressed in the study were those whose log counts per million (log2CPM) were greater than 3.0 in one or more samples. Limma-voom was carried out for differential expression analysis. Genes were designated as significantly differentially expressed (DE) when the following two conditions were met: adjusted P value of < 0.05 by the moderated t statistic and |log2FC| > 1 .
TWAS of hand OA
The TWAS of hand OA was carried out using the FUnctional Summary-based ImputatiON (FUSION) approach (http://gusevlab.org/projects/fusion/) . Given that the skeletal muscle and peripheral blood were used in previous biological studies of OA [13, 14], the FUSION precomputed functional weights of gene expression of the skeletal muscle and peripheral blood were used in our TWAS of hand OA. In brief, we took advantage of the prediction models implemented in FUSION to calculate the gene expression. Then, the calculated tissue-related expression weights were integrated with summary-level GWAS results to impute the association statistics between gene expression and the target disease. In this study, the gene expression weight panels for the peripheral blood and skeletal muscle were downloaded from the FUSION website (http://gusevlab.org/projects/fusion/).
GO and pathway enrichment analysis
The candidate genes identified by TWAS were further analyzed by the Database for Annotation, Visualization and Integrated Discovery (DAVID, https://david.ncifcrf.gov/) tool for GO and pathway enrichment analysis. Similarly, DAVID was also applied to the differentially expressed genes identified by the mRNA expression profiling of OA cartilage. To find the common GO and pathways detected both by TWAS and mRNA expression profiles, we carefully compared the DAVID analysis results of TWAS and cartilage mRNA expression profiles.
TWAS results of hand OA
The TWAS identified 177 genes with a P value < 0.05 for the skeletal muscle, including ANKRD44 (P = 0.0001), RIC3 (P = 0.0003), and AC005154.6 (P = 0.0004). TWAS identified 423 genes with a P value < 0.05 for the blood, including CRIM1 (P = 0.0002), ZNF880 (P = 0.0002), and NCKIPSD (P = 0.0003) (Supplementary Table S1). The top 10 most significant genes identified by TWAS are shown in Table 1.
GO and pathway analysis results of TWAS
DAVID analysis of TWAS identified 51 GO terms enriched in hand OA, such as protein binding (P value = 2.54 × 10− 4), cytosol (P value = 0.0012), and membrane (P value = 0.0014). We also detected 5 biological pathways for hand OA, including biosynthesis of unsaturated fatty acids (P value = 0.0034) and propanoate metabolism (P value = 0.0071) (Supplementary Table S2, S3).
Integrative analysis of TWAS and mRNA expression profiling of OA
A comparison of the genes identified by TWAS with the differentially expressed genes detected by mRNA expression profile yielded 5 common genes: DHRS3 (log2fold = − 1.85, P value = 3.31 × 10− 9), SKP2 (log2fold = 1.36, P value = 1.62 × 10− 8), IRS2 (log2fold = − 2.00, P value = 2.11 × 10− 9), TOB1 (log2fold = − 2.67, P value = 8.42 × 10− 9), PPP1R15A (log2fold = − 2.13, P value = 2.83 × 10− 8) (Table 2).
A comparison of the GO enrichment analysis results of TWAS and mRNA expression profiles yielded 6 common GO terms: protein binding (PTWAS = 2.54 × 10− 4, PmRNA = 3.42 × 10− 8), extracellular exosome (PTWAS = 0.0225, PmRNA = 1.18 × 10− 4), cytoplasm (PTWAS = 0.0183, PmRNA = 0.0048), oxidoreductase (PTWAS = 0.0312, PmRNA = 0.0259), cellular response to mechanical stimulus (PTWAS = 0.0324, PmRNA = 0.0053), and oxidation-reduction process (PTWAS = 0.0424, PmRNA = 0.0294) (Table 3).
Our study aimed to detect candidate genes closely associated with hand osteoarthritis and to explain the relationships between these genes and the disease. We conducted a TWAS and functional gene set enrichment analysis and identified multiple genes and enriched GO terms and pathways for hand OA. Further integrative analysis of TWAS and mRNA expression profiling of OA identified 5 common genes and 6 common GO terms shared by TWAS and the mRNA expression profile of OA. Our study results provide novel clues for understanding the genetic mechanism of hand OA, focusing on the roles of abnormal gene transcription in the development of OA.
TWAS identified several candidate genes for hand OA, such as IRS2, SKP2, CRIM1, and DHRS3 (Fig. 2). IRS2 (insulin receptor substrate 2) is a cytoplasmic adaptor molecule. It encodes some essential proteins related to insulin transduction and indeed plays a role in insulin resistance . Insulin-like growth factor 1 (IGF1) has important effect on bone growth and cellular differentiation. IGF1 is highly expressed in the chondrocytes and can promote the proliferation of chondrocytes. A previous study proved that IRS2 could influenced chondrocyte differentiation by targeting RUNX1 to control signaling downstream of IGF1 .
SKP2 was also detected in both the TWAS and mRNA expression profiles of OA. SKP2 gene encodes S-phase kinase-associated protein 2 (Skp2), which acts as an essential regulator of cell cycle. Previous study indicated that Skp2 also played a direct role in cellular senescence . The prevalence of OA is positively correlated with age. The relationship between cellular aging and OA has recently been well investigated [18,19,20]. In addition, Skp2 is an important molecule that mediates the ubiquitination and degradation of p27 protein. Skp2/p27 pathway is a critical mediator of ER stress-induced growth arrest. ER stress is implicated in the development of OA .
CRIM1 gene encodes human cysteine-rich motor neuron protein 1. Researchers have found that the extracellular domain of CRIM1 could bind to bone morphogenetic protein (BMP) . It has been found that chondrocyte hypertrophy could be inhibited by depressing BMP signaling pathway . In addition, previous studies suggested that BMP was involved in skeletal morphogenesis and chondrogenesis and acted as a potential chemo-protective mediators that regulated matrix metabolism .
Another identified candidate gene for OA is DHRS3. DHRS3 gene encodes dehydrogenase reductase 3 enzyme [6, 25]. This enzyme is thought to play a role in regulating all-trans-retinaldehyde in the human body . DHRS3 is detect in the biology process of vitamin A that plays an important role in bone and cartilage metabolism, and its deficiency is implicated in the pathological process of OA . In a study of rheumatoid arthritis (RA), DHRS3 was detected as one candidate genes for RA, as its expression was twice as high in the RA case group as in the control group . A previous study also found that DHRS3 mediated IL-1β-administered chondrocyte proliferation, apoptosis and ECM degradation by positively regulating GREM1 expression via competitively targeting miR-183-5p .
GO enrichment analysis was conducted to explore the functions of significant genes and how they are distributed in hand OA. In our study, oxidoreductase activity was identified by both TWAS and mRNA expression profiling. Previous research has shown that oxidoreductase activity is negatively related to collagen synthesis in osteoblasts . Studies have also indicated that osteoblast dysregulation is strongly associated with OA pathogenesis . Combined with the results of our study, these findings suggest that oxidoreductase activity may play a role in the mechanism of OA.
The oxidation-reduction process has been proven to be closely related to human aging and OA occurs more frequently in the elderly . However, the mechanism requires further research.
Biosynthesis of unsaturated fatty acids is one pathway identified as enriched in OA by TWAS. Its function in OA has been investigated in many studies. For example, oleic acid is an unsaturated fatty acid, and previous research has indicated that it could relieve OA symptoms by controlling inflammation . Eicosapentaenoic acid and docosahexaenoic acid are also unsaturated fatty acids and were both identified as substances that could reduce inflammation .
TWAS also detected the E2F1 destruction pathway as associated with hand OA. E2F1 is a transcription factor that has a critical function in controlling cell proliferation . Previous research on osteoarthritis shows that E2F1 could boost osteoclastogenesis and induce inflammation . The E2F1 destruction pathway may be a powerful therapeutic prospect for OA.
In this study, we compared the results of TWAS and mRNA expression profiling. Hand OA TWAS can detect genes with significant disease association at the DNA level, and mRNA expression profiling could provide information to explain the regulatory mechanisms of OA pathogenesis at the expression level. The combination of TWAS and mRNA expression profiling contributes to the more accurate identification of candidate genes.
Some limitations of this study should be noted. First, the GWAS summary data came from the UK Biobank cohort, in which almost all participants were of European ancestry. Therefore, it should be cautious to apply our study results to other populations. Second, to validate the TWAS results, we compared the significant genes identified by TWAS of hand OA with the differently expressed genes detected by the mRNA expression profiles of knee OA. The potential difference in the gene expression profiles between hand OA and knee OA may affect our study results, which should be interpreted with caution. Further biological studies should be conducted to confirm our findings.
In summary, we integrated the GWAS datasets of hand OA from the UK Biobank and precomputed gene expression weights of the peripheral blood and skeletal muscle to complete the TWAS. Then, we selected genes that were significantly differentially expressed in OA cases. We also compared the results of GO and pathway enrichment analyses at the gene and transcription levels. Five common genes and 6 GO terms were selected. We hope that our results will further the study of the mechanism of hand OA at the genetic and molecular levels.
Database for Annotation, Visualization and Integrated Discovery
Expression quantitative trait loci
Functional Summary-Based Imputation
Gene set enrichment analysis
Genome-wide association study
Transcriptome-wide association study
Leung GJ, Rainsford KD, Kean WF. Osteoarthritis of the hand I: aetiology and pathogenesis, risk factors, investigation and diagnosis. J Pharm Pharmacol. 2014;66(3):339–46.
Valdes K, Marik T. A systematic review of conservative interventions for osteoarthritis of the hand. J Hand Ther. 2010;23(4):334–51.
Ramonda R, Frallonardo P, Musacchio E, Vio S, Punzi L. Joint and bone assessment in hand osteoarthritis. Clin Rheumatol. 2014;33(1):11–9.
den Hollander W, Boer CG, Hart DJ, Yau MS, Ramos YFM, Metrustry S, Broer L, Deelen J, Cupples LA, Rivadeneira F, et al. Genome-wide association and functional studies identify a role for matrix Gla protein in osteoarthritis of the hand. Ann Rheum Dis. 2017;76(12):2046–53.
Taipale M, Solovieva S, Leino-Arjas P, Männikkö M. Functional polymorphisms in asporin and CILP together with joint loading predispose to hand osteoarthritis. BMC Genet. 2017;18(1):108.
Adams MK, Belyaeva OV, Wu L, Kedishvili NY. The retinaldehyde reductase activity of DHRS3 is reciprocally activated by retinol dehydrogenase 10 to control retinoid homeostasis. J Biol Chem. 2014;289(21):14868–80.
Tian J, Keller MP, Broman AT, Kendziorski C, Yandell BS, Attie AD, Broman KW. The dissection of expression quantitative trait locus hotspots. Genetics. 2016;202(4):1563–74.
Gong L, Zhang D, Lei Y, Qian Y, Tan X, Han S. Transcriptome-wide association study identifies multiple genes and pathways associated with pancreatic cancer. Cancer Med. 2018;7(11):5727–32.
Gusev A, Ko A, Shi H, Bhatia G, Chung W, Penninx BWJH, Jansen R, de Geus EJC, Boomsma DI, Wright FA, et al. Integrative approaches for large-scale transcriptome-wide association studies. Nat Genet. 2016;48(3):245–52.
Lin H, Lunetta KL, Zhao Q, Rong J, Benjamin EJ, Mendelson MM, Joehanes R, Levy D, Larson MG, Murabito JM. Transcriptome-wide association study of inflammatory biologic age. Aging (Albany NY). 2017;9(11):2288–301.
Canela-Xandri O, Rawlik K, Tenesa A. An atlas of genetic associations in UK Biobank. Nat Genet. 2018;50(11):1593–9.
Fisch KM, Gamini R, Alvarez-Garcia O, Akagi R, Saito M, Muramatsu Y, Sasho T, Koziol JA, Su AI, Lotz MK. Identification of transcription factors responsible for dysregulated networks in human osteoarthritis cartilage by global gene expression analysis. Osteoarthr Cartil. 2018;26(11):1531–8.
Cairns DM, Uchimura T, Kwon H, Lee PG, Seufert CR, Matzkin E, Zeng L. Muscle cells. Enhance resistance to pro-inflammatory cytokine-induced cartilage destruction. Biochem Biophys Res Commun. 2010;392(1):22–8.
van Vulpen LF, Schutgens RE, Coeleveld K, Alsema EC, Roosendaal G, Mastbergen SC, Lafeber FP. IL-1beta, in contrast to TNFalpha, is pivotal in blood-induced cartilage damage and is a potential target for therapy. Blood. 2015;126(19):2239–46.
Mazzoccoli G, Dagostino MP, Paroni G, Seripa D, Ciccone F, Addante F, Favuzzi G, Grandone E, Avola R, Mazza T, et al. Analysis of MTNR1B gene polymorphisms in relationship with IRS2 gene variants, epicardial fat thickness, glucose homeostasis and cognitive performance in the elderly. Chronobiol Int. 2017;34(8):1083–93.
Yang Z-Q, Zhang H-L, Duan C-C, Geng S, Wang K, Yu H-F, Yue Z-P, Guo B. IGF1 regulates RUNX1 expression via IRS1/2: implications for antler chondrocyte differentiation. Cell Cycle. 2017;16(6):522–32.
Wang G, Chan C-H, Gao Y, Lin H-K. Novel roles of Skp2 E3 ligase in cellular senescence, cancer progression, and metastasis. Chin J Cancer. 2012;31(4):169–77.
Loeser RF, Collins JA, Diekman BO. Ageing and the pathogenesis of osteoarthritis. Nat Rev Rheumatol. 2016;12(7):412–20.
McCulloch K, Litherland GJ, Rai TS. Cellular senescence in osteoarthritis pathology. Aging Cell. 2017;16(2):210–8.
Takayama K, Kawakami Y, Lee S, Greco N, Lavasani M, Mifune Y, Cummins JH, Yurube T, Kuroda R, Kurosaka M, et al. Involvement of ERCC1 in the pathogenesis of osteoarthritis through the modulation of apoptosis and cellular senescence. J Orthop Res. 2014;32(10):1326–32.
Feng Y, Li B, Li SJ, Yang X, Lv T, Shang H, Wu ZB, Zhang Y. Skp2/p27 axis regulates chondrocyte proliferation under high glucose induced endoplasmic reticulum stress. Eur Rev Med Pharmacol Sci. 2020;24(17):9129–38.
Ogasawara N, Kudo T, Sato M, Kawasaki Y, Yonezawa S, Takahashi S, Miyagi Y, Natori Y, Sugiyama A. Reduction of membrane protein CRIM1 decreases E-cadherin and increases claudin-1 and MMPs, enhancing the migration and invasion of renal carcinoma cells. Biol Pharm Bull. 2018;41(4):604–11.
Li C, Hu Q, Chen Z, Shen B, Yang J, Kang P, Zhou Z, Pei F. MicroRNA-140 suppresses human chondrocytes hypertrophy by targeting SMAD1 and controlling the bone morphogenetic protein pathway in osteoarthritis. Am J Med Sci. 2018;355(5):477–87.
Zhang J, Fu W, He M, Xie W, Lv Q, Wan G, Li G, Wang H, Lu G, Hu X, Jiang S, Li J, Lin M, Zhang Y, Kung H. MiRNA-20a promotes osteogenic differentiation of human mesenchymal stem cells by co-regulating BMP signaling. RNA Biol. 2011;8(5):829–38.
Deisenroth C, Itahana Y, Tollini L, Jin A, Zhang Y. p53-inducible DHRS3 is an endoplasmic reticulum protein associated with lipid droplet accumulation. J Biol Chem. 2011;286(32):28343–56.
Zheng X-y, Liang J, Li Y-S, Tu M. Role of fat-soluble vitamins in osteoarthritis management. JCR. 2018;24(3):132–7.
Wang H, Guo J, Jiang J, Wu W, Chang X, Zhou H, Li Z, Zhao J. New genes associated with rheumatoid arthritis identified by gene expression profiling. Int J Immunogenet. 2017;44(3):107–13.
Jiang R, Gao H, Cong F, Zhang W, Song T, Yu Z. Circ_DHRS3 positively regulates GREM1. expression by competitively targeting miR-183-5p to modulate IL-1β-administered chondrocyte proliferation, apoptosis and ECM degradation. Int Immunopharmacol. 2021;91:107293.
Zhang X, Zhao G, Zhang Y, Wang J, Wang Y, Cheng L, Sun M, Rui Y. Activation of JNK signaling in osteoblasts is inversely correlated with collagen synthesis in age-related osteoporosis. Biochem Biophys Res Commun. 2018;504(4):771–6.
Maruotti N, Corrado A, Cantatore FP. Osteoblast role in osteoarthritis pathogenesis. J Cell Physiol. 2017;232(11):2957–63.
Go Y-M, Jones DP. Redox theory of aging: implications for health and disease. Clin Sci. 2017;131(14):1669–88.
Felson DT, Bischoff-Ferrari HA. Dietary fatty acids for the treatment of OA, including fish oil. Ann Rheum Dis. 2016;75(1):1–2.
Murakami K. Potential of specialized pro-resolving lipid mediators against rheumatic diseases. Nihon Rinsho Men'eki Gakkai kaishi. 2016;39(3):155–63.
Harper JW, Elledge SJ. Skipping into the E2F1-destruction pathway. Nat Cell Biol. 1999;1(1):E5–7.
Murata K, Fang C, Terao C, Giannopoulou EG, Lee YJ, Lee MJ, Mun S-H, Bae S, Qiao Y, Yuan R, et al. Hypoxia-sensitive COMMD1 integrates signaling and cellular metabolism in human macrophages and suppresses osteoclastogenesis. Immunity. 2017;47(1):66–79.
This study was financially supported by grants from the National Natural Science Foundation of China (81672219, 81802210 and 81974347). Financial support had no impact on the outcomes of this study.
Ethics approval and consent to participate
Consent for publication
The authors declare that they have no competing interests.
Springer Nature remains neutral with regard to jurisdictional claims in published maps and institutional affiliations.
. TWAS results of hand osteoarthritis.
. GO enrichment analysis results of candidate genes detected by TWAS of hand osteoarthritis.
Pathway enrichment analysis results of candidate genes identified by the TWAS of hand osteoarthritis.
About this article
Cite this article
Xu, J., Zeng, Y., Si, H. et al. Integrating transcriptome-wide association study and mRNA expression profile identified candidate genes related to hand osteoarthritis. Arthritis Res Ther 23, 81 (2021). https://doi.org/10.1186/s13075-021-02458-2