Collection of nsSNP data and PTM site data

We collected total 517,466 human nsSNPs from the NCBI dbSNP database Build 135 [15]. From 517,466 nsSNPs, 509,183 SNPs were missense SNPs and 14,878 SNPs were stop-gain (nonsense) SNPs. We linked the location of nsSNPs on protein sequences from Ref-Seq Build 37.3. The flanking peptide of nsSNP sites from -7 to +7 retrieved by using the SNP-MapLinkProtein table in dbSNP and defined as the SNP sequence.

Human PTM site data were retrieved from dbPTM3 [16], which includes 194,886 experimentally verified PTMs from MS/MS analysis based research articles and 10 external PTM-related resources. The dbPTM also provides 404,501 computationally predicted PTMs by the KinasePhos-like method based on hidden Markov models (HMMs). It includes 20 types of PTM sites predicted with the threshold of specificity 100% for minimizing false positive. Peptides of flanking PTM sites from -7 to +7, called the PTM site sequence, were retrieved by matching the PTM site with the protein sequence from UniProt [17].

Detection of PTM-SNP

PTM-SNPs were detected by aligning the SNP sequence with the PTM site sequence because not all sequences in UniProt and RefSeq are perfectly identical (end-to-end identical). Each residue in the SNP sequence is matched with the PTM site sequence by sliding one-by-one (Figure 1). The PTM-SNP is defined as the nsSNP's position of the best-aligned SNP sequence with the PTM site sequence. As Figure 1 illustrates, an SNP sequence of nsSNP, rs1803573, is aligned well with the +3 residue of the PTM site sequence. In few cases, multiple hits with diverse locations were detected due to sequence repeats. The European Molecular Biology Open Software Suite (EMBOSS) [18] Needle, a global pairwise sequence alignment tool, was applied to handle these exceptional cases.

We also extracted and integrated PTM-SNPs from the public resource, PhosSNP 1.0. PhosSNP includes phosphorylation-associated nsSNP from computationally translating mRNA sequences. However, about 18% of SNPs in PhosSNP have not yet been defined as nsSNP by dbSNP Build 135, therefor, were removed in the integrating process.

Disease-associated PTM-SNPs

From GWASs using SNP arrays and NGS techniques, huge numbers of statistically significant disease-associated SNPs have been identified recently. National Human Genome Research Institute (NHGRI) GWAS catalog have been provided 8,771 trait-associated SNPs with p-values < 1.0 × 10^-5 (Apr 02, 2013) [19]. The genetic association database (GAD) [20] contains 29,578 disease-associated SNPs and 23,671 disease-associated SNPs are also available from Genotypes and Phenotypes (dbGaP) [21]. To find a correlation between PTM-SNPs and diseases, we collected total 52,731 distinct disease-associated SNPs from NHGRI GWAS catalog, GAD and dbGaP after then matched these SNPs with our PTM-SNPs.

Furthermore, we analyzed type 2 diabetes (T2D) associated SNPs that were identified in our previous study [22] by using Wellcome Trust Case Control Consortium (WTCCC) datasets [23]. T2D-associated SNPs with p-values < 1.0 × 10^-5 were identified from quality controlled (QC) 409,656 SNPs, based on Cochran-Armitage trend test statistics using PLINK 1.07 [24]. QC is applied as a sample missing genotype rate of > 3%, a SNP missing genotype rate of > 1%, Hardy-Weinberg Equilibrium (HWE) p-value ≤ 10^-4, and minor allele frequency (MAF) < 1%. Linkage disequilibrium pruning was leaved out to preserve the PTM-SNP candidates. T2D PTM-SNPs were detected by mapping T2D-associated SNPs and PTM-SNPs. Significantly associated T2D PTM-SNPs were further studied to interpret disease mechanisms with various biological information resources. Known disease genes are mainly collected from OMIM [25], GAD, and KEGG [26]. Drug targets are collected from DrugBank [27], PharmGKB Drug [28], and KEGG Drug.

Data Statistics

Correlation between PTM-SNP and Diseases

Not only missense and stop-gain SNPs, we also compared different functional SNPs in dbSNP. Functional classes that contain more than 100 disease-associated SNPs were considered in comparison. In addition, 414,510 SNPs on miRNA-mRNA binding site were collected from MirSNP [29]. The result shows that SNPs near the gene, such as NearGene-3, MirSNP, NearGene-5 and UTR-3 were higher disease-associated proportion than SNPs in intron because they might affect protein functions by indirect way like regulating protein expressions. On the other hand SNPs around PTM site can give a potent influence on the protein and, as a result, the ratio of disease-associated PTM-SNP was superior to other functional classes. Therefore, we can conclude that PTM-SNPs are highly associated with diseases compared with both other nsSNPs and other functional sites.

Case Study for Type 2 Diabetes associated PTM-SNP

One of T2D-associated PTM-SNPs, rs1801282, is located on PPARG which is a popular T2D gene and a T2D drug target for a thiazolidinedione and changed the Pro12 to the alanine. The SNP is reported T2D from previous study [30], however, the detailed correlation between the SNP and T2D have not been clearly known yet. According to PTM-SNP, rs1801282 may enhance the phosphorylation on Ser8 by CSNK2A1, prohibit the phosphorylation on Ser14 by MTOR or MAPK10, and also may affect the phosphorylation on Ser16. Not only PPARG but CSNK2A1, MAPK10, and MTOR are known as T2D genes and participate on parts of T2D associated pathways such as type II diabetes mellitus, insulin signaling pathway, and Wnt signaling pathway in KEGG.

The associations with T2D, body mass index (BMI), and weight of a Val66Met PTM-SNP (rs6265) located on BDNF are discovered from NHGRI GWAS catalogs and WTCCC T2D dataset analysis. The association between obesity and diabetes is discovered with rs6265 due to GWASs on BMI, weight, and T2D identified G allele as a risk allele and A allele as a non-risk allele. Adjacent site Thr62 of BDNF is predicted as a phosphorylation site and kinases for Thr62 are predicted as ARAF, RAF1, MAP3K2, and MAP3K4 by GPS2.0 and filtrated by STRING. RAF1 is known as a T2D gene and a member of insulin signaling pathway.

PTMs play crucial roles in biological mechanisms, and PTM-SNPs may influence PTM network regulation and have effects on associated diseases. From disease-associated PTM-SNPs, we can find disease-associated PTM networks regulated by PTM-SNPs that may provide insight into the mechanisms of disease. Our case studies show the possibility of using PTM-SNPs, combined with various types of biological information including known disease genes, drug targets, protein-protein interaction networks, and pathway networks, to identify mechanisms underlying disease. Novel disease-candidate genes can be found from disease-associated PTM networks assembled from the relations of PTM enzymes and PTM substrate proteins.