|
0
|
1 args <- commandArgs(trailingOnly = TRUE)
|
|
|
2
|
|
|
3 summ.file = args[1]
|
|
20
|
4 sequences.file = args[2]
|
|
|
5 aa.file = args[3]
|
|
|
6 junction.file = args[4]
|
|
28
|
7 gapped.aa.file = args[5]
|
|
|
8 out.file = args[6]
|
|
0
|
9
|
|
|
10 summ = read.table(summ.file, sep="\t", header=T, quote="", fill=T)
|
|
20
|
11 sequences = read.table(sequences.file, sep="\t", header=T, quote="", fill=T)
|
|
0
|
12 aa = read.table(aa.file, sep="\t", header=T, quote="", fill=T)
|
|
28
|
13 gapped.aa = read.table(gapped.aa.file, sep="\t", header=T, quote="", fill=T)
|
|
0
|
14 junction = read.table(junction.file, sep="\t", header=T, quote="", fill=T)
|
|
|
15
|
|
39
|
16 fix_column_names = function(df){
|
|
|
17 if("V.DOMAIN.Functionality" %in% names(df)){
|
|
|
18 names(df)[names(df) == "V.DOMAIN.Functionality"] = "Functionality"
|
|
|
19 print("found V.DOMAIN.Functionality, changed")
|
|
|
20 }
|
|
|
21 if("V.DOMAIN.Functionality.comment" %in% names(df)){
|
|
|
22 names(df)[names(df) == "V.DOMAIN.Functionality.comment"] = "Functionality.comment"
|
|
|
23 print("found V.DOMAIN.Functionality.comment, changed")
|
|
|
24 }
|
|
|
25 return(df)
|
|
|
26 }
|
|
|
27
|
|
|
28 summ = fix_column_names(summ)
|
|
|
29 sequences = fix_column_names(sequences)
|
|
|
30 aa = fix_column_names(aa)
|
|
|
31 gapped.aa = fix_column_names(gapped.aa)
|
|
|
32 junction = fix_column_names(junction)
|
|
|
33
|
|
0
|
34 old_summary_columns=c('Sequence.ID','JUNCTION.frame','V.GENE.and.allele','D.GENE.and.allele','J.GENE.and.allele','CDR1.IMGT.length','CDR2.IMGT.length','CDR3.IMGT.length','Orientation')
|
|
|
35 old_sequence_columns=c('CDR1.IMGT','CDR2.IMGT','CDR3.IMGT')
|
|
|
36 old_junction_columns=c('JUNCTION')
|
|
|
37
|
|
|
38 added_summary_columns=c('Functionality','V.REGION.identity..','V.REGION.identity.nt','D.REGION.reading.frame','AA.JUNCTION','Functionality.comment','Sequence')
|
|
|
39 added_sequence_columns=c('FR1.IMGT','FR2.IMGT','FR3.IMGT','CDR3.IMGT','JUNCTION','J.REGION','FR4.IMGT')
|
|
|
40
|
|
|
41 added_junction_columns=c('P3.V.nt.nb','N.REGION.nt.nb','N1.REGION.nt.nb','P5.D.nt.nb','P3.D.nt.nb','N2.REGION.nt.nb','P5.J.nt.nb','X3.V.REGION.trimmed.nt.nb','X5.D.REGION.trimmed.nt.nb','X3.D.REGION.trimmed.nt.nb','X5.J.REGION.trimmed.nt.nb','N.REGION','N1.REGION','N2.REGION')
|
|
|
42 added_junction_columns=c(added_junction_columns, 'P5.D1.nt.nb', 'P3.D1.nt.nb', 'N2.REGION.nt.nb', 'P5.D2.nt.nb', 'P3.D2.nt.nb', 'N3.REGION.nt.nb', 'P5.D3.nt.nb', 'P3.D2.nt.nb', 'N4.REGION.nt.nb', 'X5.D1.REGION.trimmed.nt.nb', 'X3.D1.REGION.trimmed.nt.nb', 'X5.D2.REGION.trimmed.nt.nb', 'X3.D2.REGION.trimmed.nt.nb', 'X5.D3.REGION.trimmed.nt.nb', 'X3.D3.REGION.trimmed.nt.nb', 'D.REGION.nt.nb', 'D1.REGION.nt.nb', 'D2.REGION.nt.nb', 'D3.REGION.nt.nb')
|
|
|
43
|
|
|
44 out=summ[,c("Sequence.ID","JUNCTION.frame","V.GENE.and.allele","D.GENE.and.allele","J.GENE.and.allele")]
|
|
|
45
|
|
|
46 out[,"CDR1.Seq"] = aa[,"CDR1.IMGT"]
|
|
|
47 out[,"CDR1.Length"] = summ[,"CDR1.IMGT.length"]
|
|
|
48
|
|
|
49 out[,"CDR2.Seq"] = aa[,"CDR2.IMGT"]
|
|
|
50 out[,"CDR2.Length"] = summ[,"CDR2.IMGT.length"]
|
|
|
51
|
|
31
|
52 out[,"CDR3.Seq"] = gapped.aa[,"CDR3.IMGT"]
|
|
0
|
53 out[,"CDR3.Length"] = summ[,"CDR3.IMGT.length"]
|
|
|
54
|
|
31
|
55 out[,"CDR3.IMGT"] = out[,"CDR3.Seq"]
|
|
|
56
|
|
|
57 out[,"CDR3.Seq.DNA"] = sequences[,"CDR3.IMGT"]
|
|
28
|
58 out[,"CDR3.Length.DNA"] = nchar(as.character(out[,"CDR3.Seq.DNA"]))
|
|
0
|
59 out[,"Strand"] = summ[,"Orientation"]
|
|
|
60 out[,"CDR3.Found.How"] = "a"
|
|
|
61
|
|
|
62 out[,added_summary_columns] = summ[,added_summary_columns]
|
|
|
63
|
|
|
64 out[,added_sequence_columns] = aa[,added_sequence_columns]
|
|
|
65
|
|
|
66 out[,added_junction_columns] = junction[,added_junction_columns]
|
|
|
67
|
|
|
68 out[,"Top V Gene"] = gsub(".* ", "", gsub("\\*.*", "", summ[,"V.GENE.and.allele"]))
|
|
|
69 out[,"Top D Gene"] = gsub(".* ", "", gsub("\\*.*", "", summ[,"D.GENE.and.allele"]))
|
|
|
70 out[,"Top J Gene"] = gsub(".* ", "", gsub("\\*.*", "", summ[,"J.GENE.and.allele"]))
|
|
|
71
|
|
|
72 out = out[!grepl("Less than", summ[,"V.GENE.and.allele"]),]
|
|
|
73 out = out[!grepl("Less than", summ[,"D.GENE.and.allele"]),]
|
|
|
74 out = out[!grepl("Less than", summ[,"J.GENE.and.allele"]),]
|
|
|
75
|
|
|
76 out = out[,c('Sequence.ID','JUNCTION.frame','Top V Gene','Top D Gene','Top J Gene','CDR1.Seq','CDR1.Length','CDR2.Seq','CDR2.Length','CDR3.Seq','CDR3.Length','CDR3.Seq.DNA','CDR3.Length.DNA','Strand','CDR3.Found.How','Functionality','V.REGION.identity..','V.REGION.identity.nt','D.REGION.reading.frame','AA.JUNCTION','Functionality.comment','Sequence','FR1.IMGT','FR2.IMGT','FR3.IMGT','CDR3.IMGT','JUNCTION','J.REGION','FR4.IMGT','P3.V.nt.nb','N.REGION.nt.nb','N1.REGION.nt.nb','P5.D.nt.nb','P3.D.nt.nb','N2.REGION.nt.nb','P5.J.nt.nb','X3.V.REGION.trimmed.nt.nb','X5.D.REGION.trimmed.nt.nb','X3.D.REGION.trimmed.nt.nb','X5.J.REGION.trimmed.nt.nb','N.REGION','N1.REGION','N2.REGION', 'P5.D1.nt.nb', 'P3.D1.nt.nb', 'N2.REGION.nt.nb', 'P5.D2.nt.nb', 'P3.D2.nt.nb', 'N3.REGION.nt.nb', 'P5.D3.nt.nb', 'P3.D2.nt.nb', 'N4.REGION.nt.nb', 'X5.D1.REGION.trimmed.nt.nb', 'X3.D1.REGION.trimmed.nt.nb', 'X5.D2.REGION.trimmed.nt.nb', 'X3.D2.REGION.trimmed.nt.nb', 'X5.D3.REGION.trimmed.nt.nb', 'X3.D3.REGION.trimmed.nt.nb', 'D.REGION.nt.nb', 'D1.REGION.nt.nb', 'D2.REGION.nt.nb', 'D3.REGION.nt.nb')]
|
|
|
77
|
|
|
78 names(out) = c('ID','VDJ Frame','Top V Gene','Top D Gene','Top J Gene','CDR1 Seq','CDR1 Length','CDR2 Seq','CDR2 Length','CDR3 Seq','CDR3 Length','CDR3 Seq DNA','CDR3 Length DNA','Strand','CDR3 Found How','Functionality','V-REGION identity %','V-REGION identity nt','D-REGION reading frame','AA JUNCTION','Functionality comment','Sequence','FR1-IMGT','FR2-IMGT','FR3-IMGT','CDR3-IMGT','JUNCTION','J-REGION','FR4-IMGT','P3V-nt nb','N-REGION-nt nb','N1-REGION-nt nb','P5D-nt nb','P3D-nt nb','N2-REGION-nt nb','P5J-nt nb','3V-REGION trimmed-nt nb','5D-REGION trimmed-nt nb','3D-REGION trimmed-nt nb','5J-REGION trimmed-nt nb','N-REGION','N1-REGION','N2-REGION', 'P5.D1.nt.nb', 'P3.D1.nt.nb', 'N2.REGION.nt.nb', 'P5.D2.nt.nb', 'P3.D2.nt.nb', 'N3.REGION.nt.nb', 'P5.D3.nt.nb', 'P3.D2.nt.nb', 'N4.REGION.nt.nb', 'X5.D1.REGION.trimmed.nt.nb', 'X3.D1.REGION.trimmed.nt.nb', 'X5.D2.REGION.trimmed.nt.nb', 'X3.D2.REGION.trimmed.nt.nb', 'X5.D3.REGION.trimmed.nt.nb', 'X3.D3.REGION.trimmed.nt.nb', 'D.REGION.nt.nb', 'D1.REGION.nt.nb', 'D2.REGION.nt.nb', 'D3.REGION.nt.nb')
|
|
|
79
|
|
|
80 out[,"VDJ Frame"] = as.character(out[,"VDJ Frame"])
|
|
|
81
|
|
|
82 fltr = out[,"VDJ Frame"] == "in-frame"
|
|
|
83 if(any(fltr)){
|
|
|
84 out[fltr, "VDJ Frame"] = "In-frame"
|
|
|
85 }
|
|
|
86
|
|
|
87 fltr = out[,"VDJ Frame"] == "null"
|
|
|
88 if(any(fltr)){
|
|
|
89 out[fltr, "VDJ Frame"] = "Out-of-frame"
|
|
|
90 }
|
|
|
91
|
|
|
92 fltr = out[,"VDJ Frame"] == "out-of-frame"
|
|
|
93 if(any(fltr)){
|
|
|
94 out[fltr, "VDJ Frame"] = "Out-of-frame"
|
|
|
95 }
|
|
|
96
|
|
|
97 fltr = out[,"VDJ Frame"] == ""
|
|
|
98 if(any(fltr)){
|
|
|
99 out[fltr, "VDJ Frame"] = "Out-of-frame"
|
|
|
100 }
|
|
|
101
|
|
|
102 for(col in c('Top V Gene','Top D Gene','Top J Gene')){
|
|
|
103 out[,col] = as.character(out[,col])
|
|
|
104 fltr = out[,col] == ""
|
|
|
105 fltr[is.na(fltr)] = T
|
|
|
106 if(any(fltr)){
|
|
|
107 out[fltr,col] = "NA"
|
|
|
108 }
|
|
|
109 }
|
|
|
110
|
|
|
111 write.table(out, out.file, sep="\t", quote=F, row.names=F, col.names=T)
|
