Genetic and epigenetic interplay fine-tunes expression of TGFB1, encoding TGF-Β1, contributing to osteoarthritis risk

Purpose: The application of genome wide association studies (GWAS) to osteoarthritis (OA) genetics research has facilitated the discovery of almost 100 independent single nucleotide polymorphisms (SNPs) that significantly associate with disease risk. However, it has proven challenging to identify risk genes and pathways marked by the SNPs which result in cartilage degradation. Without this understanding, the translation of GWAS discoveries into novel therapies continues to evade researchers. The overwhelming majority of OA risk SNPs fall in non-coding regions of the genome and are believed to operate by impacting upon the expression of target genes, at a detriment to cartilage health. Genotype at SNPs can correlate with tissue-specific DNA methylation (DNAm) of CG dinucleotides (CpGs), known as methylation quantitative trait loci (mQTLs). This DNAm can act as an intermediary through which SNPs regulate gene expression, however it has previously proven difficult to establish causal relationships between the two effects. One OA risk locus, at chr19q13.2, is marked by intergenic SNP rs75621460 (G>A; minor allele frequency (MAF), 0.03). Statistical fine-mapping assigned >99% probability of rs75621460 being the single causal variant at the locus. The SNP lies in a putative enhancer region downstream of TGFB1, encoding transforming growth factor beta 1 (TGF-β1). In this study we tested the region for regulatory activity, investigated the rs75621460 gene targets, and examined whether local DNAm, along with SNP genotype, could impact gene expression. Methods: An open chromatin region containing either the G- or A-allele at rs75621460 was cloned into a luciferase reporter vector. Constructs were expressed in Tc28a2 chondrocytes and luciferase activity was measured. The SNP region was deleted in the cells using CRISPR-Cas9. RNA was extracted and expression levels of the two genes flanking the SNP, TGFB1 and CCDC97, were measured by RT-qPCR. We analysed transcription factor binding to the SNP using an electrophoretic mobility shift assay (EMSA). Nuclear proteins were extracted from chondrocytes and incubated with 31bp fluorescent DNA probes which included either the G- or A-allele at rs75621460. Supershift assays were performed using antibodies for transcription factors predicted to differentially bind at the region.DNA was extracted from the cartilage, synovium, or whole blood of 319 patients undergoing hip or knee arthroplasty. We genotyped DNA at rs75621460 and quantified DNAm at six flanking CpGs by pyrosequencing. We used a catalytically dead Cas9 (dCas9) protein fused to a DNA methyltransferase (DNMT3a) to selectively methylate the six CpGs in Tc28a2 cells. Results: Luciferase analysis of the region confirmed enhancer activity. The G-allele construct conferred a 2.4-fold increase in luciferase activity (P=0.002), whilst the OA risk A-allele showed lower enhancer activity (1.7-fold increase). EMSA analysis confirmed that protein binding to the G-allele was more abundant than to the A-allele. Furthermore, protein complexes were detected that bound exclusively to either G or A. Supershift analysis identified that the transcription factor SP1 binds to both alleles in cartilage, with a greater affinity for G. CRISPR-Cas9 deletion of the SNP significantly reduced TGFB1 (P=0.01) but not CCDC97 (P=0.12) confirming TGFB1 as the gene target of the enhancer in chondrocytes.We analysed DNAm in patient samples to test for mQTLs at six CpGs flanking the SNP. Due to the low MAF at rs75621460, patient samples were screened to identify sufficient heterozygote DNA for analysis. In cartilage, significant mQTLs were identified at all CpGs (P<8x10 -5). At CpG1 (upstream of the SNP), the A-allele conferred a modest increase in median DNAm (0.8%). However, at downstream sites (CpGs 2-6), OA risk genotype had a much greater impact upon DNAm: CpG2, 14.9% median increase; CpG3, 5.8%; CpG4, 4.1%; CpG5, 5.0%; CpG6, 5.3%. Amongst heterozygous patients (GA), mean DNAm was higher in knee (n=9) samples than in hip (n=7) at all 6 CpGs.We analysed DNA from 61 knee synovium samples to test for mQTLs in a distinct tissue of the articulating joint. Significant mQTLs (P<0.0001) were identified at all CpGs. We investigated whole blood samples to test for a systemic effect but found no significant genotype-DNAm associations (P>0.14). In heterozygous cartilage and synovium samples, we correlated DNAm with TGFB1 expression to look for methylation-expression QTLs (meQTLs). In cartilage, the degree of correlation was dependent upon joint site; stronger meQTLs were measured in knee (r2=0.47 - 0.99) than in hip (r2=0.01 - 0.65). MeQTLs were observed in both cartilage and synovium at the five downstream CpGs. Notably, the impact of increasing DNAm was antagonistic between the two tissues, correlating with increased expression of TGFB1 in cartilage, yet with decreased expression in synovium.Finally, we investigated whether DNAm at the six CpGs could functionally impact upon enhancer activity. Methylation of the region in a gene reporter assay significantly reduced enhancer activity in both the G- (P=0.004) and A-allele constructs (P=0.019). We used a DNMT3a-dCas9 fusion protein for targeted editing of DNAm at the six CpGs. DNMT3a-dCas9 was expressed in Tc28a2 chondrocytes along with one of four guide RNAs (gRNAs) targeting the region. Selective increases in DNAm at CpGs 3-5 resulted in a significant decrease (P<0.01) in expression of TGFB1 (0.77-0.80-fold). Our data demonstrate that DNAm can impact binding of proteins 20-43bp downstream of the SNP, modulating TGFB1 expression. Conclusions: In this study we have confirmed that rs75621460 falls within a TGFB1 enhancer. The OA risk A-allele is associated with decreased enhancer activity, along with an increase in DNAm. We have identified directly opposing meQTLs in cartilage and synovium resulting from a pleiotropic effect at a single functional variant. Moreover, the use of modern epigenome editing tools in this study has allowed us to demonstrate a causal relationship between DNAm and gene expression. We propose a cartilage-specific mechanism through which TGFB1 is regulated in vivo. We hypothesise that a protein complex with strong transcriptional activity binds to the G-allele at rs75621460, the substitution of which alters the consensus sequence for protein binding. We hypothesise that a distinct protein complex that confers relatively weak transcriptional activation, and which also induces methylation of DNA, binds to the A-allele. This leads to a genetic deficit in chondrocyte TGFB1 expression, contributing to a breakdown in cartilage integrity over time.