Abstract
Objective
The aim of this study was to identify the susceptibility genes responsible for lumbar spondylosis (LS) in Korean patients.
Methods
Data from 1427 subjects were made available for radiographic grading and genome wide association studies (GWAS) analysis. Lateral lumbar spine radiographs were obtained and the various degrees of degenerative change were semi-quantitatively scored. A pilot GWAS was performed using the AffymetrixGenome-Wide Human single-nucleotide polymorphisms (SNPs), 500K array. A total of 352228 SNPs were analyzed and the association between the SNPs and case-control status was analyzed by stepwise logistic regression analyses.
Results
The top 100 SNPs with a cutoff p-value of less than 3.7×10−4 were selected for joint space narrowing, while a cutoff p-value of 6.0×10−4 was applied to osteophytes and the Kellgren-Lawrence (K-L) osteoarthritis grade. The SNPs with the strongest effect on disc space narrowing, osteophytes, and K-L grade were serine incorporator 1 (rs155467, odds ratio [OR]=17.58, p=1.6×10−4), stromal interaction molecule 2 (STIM1, rs210781, OR=5.53, p=5×10−4), and transient receptor potential cation channel, subfamily C (rs11224760, OR=3.99, p=4.8×10−4), respectively. Leucine-rich repeat-containing G protein-coupled receptor 4 was significantly associated with both disc space narrowing and osteophytes (rs1979400, OR=2.01, p=1.1×10−4 for disc space narrowing, OR=1.79, p=3×10−4 for osteophytes), while zinc finger and BTB domain containing 7C was significantly and negatively associated with both osteophytes and a K-L grade >2 (rs12457004,OR=0.25, p=5.8×10−4 and OR=0.27, p=5.3×10−4, respectively).
Low back pain (LBP) is an important public health problem and is associated with substantial societal costs in industrialized countries. Degenerative disease of the lumbar spine (lumbar spondylosis, LS), which is characterized radiologically by the presence of osteophytes, endplate sclerosis, and disc space narrowing13), is thought to be related to LBP. The prevalence of LS in Asian countries has scarcely been reported. However, ethnic differences in the prevalence of osteophytosis and disc degeneration have previously been noted22). Although it is widely known that the prevalence of LS increases with age, there exists a paucity of data on the risk factors for LS. LS shares many features with osteoarthritis (OA) of the peripheral joint, the most common form of arthritis affecting the elderly. Although the pathogenesis of OA has been associated with aging and environmental factors such as work activity, twin studies have revealed a significant genetic contribution, with a heritability estimate of 65% for hand OA20). Recent reports have demonstrated the significant role that genetic contribution also plays in LS. A twin study showed that the heritability of LS was 74% for the lumbar spine and 73% for the cervical spine2). Other studies also revealed significant heritability on the basis of magnetic resonance imaging (MRI) of the spine9,16). The search for specific genes associated with OA and LS has since ensued. However, due to variations in the phenotype, as well as the severity of the disease and patient characteristics such as gender and age, the definitions of the disease have varied from study to study, leading to tremendous difficulty in identifying a causative genetic factor.
A large number of genetic association studies of OA have been conducted over the last decade, and a number of significant associations have been verified by a systematic review. For example, aggrecan 1 (AGC1) is associated with hand OA, asporin (ASPN) with hip and knee OA, and vitamin D receptor (VDR) with knee and spine OA15). However, the majority of the genes that were reported to have significant associations have never been replicated, making it difficult to discern the significance of these isolated significant associations in the pathogenesis of OA.
Complex conditions such as OA and LS are caused by numerous genetic and environmental factors, any one of which can have a relatively minor effect. Genome-wide association studies (GWAS) have provided a breakthrough in deciphering the genetic influence in complex diseases, by examining hundreds of thousands of single-nucleotide polymorphisms (SNPs), thus enabling identification of SNPs implicated in a large number of robustly replicated loci of common traits.
Although GWAS has begun to unravel the genetic influence of peripheral joint OA25), there has been little attempt as yet to do the same for LS. In addition, the results of an LS GWAS among an Asian population have not yet been reported. In this study, we sought to identify the susceptibility genes for LS in a GWAS among a community-based population of Korean adults.
For this study, we selected a rural farming community (Ansung) in South Korea from among the populations incorporated into the ongoing prospective Korean Genome and Epidemiology Study (KoGES). The study methods have been described previously4). Briefly, the eligibility criteria included an age of 40–79 years, residence within the borders of the survey area for at least 6 months before testing, and the mental and physical ability to participate. Data relating to LBP were available for 4181 subjects for the years 2006 and 2007. Baseline demographic information was collected using a standard questionnaire during a face-to-face interview and included educational attainment, occupation, exercise, and co-morbidities.
Of the total number of eligible subjects, 2000 were randomly selected for spine radiography. No difference was found in the prevalence of LBP between those who underwent radiography and those who did not. After excluding 488 subjects who were unable to be evaluated due to poor film quality, and 85 patients whose genomic DNA was not obtained, data for 1427 subjects were available for radiograph grading and GWAS analysis. The study protocol was approved by the Ethics Committee of the KoGES, and written informed consent was obtained from each participant.
Lateral lumbar spine radiographs were taken according to a standard protocol with the film centered on the second lumbar vertebra. Lumbar radiographs were evaluated by a single observer. Each vertebral level from L1/2 to L4/5 was reviewed for the presence of radiographic features relating to degenerative change. Semi-quantitative scores were given for the following features using a reference atlas8) : presence and severity of anterior osteophytes (grade 0=none; grade 1=barely visible; grade 2=definite; grade 3=large), endplate sclerosis (grade 0=none; grade 1=present), and disc space narrowing (grade 0=none; grade 1=barely visible; grade 2=definite; grade 3=severe, bone to bone).
Additionally, the Kellgren-Lawrence (K-L) grading system was used for each vertebral level (grade 0=normal disc with no osteophytes; grade 1=slight anterior wear and osteophyte formation; grade 2=definite anterior wear and mild disc space narrowing with osteophyte formation; grade 3=moderate disc space narrowing with osteophytes and sclerosis; grade 4=large osteophytes, marked disc space narrowing, and sclerosis of vertebral end plates). LS was defined using one of the following criteria : 1) joint space narrowing ≥grade 2; 2) osteophytes ≥grade 2; 3) K-L grade ≥2. Radiographs were read by a single reader who was an academically-based rheumatologist. Intra-observer reproducibility was assessed by re-evaluating 50 films within 1 week of the first reading. The reproducibility of intra-reader assessments was high (for osteophyte grading, κ=0.89–0.93; for endplate sclerosis, κ=0.71–84; for joint space narrowing, κ=0.81–0.89; and for K-L grading, κ=0.69–0.80, for various vertebral levels). Films allocated different grades at each of the two readings were adjudicated by consensus between the original reader and a second reader.
Genomic DNA was isolated from peripheral blood mononuclear cells. We performed a pilot GWAS, typing cases and controls on a single platform using the Affymetrix Genome-Wide Human SNP 500K array chip (Affymetrix, Inc., Santa Clara, CA, USA). Genotype calls were determined by Bayesian robust linear modeling using the Mahalanobis distance algorithm14). We sequentially discarded 38364 markers with a Hardy-Weinberg equilibrium p-value <10−6, 17926 with genotype call rates below 95%, and 92050 with a minor allele frequency (MAF) of 0.01. This left 352228 SNPs available for subsequent analysis
The evaluation of the association between case-control status and each individual SNP was based on the odds ratio (OR) and p-values. Covariates used for multivariable-adjustment were age, sex, education, and body mass index. To optimize the joint effect of the SNPs, we conducted stepwise logistic regression analyses and also analyzed the main effects on the selected variables using a generalized linear model implemented with the R statistical software (ver. 2.14.2; R Core Development Team, Austria, Vienna).
Table 1 displays the clinical characteristics of the study population, with cases of spine OA defined as disc space narrowing, osteophytes, or K-L grade all ≥2. Subjects with spine OA were older, more likely to be male, and had a lower level of education compared to the control subjects.
We used logistic regression analysis to identify statistically significant associations with LS. For multiple comparisons, a Bonferroni-corrected p-value of 1.5×10−7 (0.05/320942) was used. Since none of the SNPs reached the very conservative Bonferroni corrected value, the top 100 SNPs with a cutoff p-value less than 3.7×10−4 were selected for use in this GWAS for joint space narrowing while a cutoff p-value of 6.0×10−4 was applied for osteophytes and K-L grade. Tables 2–4 present the detailed characteristics of the SNPs, including gene name, rs identification number, position, MAF of the cases and controls, OR, and p-value. Each feature of LS was found to be associated with distinct SNPs. The SNPs that had the strongest effect on disc space narrowing included serine incorporator 1 (rs155467, OR=17.58, p=1.6×10−4), heat shock transcription factor 2 (rs563084, OR=14.48, p=1.4×10−4), cysteine-rich hydrophobic domain (rs1568512, OR=11.57, p=4×10−5), akirin 2 (rs2787938, OR=9.86, p=8×10−5), and fibroblast growth factor receptor 2 (rs11200052, OR=9.72, p=9×10−5), while those that had the most significant effect on osteophytes included stromal interaction molecule 2 (STIM1, rs210781, OR=5.53, p=5×10−4), protein kinase C and casein kinase substrate in neurons 2 (PACSIN2, rs738379, OR=5.37, p=3.6×10−4), and ubiquinol-cytochrome c reductase complex chaperone (rs6060373, OR=3.05, p=8×10−5). Transient receptor potential cation channel, subfamily C (rs11224760, OR=3.99, p=4.8×10−4), hypocretin (orexin) receptor 1 (rs3753613, OR=2.62, p=2×10−4), and coiled-coil domain containing 60 (rs10849640, OR=2.59, p=3.8×10−4) had the strongest effect on K-L grades >2. Leucine-rich repeat-containing G protein-coupled receptor 4 was significantly associated with both disc space narrowing and osteophytes (rs1979400, OR=2.01, p=1.1×10−4 for disc space narrowing, OR=1.79, p=3×10−4 for osteophytes), while zinc finger and BTB domain containing 7C was significantly and negatively associated with both osteophytes and a K-L grade >2 (rs12457004, OR=0.25, p=5.8×10−4 for osteophytes, OR=0.27, p=5.3×10−4 for K-L grade). No SNP was found to be significantly associated with LS defined using all three criteria.
In this study, we investigated the association of clinical and genetic factors for LS defined by three features found on simple radiographs. The previously identified genetic risk factors for LS using a candidate approach included genes encoding extracellular matrix proteins expressed in the nucleus pulposus (inner structure) and annulus fibrosus (outer layer) of the disc, such as type IX collagen (COL9A2 and COL9A3), AGC1, and cartilage intermediate layer protein (CILP)15). These findings suggest that LS is caused by changes in the structural integrity of the intervertebral disc. Given that OA and LS are both degenerative diseases of the skeletal joints, and because articular cartilage and intervertebral discs share similar patterns of gene expression, it is expected that subjects with OA and LS have a similar genetic susceptibility. A systematic review of genetic associations for peripheral joint OA and degenerative disease of the spine revealed the difficulties related to complex phenotypes of these diseases. While many logical and reporting problems, including missing population details, multiple testing, and an over-reliance on subgroup analysis, were found, cases in which significant associations were replicated in independent studies were also identified15). Some of the genes thus identified are of functional importance. For example, ASPN, a member of the small leucine-rich proteoglycan family, inhibits in vitro chondrogenesis and the expression of COL2A1 and AGC1 through the inhibition of TGF-b signaling19).
A recent systematic search of the literature, including 52 studies, identified ASPN (D-repeat), COL11A1 (rs1676486), growth differentiation factor 5 (GDF5) (rs143383), SKT (rs16924573), THBS2 (rs9406328), and MMP9 (rs17576) as genes associated with LS in humans defined by MRI with a moderate level of evidence5). However, the phenotype definition of lumbar disc degeneration was highly variable among the studies, including a decrease in disc signal intensity or disc height, disc bulges, disc herniations without specification of the symptoms, Modic changes, osteophytes, and lumbar spinal stenosis. In addition, the phenotype of disc degeneration varied between the initial and replication studies. As a result, the replications were inconsistent, and most of the associations were presented with a weak level of evidence.
We used simple x-ray for diagnose and grade lumbar spondylosis. The K-L grade is most popular grading system with classification into five grade scales (0–4) where K-L grade ≥2 is the conventional standard of the diagnosis8). Epidemiological studies showed that K-L grade was associated with the degree of low back pain in elderly subjects4,11).
The variation in the identified genetic associations may reflect ethnic diversity as well as phenotypic heterogeneity. While ASPN was identified as a candidate gene for lumbar disc degeneration with a moderate level of epidemiological evidence among Asians, GDF5 was instead identified among the Northern European population24). In addition, flaws inherent in a candidate gene approach, such as isolated analyses of disparate potential associations, may merely add to the growing repertoire of weak evidence. The genome-wide association study is based on the ‘common disease-common variation’ theory, which proposes that multiple common polymorphisms with a MAF of 0.5–1%, and with small effect sizes, might be predisposed to common disorders1). Because more complex modes of inheritance involving multiple genes rather than a simple monogenic Mendelian pattern appear to be in effect for LS, the discovery of its associated genes is more likely with a GWAS than by an SNP analysis. A meta-analysis of 4 GWAS that addressed LS among 4683 individuals with European ancestries was reported by Williams et al.23) LS was defined as disc space narrowing and the formation of osteophytes, as in our study. Among the four markers identified (rs17034687, rs2187689, rs7767277, and rs926849), the rs926849 marker located in the intronic region of Parkinson protein 2, the E3 ubiquitin protein ligase (PARK2) gene, on chromosome 6 remained strongest after adjusting for age and gender. The rs2187689 and rs7767277 markers were in strong linkage disequilibrium with proteasome subunit beta type 9 (large multifunctional peptidase 2) (PSMB9). Although we used the same 500K Affymetrix kit, the genes identified in our study did not overlap with those reported by Williams et al.24). Again, the discrepancy may arise from ethnic differences between the study populations and the difference in phenotype definition.
The SNP of leucine-rich repeat-containing G protein-coupled receptor 4 (LGR4) was significantly associated with both disc space narrowing and osteophytes in our study. LGR4 encodes a receptor for R-spondins (RSPOs), which play a pleiotropic role in normal development, and the development of cancer, as well as in the survival of adult stem cells through potentiation of Wnt signaling7). Knockout of LGR4 or RSPOs in mice presents with severe developmental abnormalities, causing neonatal/embryonic lethality10). Recent studies showed that the RSPO-LGR4 axis elevates the levels of Wnt receptors through direct inhibition of the ubiquitination of Wnt receptors or though interaction with intracellular signaling proteins to potentiate the Wnt pathways3). Increased Wnt/β-catenin signaling is observed in the nucleus pulposus of dogs that suffer from premature intervertebral disc degeneration18), suggesting the role of Wnt/β-catenin signaling, and possibly the RSPO-LGR4 axis, in the pathogenesis of LS.
The SNP of the heat shock transcription factor 2 gene (HSTF2) was one of the SNPs most significantly associated with LS as defined by disc space narrowing. HSTF2 specifically binds to the heat-shock promoter element and activates transcription for heat-shock response genes under conditions of heat or other stress. Previous studies showed that the formation of heterocomplexes between HSTF1 and HSTF2 leads to enhanced activity, which activates the hsp70 promoter, and that HSTF2 was able to modulate the HSTF1-mediated expression of major heat shock protein genes12). A previous study showed that HSTF expression was more frequent in clustered cells in both the annulus fibrosus and nucleus pulposus of herniated discs17). Because the K-L grade is accounted for more by osteophytes than by joint space narrowing, it is not surprising that genes associated with LS, as defined by the K-L grade, are different from those defined by disc space narrowing. The strongest association with K-L grade was observed with the SNP for transient receptor potential cation channel (TRPC), subfamily C, which encodes a receptor-activated calcium channel. In a study using OA chondrocytes, a correlation between the appearance of TRPC6 and the state of de-differentiation of chondrocytes was identified6). Because the loss of a differentiated chondrocyte phenotype is one of the hallmarks of OA, as well as LS, the role of TRPC6 in the process of disc degeneration is plausible. On the other hand, TRPC6 was expressed in bone-derived cells and inhibitors of the TRP channel inhibited the effects of bradykinin-induced Ca2+-influx, suggesting a role in bone metabolism21).
This is the first reported GWAS in an Asian population. In terms of phenotyping, we used disc space narrowing as well as osteophytosis and K-L grading to cover all aspects of degenerative change observable in LS. Limitations included a sample size too small to detect the true effect of SNPs. Although we did not identify the same loci reported in previous studies using a candidate gene approach, differences in ethnicity and phenotype definition may have affected the results. Further studies using larger Asian samples are warranted to establish if the gene variants identified in this study are associated with an increase in the risk of LS in Asian populations.
Acknowledgements
This research was supported by a grant of the Korea Health Technology R&D Project through the Korea Health Industry Development Institute (KHIDI), funded by the Ministry of Health & Welfare, Republic of Korea (grant number : HI14C2248).
References
2. Battié MC, Videman T, Gibbons LE, Fisher LD, Mannine H, Gill K. 1995 Volvo Award in clinical sciences. Determinants of lumbar disc degeneration A study relating lifetime exposures and magnetic resonance imaging findings in identical twins. Spine. 20:2601–2612. 1995.
3. Carmon KS, Gong X, Yi J, Thomas A, Liu Q. RSPO-LGR4 functions via IQGAP1 to potentiate Wnt signaling. Proc Natl Acad Sci U S A. 111:E1221–E1229. 2014.
4. Cho NH, Jung YO, Lim SH, Chung CK, Kim HA. The prevalence and risk factors of low back pain in rural community residents of Korea. Spine (Phila Pa 1976). 37:2001–2010. 2012.
5. Eskola PJ, Lemmelä S, Kjaer P, Solovieva S, Männikkö M, Tommerup N, et al. Genetic association studies in lumbar disc degeneration: a systematic review. PLoS One. 7:e49995. 2012.
6. Gavenis K, Schumacher C, Schneider U, Eisfeld J, Mollenhauer J, Schmidt-Rohlfing B. Expression of ion channels of the TRP family in articular chondrocytes from osteoarthritic patients: changes between native and in vitro propagated chondrocytes. Mol Cell Biochem. 321:135–43. 2009.
7. Glinka A, Dolde C, Kirsch N, Huang YL, Kazanskaya O, Inqelfinger D, et al. LGR4 and LGR5 are R-spondin receptors mediating Wnt/β-catenin and Wnt/PCP signalling. EMBO Rep. 12:1055–1061. 2011.
8. Lane NE, Nevitt MC, Genant HK, Hochberg MC. Reliability of new indices of radiographic osteoarthritis of the hand and hip and lumbar disc degeneration. J Rheumatol. 20:1911–1918. 1993.
9. Matsui H, Kanamori M, Ishihara H, Yudoh K, Naruse Y, Tsuji H. Familial predisposition for lumbar degenerative disc disease. A case-control study. Spine (Phila Pa 1976). 23:1029–1034. 1998.
10. Mazerbourg S, Bouley DM, Sudo S, Klein CA, Zhang JV, Kawamura K, et al. Leucine-rich repeat-containing, G protein-coupled receptor 4 null mice exhibit intrauterine growth retardation associated with embryonic and perinatal lethality. Mol Endocrinol. 18:2241–2254. 2004.
11. Muraki S, Oka H, Akune T, Mabuchi A, En-yo Y, Yoshida M, et al. Prevalence of radiographic lumbar spondylosis and its association with low back pain in elderly subjects of population-based cohorts: the ROAD study. Ann Rheum Dis. 68:1401–1406. 2009.
12. Ostling P, Björk JK, Roos-Mattjus P, Roos-Mattjus P, Mezger V, Sistonen L. Heat shock factor 2 (HSF2) contributes to inducible expression of hsp genes through interplay with HSF1. J Biol Chem. 282:7077–7086. 2007.
13. Pye SR, Reid DM, Smith R, Adams JE, Nelson K, Silman AJ, et al. Radiographic features of lumbar disc degeneration and self-reported back pain. J Rheumatol. 31:753–758. 2004.
14. Rabbee N, Speed TP. A genotype calling algorithm for affymetrix SNP arrays. Bioinformatics. 22:7–12. 2006.
15. Ryder JJ, Garrison K, Song F, Hooper L, Skinner J, Loke Y, et al. Genetic associations in peripheral joint osteoarthritis and spinal degenerative disease: a systematic review. Ann Rheum Dis. 67:584–591. 2008.
16. Sambrook PN, MacGregor AJ, Spector TD. Genetic influences on cervical and lumbar disc degeneration: a magnetic resonance imaging study in twins. Arthritis Rheum. 42:366–372. 1999.
17. Sharp CA, Roberts S, Evans H, Brown SJ. Disc cell clusters in pathological human intervertebral discs are associated with increased stress protein immunostaining. Eur Spine J. 18:1587–1594. 2009.
18. Smolders LA, Meij BP, Onis D, Riemers FM, Bergknut N, Wubbolts R, et al. Gene expression profiling of early intervertebral disc degeneration reveals a down-regulation of canonical Wnt signaling and caveolin-1 expression: implications for development of regenerative strategies. Arthritis Res Ther. 29:R23. 2013.
19. Song YQ, Cheung KM, Ho DW, Poon SC, Chiba K, Kawaguchi Y, et al. Association of the asporin D14 allele with lumbar-disc degeneration in Asians. Am J Hum Genet. 82:744–747. 2008.
20. Spector TD, Cicuttini F, Baker J, Loughlin J, Hart D. Genetic influences on osteoarthritis in women: a twin study. BMJ. 312:940–943. 1996.
21. Suzuki Y, Kodama D, Goto S, Togari A. Involvement of TRP channels in the signal transduction of bradykinin in human osteoblasts. Biochem Biophys Res Commun. 410:317–321. 2011.
22. van Tulder MW, Assendelft WJ, Koes BW, Bouter LM. Spinal radiographic findings and nonspecific low back pain. A systematic review of observational studies. Spine (Phila Pa 1976). 22:427–434. 1997.
23. Williams FM, Bansal AT, van Meurs JB, Bell JT, Meulenbelt I, Suri P, et al. Novel genetic variants associated with lumbar disc degeneration in northern Europeans: a meta-analysis of 4600 subjects. Ann Rheum Dis. 72:1141–1148. 2013.
Table 1
Disc space narrowing | Osteophyte | K-L grade | |||||||
---|---|---|---|---|---|---|---|---|---|
|
|
|
|||||||
Case | Control | Total | Case | Control | Total | Case | Control | Total | |
Age (years) | 64.7±0.3* | 58.5±0.4 | 59.3±0.2 | 62.4±0.3* | 56.3±0.3 | 59.5±0.2 | 61.2±0.3* | 55.3±0.4 | 59.5±0.2 |
|
|||||||||
Female | 132 (57.9) | 699 (58.4) | 831 (58.3) | 369 (48.4)* | 462 (69.6) | 831 (58.2) | 536 (52.1)* | 293 (74.2) | 829 (58.2) |
|
|||||||||
BMI | 24.1±0.2 | 24.4±0.1 | 24.3±0.1 | 24.3±0.1 | 24.3±0.1 | 24.3±0.1 | 24.3±0.1 | 24.2±0.2 | 24.3±0.1 |
|
|||||||||
Smoker | 45 (19.7) | 222 (18.6) | 267 (18.7) | 168 (22.1)* | 100 (15.1) | 268 (18.8) | 215 (20.9)* | 52 (13.2) | 267 (18.8) |
|
|||||||||
Alcohol | 79 (34.7)* | 505 (42.2) | 584 (41.0) | 323 (42.3) | 261 (39.3) | 584 (40.9) | 429 (41.7) | 155 (39.2) | 584 (41.0) |
|
|||||||||
Employed | 216 (94.7) | 1124 (93.8) | 1340 (94.0) | 723 (94.8) | 617 (92.2) | 1340 (93.9) | 974 (94.7) | 364 (92.2) | 1338 (94.0) |
|
|||||||||
Education ≤6 years | 64 (28.1)* | 196 (16.4) | 260 (18.2) | 159 (20.8)* | 101 (15.2) | 260 (18.2) | 207 (20.1)* | 51 (12.9) | 258 (18.1) |