|
0
|
1 args <- commandArgs(trailingOnly = TRUE)
|
|
|
2
|
|
|
3 inFile = args[1]
|
|
|
4 outDir = args[2]
|
|
3
|
5 logfile = args[3]
|
|
|
6 min_freq = as.numeric(args[4])
|
|
|
7 min_cells = as.numeric(args[5])
|
|
|
8
|
|
|
9 cat("<html><table><tr><td>Starting analysis</td></tr>", file=logfile, append=F)
|
|
0
|
10
|
|
4
|
11 library(ggplot2)
|
|
|
12 library(reshape2)
|
|
|
13 library(data.table)
|
|
|
14 library(grid)
|
|
|
15 library(parallel)
|
|
0
|
16 #require(xtable)
|
|
3
|
17 cat("<tr><td>Reading input</td></tr>", file=logfile, append=T)
|
|
13
|
18 dat = read.table(inFile, header=T, sep="\t", dec=".", fill=T, stringsAsFactors=F)
|
|
9
|
19 dat = dat[!is.na(dat$Patient),]
|
|
13
|
20 dat$Related_to_leukemia_clone = F
|
|
9
|
21
|
|
0
|
22 setwd(outDir)
|
|
3
|
23 cat("<tr><td>Selecting first V/J Genes</td></tr>", file=logfile, append=T)
|
|
2
|
24 dat$V_Segment_Major_Gene = as.factor(as.character(lapply(strsplit(as.character(dat$V_Segment_Major_Gene), "; "), "[[", 1)))
|
|
|
25 dat$J_Segment_Major_Gene = as.factor(as.character(lapply(strsplit(as.character(dat$J_Segment_Major_Gene), "; "), "[[", 1)))
|
|
|
26
|
|
3
|
27 cat("<tr><td>Calculating Frequency</td></tr>", file=logfile, append=T)
|
|
12
|
28
|
|
13
|
29 dat$Frequency = ((10^dat$Log10_Frequency)*100)
|
|
2
|
30
|
|
3
|
31 dat = dat[dat$Frequency >= min_freq,]
|
|
|
32
|
|
13
|
33 triplets = dat[grepl("VanDongen_cALL_14696", dat$Patient) | grepl("(16278)|(26402)|(26759)", dat$Sample),]
|
|
|
34
|
|
|
35 cat("<tr><td>Normalizing to lowest cell count within locus</td></tr>", file=logfile, append=T)
|
|
|
36
|
|
|
37 dat$locus_V = substring(dat$V_Segment_Major_Gene, 0, 4)
|
|
|
38 dat$locus_J = substring(dat$J_Segment_Major_Gene, 0, 4)
|
|
|
39 min_cell_count = data.frame(data.table(dat)[, list(min_cell_count=min(.SD$Cell_Count)), by=c("Patient", "locus_V", "locus_J")])
|
|
|
40
|
|
|
41 dat$min_cell_paste = paste(dat$Patient, dat$locus_V, dat$locus_J)
|
|
|
42 min_cell_count$min_cell_paste = paste(min_cell_count$Patient, min_cell_count$locus_V, min_cell_count$locus_J)
|
|
|
43
|
|
|
44 min_cell_count = min_cell_count[,c("min_cell_paste", "min_cell_count")]
|
|
|
45
|
|
|
46 dat = merge(dat, min_cell_count, by="min_cell_paste")
|
|
|
47
|
|
|
48 dat$normalized_read_count = round(dat$Clone_Molecule_Count_From_Spikes / dat$Cell_Count * dat$min_cell_count / 2, digits=2) #??????????????????????????????????? wel of geen / 2
|
|
|
49
|
|
3
|
50 dat = dat[dat$normalized_read_count >= min_cells,]
|
|
13
|
51
|
|
|
52 dat$paste = paste(dat$Sample, dat$Clone_Sequence)
|
|
9
|
53
|
|
22
|
54 patients = split(dat, dat$Patient, drop=T)
|
|
9
|
55 intervalReads = rev(c(0,10,25,50,100,250,500,750,1000,10000))
|
|
6
|
56 intervalFreq = rev(c(0,0.01,0.05,0.1,0.5,1,5))
|
|
0
|
57 V_Segments = c(".*", "IGHV", "IGHD", "IGKV", "IGKV", "IgKINTR", "TRGV", "TRDV", "TRDD" , "TRBV")
|
|
|
58 J_Segments = c(".*", ".*", ".*", "IGKJ", "KDE", ".*", ".*", ".*", ".*", ".*")
|
|
|
59 Titles = c("Total", "IGH-Vh-Jh", "IGH-Dh-Jh", "Vk-Jk", "Vk-Kde" , "Intron-Kde", "TCRG", "TCRD-Vd-Dd", "TCRD-Dd-Dd", "TCRB-Vb-Jb")
|
|
|
60 Titles = factor(Titles, levels=Titles)
|
|
|
61 TitlesOrder = data.frame("Title"=Titles, "TitlesOrder"=1:length(Titles))
|
|
|
62
|
|
|
63 patientCountOnColumn <- function(x, product, interval, on, appendtxt=F){
|
|
28
|
64 if (!is.data.frame(x) & is.list(x)){
|
|
|
65 x = x[[1]]
|
|
|
66 }
|
|
0
|
67 x$Sample = factor(x$Sample, levels=unique(x$Sample))
|
|
|
68 onShort = "reads"
|
|
|
69 if(on == "Frequency"){
|
|
|
70 onShort = "freq"
|
|
|
71 }
|
|
18
|
72 onx = paste(on, ".x", sep="")
|
|
|
73 ony = paste(on, ".y", sep="")
|
|
0
|
74 splt = split(x, x$Sample, drop=T)
|
|
4
|
75 type="pair"
|
|
2
|
76 if(length(splt) == 1){
|
|
3
|
77 print(paste(paste(x[1,which(colnames(x) == "Patient")]), "has one sample"))
|
|
4
|
78 splt[[2]] = data.frame("Patient" = character(0), "Receptor" = character(0), "Sample" = character(0), "Cell_Count" = numeric(0), "Clone_Molecule_Count_From_Spikes" = numeric(0), "Log10_Frequency" = numeric(0), "Total_Read_Count" = numeric(0), "dsMol_per_1e6_cells" = numeric(0), "J_Segment_Major_Gene" = character(0), "V_Segment_Major_Gene" = character(0), "Clone_Sequence" = character(0), "CDR3_Sense_Sequence" = character(0), "Related_to_leukemia_clone" = logical(0), "Frequency"= numeric(0), "normalized_read_count" = numeric(0), "paste" = character(0))
|
|
|
79 type="single"
|
|
2
|
80 }
|
|
0
|
81 patient1 = splt[[1]]
|
|
|
82 patient2 = splt[[2]]
|
|
|
83
|
|
|
84 threshholdIndex = which(colnames(product) == "interval")
|
|
|
85 V_SegmentIndex = which(colnames(product) == "V_Segments")
|
|
|
86 J_SegmentIndex = which(colnames(product) == "J_Segments")
|
|
|
87 titleIndex = which(colnames(product) == "Titles")
|
|
|
88 sampleIndex = which(colnames(x) == "Sample")
|
|
|
89 patientIndex = which(colnames(x) == "Patient")
|
|
|
90 oneSample = paste(patient1[1,sampleIndex], sep="")
|
|
|
91 twoSample = paste(patient2[1,sampleIndex], sep="")
|
|
|
92 patient = paste(x[1,patientIndex])
|
|
3
|
93
|
|
0
|
94 switched = F
|
|
|
95 if(length(grep(".*_Right$", twoSample)) == 1 || length(grep(".*_Dx_BM$", twoSample)) == 1 || length(grep(".*_Dx$", twoSample)) == 1 ){
|
|
|
96 tmp = twoSample
|
|
|
97 twoSample = oneSample
|
|
|
98 oneSample = tmp
|
|
|
99 tmp = patient1
|
|
|
100 patient1 = patient2
|
|
|
101 patient2 = tmp
|
|
|
102 switched = T
|
|
|
103 }
|
|
|
104 if(appendtxt){
|
|
4
|
105 cat(paste(patient, oneSample, twoSample, type, sep="\t"), file="patients.txt", append=T, sep="", fill=3)
|
|
0
|
106 }
|
|
3
|
107 cat(paste("<tr><td>", patient, "</td></tr>", sep=""), file=logfile, append=T)
|
|
9
|
108
|
|
12
|
109 #patient1$merge = paste(patient1$V_Segment_Major_Gene, patient1$J_Segment_Major_Gene, patient1$CDR3_Sense_Sequence)
|
|
|
110 #patient2$merge = paste(patient2$V_Segment_Major_Gene, patient2$J_Segment_Major_Gene, patient2$CDR3_Sense_Sequence)
|
|
13
|
111 patient1$merge = paste(patient1$Clone_Sequence)
|
|
|
112 patient2$merge = paste(patient2$Clone_Sequence)
|
|
9
|
113
|
|
12
|
114 #patientMerge = merge(patient1, patient2, by.x="merge", by.y="merge")
|
|
9
|
115 patientMerge = merge(patient1, patient2, by.x="merge", by.y="merge")
|
|
18
|
116 patientMerge$thresholdValue = pmax(patientMerge[,onx], patientMerge[,ony])
|
|
0
|
117 res1 = vector()
|
|
|
118 res2 = vector()
|
|
|
119 resBoth = vector()
|
|
|
120 read1Count = vector()
|
|
|
121 read2Count = vector()
|
|
2
|
122 locussum1 = vector()
|
|
|
123 locussum2 = vector()
|
|
9
|
124
|
|
|
125 print(patient)
|
|
0
|
126 #for(iter in 1){
|
|
|
127 for(iter in 1:length(product[,1])){
|
|
|
128 threshhold = product[iter,threshholdIndex]
|
|
|
129 V_Segment = paste(".*", as.character(product[iter,V_SegmentIndex]), ".*", sep="")
|
|
|
130 J_Segment = paste(".*", as.character(product[iter,J_SegmentIndex]), ".*", sep="")
|
|
18
|
131 #both = (grepl(V_Segment, patientMerge$V_Segment_Major_Gene.x) & grepl(J_Segment, patientMerge$J_Segment_Major_Gene.x) & patientMerge[,onx] > threshhold & patientMerge[,ony] > threshhold) #both higher than threshold
|
|
|
132 both = (grepl(V_Segment, patientMerge$V_Segment_Major_Gene.x) & grepl(J_Segment, patientMerge$J_Segment_Major_Gene.x) & patientMerge$thresholdValue > threshhold) #highest of both higher than threshold
|
|
19
|
133 one = (grepl(V_Segment, patient1$V_Segment_Major_Gene) & grepl(J_Segment, patient1$J_Segment_Major_Gene) & patient1[,on] > threshhold & !(patient1$Clone_Sequence %in% patientMerge[both,]$merge))
|
|
|
134 two = (grepl(V_Segment, patient2$V_Segment_Major_Gene) & grepl(J_Segment, patient2$J_Segment_Major_Gene) & patient2[,on] > threshhold & !(patient2$Clone_Sequence %in% patientMerge[both,]$merge))
|
|
14
|
135 read1Count = append(read1Count, sum(patient1[one,]$normalized_read_count))
|
|
|
136 read2Count = append(read2Count, sum(patient2[two,]$normalized_read_count))
|
|
0
|
137 res1 = append(res1, sum(one))
|
|
2
|
138 res2 = append(res2, sum(two))
|
|
0
|
139 resBoth = append(resBoth, sum(both))
|
|
2
|
140 locussum1 = append(locussum1, sum(patient1[(grepl(V_Segment, patient1$V_Segment_Major_Gene) & grepl(J_Segment, patient1$J_Segment_Major_Gene)),]$normalized_read_count))
|
|
|
141 locussum2 = append(locussum2, sum(patient2[(grepl(V_Segment, patient2$V_Segment_Major_Gene) & grepl(J_Segment, patient2$J_Segment_Major_Gene)),]$normalized_read_count))
|
|
0
|
142 #threshhold = 0
|
|
|
143 if(threshhold != 0){
|
|
|
144 if(sum(one) > 0){
|
|
15
|
145 dfOne = patient1[one,c("V_Segment_Major_Gene", "J_Segment_Major_Gene", "normalized_read_count", "Frequency", "Clone_Sequence", "Related_to_leukemia_clone")]
|
|
|
146 colnames(dfOne) = c("Proximal segment", "Distal segment", "normalized_read_count", "Frequency", "Clone Sequence", "Related_to_leukemia_clone")
|
|
0
|
147 filenameOne = paste(oneSample, "_", product[iter, titleIndex], "_", threshhold, sep="")
|
|
|
148 write.table(dfOne, file=paste(filenameOne, ".txt", sep=""), quote=F, sep="\t", dec=",", row.names=F, col.names=T)
|
|
|
149 }
|
|
|
150 if(sum(two) > 0){
|
|
15
|
151 dfTwo = patient2[two,c("V_Segment_Major_Gene", "J_Segment_Major_Gene", "normalized_read_count", "Frequency", "Clone_Sequence", "Related_to_leukemia_clone")]
|
|
|
152 colnames(dfTwo) = c("Proximal segment", "Distal segment", "normalized_read_count", "Frequency", "Clone Sequence", "Related_to_leukemia_clone")
|
|
0
|
153 filenameTwo = paste(twoSample, "_", product[iter, titleIndex], "_", threshhold, sep="")
|
|
|
154 write.table(dfTwo, file=paste(filenameTwo, ".txt", sep=""), quote=F, sep="\t", dec=",", row.names=F, col.names=T)
|
|
|
155 }
|
|
|
156 }
|
|
|
157 if(sum(both) > 0){
|
|
15
|
158 dfBoth = patientMerge[both,c("V_Segment_Major_Gene.x", "J_Segment_Major_Gene.x", "normalized_read_count.x", "Frequency.x", "Related_to_leukemia_clone.x", "Clone_Sequence.x", "V_Segment_Major_Gene.y", "J_Segment_Major_Gene.y", "normalized_read_count.y", "Frequency.y", "Related_to_leukemia_clone.y")]
|
|
|
159 colnames(dfBoth) = c(paste("Proximal segment", oneSample), paste("Distal segment", oneSample), paste("Normalized_Read_Count", oneSample), paste("Frequency", oneSample), paste("Related_to_leukemia_clone", oneSample),"Clone Sequence", paste("Proximal segment", twoSample), paste("Distal segment", twoSample), paste("Normalized_Read_Count", twoSample), paste("Frequency", twoSample), paste("Related_to_leukemia_clone", twoSample))
|
|
0
|
160 filenameBoth = paste(oneSample, "_", twoSample, "_", product[iter, titleIndex], "_", threshhold, sep="")
|
|
|
161 write.table(dfBoth, file=paste(filenameBoth, ".txt", sep=""), quote=F, sep="\t", dec=",", row.names=F, col.names=T)
|
|
|
162 }
|
|
|
163 }
|
|
2
|
164 patientResult = data.frame("Locus"=product$Titles, "J_Segment"=product$J_Segments, "V_Segment"=product$V_Segments, "cut_off_value"=paste(">", product$interval, sep=""), "Both"=resBoth, "tmp1"=res1, "read_count1" = round(read1Count), "tmp2"=res2, "read_count2"= round(read2Count), "Sum"=res1 + res2 + resBoth, "percentage" = round((resBoth/(res1 + res2 + resBoth)) * 100, digits=2), "Locus_sum1"=locussum1, "Locus_sum2"=locussum2)
|
|
0
|
165 if(sum(is.na(patientResult$percentage)) > 0){
|
|
|
166 patientResult[is.na(patientResult$percentage),]$percentage = 0
|
|
|
167 }
|
|
|
168 colnames(patientResult)[6] = oneSample
|
|
|
169 colnames(patientResult)[8] = twoSample
|
|
|
170 colnamesBak = colnames(patientResult)
|
|
2
|
171 colnames(patientResult) = c("Ig/TCR gene rearrangement type", "Distal Gene segment", "Proximal gene segment", "cut_off_value", paste("Number of sequences ", patient, "_Both", sep=""), paste("Number of sequences", oneSample, sep=""), paste("Normalized Read Count", oneSample), paste("Number of sequences", twoSample, sep=""), paste("Normalized Read Count", twoSample), paste("Sum number of sequences", patient), paste("Percentage of sequences ", patient, "_Both", sep=""), paste("Locus Sum", oneSample), paste("Locus Sum", twoSample))
|
|
0
|
172 write.table(patientResult, file=paste(patient, "_", onShort, ".txt", sep=""), quote=F, sep="\t", dec=",", row.names=F, col.names=T)
|
|
|
173 colnames(patientResult) = colnamesBak
|
|
|
174
|
|
|
175 patientResult$Locus = factor(patientResult$Locus, Titles)
|
|
|
176 patientResult$cut_off_value = factor(patientResult$cut_off_value, paste(">", interval, sep=""))
|
|
|
177
|
|
|
178 plt = ggplot(patientResult[,c("Locus", "cut_off_value", "Both")])
|
|
|
179 plt = plt + geom_bar( aes( x=factor(cut_off_value), y=Both), stat='identity', position="dodge", fill="#79c36a")
|
|
|
180 plt = plt + facet_grid(.~Locus) + theme(axis.text.x = element_text(angle = 45, hjust = 1))
|
|
|
181 plt = plt + geom_text(aes(ymax=max(Both), x=cut_off_value,y=Both,label=Both), angle=90, hjust=0)
|
|
|
182 plt = plt + xlab("Reads per locus") + ylab("Count") + ggtitle("Number of clones in both")
|
|
|
183 plt = plt + theme(plot.margin = unit(c(1,8.8,0.5,1.5), "lines"))
|
|
|
184 png(paste(patient, "_", onShort, ".png", sep=""), width=1920, height=1080)
|
|
|
185 print(plt)
|
|
|
186 dev.off()
|
|
|
187 #(t,r,b,l)
|
|
|
188 plt = ggplot(patientResult[,c("Locus", "cut_off_value", "percentage")])
|
|
|
189 plt = plt + geom_bar( aes( x=factor(cut_off_value), y=percentage), stat='identity', position="dodge", fill="#79c36a")
|
|
|
190 plt = plt + facet_grid(.~Locus) + theme(axis.text.x = element_text(angle = 45, hjust = 1))
|
|
|
191 plt = plt + geom_text(aes(ymax=max(percentage), x=cut_off_value,y=percentage,label=percentage), angle=90, hjust=0)
|
|
|
192 plt = plt + xlab("Reads per locus") + ylab("Count") + ggtitle("% clones in both left and right")
|
|
|
193 plt = plt + theme(plot.margin = unit(c(1,8.8,0.5,1.5), "lines"))
|
|
|
194 png(paste(patient, "_percent_", onShort, ".png", sep=""), width=1920, height=1080)
|
|
|
195 print(plt)
|
|
|
196 dev.off()
|
|
|
197
|
|
|
198 patientResult = melt(patientResult[,c('Locus','cut_off_value', oneSample, twoSample)] ,id.vars=1:2)
|
|
|
199 patientResult$relativeValue = patientResult$value * 10
|
|
|
200 patientResult[patientResult$relativeValue == 0,]$relativeValue = 1
|
|
|
201 plt = ggplot(patientResult)
|
|
|
202 plt = plt + geom_bar( aes( x=factor(cut_off_value), y=relativeValue, fill=variable), stat='identity', position="dodge")
|
|
|
203 plt = plt + facet_grid(.~Locus) + theme(axis.text.x = element_text(angle = 45, hjust = 1))
|
|
|
204 plt = plt + scale_y_continuous(trans="log", breaks=10^c(0:10), labels=c(0, 10^c(0:9)))
|
|
|
205 plt = plt + geom_text(data=patientResult[patientResult$variable == oneSample,], aes(ymax=max(value), x=cut_off_value,y=relativeValue,label=value), angle=90, position=position_dodge(width=0.9), hjust=0, vjust=-0.2)
|
|
|
206 plt = plt + geom_text(data=patientResult[patientResult$variable == twoSample,], aes(ymax=max(value), x=cut_off_value,y=relativeValue,label=value), angle=90, position=position_dodge(width=0.9), hjust=0, vjust=0.8)
|
|
|
207 plt = plt + xlab("Reads per locus") + ylab("Count") + ggtitle(paste("Number of clones in only ", oneSample, " and only ", twoSample, sep=""))
|
|
|
208 png(paste(patient, "_", onShort, "_both.png", sep=""), width=1920, height=1080)
|
|
|
209 print(plt)
|
|
|
210 dev.off()
|
|
|
211 }
|
|
|
212
|
|
3
|
213 cat("<tr><td>Starting Frequency analysis</td></tr>", file=logfile, append=T)
|
|
|
214
|
|
0
|
215 interval = intervalFreq
|
|
|
216 intervalOrder = data.frame("interval"=paste(">", interval, sep=""), "intervalOrder"=1:length(interval))
|
|
4
|
217 product = data.frame("Titles"=rep(Titles, each=length(interval)), "interval"=rep(interval, times=10), "V_Segments"=rep(V_Segments, each=length(interval)), "J_Segments"=rep(J_Segments, each=length(interval)))
|
|
|
218 mclapply(patients, FUN=patientCountOnColumn, product = product, interval=interval, on="Frequency", appendtxt=T)
|
|
0
|
219
|
|
3
|
220 cat("<tr><td>Starting Cell Count analysis</td></tr>", file=logfile, append=T)
|
|
|
221
|
|
0
|
222 interval = intervalReads
|
|
|
223 intervalOrder = data.frame("interval"=paste(">", interval, sep=""), "intervalOrder"=1:length(interval))
|
|
4
|
224 product = data.frame("Titles"=rep(Titles, each=length(interval)), "interval"=rep(interval, times=10), "V_Segments"=rep(V_Segments, each=length(interval)), "J_Segments"=rep(J_Segments, each=length(interval)))
|
|
9
|
225 mclapply(patients, FUN=patientCountOnColumn, product = product, interval=interval, on="normalized_read_count")
|
|
0
|
226
|
|
3
|
227 cat("</table></html>", file=logfile, append=T)
|
|
|
228
|
|
7
|
229
|
|
13
|
230
|
|
7
|
231 tripletAnalysis <- function(patient1, label1, patient2, label2, patient3, label3, product, interval, on, appendTriplets= FALSE){
|
|
|
232 onShort = "reads"
|
|
|
233 if(on == "Frequency"){
|
|
|
234 onShort = "freq"
|
|
|
235 }
|
|
18
|
236 onx = paste(on, ".x", sep="")
|
|
|
237 ony = paste(on, ".y", sep="")
|
|
|
238 onz = paste(on, ".z", sep="")
|
|
7
|
239 type="triplet"
|
|
|
240
|
|
|
241 threshholdIndex = which(colnames(product) == "interval")
|
|
|
242 V_SegmentIndex = which(colnames(product) == "V_Segments")
|
|
|
243 J_SegmentIndex = which(colnames(product) == "J_Segments")
|
|
|
244 titleIndex = which(colnames(product) == "Titles")
|
|
|
245 sampleIndex = which(colnames(patient1) == "Sample")
|
|
|
246 patientIndex = which(colnames(patient1) == "Patient")
|
|
|
247 oneSample = paste(patient1[1,sampleIndex], sep="")
|
|
|
248 twoSample = paste(patient2[1,sampleIndex], sep="")
|
|
|
249 threeSample = paste(patient3[1,sampleIndex], sep="")
|
|
|
250
|
|
12
|
251 #patient1$merge = paste(patient1$V_Segment_Major_Gene, patient1$J_Segment_Major_Gene, patient1$CDR3_Sense_Sequence)
|
|
|
252 #patient2$merge = paste(patient2$V_Segment_Major_Gene, patient2$J_Segment_Major_Gene, patient2$CDR3_Sense_Sequence)
|
|
|
253 #patient3$merge = paste(patient3$V_Segment_Major_Gene, patient3$J_Segment_Major_Gene, patient3$CDR3_Sense_Sequence)
|
|
|
254
|
|
15
|
255 patient1$merge = paste(patient1$Clone_Sequence)
|
|
|
256 patient2$merge = paste(patient2$Clone_Sequence)
|
|
|
257 patient3$merge = paste(patient3$Clone_Sequence)
|
|
9
|
258
|
|
|
259 patientMerge = merge(patient1, patient2, by="merge")
|
|
|
260 patientMerge = merge(patientMerge, patient3, by="merge")
|
|
28
|
261 colnames(patientMerge)[which(!grepl("(\\.x$)|(\\.y$)|(merge)", names(patientMerge)))] = paste(colnames(patientMerge)[which(!grepl("(\\.x$)|(\\.y$)|(merge)", names(patientMerge), perl=T))], ".z", sep="")
|
|
18
|
262 patientMerge$thresholdValue = pmax(patientMerge[,onx], patientMerge[,ony], patientMerge[,onz])
|
|
20
|
263 patientMerge12 = merge(patient1, patient2, by="merge")
|
|
|
264 patientMerge12$thresholdValue = pmax(patientMerge12[,onx], patientMerge12[,ony])
|
|
|
265 patientMerge13 = merge(patient1, patient3, by="merge")
|
|
|
266 patientMerge13$thresholdValue = pmax(patientMerge13[,onx], patientMerge13[,ony])
|
|
|
267 patientMerge23 = merge(patient2, patient3, by="merge")
|
|
|
268 patientMerge23$thresholdValue = pmax(patientMerge23[,onx], patientMerge23[,ony])
|
|
|
269
|
|
28
|
270
|
|
24
|
271 scatterplot_data_columns = c("Clone_Sequence", "Frequency", "normalized_read_count", "V_Segment_Major_Gene", "J_Segment_Major_Gene")
|
|
20
|
272 scatterplot_data = rbind(patient1[,scatterplot_data_columns], patient2[,scatterplot_data_columns], patient3[,scatterplot_data_columns])
|
|
|
273 scatterplot_data = scatterplot_data[!duplicated(scatterplot_data$Clone_Sequence),]
|
|
27
|
274 scatterplot_data$type = factor(x="In one", levels=c("In one", "In two", "In three", "In multiple"))
|
|
20
|
275
|
|
7
|
276 res1 = vector()
|
|
|
277 res2 = vector()
|
|
|
278 res3 = vector()
|
|
20
|
279 res12 = vector()
|
|
|
280 res13 = vector()
|
|
|
281 res23 = vector()
|
|
7
|
282 resAll = vector()
|
|
|
283 read1Count = vector()
|
|
|
284 read2Count = vector()
|
|
|
285 read3Count = vector()
|
|
|
286
|
|
|
287 if(appendTriplets){
|
|
9
|
288 cat(paste(label1, label2, label3, sep="\t"), file="triplets.txt", append=T, sep="", fill=3)
|
|
7
|
289 }
|
|
|
290 for(iter in 1:length(product[,1])){
|
|
|
291 threshhold = product[iter,threshholdIndex]
|
|
|
292 V_Segment = paste(".*", as.character(product[iter,V_SegmentIndex]), ".*", sep="")
|
|
|
293 J_Segment = paste(".*", as.character(product[iter,J_SegmentIndex]), ".*", sep="")
|
|
18
|
294 #all = (grepl(V_Segment, patientMerge$V_Segment_Major_Gene.x) & grepl(J_Segment, patientMerge$J_Segment_Major_Gene.x) & patientMerge[,onx] > threshhold & patientMerge[,ony] > threshhold & patientMerge[,onz] > threshhold)
|
|
|
295 all = (grepl(V_Segment, patientMerge$V_Segment_Major_Gene.x) & grepl(J_Segment, patientMerge$J_Segment_Major_Gene.x) & patientMerge$thresholdValue > threshhold)
|
|
7
|
296
|
|
25
|
297 one_two = (grepl(V_Segment, patientMerge12$V_Segment_Major_Gene.x) & grepl(J_Segment, patientMerge12$J_Segment_Major_Gene.x) & patientMerge12$thresholdValue > threshhold & !(patientMerge12$Clone_Sequence.x %in% patientMerge[all,]$merge))
|
|
|
298 one_three = (grepl(V_Segment, patientMerge13$V_Segment_Major_Gene.x) & grepl(J_Segment, patientMerge13$J_Segment_Major_Gene.x) & patientMerge13$thresholdValue > threshhold & !(patientMerge13$Clone_Sequence.x %in% patientMerge[all,]$merge))
|
|
|
299 two_three = (grepl(V_Segment, patientMerge23$V_Segment_Major_Gene.x) & grepl(J_Segment, patientMerge23$J_Segment_Major_Gene.x) & patientMerge23$thresholdValue > threshhold & !(patientMerge23$Clone_Sequence.x %in% patientMerge[all,]$merge))
|
|
24
|
300
|
|
25
|
301 one = (grepl(V_Segment, patient1$V_Segment_Major_Gene) & grepl(J_Segment, patient1$J_Segment_Major_Gene) & patient1[,on] > threshhold & !(patient1$Clone_Sequence %in% patientMerge[all,]$merge) & !(patient1$Clone_Sequence %in% patientMerge12[one_two,]$merge) & !(patient1$Clone_Sequence %in% patientMerge13[one_three,]$merge))
|
|
|
302 two = (grepl(V_Segment, patient2$V_Segment_Major_Gene) & grepl(J_Segment, patient2$J_Segment_Major_Gene) & patient2[,on] > threshhold & !(patient2$Clone_Sequence %in% patientMerge[all,]$merge) & !(patient2$Clone_Sequence %in% patientMerge12[one_two,]$merge) & !(patient2$Clone_Sequence %in% patientMerge23[two_three,]$merge))
|
|
|
303 three = (grepl(V_Segment, patient3$V_Segment_Major_Gene) & grepl(J_Segment, patient3$J_Segment_Major_Gene) & patient3[,on] > threshhold & !(patient3$Clone_Sequence %in% patientMerge[all,]$merge) & !(patient3$Clone_Sequence %in% patientMerge13[one_three,]$merge) & !(patient3$Clone_Sequence %in% patientMerge23[two_three,]$merge))
|
|
20
|
304
|
|
18
|
305 read1Count = append(read1Count, sum(patient1[one,]$normalized_read_count) + sum(patientMerge[all,]$normalized_read_count.x))
|
|
|
306 read2Count = append(read2Count, sum(patient2[two,]$normalized_read_count) + sum(patientMerge[all,]$normalized_read_count.y))
|
|
|
307 read3Count = append(read3Count, sum(patient3[three,]$normalized_read_count) + sum(patientMerge[all,]$normalized_read_count.z))
|
|
7
|
308 res1 = append(res1, sum(one))
|
|
|
309 res2 = append(res2, sum(two))
|
|
|
310 res3 = append(res3, sum(three))
|
|
|
311 resAll = append(resAll, sum(all))
|
|
20
|
312 res12 = append(res12, sum(one_two))
|
|
|
313 res13 = append(res13, sum(one_three))
|
|
|
314 res23 = append(res23, sum(two_three))
|
|
7
|
315 #threshhold = 0
|
|
|
316 if(threshhold != 0){
|
|
|
317 if(sum(one) > 0){
|
|
20
|
318 dfOne = patient1[one,c("V_Segment_Major_Gene", "J_Segment_Major_Gene", "normalized_read_count", "Frequency", "Clone_Sequence", "Related_to_leukemia_clone")]
|
|
|
319 colnames(dfOne) = c("Proximal segment", "Distal segment", "normalized_read_count", "Frequency", "Clone_Sequence", "Related_to_leukemia_clone")
|
|
7
|
320 filenameOne = paste(label1, "_", product[iter, titleIndex], "_", threshhold, sep="")
|
|
|
321 write.table(dfOne, file=paste(filenameOne, ".txt", sep=""), quote=F, sep="\t", dec=",", row.names=F, col.names=T)
|
|
|
322 }
|
|
|
323 if(sum(two) > 0){
|
|
20
|
324 dfTwo = patient2[two,c("V_Segment_Major_Gene", "J_Segment_Major_Gene", "normalized_read_count", "Frequency", "Clone_Sequence", "Related_to_leukemia_clone")]
|
|
|
325 colnames(dfTwo) = c("Proximal segment", "Distal segment", "normalized_read_count", "Frequency", "Clone_Sequence", "Related_to_leukemia_clone")
|
|
7
|
326 filenameTwo = paste(label2, "_", product[iter, titleIndex], "_", threshhold, sep="")
|
|
|
327 write.table(dfTwo, file=paste(filenameTwo, ".txt", sep=""), quote=F, sep="\t", dec=",", row.names=F, col.names=T)
|
|
|
328 }
|
|
|
329 if(sum(three) > 0){
|
|
20
|
330 dfThree = patient3[three,c("V_Segment_Major_Gene", "J_Segment_Major_Gene", "normalized_read_count", "Frequency", "Clone_Sequence", "Related_to_leukemia_clone")]
|
|
|
331 colnames(dfThree) = c("Proximal segment", "Distal segment", "normalized_read_count", "Frequency", "Clone_Sequence", "Related_to_leukemia_clone")
|
|
7
|
332 filenameThree = paste(label3, "_", product[iter, titleIndex], "_", threshhold, sep="")
|
|
|
333 write.table(dfThree, file=paste(filenameThree, ".txt", sep=""), quote=F, sep="\t", dec=",", row.names=F, col.names=T)
|
|
|
334 }
|
|
20
|
335 if(sum(one_two) > 0){
|
|
|
336 dfOne_two = patientMerge12[one_two,c("V_Segment_Major_Gene.x", "J_Segment_Major_Gene.x", "normalized_read_count.x", "Frequency.x", "Related_to_leukemia_clone.x", "Clone_Sequence.x", "V_Segment_Major_Gene.y", "J_Segment_Major_Gene.y", "normalized_read_count.y", "Frequency.y", "Related_to_leukemia_clone.y")]
|
|
|
337 colnames(dfOne_two) = c(paste("Proximal segment", oneSample), paste("Distal segment", oneSample), paste("Normalized_Read_Count", oneSample), paste("Frequency", oneSample), paste("Related_to_leukemia_clone", oneSample),"Clone_Sequence", paste("Proximal segment", twoSample), paste("Distal segment", twoSample), paste("Normalized_Read_Count", twoSample), paste("Frequency", twoSample), paste("Related_to_leukemia_clone", twoSample))
|
|
|
338 filenameOne_two = paste(label1, "_", label2, "_", product[iter, titleIndex], "_", threshhold, onShort, sep="")
|
|
|
339 write.table(dfOne_two, file=paste(filenameOne_two, ".txt", sep=""), quote=F, sep="\t", dec=",", row.names=F, col.names=T)
|
|
|
340 }
|
|
|
341 if(sum(one_three) > 0){
|
|
|
342 dfOne_three = patientMerge13[one_three,c("V_Segment_Major_Gene.x", "J_Segment_Major_Gene.x", "normalized_read_count.x", "Frequency.x", "Related_to_leukemia_clone.x", "Clone_Sequence.x", "V_Segment_Major_Gene.y", "J_Segment_Major_Gene.y", "normalized_read_count.y", "Frequency.y", "Related_to_leukemia_clone.y")]
|
|
|
343 colnames(dfOne_three) = c(paste("Proximal segment", oneSample), paste("Distal segment", oneSample), paste("Normalized_Read_Count", oneSample), paste("Frequency", oneSample), paste("Related_to_leukemia_clone", oneSample),"Clone_Sequence", paste("Proximal segment", threeSample), paste("Distal segment", threeSample), paste("Normalized_Read_Count", threeSample), paste("Frequency", threeSample), paste("Related_to_leukemia_clone", threeSample))
|
|
|
344 filenameOne_three = paste(label1, "_", label3, "_", product[iter, titleIndex], "_", threshhold, onShort, sep="")
|
|
|
345 write.table(dfOne_three, file=paste(filenameOne_three, ".txt", sep=""), quote=F, sep="\t", dec=",", row.names=F, col.names=T)
|
|
|
346 }
|
|
|
347 if(sum(two_three) > 0){
|
|
|
348 dfTwo_three = patientMerge23[two_three,c("V_Segment_Major_Gene.x", "J_Segment_Major_Gene.x", "normalized_read_count.x", "Frequency.x", "Related_to_leukemia_clone.x", "Clone_Sequence.x", "V_Segment_Major_Gene.y", "J_Segment_Major_Gene.y", "normalized_read_count.y", "Frequency.y", "Related_to_leukemia_clone.y")]
|
|
|
349 colnames(dfTwo_three) = c(paste("Proximal segment", twoSample), paste("Distal segment", twoSample), paste("Normalized_Read_Count", twoSample), paste("Frequency", twoSample), paste("Related_to_leukemia_clone", twoSample),"Clone_Sequence", paste("Proximal segment", threeSample), paste("Distal segment", threeSample), paste("Normalized_Read_Count", threeSample), paste("Frequency", threeSample), paste("Related_to_leukemia_clone", threeSample))
|
|
|
350 filenameTwo_three = paste(label2, "_", label3, "_", product[iter, titleIndex], "_", threshhold, onShort, sep="")
|
|
|
351 write.table(dfTwo_three, file=paste(filenameTwo_three, ".txt", sep=""), quote=F, sep="\t", dec=",", row.names=F, col.names=T)
|
|
|
352 }
|
|
|
353 } else { #scatterplot data
|
|
24
|
354 scatterplot_locus_data = scatterplot_data[grepl(V_Segment, scatterplot_data$V_Segment_Major_Gene) & grepl(J_Segment, scatterplot_data$J_Segment_Major_Gene),]
|
|
27
|
355 in_two = (scatterplot_locus_data$Clone_Sequence %in% patientMerge12[one_two,]$Clone_Sequence.x) | (scatterplot_locus_data$Clone_Sequence %in% patientMerge13[one_three,]$Clone_Sequence.x) | (scatterplot_locus_data$Clone_Sequence %in% patientMerge23[two_three,]$Clone_Sequence.x)
|
|
|
356 if(sum(in_two) > 0){
|
|
|
357 scatterplot_locus_data[in_two,]$type = "In two"
|
|
|
358 }
|
|
|
359 in_three = (scatterplot_locus_data$Clone_Sequence %in% patientMerge[all,]$Clone_Sequence.x)
|
|
|
360 if(sum(in_three)> 0){
|
|
|
361 scatterplot_locus_data[in_three,]$type = "In three"
|
|
|
362 }
|
|
|
363 not_in_one = scatterplot_locus_data$type != "In one"
|
|
|
364 if(sum(not_in_one) > 0){
|
|
|
365 scatterplot_locus_data[not_in_one,]$type = "In multiple"
|
|
|
366 }
|
|
20
|
367 p = NULL
|
|
|
368 if(nrow(scatterplot_locus_data) != 0){
|
|
|
369 if(on == "normalized_read_count"){
|
|
|
370 scales = 10^(0:ceiling(log10(max(scatterplot_locus_data$normalized_read_count))))
|
|
|
371 p = ggplot(scatterplot_locus_data, aes(type, normalized_read_count)) + scale_y_log10(breaks=scales,labels=scales)
|
|
|
372 } else {
|
|
|
373 p = ggplot(scatterplot_locus_data, aes(type, Frequency))
|
|
|
374 }
|
|
|
375 p = p + geom_point(aes(colour=type), position="jitter")
|
|
25
|
376 p = p + xlab("In one or in multiple samples") + ylab(onShort) + ggtitle(paste(label1, label2, label3, onShort, product[iter, titleIndex]))
|
|
20
|
377 } else {
|
|
|
378 p = ggplot(NULL, aes(x=c("In one", "In multiple"),y=0)) + geom_blank(NULL) + xlab("In two or in three of the samples") + ylab(onShort) + ggtitle(paste(label1, label2, label3, onShort, product[iter, titleIndex]))
|
|
|
379 }
|
|
|
380 png(paste(label1, "_", label2, "_", label3, "_", onShort, "_", product[iter, titleIndex],"_scatter.png", sep=""))
|
|
|
381 print(p)
|
|
|
382 dev.off()
|
|
|
383 }
|
|
7
|
384 if(sum(all) > 0){
|
|
20
|
385 dfAll = patientMerge[all,c("V_Segment_Major_Gene.x", "J_Segment_Major_Gene.x", "normalized_read_count.x", "Frequency.x", "Related_to_leukemia_clone.x", "Clone_Sequence.x", "V_Segment_Major_Gene.y", "J_Segment_Major_Gene.y", "normalized_read_count.y", "Frequency.y", "Related_to_leukemia_clone.y", "V_Segment_Major_Gene.z", "J_Segment_Major_Gene.z", "normalized_read_count.z", "Frequency.z", "Related_to_leukemia_clone.z")]
|
|
|
386 colnames(dfAll) = c(paste("Proximal segment", oneSample), paste("Distal segment", oneSample), paste("Normalized_Read_Count", oneSample), paste("Frequency", oneSample), paste("Related_to_leukemia_clone", oneSample),"Clone_Sequence", paste("Proximal segment", twoSample), paste("Distal segment", twoSample), paste("Normalized_Read_Count", twoSample), paste("Frequency", twoSample), paste("Related_to_leukemia_clone", twoSample), paste("Proximal segment", threeSample), paste("Distal segment", threeSample), paste("Normalized_Read_Count", threeSample), paste("Frequency", threeSample), paste("Related_to_leukemia_clone", threeSample))
|
|
7
|
387 filenameAll = paste(label1, "_", label2, "_", label3, "_", product[iter, titleIndex], "_", threshhold, sep="")
|
|
|
388 write.table(dfAll, file=paste(filenameAll, ".txt", sep=""), quote=F, sep="\t", dec=",", row.names=F, col.names=T)
|
|
|
389 }
|
|
|
390 }
|
|
20
|
391 #patientResult = data.frame("Locus"=product$Titles, "J_Segment"=product$J_Segments, "V_Segment"=product$V_Segments, "cut_off_value"=paste(">", product$interval, sep=""), "All"=resAll, "tmp1"=res1, "read_count1" = round(read1Count), "tmp2"=res2, "read_count2"= round(read2Count), "tmp3"=res3, "read_count3"=round(read3Count))
|
|
|
392 patientResult = data.frame("Locus"=product$Titles, "J_Segment"=product$J_Segments, "V_Segment"=product$V_Segments, "cut_off_value"=paste(">", product$interval, sep=""), "All"=resAll, "tmp1"=res1, "tmp2"=res2, "tmp3"=res3, "tmp12"=res12, "tmp13"=res13, "tmp23"=res23)
|
|
7
|
393 colnames(patientResult)[6] = oneSample
|
|
20
|
394 colnames(patientResult)[7] = twoSample
|
|
|
395 colnames(patientResult)[8] = threeSample
|
|
|
396 colnames(patientResult)[9] = paste(oneSample, twoSample, sep="_")
|
|
|
397 colnames(patientResult)[10] = paste(oneSample, twoSample, sep="_")
|
|
|
398 colnames(patientResult)[11] = paste(oneSample, twoSample, sep="_")
|
|
7
|
399
|
|
|
400 colnamesBak = colnames(patientResult)
|
|
20
|
401 colnames(patientResult) = c("Ig/TCR gene rearrangement type", "Distal Gene segment", "Proximal gene segment", "cut_off_value", "Number of sequences All", paste("Number of sequences", oneSample), paste("Number of sequences", twoSample), paste("Number of sequences", threeSample), paste("Number of sequences", oneSample, twoSample), paste("Number of sequences", oneSample, threeSample), paste("Number of sequences", twoSample, threeSample))
|
|
7
|
402 write.table(patientResult, file=paste(label1, "_", label2, "_", label3, "_", onShort, ".txt", sep=""), quote=F, sep="\t", dec=",", row.names=F, col.names=T)
|
|
|
403 colnames(patientResult) = colnamesBak
|
|
|
404
|
|
|
405 patientResult$Locus = factor(patientResult$Locus, Titles)
|
|
|
406 patientResult$cut_off_value = factor(patientResult$cut_off_value, paste(">", interval, sep=""))
|
|
|
407
|
|
|
408 plt = ggplot(patientResult[,c("Locus", "cut_off_value", "All")])
|
|
|
409 plt = plt + geom_bar( aes( x=factor(cut_off_value), y=All), stat='identity', position="dodge", fill="#79c36a")
|
|
|
410 plt = plt + facet_grid(.~Locus) + theme(axis.text.x = element_text(angle = 45, hjust = 1))
|
|
|
411 plt = plt + geom_text(aes(ymax=max(All), x=cut_off_value,y=All,label=All), angle=90, hjust=0)
|
|
|
412 plt = plt + xlab("Reads per locus") + ylab("Count") + ggtitle("Number of clones in All")
|
|
|
413 plt = plt + theme(plot.margin = unit(c(1,8.8,0.5,1.5), "lines"))
|
|
|
414 png(paste(label1, "_", label2, "_", label3, "_", onShort, "_total_all.png", sep=""), width=1920, height=1080)
|
|
|
415 print(plt)
|
|
|
416 dev.off()
|
|
|
417
|
|
|
418 fontSize = 4
|
|
|
419
|
|
|
420 bak = patientResult
|
|
|
421 patientResult = melt(patientResult[,c('Locus','cut_off_value', oneSample, twoSample, threeSample)] ,id.vars=1:2)
|
|
|
422 patientResult$relativeValue = patientResult$value * 10
|
|
|
423 patientResult[patientResult$relativeValue == 0,]$relativeValue = 1
|
|
|
424 plt = ggplot(patientResult)
|
|
|
425 plt = plt + geom_bar( aes( x=factor(cut_off_value), y=relativeValue, fill=variable), stat='identity', position="dodge")
|
|
|
426 plt = plt + facet_grid(.~Locus) + theme(axis.text.x = element_text(angle = 45, hjust = 1))
|
|
|
427 plt = plt + scale_y_continuous(trans="log", breaks=10^c(0:10), labels=c(0, 10^c(0:9)))
|
|
|
428 plt = plt + geom_text(data=patientResult[patientResult$variable == oneSample,], aes(ymax=max(value), x=cut_off_value,y=relativeValue,label=value), angle=90, position=position_dodge(width=0.9), hjust=0, vjust=-0.7, size=fontSize)
|
|
|
429 plt = plt + geom_text(data=patientResult[patientResult$variable == twoSample,], aes(ymax=max(value), x=cut_off_value,y=relativeValue,label=value), angle=90, position=position_dodge(width=0.9), hjust=0, vjust=0.4, size=fontSize)
|
|
|
430 plt = plt + geom_text(data=patientResult[patientResult$variable == threeSample,], aes(ymax=max(value), x=cut_off_value,y=relativeValue,label=value), angle=90, position=position_dodge(width=0.9), hjust=0, vjust=1.5, size=fontSize)
|
|
|
431 plt = plt + xlab("Reads per locus") + ylab("Count") + ggtitle("Number of clones in only one sample")
|
|
|
432 png(paste(label1, "_", label2, "_", label3, "_", onShort, "_indiv_all.png", sep=""), width=1920, height=1080)
|
|
|
433 print(plt)
|
|
|
434 dev.off()
|
|
|
435 }
|
|
|
436
|
|
9
|
437 triplets$uniqueID = "ID"
|
|
|
438
|
|
|
439 triplets[grepl("16278_Left", triplets$Sample),]$uniqueID = "16278_26402_26759_Left"
|
|
|
440 triplets[grepl("26402_Left", triplets$Sample),]$uniqueID = "16278_26402_26759_Left"
|
|
|
441 triplets[grepl("26759_Left", triplets$Sample),]$uniqueID = "16278_26402_26759_Left"
|
|
|
442
|
|
|
443 triplets[grepl("16278_Right", triplets$Sample),]$uniqueID = "16278_26402_26759_Right"
|
|
|
444 triplets[grepl("26402_Right", triplets$Sample),]$uniqueID = "16278_26402_26759_Right"
|
|
|
445 triplets[grepl("26759_Right", triplets$Sample),]$uniqueID = "16278_26402_26759_Right"
|
|
|
446
|
|
|
447 triplets[grepl("14696", triplets$Patient),]$uniqueID = "14696"
|
|
|
448
|
|
13
|
449 triplets$locus_V = substring(triplets$V_Segment_Major_Gene, 0, 4)
|
|
|
450 triplets$locus_J = substring(triplets$J_Segment_Major_Gene, 0, 4)
|
|
|
451 min_cell_count = data.frame(data.table(triplets)[, list(min_cell_count=min(.SD$Cell_Count)), by=c("uniqueID", "locus_V", "locus_J")])
|
|
|
452
|
|
|
453 triplets$min_cell_paste = paste(triplets$uniqueID, triplets$locus_V, triplets$locus_J)
|
|
|
454 min_cell_count$min_cell_paste = paste(min_cell_count$uniqueID, min_cell_count$locus_V, min_cell_count$locus_J)
|
|
|
455
|
|
|
456 min_cell_count = min_cell_count[,c("min_cell_paste", "min_cell_count")]
|
|
9
|
457
|
|
13
|
458 triplets = merge(triplets, min_cell_count, by="min_cell_paste")
|
|
|
459
|
|
|
460 triplets$normalized_read_count = round(triplets$Clone_Molecule_Count_From_Spikes / triplets$Cell_Count * triplets$min_cell_count / 2, digits=2) #??????????????????????????????????? wel of geen / 2
|
|
|
461
|
|
|
462 triplets = triplets[triplets$normalized_read_count >= min_cells,]
|
|
|
463
|
|
|
464 column_drops = c("locus_V", "locus_J", "min_cell_count", "min_cell_paste")
|
|
|
465
|
|
|
466 triplets = triplets[,!(colnames(triplets) %in% column_drops)]
|
|
9
|
467
|
|
7
|
468 interval = intervalReads
|
|
|
469 intervalOrder = data.frame("interval"=paste(">", interval, sep=""), "intervalOrder"=1:length(interval))
|
|
|
470 product = data.frame("Titles"=rep(Titles, each=length(interval)), "interval"=rep(interval, times=10), "V_Segments"=rep(V_Segments, each=length(interval)), "J_Segments"=rep(J_Segments, each=length(interval)))
|
|
|
471
|
|
9
|
472 one = triplets[triplets$Sample == "14696_reg_BM",]
|
|
|
473 two = triplets[triplets$Sample == "24536_reg_BM",]
|
|
|
474 three = triplets[triplets$Sample == "24062_reg_BM",]
|
|
8
|
475 tripletAnalysis(one, "14696_1", two, "14696_2", three, "14696_3", product=product, interval=interval, on="normalized_read_count", T)
|
|
7
|
476
|
|
9
|
477 one = triplets[triplets$Sample == "16278_Left",]
|
|
|
478 two = triplets[triplets$Sample == "26402_Left",]
|
|
|
479 three = triplets[triplets$Sample == "26759_Left",]
|
|
8
|
480 tripletAnalysis(one, "16278_Left", two, "26402_Left", three, "26759_Left", product=product, interval=interval, on="normalized_read_count", T)
|
|
7
|
481
|
|
9
|
482 one = triplets[triplets$Sample == "16278_Right",]
|
|
|
483 two = triplets[triplets$Sample == "26402_Right",]
|
|
|
484 three = triplets[triplets$Sample == "26759_Right",]
|
|
8
|
485 tripletAnalysis(one, "16278_Right", two, "26402_Right", three, "26759_Right", product=product, interval=interval, on="normalized_read_count", T)
|
|
7
|
486
|
|
|
487
|
|
|
488 interval = intervalFreq
|
|
|
489 intervalOrder = data.frame("interval"=paste(">", interval, sep=""), "intervalOrder"=1:length(interval))
|
|
|
490 product = data.frame("Titles"=rep(Titles, each=length(interval)), "interval"=rep(interval, times=10), "V_Segments"=rep(V_Segments, each=length(interval)), "J_Segments"=rep(J_Segments, each=length(interval)))
|
|
|
491
|
|
9
|
492 one = triplets[triplets$Sample == "14696_reg_BM",]
|
|
|
493 two = triplets[triplets$Sample == "24536_reg_BM",]
|
|
|
494 three = triplets[triplets$Sample == "24062_reg_BM",]
|
|
8
|
495 tripletAnalysis(one, "14696_1", two, "14696_2", three, "14696_3", product=product, interval=interval, on="Frequency", F)
|
|
7
|
496
|
|
9
|
497 one = triplets[triplets$Sample == "16278_Left",]
|
|
|
498 two = triplets[triplets$Sample == "26402_Left",]
|
|
|
499 three = triplets[triplets$Sample == "26759_Left",]
|
|
8
|
500 tripletAnalysis(one, "16278_Left", two, "26402_Left", three, "26759_Left", product=product, interval=interval, on="Frequency", F)
|
|
7
|
501
|
|
9
|
502 one = triplets[triplets$Sample == "16278_Right",]
|
|
|
503 two = triplets[triplets$Sample == "26402_Right",]
|
|
|
504 three = triplets[triplets$Sample == "26759_Right",]
|
|
8
|
505 tripletAnalysis(one, "16278_Right", two, "26402_Right", three, "26759_Right", product=product, interval=interval, on="Frequency", F)
|
