76
|
1 library(reshape2)
|
|
2
|
|
3 args <- commandArgs(trailingOnly = TRUE)
|
|
4
|
|
5 gene.matches = args[1]
|
|
6 sequence.file = args[2]
|
|
7 outputdir = args[3]
|
81
|
8 NToverview.file = paste(outputdir, "ntoverview.txt", sep="/")
|
|
9 NTsum.file = paste(outputdir, "ntsum.txt", sep="/")
|
76
|
10 main.html = "index.html"
|
|
11
|
|
12 setwd(outputdir)
|
|
13
|
|
14 genes = read.table(gene.matches, header=T, sep="\t", fill=T)
|
|
15 sequences = read.table(sequence.file, header=T, sep="\t", fill=T, stringsAsFactors=F, quote="")
|
|
16
|
|
17 dat = merge(sequences, genes, by="Sequence.ID")
|
|
18
|
|
19 dat$seq_conc = paste(dat$CDR1.IMGT, dat$CDR2.IMGT, dat$CDR3.IMGT, dat$FR2.IMGT, dat$FR3.IMGT)
|
|
20
|
78
|
21 IDs = dat[,c("Sequence.ID", "seq_conc", "best_match", "Functionality")]
|
76
|
22 IDs$best_match = as.character(IDs$best_match)
|
|
23
|
|
24 #dat = data.frame(data.table(dat)[, list(freq=.N), by=c("best_match", "seq_conc")])
|
|
25
|
81
|
26 dat = data.frame(table(dat$seq_conc))
|
76
|
27
|
|
28 dat = dat[dat$Freq > 1,]
|
|
29
|
81
|
30 names(dat) = c("seq_conc", "Freq")
|
76
|
31
|
|
32 dat$seq_conc = factor(dat$seq_conc)
|
|
33
|
|
34 dat = dat[order(nchar(as.character(dat$seq_conc))),]
|
|
35
|
|
36 #writing html from R...
|
|
37 td = function(val) { paste("<td>", val, "</td>", sep="") }
|
|
38 tr = function(val) { capture.output(cat("<tr>", td(val), "</tr>", sep="")) }
|
|
39 make.link = function(id, clss, val) { paste("<a href='", clss, "_", id, ".html'>", val, "</a>", sep="") }
|
|
40 tbl = function(df) { res = "<table border='1'>"; for(i in 1:nrow(df)){ res = paste(res, tr(df[i,]), sep=""); }; res = paste(res, "</table>"); }
|
|
41
|
|
42 cat("<table border='1'>", file=main.html, append=F)
|
77
|
43 cat("<caption>CDR1+CDR2+CDR3+FR2+FR3 sequences that show up more than once</caption>", file=main.html, append=T)
|
79
|
44 cat("<tr><th>Sequence</th><th>Functionality</th><th>ca1</th><th>ca2</th><th>cg1</th><th>cg2</th><th>cg3</th><th>cg4</th><th>cm</th></tr>", file=main.html, append=T)
|
76
|
45
|
|
46 for(i in 1:nrow(dat)){
|
|
47 ca1 = IDs[IDs$seq_conc == dat[i,c("seq_conc")] & IDs$best_match == "ca1",]
|
|
48 ca2 = IDs[IDs$seq_conc == dat[i,c("seq_conc")] & IDs$best_match == "ca2",]
|
|
49
|
|
50 cg1 = IDs[IDs$seq_conc == dat[i,c("seq_conc")] & IDs$best_match == "cg1",]
|
|
51 cg2 = IDs[IDs$seq_conc == dat[i,c("seq_conc")] & IDs$best_match == "cg2",]
|
|
52 cg3 = IDs[IDs$seq_conc == dat[i,c("seq_conc")] & IDs$best_match == "cg3",]
|
|
53 cg4 = IDs[IDs$seq_conc == dat[i,c("seq_conc")] & IDs$best_match == "cg4",]
|
|
54
|
|
55 cm = IDs[IDs$seq_conc == dat[i,c("seq_conc")] & IDs$best_match == "cm",]
|
|
56
|
|
57 id = as.numeric(dat[i,"seq_conc"])
|
79
|
58
|
|
59 functionality = dat[i,"Functionality"]
|
76
|
60
|
|
61 if(nrow(ca1) > 0){
|
|
62 cat(tbl(ca1), file=paste("ca1_", id, ".html", sep=""))
|
|
63 }
|
|
64
|
|
65 if(nrow(ca2) > 0){
|
|
66 cat(tbl(ca2), file=paste("ca2_", id, ".html", sep=""))
|
|
67 }
|
|
68
|
|
69 if(nrow(cg1) > 0){
|
|
70 cat(tbl(cg1), file=paste("cg1_", id, ".html", sep=""))
|
|
71 }
|
|
72
|
|
73 if(nrow(cg2) > 0){
|
|
74 cat(tbl(cg2), file=paste("cg2_", id, ".html", sep=""))
|
|
75 }
|
|
76
|
|
77 if(nrow(cg3) > 0){
|
|
78 cat(tbl(cg3), file=paste("cg3_", id, ".html", sep=""))
|
|
79 }
|
|
80
|
|
81 if(nrow(cg4) > 0){
|
|
82 cat(tbl(cg4), file=paste("cg4_", id, ".html", sep=""))
|
|
83 }
|
|
84
|
|
85 if(nrow(cm) > 0){
|
|
86 cat(tbl(cm), file=paste("cm_", id, ".html", sep=""))
|
|
87 }
|
|
88
|
|
89 ca1.html = make.link(id, "ca1", nrow(ca1))
|
|
90 ca2.html = make.link(id, "ca2", nrow(ca2))
|
|
91
|
|
92 cg1.html = make.link(id, "cg1", nrow(cg1))
|
|
93 cg2.html = make.link(id, "cg2", nrow(cg2))
|
|
94 cg3.html = make.link(id, "cg3", nrow(cg3))
|
|
95 cg4.html = make.link(id, "cg4", nrow(cg4))
|
|
96
|
|
97 cm.html = make.link(id, "cm", nrow(cm))
|
|
98
|
79
|
99 rw = c(as.character(dat[i,"seq_conc"]), as.character(functionality), ca1.html, ca2.html, cg1.html, cg2.html, cg3.html, cg4.html, cm.html)
|
76
|
100
|
|
101 cat(tr(rw), file=main.html, append=T)
|
|
102 }
|
|
103
|
|
104 cat("</table>", file=main.html, append=T)
|
81
|
105
|
|
106
|
|
107 #ACGT overview
|
|
108
|
|
109
|
|
110
|
|
111 NToverview = genes[,c("Sequence.ID", "best_match")]
|
|
112 sequences$seq = paste(sequences$CDR2.IMGT, sequences$CDR2.IMGT, sequences$FR2.IMGT, sequences$FR3.IMGT, sep="_")
|
|
113
|
|
114 NToverview = merge(NToverview, sequences[,c("Sequence.ID", "seq")], by="Sequence.ID")
|
|
115
|
|
116 NToverview$A = nchar(gsub("[^Aa]", "", NToverview$seq))
|
|
117 NToverview$C = nchar(gsub("[^Cc]", "", NToverview$seq))
|
|
118 NToverview$G = nchar(gsub("[^Gg]", "", NToverview$seq))
|
|
119 NToverview$T = nchar(gsub("[^Tt]", "", NToverview$seq))
|
|
120
|
|
121 NTsum = data.frame(Sequence.ID="-", best_match="Sum", seq="-", A = sum(NToverview$A), C = sum(NToverview$C), G = sum(NToverview$G), T = sum(NToverview$T))
|
|
122
|
|
123 print(names(NToverview))
|
|
124 print(names(NTsum))
|
|
125
|
|
126 NToverview = rbind(NToverview, NTsum)
|
|
127
|
|
128 write.table(NToverview, NToverview.file, quote=F, sep="\t", row.names=F, col.names=T)
|
|
129 #write.table(NTsum, NTsum.file, quote=F, sep="\t", row.names=F, col.names=T)
|
|
130
|
|
131
|
|
132
|
|
133
|
|
134
|
|
135
|
|
136
|
|
137
|
|
138
|
|
139
|
|
140
|
|
141
|
|
142
|
|
143
|
|
144
|
|
145
|
|
146
|
|
147
|
|
148
|
|
149
|
|
150
|
|
151
|
|
152
|
|
153
|
|
154
|
|
155
|
|
156
|
|
157
|
|
158
|
|
159
|