Gene-gene interactions between TGF-β/Smad3 signalling pathway polymorphisms affect susceptibility to knee osteoarthritis.

OBJECTIVE: Transforming growth factor/Smad family member 3 (TGF)-β/Smad3 signalling is essential for maintaining articular cartilage. A relationship between the genetic variants of TGF-β itself, TGF-β signalling and binding molecules, and osteoarthritis (OA) has been reported. Although variants of candidate genes have become prime targets for genetic analysis, their detailed interplay has not been documented. Our goal was to establish whether single nucleotide polymorphisms (SNPs) of TGF-β1, TGF-βRI, Smad3 and tissue inhibitor of metalloproteinases 3 (TIMP3), and their interactions, are associated with knee OA.
DESIGN: We performed a case-control association study and genotyped 518 knee patients with OA and 468 healthy controls. All participants were genotyped for TGF-β1 (rs1800469C/T), TGF-βRI (rs1590A/G), Smad3 (rs12901499A/G and rs6494629T/C), and TIMP3 (rs715572G/A and rs1962223G/C) polymorphisms by polymerase chain reaction-restriction fragment length polymorphism analysis. Multifactor dimensionality reduction (MDR) was used to identify gene-gene interactions.
RESULTS: Significant associations were observed for TIMP3 rs715572G/A polymorphisms in knee patients with OA and healthy individuals. The GA heterozygote in TIMP3 (rs715572G/A) was significantly associated with OA (p=0.007). Patient stratification using the Kellgren-Lawrence grading scale showed significant differences in TIMP3 rs715572G/A genotypes between grade 4 knee OA and controls. By MDR analysis, a two-locus model (Smad3 rs6494629T/C and TIMP3 rs715572G/A) of gene-gene interaction was the best for predicting knee OA risk, and its maximum testing accuracy was 57.55% and maximum cross-validation consistency was 10/10.
CONCLUSIONS: TIMP3 rs715572G/A is a candidate protective gene for severe knee OA. Gene-gene interactions between Smad3 rs6494629T/C and TIMP3 rs715572G/A polymorphisms may play more important protective roles in knee OA. Published by the BMJ Publishing Group Limited. For permission to use (where not already granted under a licence) please go to http://group.bmj.com/group/rights-licensing/permissions.

This study is the first population-based study to evaluate the interactions between single nucleotide polymorphism variants of the transforming growth factor/Smad family member 3 (TGF-β)/Smad3 signalling pathway for knee osteoarthritis (OA).

Our results indicate tissue inhibitor of metalloproteinases 3 (TIMP3) rs715572G/A is associated with more severe knee OA.

Our study highlights the importance of the effect of interactions between Smad3 rs6494629T/C and TIMP3 rs715572G/A polymorphisms for knee OA, which would be likely to be missed if genes are individually examined without considering potential related pathways.

In future research, the mechanisms of interactions between Smad3 rs6494629T/C and TIMP3 rs715572G/A polymorphisms and their effects on knee OA need to be established.

Introduction

Osteoarthritis (OA) is the most common form of arthritis and is a leading cause of disability in the elderly. An increasing body of evidence suggests that ageing, genetic predisposition, obesity, inflammation and excessive mechanical loading predispose to OA development.1 The mechanisms by which these risk factors predispose to the development of OA are beginning to be explored and understood. Irrespective of the initiating event, OA results from an imbalance in catabolic and anabolic processes, which leads to progressive cartilage damage and destruction.2

The heritable component of OA is estimated to be around 40–65%. Candidate gene studies and, more recently, genome-wide association studies, are beginning to help identify key genetic factors that may influence susceptibility to onset and progression of OA.3–5 Candidate gene studies and subsequent large-scale studies and meta-analyses suggest that polymorphisms ASPN and GDF5 are associated with OA.6–8 The gene for GDF5 codes for growth differentiation factor 5 is a member of the TGF-β superfamily and has important roles in skeletal and joint development with mutations resulting in a range of skeletal abnormalities.9 10 Biological studies indicate that the rs143383 single nucleotide polymorphisms (SNP) in GDF5 results in reduced GDF5 transcription in joint tissues, which in turn may be important in OA development.11 ASPN in turn encodes for asporin, a member of the sub family of small leucine-rich proteoglycans. Functionally, asporin binds to transforming growth factor-β (TGF-β), preventing its binding to the TGF-β type II receptor and inhibiting TGF-β-induced expression of anabolic cartilage molecules including aggrecan and type II collagen.12 The effect on TGF-β activity is allele-specific, with the D14 allele, which is associated with OA, causing a greater inhibition of TGF-β activity than other alleles.13

TGF-β is a pleiotropic cytokine/growth factor with important anabolic effects on chondrocytes14 and, as such ,TGF-β signalling, especially via the Smad family member 3 (Smad3), which plays a pivotal role in the homeostasis of synovial joints.15 In the classical TGF-β/Smad signalling pathway, phosphorylated Smad3 forms a complex with Smad4; this complex then translocates to the nucleus to regulate gene expression and promote an anabolic phenotype in cartilage.16 This includes TGF-β-induced production of a tissue inhibitor of metalloproteinases 3 (TIMP3) via the PI3K/Akt signalling pathway.17 By inhibiting activity of matrix metalloproteinases, a disintegrin and metalloproteinase with thrombospondin motifs 4/5 (ADAMTS-4/5) and tumour necrosis factor (TNF-α) converting enzyme (TACE/ADAM-17) TIMP3 acts to reduce joint inflammation and cartilage matrix resorption.18

A relationship between the genetic variants of TGF-β itself, TGF-β signalling and binding molecules, and OA, has been reported in humans.19 Polymorphic variants of TGF-β1, TGF-βRI, Smad3 and TIMP3 may be functionally expressed, suggesting that SNPs are among the factors associated with susceptibility to OA. The genetic aetiology of OA is likely to involve interactions between multiple genetic variants of molecules within important chondroprotective pathways such as the TGF-β/Smad3 axis. The current study was therefore undertaken to assess whether interactions between multiple SNP variants of TGF-β1, TGF-βRI, Smad3 and TIMP3, were associated with knee OA.

Materials and methods

Subjects

This case–control study included 518 knee patients with OA (328 females and 190 males; age 72.98±7.57 years) received at Tri-Service General Hospital, Taipei, Taiwan. Disease severity in the knee OA population was assessed using the Kellgren–Lawrence (K–L) grading scale. All patients had a K–L grade ≥2. Knee joint diseases of other aetiologies such as inflammatory arthritis, post-traumatic or postseptic arthritis, skeletal dysplasia or developmental dysplasia, were excluded. The study also included 468 healthy control participants (261 females and 207 males; mean age 69.59±9.30 years) with no symptoms of joint disease (pain, swelling, tenderness or restriction of movement) in whom standard X-rays of knee joints confirmed the absence of radiographical knee OA. All clinical and biological samples were collected, and DNA was genotyped after obtaining the approval of this committee. After full explanation of the study, written informed consent was obtained from all participants.

Radiographic assessment

All participants underwent weight-bearing anteroposterior radiographs to assess the structural changes of the affected knee. Radiographic severity was assessed according to the Kellgren-Lawrence (K-L) grading system: grade 1, doubtful narrowing of joint space and possible osteophytic lipping; grade 2, definite osteophytes and possible narrowing of joint space; grade 3, moderate multiple osteophytes, definite narrowing of joints space, some sclerosis and possible deformity of bone contour; grade 4, large osteophytes, marked narrowing of joint space, severe sclerosis and definite deformity of bone contour. An experienced observer, who was blinded to the source of participants, scored the grading of radiographs. All participants whose K-L grade was less than 2 were included in this study as normal controls.

SNP selection and genotyping

We selected TGF-β1, TGF-βRI, Smad3 and TIMP3 as candidate genes based on the published literature.20–23 To select the most representative SNPs by capturing the majority of genetic variations, SNP genotype information was downloaded from the HapMap database (http://www.hapmap.ncbi.nlm.nih.gov/) and the National Center for Biotechnology Information dbSNP database (http://www.ncbi.nlm.nih.gov/snp). Tag SNPs were selected for TGF-β1, TGF-βRI, Smad3 and TIMP3, using the criterion of minor allele frequency (MAF) >10%. We also examined SNPs in regulatory regions and those reported by other investigators. Genomic DNA was extracted from the peripheral blood of patients and controls using the QIAamp DNA Blood Mini Kit (QIAGEN Inc, Hilden, Germany) and stored at −20°C until genotyping. TGF-β1 (rs1800469C/T), TGF-βRI (rs1590A/G), Smad3 (rs12901499A/G and rs6494629T/C) and TIMP3 (rs715572G/A and rs1962223G/C) polymorphisms were screened by polymerase chain reaction–restriction fragment length polymorphism (PCR–RFLP) analysis. The primer design was based on published sequences24 or designed using the Primer Z software (http://genepipe.ngc.sinica.edu.tw/primerz/beginDesign.do). PCR cycling conditions were: an initial denaturation at 95°C for 5 min, followed by 35 denaturation cycles at 95°C for 30 s, annealing at 55°C for 30 s, extension at 72°C for 30 s and final extension at 72°C for 7 min. The PCR products were digested with appropriate restriction endonucleases (New England Biolabs, Inc, Ipswich, USA). Resulting fragments were separated in 2.5% agarose gel containing 0.5 μg/mL ethidium bromide by electrophoresis at 100V and visualised under UV light. Genotyping was performed by laboratory personnel blinded to the case status, and 10% of the samples were randomly selected for repeated testing to validate genotyping procedures. Two authors independently reviewed the genotyping results, data entry and statistical analyses. Online supplementary I summarises SNP description and RFLP condition.

Statistical methods

Demographics were evaluated by Student's t test for continuous variables and expressed as mean±SD. The Hardy–Weinberg equilibrium (HWE) test was assessed by a goodness-of-fit χ2 test and was performed to examine possible genotyping error for each SNP among the controls. Genotypes and allelic frequencies were compared between knee patients with OA and healthy controls using the χ2 or Fisher's exact test, when appropriate. Logistic regression was used to estimate crude and adjusted (age, gender and body mass index) ORs and 95% CIs as a measure of the association with knee OA risk.

The level of significance was determined by Bonferroni's method for correcting multiple testing errors. Under the selected six SNPs, a p value of less than 0.0083 (0.05 divided by 6) was considered statistically significant. Statistical analysis was performed with SPSS for Windows, V.18.0 (SPSS, Chicago, Illinois, USA). To investigate the effect of gene–gene interaction on OA, multifactor dimensionality reduction (MDR) (V.2.0 β) and MDR–permutation testing software applications (V.1.0 β) were employed. In addition, the logistic regression model was performed to confirm the results of gene–gene interaction analyses.

Results

Basic characteristics of the study population

The demographic and clinical characteristics of knee OA cases (n=518) and the controls (n=468) are shown in table 1. Overall, patients with OA were significantly older than control individuals, and were more likely to be obese.

Table 1

Characteristics of the study population in case and control participants

	Case	Control	p Value
Number	518	468
Age	72.98±7.57	69.59±9.30	<0.001
Gender
Male	190 (36.7%)	207 (44.2%)	0.016
Female	328 (63.3%)	261 (55.8%)
BMI	25.81±3.33	24.40±3.72	<0.001
K–L grade
0	0	246
1	0	222
2	194	0
3	104	0
4	220	0

BMI, body mass index; K–L, Kellgren–Lawrence.

TGF-β, TGF-βRI, Smad3 and TIMP3 allele and genotype frequencies

TGF-β1 (rs1800469C/T), TGF-βRI (rs1590A/G), Smad3 (rs12901499A/G and rs6494629T/C) and TIMP3 (rs715572G/A and rs1962223G/C) genotype distributions were compatible with the HWE in knee OA cases and controls (p>0.05). This indicates that the study participants were representative of the study field. The genotype and allele frequencies of six SNPs in knee patients with OA and healthy controls are presented in table 2. There were no significant differences between the genotype or allele frequencies of TGF-β1 (rs1800469C/T), TGF-βRI (rs1590A/G), Smad3 (rs12901499A/G and rs6494629T/C) and TIMP3 (rs1962223G/C) polymorphisms in the patient and control groups. SNPs in the dominant and recessive modes showed no significant differences (data not shown). The genotypic distributions of rs715572G/A in TIMP3 significantly differed between knee OA cases and healthy controls (p<0.05). When the TIMP3 rs715572GG genotype was used as the reference group, the TIMP3 rs715572GA heterozygotes appeared to have a lower risk for knee OA (adjusted OR=0.64, 95% CI=0.46 to 0.88; p=0.007). After the correction for multiple comparisons, the TIMP3 rs715572GA genotype still appeared to have a lower risk for knee OA.

Table 2

Analyses of the association of six SNPs with knee OA

SNP	Case	Control	Crude OR (95% CI)	Adjusted OR (95% CI)*	p Value
TGF-β1 rs1800469C/T
T/T	166	167	1	1
T/C	238	212	1.03 (0.96 to 1.11)	1.19 (0.88 to 1.61)	0.254
C/C	114	89	1.07 (0.98 to 1.16)	1.27 (0.87 to 1.84)	0.214
C allele	0.45	0.42	1.16 (0.96 to 1.37)	1.14 (0.95 to 1.38)	0.167
TGF-βRI rs1590A/G
A/A	199	167	1	1
C/A	227	208	0.98 (0.91 to 1.05)	0.95 (0.70 to 1.28)	0.717
C/C	92	93	0.96 (0.87 to 1.04)	0.91 (0.63 to 1.33)	0.636
C allele	0.40	0.42	0.91 (0.76 to 1.08)	0.95 (0.78 to 1.51)	0.602
Smad3 rs12901499A/G
A/A	142	116	1	1
G/A	274	228	0.97 (0.90 to 1.05)	0.86 (0.62 to 1.19)	0.360
G/G	129	124	0.96 (0.88 to 1.05)	0.87 (0.60 to 1.26)	0.447
G allele	0.49	0.51	0.92 (0.77 to 1.8)	0.93 (0.77 to 1.12)	0.434
Smad3 rs6494629T/C
T/T	241	184	1	1
C/T	215	214	0.94 (0.88 to 1.00)	0.79 (0.59 to 1.05)	0.107
C/C	62	70	0.91 (0.82 to 1.00)	0.79 (0.52 to 1.20)	0.262
C allele	0.33	0.38	0.80 (0.66 to 0.96)	0.86 (0.70 to 1.04)	0.119
TIMP3 rs715572G/A
G/G	157	100	1	1
G/A	236	242	0.89 (0.83 to 0.96)	0.64 (0.46 to 0.88)	0.007†
A/A	125	126	0.89 (0.82 to 0.97)	0.71 (0.49 to 1.03)	0.071
A allele	0.47	0.53	0.79 (0.66 to 0.94)	0.84 (0.69 to 1.01)	0.065
TIMP3 rs1962223G/C
G/G	173	155	1	1
C/G	259	240	0.99 (0.93 to 1.06)	0.96 (0.71 to 1.58)	0.792
C/C	86	73	1.01 (0.92 to 1.11)	1.05 (0.71 to 1.29)	0.790
C allele	0.42	0.41	1.02 (0.85 to 1.22)	1.01 (0.84 to 1.23)	0.840

*Adjusted for age, gender and BMI.

†p Values were based on Bonferroni's method.

BMI, body mass index; OA, osteoarthritis; Smad3, Smad family member 3; SNPs, single nucleotide polymorphisms; TGF, transforming growth factor; TIMP3, tissue inhibitor of metalloproteinases 3.

Stratification analysis according to disease severity

We conducted an analysis of associations between the TIMP3 rs715572G/A genotypes and knee OA risk after stratifying the patients using the K–L grading scale. The results revealed significant differences between patients with grade 4 knee OA after the correction for multiple comparisons (GA/GG, adjusted OR=0.53, 95% CI=0.35 to 0.80), and controls (table 3).

Table 3

Stratified analysis of associations between TIMP3 rs715572G/A genotypes and knee OA risk using the K–L grading scale

Genotype	Model	K–L grading scale*
		K–L	OR (95% CI)	Adjusted OR (95% CI)†
GG	AA/GG	2	0.78 (0.49 to 1.24)	0.86 (0.53 to 1.41)
GA	AA/GG	3	0.56 (0.31 to 1.03)	0.60 (0.32 to 1.12)
AA	AA/GG	4	0.56 (0.36 to 0.87)	0.57 (0.35 to 0.92)
	GA/GG	2	0.72 (0.48 to 1.09)	0.68 (0.44 to 1.05)
	GA/GG	3	0.68 (0.41 to 1.13)	0.60 (0.36 to 1.02)
	GA/GG	4	0.53 (0.36 to 0.77)‡	0.53 (0.35 to 0.80)‡

*Grade 0, 1 as a reference category.

†Adjusted for age, gender and BMI.

‡p Values were based on Bonferroni's method.

BMI, body mass index; K–L, Kellgren–Lawrence; OA, osteoarthritis; TIMP3, tissue inhibitor of metalloproteinases 3.

Evaluation of gene–gene interactions: MDR

Table 4 summarises the results of exhaustive MDR analysis evaluating all possible combinations of the studied polymorphisms. According to MDR analysis, the best MDR model included Smad3 rs6494629T/C and TIMP3 rs715572G/A polymorphisms. This model had a maximum testing accuracy of 0.5755 and a maximum cross-validation consistency of 10/10. The model was significant at the 0.010 level, which indicates that a model as good or better was observed only once in 1000 permutations; thus, this was unlikely under the null hypothesis of no association. The significance was also confirmed by a logistic regression model (p for interaction=0.021 for the interaction term, data not shown). Figure 1 depicts the interaction map of all genes, based on entropy measures between individual variables. A strong interaction effect was observed for Smad3 rs6494629T/C and TIMP3 rs715572G/A, which had information gain values of 0.60%.

Table 4

Results of MDR analysis

Locus number	Model	Training Bal ACC	Testing Bal ACC	Cross-validation consistency	p Value*
1	rs715572G/A	0.5447	0.5360	9/10	0.3900
2	rs6494629T/C, rs715572G/A	0.5780	0.5755	10/10	0.0100
3	rs6494629T/C, rs715572 G/A, rs1962223G/C	0.5980	0.5403	5/10	0.3090
4	rs1800469C/T, rs6494629T/C, rs715572G/A, rs1962223G/C	0.6331	0.5105	8/10	0.8360
5	rs1800469C/T, rs1590A/C, rs6494629T/C, rs715572G/A, rs1962223G/C	0.6939	0.5092	5/10	0.9110
6	rs1800469C/T, rs1590A/C, rs6494629T/C, rs12901499A/G, rs715572G/A, rs1962223G/C	0.7723	0.5080	10/10	0.8620

*p Values were based on 1000 permutations.

MDR, multifactor dimensionality reduction; Testing Bal ACC, testing-balanced accuracy.

Figure 1

Interaction map for osteoarthritis risk. Values inside nodes indicate information gain (IG) of individual attributes or main effects, whereas values between nodes show IG of pairwise combinations of attributes or interaction effects. Positive entropy (plotted in red or orange) indicates interaction, while negative entropy (plotted in green) indicates redundancy. Smad3, Smad family member 3; TGF, transforming growth factor; TIMP3, tissue inhibitor of metalloproteinases 3.

Discussion

We investigated TGF-β1, TGF-βRI, Smad3 and TIMP3 polymorphisms in knee patients with OA and identified a significant association between knee OA and TIMP3 rs715572G/A. We also presented statistical evidence of significant interaction between Smad3 rs6494629T/C and TIMP3 rs715572G/A affecting knee OA risk. This interaction was also echoed in the logistic regression approach.

TGF-β is an important anabolic and anticatabolic factor in the maintenance of articular cartilage. Previous gene-association studies have reported that TGF-β is independently associated with knee OA25 and spinal osteophytosis.26 In knee OA, TGFβ1 rs2278422 and rs8179181 were found to have a possible role in susceptibility to knee OA in a British Caucasian population,25 whereas a variation on Leu10Pro or SNP rs1982073 was implicated with spinal osteophytosis in Japanese women.26 An interaction between TGFβ1 rs1800469C/T polymorphism and obesity with risk of hip OA has been identified.27 Although no association with this SNP was seen in the knee OA group, there was an association between obesity and risk of OA with the rs2278422 TGFβ1 polymorphism.27 In our current study, the frequencies of TGF-β1 rs1800469C/T and TGF-βRI rs1590A/G genotypes and alleles, did not differ between knee patients with OA and control groups, consistent with previous findings. In our study population, the patient group was, on average, mildly overweight (BMI=25.81±3.33) rather than obese. It is possible that TGF-β1 and TGF-βRI polymorphisms may only have an effect on development of OA in specific joints, and then only when appropriate additional conditions such as obesity are present.

Smad3 is an intracellular molecule that links the extracellular TGF-β signal with changes in gene transcription. A number of studies have extensively investigated the role of Smad3 protein in OA. A reduction in Smad3 activity results in OA phenotype in some model systems.28–30 Smad3 variants have been recently reported as being associated with OA in European populations, supporting results from animal studies suggesting an important role for this molecule in OA pathogenesis.31 Valdes et al31 showed that four SNPs (rs266335G/A, rs12901499A/G, rs6494629T/C and rs2289263A/C) in Smad3 were found to be significantly associated with knee OA, but only one of them, rs12901499A/G, was associated also with hip OA. A recent study on the role of SMAD3 in graft-versus-host disease suggested that inter-individual differences in SMAD3 expression levels could not be attributed to in-cis genetic interactions in a panel of 22 SNPs tested.32 In the northeastern Chinese population, the Smad3 rs12901499A/G appears to be involved in OA pathogenesis.33 However, no associations were found between knee OA and Smad3 polymorphisms (rs12901499A/G and rs6494629T/C) in this study. These inconsistencies or contradictory findings in different studies may be due to factors such as the size of the sample set and ethnic factors. Small sample size is a common factor leading to different findings. Therefore, more association studies with larger numbers of participants are needed to confirm the association between Smad3 SNPs and knee OA.

Polymorphisms of the tissue inhibitor of TIMP3 have been associated with a range of conditions including resistance to high-altitude pulmonary oedema,34 and susceptibility and survival of patients with breast carcinoma,35 and adenocarcinoma of the gastro-oesophageal junction.36 As far as we are aware, no previous study has yet shown an association with this gene and OA. TIMP3 is potentially chondroprotective. It is closely associated with chondrocytes in articular cartilage, and expression by chondrocytes in vitro is increased following exposure to TGFβ1.18 Also, TIMP-3 deficiency in mice results in cartilage degradation similar to changes seen in patients with OA, indicating TIMP-3 may play a pathophysiologic role in the development of OA.37 In the current study, TIMP3 rs715572G/A was associated with more severe knee OA. The genetics of OA are complex and considered to involve interactions between multiple genetic variants.

The magnitude of the effect of any single polymorphism is likely to be missed if genes are individually examined without considering potential interactions with other genes, especially those in related pathways. As such, our findings suggesting that Smad3 rs6494629T/C and TIMP3 rs715572G/A may cooperate in the determination of individual knee OA susceptibility profiles, is relevant. The evaluation of gene–gene interactions not only increases the detection power but also helps in understanding the genetics of the biological and biochemical pathways underlying the disease. Additional studies are needed to establish the mechanisms of interactions between Smad3 rs6494629T/C and TIMP3 rs715572G/A polymorphisms, and their effects on knee OA predisposition.

Our results suggest that a TIMP3 polymorphism is associated with severe knee OA in a Chinese Han population. The effect of interactions between Smad3 rs6494629T/C and TIMP3 rs715572G/A polymorphisms may be more important in knee OA. Further studies have to replicate our findings and investigate whether environmental factors act on different SNPs, whereas functional studies have to investigate the exact biological mechanism of these gene–gene interactions.