annotate shm_csr.r @ 9:372ccdcf0b2d draft

Uploaded
author davidvanzessen
date Tue, 08 Nov 2016 07:32:54 -0500
parents 3968d04b5724
children 81453585dfc3
Ignore whitespace changes - Everywhere: Within whitespace: At end of lines:
rev   line source
1
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
1 library(data.table)
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
2 library(ggplot2)
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
3 library(reshape2)
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
4
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
5 args <- commandArgs(trailingOnly = TRUE)
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
6
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
7 input = args[1]
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
8 genes = unlist(strsplit(args[2], ","))
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
9 outputdir = args[3]
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
10 empty.region.filter = args[4]
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
11 setwd(outputdir)
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
12
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
13 dat = read.table(input, header=T, sep="\t", fill=T, stringsAsFactors=F)
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
14
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
15 if(length(dat$Sequence.ID) == 0){
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
16 setwd(outputdir)
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
17 result = data.frame(x = rep(0, 5), y = rep(0, 5), z = rep(NA, 5))
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
18 row.names(result) = c("Number of Mutations (%)", "Transition (%)", "Transversions (%)", "Transitions at G C (%)", "Targeting of C G (%)")
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
19 write.table(x=result, file="mutations.txt", sep=",",quote=F,row.names=T,col.names=F)
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
20 transitionTable = data.frame(A=rep(0, 4),C=rep(0, 4),G=rep(0, 4),T=rep(0, 4))
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
21 row.names(transitionTable) = c("A", "C", "G", "T")
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
22 transitionTable["A","A"] = NA
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
23 transitionTable["C","C"] = NA
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
24 transitionTable["G","G"] = NA
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
25 transitionTable["T","T"] = NA
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
26
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
27 write.table(x=transitionTable, file="transitions.txt", sep=",",quote=F,row.names=T,col.names=NA)
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
28 cat("0", file="n.txt")
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
29 stop("No data")
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
30 }
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
31
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
32 cleanup_columns = c("FR1.IMGT.c.a",
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
33 "FR2.IMGT.g.t",
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
34 "CDR1.IMGT.Nb.of.nucleotides",
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
35 "CDR2.IMGT.t.a",
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
36 "FR1.IMGT.c.g",
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
37 "CDR1.IMGT.c.t",
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
38 "FR2.IMGT.a.c",
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
39 "FR2.IMGT.Nb.of.mutations",
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
40 "FR2.IMGT.g.c",
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
41 "FR2.IMGT.a.g",
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
42 "FR3.IMGT.t.a",
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
43 "FR3.IMGT.t.c",
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
44 "FR2.IMGT.g.a",
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
45 "FR3.IMGT.c.g",
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
46 "FR1.IMGT.Nb.of.mutations",
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
47 "CDR1.IMGT.g.a",
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
48 "CDR1.IMGT.t.g",
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
49 "CDR1.IMGT.g.c",
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
50 "CDR2.IMGT.Nb.of.nucleotides",
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
51 "FR2.IMGT.a.t",
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
52 "CDR1.IMGT.Nb.of.mutations",
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
53 "CDR3.IMGT.Nb.of.nucleotides",
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
54 "CDR1.IMGT.a.g",
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
55 "FR3.IMGT.a.c",
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
56 "FR1.IMGT.g.a",
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
57 "FR3.IMGT.a.g",
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
58 "FR1.IMGT.a.t",
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
59 "CDR2.IMGT.a.g",
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
60 "CDR2.IMGT.Nb.of.mutations",
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
61 "CDR2.IMGT.g.t",
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
62 "CDR2.IMGT.a.c",
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
63 "CDR1.IMGT.t.c",
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
64 "FR3.IMGT.g.c",
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
65 "FR1.IMGT.g.t",
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
66 "FR3.IMGT.g.t",
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
67 "CDR1.IMGT.a.t",
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
68 "FR1.IMGT.a.g",
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
69 "FR3.IMGT.a.t",
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
70 "FR3.IMGT.Nb.of.nucleotides",
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
71 "FR2.IMGT.t.c",
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
72 "CDR2.IMGT.g.a",
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
73 "FR2.IMGT.t.a",
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
74 "CDR1.IMGT.t.a",
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
75 "FR2.IMGT.t.g",
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
76 "FR3.IMGT.t.g",
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
77 "FR2.IMGT.Nb.of.nucleotides",
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
78 "FR1.IMGT.t.a",
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
79 "FR1.IMGT.t.g",
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
80 "FR3.IMGT.c.t",
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
81 "FR1.IMGT.t.c",
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
82 "CDR2.IMGT.a.t",
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
83 "FR2.IMGT.c.t",
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
84 "CDR1.IMGT.g.t",
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
85 "CDR2.IMGT.t.g",
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
86 "FR1.IMGT.Nb.of.nucleotides",
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
87 "CDR1.IMGT.c.g",
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
88 "CDR2.IMGT.t.c",
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
89 "FR3.IMGT.g.a",
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
90 "CDR1.IMGT.a.c",
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
91 "FR2.IMGT.c.a",
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
92 "FR3.IMGT.Nb.of.mutations",
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
93 "FR2.IMGT.c.g",
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
94 "CDR2.IMGT.g.c",
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
95 "FR1.IMGT.g.c",
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
96 "CDR2.IMGT.c.t",
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
97 "FR3.IMGT.c.a",
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
98 "CDR1.IMGT.c.a",
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
99 "CDR2.IMGT.c.g",
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
100 "CDR2.IMGT.c.a",
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
101 "FR1.IMGT.c.t",
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
102 "FR1.IMGT.Nb.of.silent.mutations",
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
103 "FR2.IMGT.Nb.of.silent.mutations",
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
104 "FR3.IMGT.Nb.of.silent.mutations",
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
105 "FR1.IMGT.Nb.of.nonsilent.mutations",
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
106 "FR2.IMGT.Nb.of.nonsilent.mutations",
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
107 "FR3.IMGT.Nb.of.nonsilent.mutations")
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
108
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
109 print("Cleaning up columns")
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
110
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
111 for(col in cleanup_columns){
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
112 dat[,col] = gsub("\\(.*\\)", "", dat[,col])
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
113 #dat[dat[,col] == "",] = "0"
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
114 dat[,col] = as.numeric(dat[,col])
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
115 dat[is.na(dat[,col]),col] = 0
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
116 }
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
117
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
118 regions = c("FR1", "CDR1", "FR2", "CDR2", "FR3")
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
119 if(empty.region.filter == "FR1") {
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
120 regions = c("CDR1", "FR2", "CDR2", "FR3")
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
121 } else if (empty.region.filter == "CDR1") {
7
ad9be244b104 Uploaded
davidvanzessen
parents: 5
diff changeset
122 regions = c("FR2", "CDR2", "FR3")
1
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
123 } else if (empty.region.filter == "FR2") {
7
ad9be244b104 Uploaded
davidvanzessen
parents: 5
diff changeset
124 regions = c("CDR2", "FR3")
1
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
125 }
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
126
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
127 sum_by_row = function(x, columns) { sum(as.numeric(x[columns]), na.rm=T) }
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
128
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
129 print("aggregating data into new columns")
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
130
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
131 VRegionMutations_columns = paste(regions, ".IMGT.Nb.of.mutations", sep="")
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
132 dat$VRegionMutations = apply(dat, FUN=sum_by_row, 1, columns=VRegionMutations_columns)
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
133
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
134 VRegionNucleotides_columns = paste(regions, ".IMGT.Nb.of.nucleotides", sep="")
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
135 dat$FR3.IMGT.Nb.of.nucleotides = nchar(dat$FR3.IMGT.seq)
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
136 dat$VRegionNucleotides = apply(dat, FUN=sum_by_row, 1, columns=VRegionNucleotides_columns)
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
137
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
138 transitionMutations_columns = paste(rep(regions, each=4), c(".IMGT.a.g", ".IMGT.g.a", ".IMGT.c.t", ".IMGT.t.c"), sep="")
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
139 dat$transitionMutations = apply(dat, FUN=sum_by_row, 1, columns=transitionMutations_columns)
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
140
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
141 transversionMutations_columns = paste(rep(regions, each=8), c(".IMGT.a.c",".IMGT.c.a",".IMGT.a.t",".IMGT.t.a",".IMGT.g.c",".IMGT.c.g",".IMGT.g.t",".IMGT.t.g"), sep="")
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
142 dat$transversionMutations = apply(dat, FUN=sum_by_row, 1, columns=transversionMutations_columns)
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
143
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
144 transitionMutationsAtGC_columns = paste(rep(regions, each=2), c(".IMGT.g.a",".IMGT.c.t"), sep="")
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
145 dat$transitionMutationsAtGC = apply(dat, FUN=sum_by_row, 1, columns=transitionMutationsAtGC_columns)
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
146
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
147 totalMutationsAtGC_columns = paste(rep(regions, each=6), c(".IMGT.c.g",".IMGT.c.t",".IMGT.c.a",".IMGT.g.c",".IMGT.g.a",".IMGT.g.t"), sep="")
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
148 #totalMutationsAtGC_columns = paste(rep(regions, each=6), c(".IMGT.g.a",".IMGT.c.t",".IMGT.c.a",".IMGT.c.g",".IMGT.g.t"), sep="")
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
149 dat$totalMutationsAtGC = apply(dat, FUN=sum_by_row, 1, columns=totalMutationsAtGC_columns)
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
150
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
151 transitionMutationsAtAT_columns = paste(rep(regions, each=2), c(".IMGT.a.g",".IMGT.t.c"), sep="")
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
152 dat$transitionMutationsAtAT = apply(dat, FUN=sum_by_row, 1, columns=transitionMutationsAtAT_columns)
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
153
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
154 totalMutationsAtAT_columns = paste(rep(regions, each=6), c(".IMGT.a.g",".IMGT.a.c",".IMGT.a.t",".IMGT.t.g",".IMGT.t.c",".IMGT.t.a"), sep="")
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
155 #totalMutationsAtAT_columns = paste(rep(regions, each=5), c(".IMGT.a.g",".IMGT.t.c",".IMGT.a.c",".IMGT.g.c",".IMGT.t.g"), sep="")
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
156 dat$totalMutationsAtAT = apply(dat, FUN=sum_by_row, 1, columns=totalMutationsAtAT_columns)
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
157
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
158 FRRegions = regions[grepl("FR", regions)]
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
159 CDRRegions = regions[grepl("CDR", regions)]
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
160
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
161 FR_silentMutations_columns = paste(FRRegions, ".IMGT.Nb.of.silent.mutations", sep="")
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
162 dat$silentMutationsFR = apply(dat, FUN=sum_by_row, 1, columns=FR_silentMutations_columns)
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
163
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
164 CDR_silentMutations_columns = paste(CDRRegions, ".IMGT.Nb.of.silent.mutations", sep="")
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
165 dat$silentMutationsCDR = apply(dat, FUN=sum_by_row, 1, columns=CDR_silentMutations_columns)
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
166
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
167 FR_nonSilentMutations_columns = paste(FRRegions, ".IMGT.Nb.of.nonsilent.mutations", sep="")
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
168 dat$nonSilentMutationsFR = apply(dat, FUN=sum_by_row, 1, columns=FR_nonSilentMutations_columns)
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
169
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
170 CDR_nonSilentMutations_columns = paste(CDRRegions, ".IMGT.Nb.of.nonsilent.mutations", sep="")
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
171 dat$nonSilentMutationsCDR = apply(dat, FUN=sum_by_row, 1, columns=CDR_nonSilentMutations_columns)
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
172
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
173 mutation.sum.columns = c("Sequence.ID", "VRegionMutations", "VRegionNucleotides", "transitionMutations", "transversionMutations", "transitionMutationsAtGC", "transitionMutationsAtAT", "silentMutationsFR", "nonSilentMutationsFR", "silentMutationsCDR", "nonSilentMutationsCDR")
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
174 write.table(dat[,mutation.sum.columns], "mutation_by_id.txt", sep="\t",quote=F,row.names=F,col.names=T)
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
175
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
176 setwd(outputdir)
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
177
8
3968d04b5724 Uploaded
davidvanzessen
parents: 7
diff changeset
178 write.table(dat, input, sep="\t",quote=F,row.names=F,col.names=T)
3968d04b5724 Uploaded
davidvanzessen
parents: 7
diff changeset
179
1
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
180 base.order = data.frame(base=c("A", "T", "C", "G"), order=1:4)
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
181
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
182 calculate_result = function(i, gene, dat, matrx, f, fname, name){
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
183 tmp = dat[grepl(paste("^", gene, ".*", sep=""), dat$best_match),]
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
184
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
185 j = i - 1
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
186 x = (j * 3) + 1
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
187 y = (j * 3) + 2
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
188 z = (j * 3) + 3
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
189
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
190 if(nrow(tmp) > 0){
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
191 if(fname == "sum"){
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
192 matrx[1,x] = round(f(tmp$VRegionMutations, na.rm=T), digits=1)
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
193 matrx[1,y] = round(f(tmp$VRegionNucleotides, na.rm=T), digits=1)
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
194 matrx[1,z] = round(f(matrx[1,x] / matrx[1,y]) * 100, digits=1)
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
195 } else {
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
196 matrx[1,x] = round(f(tmp$VRegionMutations, na.rm=T), digits=1)
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
197 matrx[1,y] = round(f(tmp$VRegionNucleotides, na.rm=T), digits=1)
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
198 matrx[1,z] = round(f(tmp$VRegionMutations / tmp$VRegionNucleotides) * 100, digits=1)
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
199 }
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
200
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
201 matrx[2,x] = round(f(tmp$transitionMutations, na.rm=T), digits=1)
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
202 matrx[2,y] = round(f(tmp$VRegionMutations, na.rm=T), digits=1)
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
203 matrx[2,z] = round(matrx[2,x] / matrx[2,y] * 100, digits=1)
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
204
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
205 matrx[3,x] = round(f(tmp$transversionMutations, na.rm=T), digits=1)
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
206 matrx[3,y] = round(f(tmp$VRegionMutations, na.rm=T), digits=1)
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
207 matrx[3,z] = round(matrx[3,x] / matrx[3,y] * 100, digits=1)
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
208
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
209 matrx[4,x] = round(f(tmp$transitionMutationsAtGC, na.rm=T), digits=1)
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
210 matrx[4,y] = round(f(tmp$totalMutationsAtGC, na.rm=T), digits=1)
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
211 matrx[4,z] = round(matrx[4,x] / matrx[4,y] * 100, digits=1)
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
212
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
213 matrx[5,x] = round(f(tmp$totalMutationsAtGC, na.rm=T), digits=1)
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
214 matrx[5,y] = round(f(tmp$VRegionMutations, na.rm=T), digits=1)
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
215 matrx[5,z] = round(matrx[5,x] / matrx[5,y] * 100, digits=1)
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
216
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
217 matrx[6,x] = round(f(tmp$transitionMutationsAtAT, na.rm=T), digits=1)
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
218 matrx[6,y] = round(f(tmp$totalMutationsAtAT, na.rm=T), digits=1)
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
219 matrx[6,z] = round(matrx[6,x] / matrx[6,y] * 100, digits=1)
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
220
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
221 matrx[7,x] = round(f(tmp$totalMutationsAtAT, na.rm=T), digits=1)
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
222 matrx[7,y] = round(f(tmp$VRegionMutations, na.rm=T), digits=1)
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
223 matrx[7,z] = round(matrx[7,x] / matrx[7,y] * 100, digits=1)
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
224
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
225 matrx[8,x] = round(f(tmp$nonSilentMutationsFR, na.rm=T), digits=1)
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
226 matrx[8,y] = round(f(tmp$silentMutationsFR, na.rm=T), digits=1)
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
227 matrx[8,z] = round(matrx[8,x] / matrx[8,y], digits=1)
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
228
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
229 matrx[9,x] = round(f(tmp$nonSilentMutationsCDR, na.rm=T), digits=1)
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
230 matrx[9,y] = round(f(tmp$silentMutationsCDR, na.rm=T), digits=1)
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
231 matrx[9,z] = round(matrx[9,x] / matrx[9,y], digits=1)
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
232
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
233 if(fname == "sum"){
7
ad9be244b104 Uploaded
davidvanzessen
parents: 5
diff changeset
234
ad9be244b104 Uploaded
davidvanzessen
parents: 5
diff changeset
235 regions.fr = regions[grepl("FR", regions)]
ad9be244b104 Uploaded
davidvanzessen
parents: 5
diff changeset
236 regions.fr = paste(regions.fr, ".IMGT.Nb.of.nucleotides", sep="")
ad9be244b104 Uploaded
davidvanzessen
parents: 5
diff changeset
237 regions.cdr = regions[grepl("CDR", regions)]
ad9be244b104 Uploaded
davidvanzessen
parents: 5
diff changeset
238 regions.cdr = paste(regions.cdr, ".IMGT.Nb.of.nucleotides", sep="")
ad9be244b104 Uploaded
davidvanzessen
parents: 5
diff changeset
239
ad9be244b104 Uploaded
davidvanzessen
parents: 5
diff changeset
240 if(length(regions.fr) > 1){ #in case there is only on FR region (rowSums needs >1 column)
ad9be244b104 Uploaded
davidvanzessen
parents: 5
diff changeset
241 matrx[10,x] = round(f(rowSums(tmp[,regions.fr], na.rm=T)), digits=1)
ad9be244b104 Uploaded
davidvanzessen
parents: 5
diff changeset
242 } else {
ad9be244b104 Uploaded
davidvanzessen
parents: 5
diff changeset
243 matrx[10,x] = round(f(tmp[,regions.fr], na.rm=T), digits=1)
ad9be244b104 Uploaded
davidvanzessen
parents: 5
diff changeset
244 }
1
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
245 matrx[10,y] = round(f(tmp$VRegionNucleotides, na.rm=T), digits=1)
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
246 matrx[10,z] = round(matrx[10,x] / matrx[10,y] * 100, digits=1)
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
247
7
ad9be244b104 Uploaded
davidvanzessen
parents: 5
diff changeset
248 if(length(regions.cdr) > 1){ #in case there is only on CDR region
ad9be244b104 Uploaded
davidvanzessen
parents: 5
diff changeset
249 matrx[11,x] = round(f(rowSums(tmp[,regions.cdr], na.rm=T)), digits=1)
ad9be244b104 Uploaded
davidvanzessen
parents: 5
diff changeset
250 } else {
ad9be244b104 Uploaded
davidvanzessen
parents: 5
diff changeset
251 matrx[11,x] = round(f(tmp[,regions.cdr], na.rm=T), digits=1)
ad9be244b104 Uploaded
davidvanzessen
parents: 5
diff changeset
252 }
1
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
253 matrx[11,y] = round(f(tmp$VRegionNucleotides, na.rm=T), digits=1)
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
254 matrx[11,z] = round(matrx[11,x] / matrx[11,y] * 100, digits=1)
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
255 }
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
256 }
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
257
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
258 transitionTable = data.frame(A=zeros,C=zeros,G=zeros,T=zeros)
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
259 row.names(transitionTable) = c("A", "C", "G", "T")
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
260 transitionTable["A","A"] = NA
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
261 transitionTable["C","C"] = NA
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
262 transitionTable["G","G"] = NA
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
263 transitionTable["T","T"] = NA
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
264
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
265 if(nrow(tmp) > 0){
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
266 for(nt1 in nts){
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
267 for(nt2 in nts){
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
268 if(nt1 == nt2){
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
269 next
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
270 }
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
271 NT1 = LETTERS[letters == nt1]
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
272 NT2 = LETTERS[letters == nt2]
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
273 FR1 = paste("FR1.IMGT.", nt1, ".", nt2, sep="")
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
274 CDR1 = paste("CDR1.IMGT.", nt1, ".", nt2, sep="")
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
275 FR2 = paste("FR2.IMGT.", nt1, ".", nt2, sep="")
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
276 CDR2 = paste("CDR2.IMGT.", nt1, ".", nt2, sep="")
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
277 FR3 = paste("FR3.IMGT.", nt1, ".", nt2, sep="")
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
278 if (empty.region.filter == "leader"){
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
279 transitionTable[NT1,NT2] = sum(tmp[,c(FR1, CDR1, FR2, CDR2, FR3)])
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
280 } else if (empty.region.filter == "FR1") {
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
281 transitionTable[NT1,NT2] = sum(tmp[,c(CDR1, FR2, CDR2, FR3)])
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
282 } else if (empty.region.filter == "CDR1") {
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
283 transitionTable[NT1,NT2] = sum(tmp[,c(FR2, CDR2, FR3)])
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
284 } else if (empty.region.filter == "FR2") {
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
285 transitionTable[NT1,NT2] = sum(tmp[,c(CDR2, FR3)])
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
286 }
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
287 }
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
288 }
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
289 transition = transitionTable
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
290 transition$id = names(transition)
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
291
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
292 transition2 = melt(transition, id.vars="id")
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
293
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
294 transition2 = merge(transition2, base.order, by.x="id", by.y="base")
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
295
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
296 transition2 = merge(transition2, base.order, by.x="variable", by.y="base")
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
297
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
298 transition2[is.na(transition2$value),]$value = 0
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
299
4
477e95b098fd Uploaded
davidvanzessen
parents: 3
diff changeset
300 if(any(transition2$value != 0)){ #having rows of data but a transition table filled with 0 is bad
1
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
301 print("Plotting stacked transition")
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
302 png(filename=paste("transitions_stacked_", name, ".png", sep=""))
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
303 p = ggplot(transition2, aes(factor(reorder(id, order.x)), y=value, fill=factor(reorder(variable, order.y)))) + geom_bar(position="fill", stat="identity", colour="black") #stacked bar
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
304 p = p + xlab("From base") + ylab("To base") + ggtitle("Mutations frequency from base to base") + guides(fill=guide_legend(title=NULL))
5
012a738edf5a Uploaded
davidvanzessen
parents: 4
diff changeset
305 p = p + theme(panel.background = element_rect(fill = "white", colour="black"), text = element_text(size=16, colour="black")) + scale_fill_manual(values=c("A" = "blue4", "G" = "lightblue1", "C" = "olivedrab3", "T" = "olivedrab4"))
1
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
306 #p = p + scale_colour_manual(values=c("A" = "black", "G" = "black", "C" = "black", "T" = "black"))
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
307 print(p)
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
308 dev.off()
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
309
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
310 print("Plotting heatmap transition")
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
311
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
312 png(filename=paste("transitions_heatmap_", name, ".png", sep=""))
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
313 p = ggplot(transition2, aes(factor(reorder(id, order.x)), factor(reorder(variable, order.y)))) + geom_tile(aes(fill = value)) + scale_fill_gradient(low="white", high="steelblue") #heatmap
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
314 p = p + xlab("From base") + ylab("To base") + ggtitle("Mutations frequency from base to base") + theme(panel.background = element_rect(fill = "white", colour="black"), text = element_text(size=13, colour="black"))
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
315 print(p)
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
316 dev.off()
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
317 } else {
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
318 print("No data to plot")
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
319 }
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
320 }
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
321
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
322 #print(paste("writing value file: ", name, "_", fname, "_value.txt" ,sep=""))
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
323 write.table(x=transitionTable, file=paste("transitions_", name ,"_", fname, ".txt", sep=""), sep=",",quote=F,row.names=T,col.names=NA)
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
324 write.table(x=tmp[,c("Sequence.ID", "best_match", "chunk_hit_percentage", "nt_hit_percentage", "start_locations")], file=paste("matched_", name , "_", fname, ".txt", sep=""), sep="\t",quote=F,row.names=F,col.names=T)
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
325 cat(matrx[1,x], file=paste(name, "_", fname, "_value.txt" ,sep=""))
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
326 cat(nrow(tmp), file=paste(name, "_", fname, "_n.txt" ,sep=""))
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
327 #print(paste(fname, name, nrow(tmp)))
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
328 matrx
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
329 }
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
330 nts = c("a", "c", "g", "t")
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
331 zeros=rep(0, 4)
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
332 funcs = c(median, sum, mean)
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
333 fnames = c("median", "sum", "mean")
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
334
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
335 print("Creating result tables")
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
336
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
337 for(i in 1:length(funcs)){
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
338 func = funcs[[i]]
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
339 fname = fnames[[i]]
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
340
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
341 rows = 9
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
342 if(fname == "sum"){
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
343 rows = 11
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
344 }
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
345 matrx = matrix(data = 0, ncol=((length(genes) + 1) * 3),nrow=rows)
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
346 for(i in 1:length(genes)){
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
347 print(paste("Creating table for", fname, genes[i]))
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
348 matrx = calculate_result(i, genes[i], dat, matrx, func, fname, genes[i])
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
349 }
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
350 matrx = calculate_result(i + 1, ".*", dat[!grepl("unmatched", dat$best_match),], matrx, func, fname, name="all")
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
351
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
352 result = data.frame(matrx)
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
353 if(fname == "sum"){
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
354 row.names(result) = c("Number of Mutations (%)", "Transitions (%)", "Transversions (%)", "Transitions at G C (%)", "Targeting of C G (%)", "Transitions at A T (%)", "Targeting of A T (%)", "FR R/S (ratio)", "CDR R/S (ratio)", "nt in FR", "nt in CDR")
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
355 } else {
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
356 row.names(result) = c("Number of Mutations (%)", "Transitions (%)", "Transversions (%)", "Transitions at G C (%)", "Targeting of C G (%)", "Transitions at A T (%)", "Targeting of A T (%)", "FR R/S (ratio)", "CDR R/S (ratio)")
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
357 }
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
358 write.table(x=result, file=paste("mutations_", fname, ".txt", sep=""), sep=",",quote=F,row.names=T,col.names=F)
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
359 }
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
360
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
361 print("Adding median number of mutations to sum table")
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
362 sum.table = read.table("mutations_sum.txt", sep=",", header=F)
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
363 median.table = read.table("mutations_median.txt", sep=",", header=F)
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
364
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
365 new.table = sum.table[1,]
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
366 new.table[2,] = median.table[1,]
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
367 new.table[3:12,] = sum.table[2:11,]
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
368 new.table[,1] = as.character(new.table[,1])
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
369 new.table[2,1] = "Median of Number of Mutations (%)"
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
370
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
371 #sum.table = sum.table[c("Number of Mutations (%)", "Median of Number of Mutations (%)", "Transition (%)", "Transversions (%)", "Transitions at G C (%)", "Targeting of C G (%)", "Transitions at A T (%)", "Targeting of A T (%)", "FR R/S (ratio)", "CDR R/S (ratio)", "nt in FR", "nt in CDR"),]
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
372
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
373 write.table(x=new.table, file="mutations_sum.txt", sep=",",quote=F,row.names=F,col.names=F)
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
374
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
375 print("Plotting IGA piechart")
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
376
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
377 dat = dat[!grepl("^unmatched", dat$best_match),]
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
378
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
379 #blegh
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
380
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
381 genesForPlot = dat[grepl("IGA", dat$best_match),]$best_match
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
382
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
383 if(length(genesForPlot) > 0){
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
384 genesForPlot = data.frame(table(genesForPlot))
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
385 colnames(genesForPlot) = c("Gene","Freq")
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
386 genesForPlot$label = paste(genesForPlot$Gene, "-", genesForPlot$Freq)
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
387
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
388 pc = ggplot(genesForPlot, aes(x = factor(1), y=Freq, fill=Gene))
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
389 pc = pc + geom_bar(width = 1, stat = "identity") + scale_fill_manual(labels=genesForPlot$label, values=c("IGA1" = "lightblue1", "IGA2" = "blue4"))
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
390 pc = pc + coord_polar(theta="y") + scale_y_continuous(breaks=NULL)
5
012a738edf5a Uploaded
davidvanzessen
parents: 4
diff changeset
391 pc = pc + theme(panel.background = element_rect(fill = "white", colour="black"), text = element_text(size=16, colour="black"), axis.title=element_blank(), axis.text=element_blank(), axis.ticks=element_blank())
1
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
392 pc = pc + xlab(" ") + ylab(" ") + ggtitle(paste("IGA subclasses", "( n =", sum(genesForPlot$Freq), ")"))
2
e85fec274cde Uploaded
davidvanzessen
parents: 1
diff changeset
393 write.table(genesForPlot, "IGA_pie.txt", sep="\t",quote=F,row.names=F,col.names=T)
1
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
394
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
395 png(filename="IGA.png")
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
396 print(pc)
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
397 dev.off()
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
398 }
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
399
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
400 print("Plotting IGG piechart")
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
401
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
402 genesForPlot = dat[grepl("IGG", dat$best_match),]$best_match
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
403
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
404 if(length(genesForPlot) > 0){
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
405 genesForPlot = data.frame(table(genesForPlot))
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
406 colnames(genesForPlot) = c("Gene","Freq")
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
407 genesForPlot$label = paste(genesForPlot$Gene, "-", genesForPlot$Freq)
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
408
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
409 pc = ggplot(genesForPlot, aes(x = factor(1), y=Freq, fill=Gene))
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
410 pc = pc + geom_bar(width = 1, stat = "identity") + scale_fill_manual(labels=genesForPlot$label, values=c("IGG1" = "olivedrab3", "IGG2" = "red", "IGG3" = "gold", "IGG4" = "darkred"))
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
411 pc = pc + coord_polar(theta="y") + scale_y_continuous(breaks=NULL)
5
012a738edf5a Uploaded
davidvanzessen
parents: 4
diff changeset
412 pc = pc + theme(panel.background = element_rect(fill = "white", colour="black"), text = element_text(size=16, colour="black"), axis.title=element_blank(), axis.text=element_blank(), axis.ticks=element_blank())
1
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
413 pc = pc + xlab(" ") + ylab(" ") + ggtitle(paste("IGG subclasses", "( n =", sum(genesForPlot$Freq), ")"))
2
e85fec274cde Uploaded
davidvanzessen
parents: 1
diff changeset
414 write.table(genesForPlot, "IGG_pie.txt", sep="\t",quote=F,row.names=F,col.names=T)
1
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
415
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
416 png(filename="IGG.png")
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
417 print(pc)
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
418 dev.off()
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
419 }
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
420
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
421 print("Plotting scatterplot")
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
422
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
423 dat$percentage_mutations = round(dat$VRegionMutations / dat$VRegionNucleotides * 100, 2)
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
424 dat.clss = dat
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
425
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
426 dat.clss$best_match = substr(dat.clss$best_match, 0, 3)
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
427
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
428 dat.clss = rbind(dat, dat.clss)
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
429
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
430 p = ggplot(dat.clss, aes(best_match, percentage_mutations))
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
431 p = p + geom_point(aes(colour=best_match), position="jitter") + geom_boxplot(aes(middle=mean(percentage_mutations)), alpha=0.1, outlier.shape = NA)
5
012a738edf5a Uploaded
davidvanzessen
parents: 4
diff changeset
432 p = p + xlab("Subclass") + ylab("Frequency") + ggtitle("Frequency scatter plot") + theme(panel.background = element_rect(fill = "white", colour="black"), text = element_text(size=16, colour="black"))
9
372ccdcf0b2d Uploaded
davidvanzessen
parents: 8
diff changeset
433 p = p + scale_fill_manual(values=c("IGA" = "blue4", "IGA1" = "lightblue1", "IGA2" = "blue4", "IGG" = "olivedrab3", "IGG1" = "olivedrab3", "IGG2" = "red", "IGG3" = "gold", "IGG4" = "darkred", "IGM" = "darkviolet", "IGE" = "darkorange", "all" = "blue4"))
372ccdcf0b2d Uploaded
davidvanzessen
parents: 8
diff changeset
434 p = p + scale_colour_manual(values=c("IGA" = "blue4", "IGA1" = "lightblue1", "IGA2" = "blue4", "IGG" = "olivedrab3", "IGG1" = "olivedrab3", "IGG2" = "red", "IGG3" = "gold", "IGG4" = "darkred", "IGM" = "darkviolet", "IGE" = "darkorange", "all" = "blue4"))
1
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
435
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
436 png(filename="scatter.png")
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
437 print(p)
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
438 dev.off()
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
439
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
440 write.table(dat[,c("Sequence.ID", "best_match", "VRegionMutations", "VRegionNucleotides", "percentage_mutations")], "scatter.txt", sep="\t",quote=F,row.names=F,col.names=T)
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
441
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
442 print("Plotting frequency ranges plot")
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
443
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
444 dat$best_match_class = substr(dat$best_match, 0, 3)
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
445 freq_labels = c("0", "0-2", "2-5", "5-10", "10-15", "15-20", "20")
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
446 dat$frequency_bins = cut(dat$percentage_mutations, breaks=c(-Inf, 0, 2,5,10,15,20, Inf), labels=freq_labels)
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
447
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
448 frequency_bins_sum = data.frame(data.table(dat)[, list(class_sum=sum(.N)), by=c("best_match_class")])
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
449
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
450 frequency_bins_data = data.frame(data.table(dat)[, list(frequency_count=.N), by=c("best_match_class", "frequency_bins")])
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
451
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
452 frequency_bins_data = merge(frequency_bins_data, frequency_bins_sum, by="best_match_class")
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
453
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
454 frequency_bins_data$frequency = round(frequency_bins_data$frequency_count / frequency_bins_data$class_sum * 100, 2)
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
455
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
456 p = ggplot(frequency_bins_data, aes(frequency_bins, frequency))
5
012a738edf5a Uploaded
davidvanzessen
parents: 4
diff changeset
457 p = p + geom_bar(aes(fill=best_match_class), stat="identity", position="dodge") + theme(panel.background = element_rect(fill = "white", colour="black"), text = element_text(size=16, colour="black"))
9
372ccdcf0b2d Uploaded
davidvanzessen
parents: 8
diff changeset
458 p = p + xlab("Frequency ranges") + ylab("Frequency") + ggtitle("Mutation Frequencies by class") + scale_fill_manual(values=c("IGA" = "blue4", "IGG" = "olivedrab3", "IGM" = "darkviolet", "IGE" = "darkorange", "all" = "blue4"))
1
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
459
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
460 png(filename="frequency_ranges.png")
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
461 print(p)
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
462 dev.off()
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
463
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
464 frequency_bins_data_by_class = frequency_bins_data
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
465
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
466 write.table(frequency_bins_data_by_class, "frequency_ranges_classes.txt", sep="\t",quote=F,row.names=F,col.names=T)
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
467
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
468 frequency_bins_data = data.frame(data.table(dat)[, list(frequency_count=.N), by=c("best_match", "best_match_class", "frequency_bins")])
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
469
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
470 frequency_bins_data = merge(frequency_bins_data, frequency_bins_sum, by="best_match_class")
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
471
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
472 frequency_bins_data$frequency = round(frequency_bins_data$frequency_count / frequency_bins_data$class_sum * 100, 2)
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
473
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
474 write.table(frequency_bins_data, "frequency_ranges_subclasses.txt", sep="\t",quote=F,row.names=F,col.names=T)
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
475
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
476
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
477
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
478
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
479
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
480
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
481
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
482
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
483
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
484
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
485
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
486
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
487
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
488
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
489
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
490
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
491
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
492
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
493
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
494
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
495
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
496
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
497
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
498
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
499
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
500
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
501
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
502
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
503
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
504
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
505
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
506
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
507
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
508
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
509
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
510
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
511
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
512
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
513
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
514
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
515
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
516
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
517
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
518
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
519
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
520
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
521
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
522
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
523
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
524
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
525
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
526
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
527
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
528
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
529
faae21ba5c63 Uploaded
davidvanzessen
parents: 0
diff changeset
530