comparison kinatestid_r/Kinatest-R.R @ 17:26ef4add9f7b draft

Uploaded
author jfb
date Wed, 28 Feb 2018 14:09:19 -0500
parents 15b5d4ae4480
children e16ca3e9fc49
comparison
equal deleted inserted replaced
16:15b5d4ae4480 17:26ef4add9f7b
1
2 ImportedSubstrateList<- read.csv("input1", stringsAsFactors=FALSE) 1 ImportedSubstrateList<- read.csv("input1", stringsAsFactors=FALSE)
3 NegativeSubstrateList<- read.csv("input2", stringsAsFactors=FALSE) 2 NegativeSubstrateList<- read.csv("input2", stringsAsFactors=FALSE)
4 SubstrateBackgroundFrequency<- read.csv("input3", stringsAsFactors=FALSE) 3 SubstrateBackgroundFrequency<- read.csv("input3", stringsAsFactors=FALSE, header = FALSE)
5 4
6 ScreenerFilename<-"screener" 5 ScreenerFilename<-"screener"
7 6
8 7
9 8
11 FILENAME2<-"output2" 10 FILENAME2<-"output2"
12 FILENAME3<-"output3" 11 FILENAME3<-"output3"
13 12
14 13
15 14
16 15 SubstrateBackgroundFrequency<-t(SubstrateBackgroundFrequency)
17 16 # number<-nrow(SubstrateBackgroundFrequency)-1
18 17 SubstrateBackgroundFrequency<-SubstrateBackgroundFrequency[2:nrow(SubstrateBackgroundFrequency),]
19 18 Sub<-na.omit(SubstrateBackgroundFrequency)
20 19 SubstrateBackgroundFrequency<-Sub
21 20
22 21 args = commandArgs(trailingOnly=TRUE)
23 22 TodaysKinase<-args[1]
24
25
26
27
28
29
30
31 23
32 24
33 25
34 26
35 27
592 SDtable<-data.frame(SetOfAAs,SDtable) 584 SDtable<-data.frame(SetOfAAs,SDtable)
593 585
594 PercentTable<-rbind(HeaderSD,PercentTable) 586 PercentTable<-rbind(HeaderSD,PercentTable)
595 row.names(PercentTable)<-NULL 587 row.names(PercentTable)<-NULL
596 PercentTable<-data.frame(SetOfAAs,PercentTable) 588 PercentTable<-data.frame(SetOfAAs,PercentTable)
597 numberofY<-as.numeric(SubstrateBackgroundFrequency$Number.of.Y) 589 numberofY<-as.numeric(SubstrateBackgroundFrequency[,34])
598 numberofY<-numberofY[!is.na(numberofY)] 590 numberofY<-numberofY[!is.na(numberofY)]
599 591
600 numberofPY<-as.numeric(SubstrateBackgroundFrequency$Number.of.pY) 592 numberofPY<-as.numeric(SubstrateBackgroundFrequency[,35])
601 numberofPY<-numberofPY[!is.na(numberofPY)] 593 numberofPY<-numberofPY[!is.na(numberofPY)]
602 594
603 NormalizationScore<-sum(numberofPY)/sum(numberofY) 595 NormalizationScore<-sum(numberofPY)/sum(numberofY)
604 596
605 # positions<-matrix(data = NA, nrow=20,ncol = 15) 597 # positions<-matrix(data = NA, nrow=20,ncol = 15)
1159 1151
1160 #View(SDtable) 1152 #View(SDtable)
1161 bareSDs<-SDtable[2:21,2:16] 1153 bareSDs<-SDtable[2:21,2:16]
1162 goodones<-bareSDs>2 1154 goodones<-bareSDs>2
1163 1155
1156 # Positionm7<-which(goodones[,1] %in% TRUE)
1157 # if (length(Positionm7)<1){Positionm7<-which(bareSDs[,1]==max(bareSDs[,1]))}
1158 # Positionm6<-which(goodones[,2] %in% TRUE)
1159 # if (length(Positionm6)<1){Positionm6<-which(bareSDs[,2]==max(bareSDs[,2]))}
1160 # Positionm5<-which(goodones[,3] %in% TRUE)
1161 # if (length(Positionm5)<1){Positionm5<-which(bareSDs[,3]==max(bareSDs[,3]))}
1162 # Positionm4<-which(goodones[,4] %in% TRUE)
1163 # if (length(Positionm4)<2){Positionm4<-bareSDs[,4][order(bareSDs[,4])[1:2]]}
1164 # Positionm3<-which(goodones[,5] %in% TRUE)
1165 # if (length(Positionm3)<2){Positionm3<-bareSDs[,5][order(bareSDs[,5])[1:2]]}
1166 # Positionm2<-which(goodones[,6] %in% TRUE)
1167 # if (length(Positionm2)<2){Positionm2<-bareSDs[,6][order(bareSDs[,6])[1:2]]}
1168 # Positionm1<-which(goodones[,7] %in% TRUE)
1169 # if (length(Positionm1)<2){Positionm1<-bareSDs[,7][order(bareSDs[,7])[1:2]]}
1170 #
1171 # Positiond0<-which(goodones[,8] %in% TRUE)
1172 # if (length(Positiond0)<1){Positiond0<-which(bareSDs[,8]==max(bareSDs[,8]))}
1173 #
1174 # Positionp1<-which(goodones[,9] %in% TRUE)
1175 # if (length(Positionp1)<2){Positionp1<-bareSDs[,9][order(bareSDs[,9])[1:2]]}
1176 # Positionp2<-which(goodones[,10] %in% TRUE)
1177 # if (length(Positionp2)<2){Positionp2<-bareSDs[,10][order(bareSDs[,10])[1:2]]}
1178 # Positionp3<-which(goodones[,11] %in% TRUE)
1179 # if (length(Positionp3)<2){Positionp3<-bareSDs[,11][order(bareSDs[,11])[1:2]]}
1180 # Positionp4<-which(goodones[,12] %in% TRUE)
1181 # if (length(Positionp4)<2){Positionp4<-bareSDs[,12][order(bareSDs[,12])[1:2]]}
1182 # Positionp5<-which(goodones[,13] %in% TRUE)
1183 # if (length(Positionp5)<1){Positionp5<-which(bareSDs[,13]==max(bareSDs[,13]))}
1184 # Positionp6<-which(goodones[,14] %in% TRUE)
1185 # if (length(Positionp6)<1){Positionp6<-which(bareSDs[,14]==max(bareSDs[,14]))}
1186 # Positionp7<-which(goodones[,15] %in% TRUE)
1187 # if (length(Positionp7)<1){Positionp7<-which(bareSDs[,15]==max(bareSDs[,15]))}
1188
1189
1190
1191
1192 # Positionm7<-which(goodones[,1] %in% TRUE)
1193 # if (length(Positionm7)<1){Positionm7<-which(bareSDs[,1]==max(bareSDs[,1]))}
1194 # Positionm6<-which(goodones[,2] %in% TRUE)
1195 # if (length(Positionm6)<1){Positionm6<-which(bareSDs[,2]==max(bareSDs[,2]))}
1196 # Positionm5<-which(goodones[,3] %in% TRUE)
1197 # if (length(Positionm5)<1){Positionm5<-which(bareSDs[,3]==max(bareSDs[,3]))}
1198 # Positionm4<-which(goodones[,4] %in% TRUE)
1199 # if (length(Positionm4)<1){Positionm4<-which(bareSDs[,4]==max(bareSDs[,4]))}
1200 # Positionm3<-which(goodones[,5] %in% TRUE)
1201 # if (length(Positionm3)<1){Positionm3<-which(bareSDs[,5]==max(bareSDs[,5]))}
1202 # Positionm2<-which(goodones[,6] %in% TRUE)
1203 # if (length(Positionm2)<1){Positionm2<-which(bareSDs[,6]==max(bareSDs[,6]))}
1204 # Positionm1<-which(goodones[,7] %in% TRUE)
1205 # if (length(Positionm1)<1){Positionm1<-which(bareSDs[,7]==max(bareSDs[,7]))}
1206 #
1207 # Positiond0<-which(goodones[,8] %in% TRUE)
1208 # if (length(Positiond0)<1){Positiond0<-which(bareSDs[,8]==max(bareSDs[,8]))}
1209 #
1210 # Positionp1<-which(goodones[,9] %in% TRUE)
1211 # if (length(Positionp1)<1){Positionp1<-which(bareSDs[,9]==max(bareSDs[,9]))}
1212 # Positionp2<-which(goodones[,10] %in% TRUE)
1213 # if (length(Positionp2)<1){Positionp2<-which(bareSDs[,10]==max(bareSDs[,10]))}
1214 # Positionp3<-which(goodones[,11] %in% TRUE)
1215 # if (length(Positionp3)<1){Positionp3<-which(bareSDs[,11]==max(bareSDs[,11]))}
1216 # Positionp4<-which(goodones[,12] %in% TRUE)
1217 # if (length(Positionp4)<1){Positionp4<-which(bareSDs[,12]==max(bareSDs[,12]))}
1218 # Positionp5<-which(goodones[,13] %in% TRUE)
1219 # if (length(Positionp5)<1){Positionp5<-which(bareSDs[,13]==max(bareSDs[,13]))}
1220 # Positionp6<-which(goodones[,14] %in% TRUE)
1221 # if (length(Positionp6)<1){Positionp6<-which(bareSDs[,14]==max(bareSDs[,14]))}
1222 # Positionp7<-which(goodones[,15] %in% TRUE)
1223 # if (length(Positionp7)<1){Positionp7<-which(bareSDs[,15]==max(bareSDs[,15]))}
1224
1225 match(c(bareSDs[,2][order(bareSDs[,2])[1:2]]),bareSDs[,2])
1226
1164 Positionm7<-which(goodones[,1] %in% TRUE) 1227 Positionm7<-which(goodones[,1] %in% TRUE)
1165 if (length(Positionm7)<1){Positionm7<-which(bareSDs[,1]==max(bareSDs[,1]))} 1228 if (length(Positionm7)<3){Positionm7<-match(c(bareSDs[,1][order(bareSDs[,1])[19:20]]),bareSDs[,1])}
1166 Positionm6<-which(goodones[,2] %in% TRUE) 1229 Positionm6<-which(goodones[,2] %in% TRUE)
1167 if (length(Positionm6)<1){Positionm6<-which(bareSDs[,2]==max(bareSDs[,2]))} 1230 if (length(Positionm6)<3){Positionm6<-match(c(bareSDs[,2][order(bareSDs[,2])[19:20]]),bareSDs[,2])}
1168 Positionm5<-which(goodones[,3] %in% TRUE) 1231 Positionm5<-which(goodones[,3] %in% TRUE)
1169 if (length(Positionm5)<1){Positionm5<-which(bareSDs[,3]==max(bareSDs[,3]))} 1232 if (length(Positionm5)<3){Positionm5<-match(c(bareSDs[,3][order(bareSDs[,3])[19:20]]),bareSDs[,3])}
1170 Positionm4<-which(goodones[,4] %in% TRUE) 1233 Positionm4<-which(goodones[,4] %in% TRUE)
1171 if (length(Positionm4)<1){Positionm4<-which(bareSDs[,4]==max(bareSDs[,4]))} 1234 if (length(Positionm4)<3){Positionm4<-match(c(bareSDs[,4][order(bareSDs[,4])[19:20]]),bareSDs[,4])}
1172 Positionm3<-which(goodones[,5] %in% TRUE) 1235 Positionm3<-which(goodones[,5] %in% TRUE)
1173 if (length(Positionm3)<1){Positionm3<-which(bareSDs[,5]==max(bareSDs[,5]))} 1236 if (length(Positionm3)<3){Positionm3<-match(c(bareSDs[,5][order(bareSDs[,5])[19:20]]),bareSDs[,5])}
1174 Positionm2<-which(goodones[,6] %in% TRUE) 1237 Positionm2<-which(goodones[,6] %in% TRUE)
1175 if (length(Positionm2)<1){Positionm2<-which(bareSDs[,6]==max(bareSDs[,6]))} 1238 if (length(Positionm2)<3){Positionm2<-match(c(bareSDs[,6][order(bareSDs[,6])[19:20]]),bareSDs[,6])}
1176 Positionm1<-which(goodones[,7] %in% TRUE) 1239 Positionm1<-which(goodones[,7] %in% TRUE)
1177 if (length(Positionm1)<1){Positionm1<-which(bareSDs[,7]==max(bareSDs[,7]))} 1240 if (length(Positionm1)<3){Positionm1<-match(c(bareSDs[,7][order(bareSDs[,7])[19:20]]),bareSDs[,7])}
1178 1241
1179 Positiond0<-which(goodones[,8] %in% TRUE) 1242 Positiond0<-which(goodones[,8] %in% TRUE)
1180 if (length(Positiond0)<1){Positiond0<-which(bareSDs[,8]==max(bareSDs[,8]))} 1243 #if (length(Positiond0)<3){Positiond0<-bareSDs[,8][order(bareSDs[,8])[1:2]]}
1181 1244
1182 Positionp1<-which(goodones[,9] %in% TRUE) 1245 Positionp1<-which(goodones[,9] %in% TRUE)
1183 if (length(Positionp1)<1){Positionp1<-which(bareSDs[,9]==max(bareSDs[,9]))} 1246 if (length(Positionp1)<3){Positionp1<-match(c(bareSDs[,9][order(bareSDs[,9])[19:20]]),bareSDs[,9])}
1184 Positionp2<-which(goodones[,10] %in% TRUE) 1247 Positionp2<-which(goodones[,10] %in% TRUE)
1185 if (length(Positionp2)<1){Positionp2<-which(bareSDs[,10]==max(bareSDs[,10]))} 1248 if (length(Positionp2)<3){Positionp2<-match(c(bareSDs[,10][order(bareSDs[,10])[19:20]]),bareSDs[,10])}
1186 Positionp3<-which(goodones[,11] %in% TRUE) 1249 Positionp3<-which(goodones[,11] %in% TRUE)
1187 if (length(Positionp3)<1){Positionp3<-which(bareSDs[,11]==max(bareSDs[,11]))} 1250 if (length(Positionp3)<3){Positionp3<-match(c(bareSDs[,11][order(bareSDs[,11])[19:20]]),bareSDs[,11])}
1188 Positionp4<-which(goodones[,12] %in% TRUE) 1251 Positionp4<-which(goodones[,12] %in% TRUE)
1189 if (length(Positionp4)<1){Positionp4<-which(bareSDs[,12]==max(bareSDs[,12]))} 1252 if (length(Positionp4)<3){Positionp4<-match(c(bareSDs[,12][order(bareSDs[,12])[19:20]]),bareSDs[,12])}
1190 Positionp5<-which(goodones[,13] %in% TRUE) 1253 Positionp5<-which(goodones[,13] %in% TRUE)
1191 if (length(Positionp5)<1){Positionp5<-which(bareSDs[,13]==max(bareSDs[,13]))} 1254 if (length(Positionp5)<3){Positionp5<-match(c(bareSDs[,13][order(bareSDs[,13])[19:20]]),bareSDs[,13])}
1192 Positionp6<-which(goodones[,14] %in% TRUE) 1255 Positionp6<-which(goodones[,14] %in% TRUE)
1193 if (length(Positionp6)<1){Positionp6<-which(bareSDs[,14]==max(bareSDs[,14]))} 1256 if (length(Positionp6)<3){Positionp6<-match(c(bareSDs[,14][order(bareSDs[,14])[19:20]]),bareSDs[,14])}
1194 Positionp7<-which(goodones[,15] %in% TRUE) 1257 Positionp7<-which(goodones[,15] %in% TRUE)
1195 if (length(Positionp7)<1){Positionp7<-which(bareSDs[,15]==max(bareSDs[,15]))} 1258 if (length(Positionp7)<3){Positionp7<-match(c(bareSDs[,15][order(bareSDs[,15])[19:20]]),bareSDs[,15])}
1259
1196 1260
1197 aa_props2 <- c("1"="A", "2"="C", "3"="D", "4"="E", "5"="F", "6"="G", "7"="H", "8"="I", "9"="K", "10"="L", "11"="M", "12"="N", 1261 aa_props2 <- c("1"="A", "2"="C", "3"="D", "4"="E", "5"="F", "6"="G", "7"="H", "8"="I", "9"="K", "10"="L", "11"="M", "12"="N",
1198 "13"="P", "14"="Q", "15"="R", "16"="S", "17"="T", "18"="V", "19"="W", "20"="Y") 1262 "13"="P", "14"="Q", "15"="R", "16"="S", "17"="T", "18"="V", "19"="W", "20"="Y")
1199 1263
1200 Positionm7<-sapply(Positionm7, function (x) aa_props2[x]) 1264 Positionm7<-sapply(Positionm7, function (x) aa_props2[x])
1712 1776
1713 AblNorm<-1/as.numeric(Abl[22,1]) 1777 AblNorm<-1/as.numeric(Abl[22,1])
1714 AblThresh<-as.numeric(Abl[24,1]) 1778 AblThresh<-as.numeric(Abl[24,1])
1715 AblTrueThresh<-((AblThresh*AblNorm)/(100-AblThresh)) 1779 AblTrueThresh<-((AblThresh*AblNorm)/(100-AblThresh))
1716 AblActive<-unlist(AblGeneratedScores)>AblTrueThresh 1780 AblActive<-unlist(AblGeneratedScores)>AblTrueThresh
1781 if (TodaysKinase=="ABL"){AblActive<-rep(0,times=nrow(GeneratedPeptides))}
1717 1782
1718 ArgNorm<-1/as.numeric(Arg[22,1]) 1783 ArgNorm<-1/as.numeric(Arg[22,1])
1719 ArgThresh<-as.numeric(Arg[24,1]) 1784 ArgThresh<-as.numeric(Arg[24,1])
1720 ArgTrueThresh<-((ArgThresh*ArgNorm)/(100-ArgThresh)) 1785 ArgTrueThresh<-((ArgThresh*ArgNorm)/(100-ArgThresh))
1721 ArgActive<-unlist(ArgGeneratedScores)>ArgTrueThresh 1786 ArgActive<-unlist(ArgGeneratedScores)>ArgTrueThresh
1787 if (TodaysKinase=="ARG"){ArgActive<-rep(0,times=nrow(GeneratedPeptides))}
1722 1788
1723 BtkNorm<-1/as.numeric(Btk[22,1]) 1789 BtkNorm<-1/as.numeric(Btk[22,1])
1724 BtkThresh<-as.numeric(Btk[24,1]) 1790 BtkThresh<-as.numeric(Btk[24,1])
1725 BtkTrueThresh<-((BtkThresh*BtkNorm)/(100-BtkThresh)) 1791 BtkTrueThresh<-((BtkThresh*BtkNorm)/(100-BtkThresh))
1726 BtkActive<-unlist(BtkGeneratedScores)>BtkTrueThresh 1792 BtkActive<-unlist(BtkGeneratedScores)>BtkTrueThresh
1793 if (TodaysKinase=="BTK"){BtkActive<-rep(0,times=nrow(GeneratedPeptides))}
1727 1794
1728 CskNorm<-1/as.numeric(Csk[22,1]) 1795 CskNorm<-1/as.numeric(Csk[22,1])
1729 CskThresh<-as.numeric(Csk[24,1]) 1796 CskThresh<-as.numeric(Csk[24,1])
1730 CskTrueThresh<-((CskThresh*CskNorm)/(100-CskThresh)) 1797 CskTrueThresh<-((CskThresh*CskNorm)/(100-CskThresh))
1731 CskActive<-(CskGeneratedScores)>CskTrueThresh 1798 CskActive<-(CskGeneratedScores)>CskTrueThresh
1799 if (TodaysKinase=="CSK"){CskActive<-rep(0,times=nrow(GeneratedPeptides))}
1732 1800
1733 FynNorm<-1/as.numeric(Fyn[22,1]) 1801 FynNorm<-1/as.numeric(Fyn[22,1])
1734 FynThresh<-as.numeric(Fyn[24,1]) 1802 FynThresh<-as.numeric(Fyn[24,1])
1735 FynTrueThresh<-((FynThresh*FynNorm)/(100-FynThresh)) 1803 FynTrueThresh<-((FynThresh*FynNorm)/(100-FynThresh))
1736 FynActive<-unlist(FynGeneratedScores)>FynTrueThresh 1804 FynActive<-unlist(FynGeneratedScores)>FynTrueThresh
1805 if (TodaysKinase=="FYN"){FynActive<-rep(0,times=nrow(GeneratedPeptides))}
1737 1806
1738 HckNorm<-1/as.numeric(Hck[22,1]) 1807 HckNorm<-1/as.numeric(Hck[22,1])
1739 HckThresh<-as.numeric(Hck[24,1]) 1808 HckThresh<-as.numeric(Hck[24,1])
1740 HckTrueThresh<-((HckThresh*HckNorm)/(100-HckThresh)) 1809 HckTrueThresh<-((HckThresh*HckNorm)/(100-HckThresh))
1741 HckActive<-unlist(HckGeneratedScores)>HckTrueThresh 1810 HckActive<-unlist(HckGeneratedScores)>HckTrueThresh
1811 if (TodaysKinase=="HCK"){HckActive<-rep(0,times=nrow(GeneratedPeptides))}
1742 1812
1743 JAK2Norm<-1/as.numeric(JAK2[22,1]) 1813 JAK2Norm<-1/as.numeric(JAK2[22,1])
1744 JAK2Thresh<-as.numeric(JAK2[24,1]) 1814 JAK2Thresh<-as.numeric(JAK2[24,1])
1745 JAK2TrueThresh<-((JAK2Thresh*JAK2Norm)/(100-JAK2Thresh)) 1815 JAK2TrueThresh<-((JAK2Thresh*JAK2Norm)/(100-JAK2Thresh))
1746 JAk2Active<-unlist(JAK2GeneratedScores)>JAK2TrueThresh 1816 JAk2Active<-unlist(JAK2GeneratedScores)>JAK2TrueThresh
1817 if (TodaysKinase=="JAK2"){JAk2Active<-rep(0,times=nrow(GeneratedPeptides))}
1747 1818
1748 LckNorm<-1/as.numeric(Lck[22,1]) 1819 LckNorm<-1/as.numeric(Lck[22,1])
1749 LckThresh<-as.numeric(Lck[24,1]) 1820 LckThresh<-as.numeric(Lck[24,1])
1750 LckTrueThresh<-((LckThresh*LckNorm)/(100-LckThresh)) 1821 LckTrueThresh<-((LckThresh*LckNorm)/(100-LckThresh))
1751 LckActive<-unlist(LckGeneratedScores)>LckTrueThresh 1822 LckActive<-unlist(LckGeneratedScores)>LckTrueThresh
1823 if (TodaysKinase=="LCK"){LckActive<-rep(0,times=nrow(GeneratedPeptides))}
1752 1824
1753 LynNorm<-1/as.numeric(Lyn[22,1]) 1825 LynNorm<-1/as.numeric(Lyn[22,1])
1754 LynThresh<-as.numeric(Lyn[24,1]) 1826 LynThresh<-as.numeric(Lyn[24,1])
1755 LynTrueThresh<-((LynThresh*LynNorm)/(100-LynThresh)) 1827 LynTrueThresh<-((LynThresh*LynNorm)/(100-LynThresh))
1756 LynActive<-unlist(LynGeneratedScores)>LynTrueThresh 1828 LynActive<-unlist(LynGeneratedScores)>LynTrueThresh
1829 if (TodaysKinase=="LYN"){LynActive<-rep(0,times=nrow(GeneratedPeptides))}
1757 1830
1758 Pyk2Norm<-1/as.numeric(Pyk2[22,1]) 1831 Pyk2Norm<-1/as.numeric(Pyk2[22,1])
1759 Pyk2Thresh<-as.numeric(Pyk2[24,1]) 1832 Pyk2Thresh<-as.numeric(Pyk2[24,1])
1760 Pyk2TrueThresh<-((Pyk2Thresh*Pyk2Norm)/(100-Pyk2Thresh)) 1833 Pyk2TrueThresh<-((Pyk2Thresh*Pyk2Norm)/(100-Pyk2Thresh))
1761 Pyk2Active<-unlist(Pyk2GeneratedScores)>Pyk2TrueThresh 1834 Pyk2Active<-unlist(Pyk2GeneratedScores)>Pyk2TrueThresh
1835 if (TodaysKinase=="PYK2"){Pyk2Active<-rep(0,times=nrow(GeneratedPeptides))}
1762 1836
1763 SrcNorm<-1/as.numeric(Src[22,1]) 1837 SrcNorm<-1/as.numeric(Src[22,1])
1764 SrcThresh<-as.numeric(Src[24,1]) 1838 SrcThresh<-as.numeric(Src[24,1])
1765 SrcTrueThresh<-((SrcThresh*SrcNorm)/(100-SrcThresh)) 1839 SrcTrueThresh<-((SrcThresh*SrcNorm)/(100-SrcThresh))
1766 SrcActive<-unlist(SrcGeneratedScores)>SrcTrueThresh 1840 SrcActive<-unlist(SrcGeneratedScores)>SrcTrueThresh
1841 if (TodaysKinase=="SRC"){SrcActive<-rep(0,times=nrow(GeneratedPeptides))}
1767 1842
1768 SykNorm<-1/as.numeric(Syk[22,1]) 1843 SykNorm<-1/as.numeric(Syk[22,1])
1769 SykThresh<-as.numeric(Syk[24,1]) 1844 SykThresh<-as.numeric(Syk[24,1])
1770 SykTrueThresh<-((SykThresh*SykNorm)/(100-SykThresh)) 1845 SykTrueThresh<-((SykThresh*SykNorm)/(100-SykThresh))
1771 SykActive<-unlist(SykGeneratedScores)>SykTrueThresh 1846 SykActive<-unlist(SykGeneratedScores)>SykTrueThresh
1847 if (TodaysKinase=="SYK"){SykActive<-rep(0,times=nrow(GeneratedPeptides))}
1772 1848
1773 YesNorm<-1/as.numeric(Yes[22,1]) 1849 YesNorm<-1/as.numeric(Yes[22,1])
1774 YesThresh<-as.numeric(Yes[24,1]) 1850 YesThresh<-as.numeric(Yes[24,1])
1775 YesTrueThresh<-((YesThresh*YesNorm)/(100-YesThresh)) 1851 YesTrueThresh<-((YesThresh*YesNorm)/(100-YesThresh))
1776 YesActive<-unlist(YesGeneratedScores)>YesTrueThresh 1852 YesActive<-unlist(YesGeneratedScores)>YesTrueThresh
1853 if (TodaysKinase=="YES"){YesActive<-rep(0,times=nrow(GeneratedPeptides))}
1777 1854
1778 AllActive<-AblActive+ArgActive+BtkActive+CskActive+FynActive+HckActive+JAk2Active+LckActive+LynActive+Pyk2Active+SrcActive+SykActive+YesActive 1855 AllActive<-AblActive+ArgActive+BtkActive+CskActive+FynActive+HckActive+JAk2Active+LckActive+LynActive+Pyk2Active+SrcActive+SykActive+YesActive
1779 #Btkactive+ 1856 #Btkactive+
1780 1857
1781 Scores<-ThisKinGeneratedScores 1858 Scores<-ThisKinGeneratedScores
1882 # TNninetyone<-length(NegativeWeirdScores[NegativeWeirdScores<91]) 1959 # TNninetyone<-length(NegativeWeirdScores[NegativeWeirdScores<91])
1883 # Specninetyone<-TNninetyone/100 1960 # Specninetyone<-TNninetyone/100
1884 1961
1885 #create the MCC table 1962 #create the MCC table
1886 1963
1887 threshold<-c(1:100) 1964 threshold<-c(1:100,(1:9)/10,(1:9)/100,0,-.1)
1888 threshold<-order(threshold,decreasing = TRUE) 1965 threshold<-threshold[order(threshold,decreasing = TRUE)]
1889 1966 threshold
1890 Truepositives<-c(1:100) 1967
1891 Falsenegatives<-c(1:100) 1968 Truepositives<-c(1:120)
1892 Sensitivity<-c(1:100) 1969 Falsenegatives<-c(1:120)
1893 TrueNegatives<-c(1:100) 1970 Sensitivity<-c(1:120)
1894 FalsePositives<-c(1:100) 1971 TrueNegatives<-c(1:120)
1895 Specificity<-c(1:100) 1972 FalsePositives<-c(1:120)
1896 Accuracy<-c(1:100) 1973 Specificity<-c(1:120)
1897 MCC<-c(1:100) 1974 Accuracy<-c(1:120)
1898 EER<-c(1:100) 1975 MCC<-c(1:120)
1976 EER<-c(1:120)
1899 1977
1900 #MAKE DAMN SURE THAT THE ACCESSION NUMBERS FOLLOW THE MOTIFS 1978 #MAKE DAMN SURE THAT THE ACCESSION NUMBERS FOLLOW THE MOTIFS
1901 1979
1902 for (z in 1:100) { 1980 for (z in 1:120) {
1903 thres<-101-z 1981 thres<-threshold[z]
1904 Truepositives[z]<-length(PositiveWeirdScores[PositiveWeirdScores>=(thres)]) 1982 Truepositives[z]<-length(PositiveWeirdScores[PositiveWeirdScores>=(thres)])
1905 Falsenegatives[z]<-nrow(positivesubstrates)-Truepositives[z] 1983 Falsenegatives[z]<-nrow(positivesubstrates)-Truepositives[z]
1906 Sensitivity[z]<-Truepositives[z]/(Falsenegatives[z]+Truepositives[z]) 1984 Sensitivity[z]<-Truepositives[z]/(Falsenegatives[z]+Truepositives[z])
1907 TrueNegatives[z]<-length(NegativeWeirdScores[NegativeWeirdScores<(thres)]) 1985 TrueNegatives[z]<-length(NegativeWeirdScores[NegativeWeirdScores<(thres)])
1908 # at thresh 100 this should be 0, because it is total minus true negatives 1986 # at thresh 100 this should be 0, because it is total minus true negatives
1909 FalsePositives[z]<-nrow(NegativeSubstrateList)-TrueNegatives[z] 1987 FalsePositives[z]<-nrow(NegativeSubstrateList)-TrueNegatives[z]
1910 Specificity[z]<-1-(TrueNegatives[z]/(FalsePositives[z]+TrueNegatives[z])) 1988 Specificity[z]<-1-(TrueNegatives[z]/(FalsePositives[z]+TrueNegatives[z]))
1911 Accuracy[z]<-100*(Truepositives[z]+TrueNegatives[z])/(Falsenegatives[z]+FalsePositives[z]+TrueNegatives[z]+Truepositives[z]) 1989 Accuracy[z]<-100*(Truepositives[z]+TrueNegatives[z])/(Falsenegatives[z]+FalsePositives[z]+TrueNegatives[z]+Truepositives[z])
1912 MCC[z]<-((Truepositives[z]+TrueNegatives[z])-(Falsenegatives[z]+FalsePositives[z]))/sqrt(round(round(Truepositives[z]+Falsenegatives[z])*round(TrueNegatives[z]+FalsePositives[z])*round(Truepositives[z]+FalsePositives[z])*round(TrueNegatives[z]+Falsenegatives[z]))) 1990 MCC[z]<-((Truepositives[z]*TrueNegatives[z])-(Falsenegatives[z]*FalsePositives[z]))/sqrt(round(round(Truepositives[z]+Falsenegatives[z])*round(TrueNegatives[z]+FalsePositives[z])*round(Truepositives[z]+FalsePositives[z])*round(TrueNegatives[z]+Falsenegatives[z])))
1913 EER[z]<-.01*(((1-(Sensitivity[z]))*(Truepositives[z]+Falsenegatives[z]))+(Specificity[z]*(1-(Truepositives[z]+Falsenegatives[z])))) 1991 EER[z]<-.01*(((1-(Sensitivity[z]))*(Truepositives[z]+Falsenegatives[z]))+(Specificity[z]*(1-(Truepositives[z]+Falsenegatives[z]))))
1914 } 1992 }
1915 Characterization<-cbind.data.frame(threshold,Truepositives,Falsenegatives,Sensitivity,TrueNegatives,FalsePositives,Specificity,Accuracy,MCC,EER) 1993 Characterization<-cbind.data.frame(threshold,Truepositives,Falsenegatives,Sensitivity,TrueNegatives,FalsePositives,Specificity,MCC,EER)
1916 1994
1917 positiveheader<-c(1,2,3,4,5,6,7,8,9,10,11,12,13,14,15,"RPMS","PMS") 1995 positiveheader<-c(1,2,3,4,5,6,7,8,9,10,11,12,13,14,15,"RPMS","PMS")
1918 positivewithscores<-rbind.data.frame(positiveheader,positivewithscores) 1996 positivewithscores<-rbind.data.frame(positiveheader,positivewithscores)
1919 1997
1920 negativeheader<-c("Substrate","RPMS","PMS") 1998 negativeheader<-c("Substrate","RPMS","PMS")
1929 Characterization<-rbind.data.frame(header,Characterization) 2007 Characterization<-rbind.data.frame(header,Characterization)
1930 write.table(Characterization,file = FILENAME2, col.names = FALSE,row.names = FALSE, append = TRUE,sep = ",") 2008 write.table(Characterization,file = FILENAME2, col.names = FALSE,row.names = FALSE, append = TRUE,sep = ",")
1931 2009
1932 # header<-colnames(RanksPeptides) 2010 # header<-colnames(RanksPeptides)
1933 # RanksPeptides<-rbind.data.frame(header,RanksPeptides) 2011 # RanksPeptides<-rbind.data.frame(header,RanksPeptides)
2012 write.table(x="Off Target Kinase activity (your kinase of interest should have zeros here because it is ON-target)",file = FILENAME3,append = FALSE,row.names = FALSE,col.names = TRUE,sep = ",")
1934 write.table(RanksPeptides,file = FILENAME3,append = FALSE,row.names = FALSE,col.names = TRUE,sep = ",") 2013 write.table(RanksPeptides,file = FILENAME3,append = FALSE,row.names = FALSE,col.names = TRUE,sep = ",")