Quantification of the effects of single nucleotide variants in NKX2.1 transcription factor binding sites

Transcription factors recognise and bind specific DNA sequence patterns in promoters and enhancers thereby regulating gene expression. Variations in the DNA sequence of transcription factor binding sites (TFBSs) can alter gene regulation and may disrupt development. The transcription factor NKX2.1 is a crucial regulator of thyroid, lung, and neural development. Mutations in its coding gene NKX2-1 may cause choreoathetosis and congenital hypothyroidism with or without pulmonary dysfunction (CAHTP, OMIM #610978). Most genetically solved patients carry mutations in the coding regions of NKX2-1 that affect DNA binding, while the majority of patients with CAHTP-like symptoms do not carry mutations in the NKX2-1 coding sequence. We hypothesise that variations in NKX2.1 DNA-binding at promoter and enhancer sites could cause disease as well. We employ EMSA-seq to quantify the effects of genetic variation on NKX2.1 binding strength. We used this data to train neural network models to forecast the influence of DNA variation on NKX2.1 binding and validated our models using microscale thermophoresis, X-ray crystallography, and publicly available ChIP-seq data. The neural networks were able to detect TFBSs in ChIP-seq data and can thus be used to evaluate whole genome sequencing data of CAHTP-patients in order to prioritise potential disease-causing variants in regulatory elements. Overall design: We created mutant libraries that contain single nucleotide variants or combinations thereof in the NKX2.1 binding site. To design the libraries, we took a 24 bp sequence of the NKX2.1 binding site (see below: middle sequence) in the rat thyroglobulin promoter and added constant primer regions to either end (see below: italics) for library preparation. We performed EMSA-seq to determine the effect of these variants on protein binding. In this experiment we incubated the NKX2.1 DNA binding domain with the mutant libraries and separated the bound sequences using electromobility shift assays (EMSAs). The bound DNA was purified from the gel and used in high-throughput sequencing. EMSA-seq was performed with three different mutant libraries [CORE: the 4 "CAAG" core nucleotides are randomised, FLANK: 5 nucleotides on each side of the core motif are randomised, ALL: all 14 flanking and core nucleotides are randomised]. The binding reactions either contained no protein (control for sequence bias) or with protein at different DNA:protein molar ratios (e.g. 1:0.1 is 10 times more DNA than protein). All experiments were performed in three technical replicates. Binding site sequence: 5'-TCGTCGGCAGCGTCAGATGTGTAT CACTGCCCAGTCAAGTGTTCTTGA TCACCGCTCCGACTGCAGAAAA-3'. CORE: 5'-TCGTCGGCAGCGTCAGATGTGTAT CACTGCCCAGTNNNNTGTTCTTGA TCACCGCTCCGACTGCAGAAAA-3'. FLANK: 5'-TCGTCGGCAGCGTCAGATGTGTAT CACTGCNNNNNCAAGNNNNNTTGA TCACCGCTCCGACTGCAGAAAA-3'. ALL: 5'-TCGTCGGCAGCGTCAGATGTGTAT CACTGCNNNNNNNNNNNNNNTTGA TCACCGCTCCGACTGCAGAAAA-3'