Mercurial > repos > davidvanzessen > imgt_loader

changeset 0:bd23170af386 draft default tip

Find changesets by keywords (author, files, the commit message), revision number or hash, or revset expression.
Uploaded
author davidvanzessen
date Fri, 18 Aug 2017 10:41:51 -0400
parents
children
files imgt_loader.xml imgt_loader/imgt_loader.r imgt_loader/imgt_loader.sh
diffstat 3 files changed, 242 insertions(+), 0 deletions(-) [+]
line wrap: on
line diff
--- /dev/null	Thu Jan 01 00:00:00 1970 +0000
+++ b/imgt_loader.xml	Fri Aug 18 10:41:51 2017 -0400
@@ -0,0 +1,48 @@
+<tool id="imgt_loader_igg" name="IMGT Loader" version="1.0">
+	<description> </description>
+	<command interpreter="bash">
+		imgt_loader/imgt_loader.sh $in_file $out_file "tmp"
+	</command>
+	<inputs>
+        <param name="in_file" type="data" label="Archive with files" />
+	</inputs>
+	<outputs>
+		<data format="tabular" name="out_file" label="IMGT Loader on ${in_file.name}"/>
+	</outputs>
+	<help>
+**INPUT**
+
+This tool accepts an IMGT/HIGHV-QUEST ZIP file
+
+**OUTPUT**
+
+The following data is used for ARGalaxy
+
++-----------------+----------------------------------------------+
+| Column name     | Column contents                              |
++-----------------+----------------------------------------------+
+| ID              | The Sequence ID provided by the sequencer.   |
++-----------------+----------------------------------------------+
+| VDJ Frame       | In-frame/Out-frame                           |
++-----------------+----------------------------------------------+
+| Top V Gene      | The best matching V gene found.              |
++-----------------+----------------------------------------------+
+| Top D Gene      | The best matching D gene found.              |
++-----------------+----------------------------------------------+
+| Top J Gene      | The best matching J gene found.              |
++-----------------+----------------------------------------------+
+| CDR3 Seq        | The CDR3 region.                             |
++-----------------+----------------------------------------------+
+| CDR3 Length     | The length of the CDR3 region.               |
++-----------------+----------------------------------------------+
+| CDR3 Seq DNA    | The CDR3 sequence region.                    |
++-----------------+----------------------------------------------+
+| CDR3 Length DNA | The length of the CDR3 sequence region.      |
++-----------------+----------------------------------------------+
+| Functionality   | If sequence is productive/unproductive       |
++-----------------+----------------------------------------------+
+
+
+	</help>
+
+</tool>
--- /dev/null	Thu Jan 01 00:00:00 1970 +0000
+++ b/imgt_loader/imgt_loader.r	Fri Aug 18 10:41:51 2017 -0400
@@ -0,0 +1,120 @@
+args <- commandArgs(trailingOnly = TRUE)
+options(show.error.locations = TRUE)
+
+summ.file = args[1]
+sequences.file = args[2]
+aa.file = args[3]
+junction.file = args[4]
+gapped.aa.file = args[5]
+out.file = args[6]
+
+print(summ.file)
+print(out.file)
+
+summ = read.table(summ.file, sep="\t", header=T, quote="", fill=T)
+sequences = read.table(sequences.file, sep="\t", header=T, quote="", fill=T)
+aa = read.table(aa.file, sep="\t", header=T, quote="", fill=T)
+gapped.aa = read.table(gapped.aa.file, sep="\t", header=T, quote="", fill=T)
+junction = read.table(junction.file, sep="\t", header=T, quote="", fill=T)
+
+print(paste("nrow(summ)", nrow(summ)))
+
+write.table(summ, out.file, sep="\t", quote=F, row.names=F, col.names=T)
+
+fix_column_names = function(df){
+    if("V.DOMAIN.Functionality" %in% names(df)){
+        names(df)[names(df) == "V.DOMAIN.Functionality"] = "Functionality"
+        print("found V.DOMAIN.Functionality, changed")
+    }
+    if("V.DOMAIN.Functionality.comment" %in% names(df)){
+        names(df)[names(df) == "V.DOMAIN.Functionality.comment"] = "Functionality.comment"
+        print("found V.DOMAIN.Functionality.comment, changed")
+    }
+    return(df)
+}
+
+summ = fix_column_names(summ)
+sequences = fix_column_names(sequences)
+aa = fix_column_names(aa)
+gapped.aa = fix_column_names(gapped.aa)
+junction = fix_column_names(junction)
+
+print(paste("nrow(summ)", nrow(summ)))
+
+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')
+old_sequence_columns=c('CDR1.IMGT','CDR2.IMGT','CDR3.IMGT')
+old_junction_columns=c('JUNCTION')
+
+added_summary_columns=c('Functionality','V.REGION.identity..','V.REGION.identity.nt','D.REGION.reading.frame','AA.JUNCTION','Functionality.comment','Sequence')
+added_sequence_columns=c('FR1.IMGT','FR2.IMGT','FR3.IMGT','CDR3.IMGT','JUNCTION','J.REGION','FR4.IMGT')
+
+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')
+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')
+
+out=summ[,c("Sequence.ID","JUNCTION.frame","V.GENE.and.allele","D.GENE.and.allele","J.GENE.and.allele")]
+
+print(paste("nrow(summ)", nrow(summ)))
+print(paste("nrow(aa)", nrow(aa)))
+
+out[,"CDR1.Seq"] = aa[,"CDR1.IMGT"]
+out[,"CDR1.Length"] = summ[,"CDR1.IMGT.length"]
+
+out[,"CDR2.Seq"] = aa[,"CDR2.IMGT"]
+out[,"CDR2.Length"] = summ[,"CDR2.IMGT.length"]
+
+out[,"CDR3.Seq"] = gapped.aa[,"CDR3.IMGT"]
+out[,"CDR3.Length"] = summ[,"CDR3.IMGT.length"]
+
+out[,"CDR3.IMGT"] = out[,"CDR3.Seq"]
+
+out[,"CDR3.Seq.DNA"] = sequences[,"CDR3.IMGT"]
+out[,"CDR3.Length.DNA"] = nchar(as.character(out[,"CDR3.Seq.DNA"]))
+out[,"Strand"] = summ[,"Orientation"]
+out[,"CDR3.Found.How"] = "a"
+
+out[,added_summary_columns] = summ[,added_summary_columns]
+
+out[,added_sequence_columns] = aa[,added_sequence_columns]
+
+out[,added_junction_columns] = junction[,added_junction_columns]
+
+out[,"Top V Gene"] = gsub(".* ", "", gsub("\\*.*", "", summ[,"V.GENE.and.allele"]))
+out[,"Top D Gene"] = gsub(".* ", "", gsub("\\*.*", "", summ[,"D.GENE.and.allele"]))
+out[,"Top J Gene"] = gsub(".* ", "", gsub("\\*.*", "", summ[,"J.GENE.and.allele"]))
+
+out = out[!grepl("Less than", summ[,"V.GENE.and.allele"]),]
+out = out[!grepl("Less than", summ[,"D.GENE.and.allele"]),]
+out = out[!grepl("Less than", summ[,"J.GENE.and.allele"]),]
+
+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')]
+
+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')
+
+out[,"VDJ Frame"] = as.character(out[,"VDJ Frame"])
+fltr = out[,"VDJ Frame"] == "in-frame"
+if(any(fltr, na.rm=T)){
+	out[fltr, "VDJ Frame"] = "In-frame"
+}
+fltr = out[,"VDJ Frame"] == "null"
+if(any(fltr, na.rm = T)){
+	out[fltr, "VDJ Frame"] = "Out-of-frame"
+}
+fltr = out[,"VDJ Frame"] == "out-of-frame"
+if(any(fltr, na.rm = T)){
+	out[fltr, "VDJ Frame"] = "Out-of-frame"
+}
+fltr = out[,"VDJ Frame"] == ""
+if(any(fltr, na.rm = T)){
+	out[fltr, "VDJ Frame"] = "Out-of-frame"
+}
+
+for(col in c('Top V Gene','Top D Gene','Top J Gene')){
+	out[,col] = as.character(out[,col])
+	fltr = out[,col] == ""
+	fltr[is.na(fltr)] = T
+	if(any(fltr, na.rm = T)){
+		out[fltr,col] = "NA"
+	}
+}
+
+write.table(out, out.file, sep="\t", quote=F, row.names=F, col.names=T)
--- /dev/null	Thu Jan 01 00:00:00 1970 +0000
+++ b/imgt_loader/imgt_loader.sh	Fri Aug 18 10:41:51 2017 -0400
@@ -0,0 +1,74 @@
+#!/bin/bash
+input=$1
+output=$2
+name=$3
+dir="$(cd "$(dirname "$0")" && pwd)"
+mkdir -p $PWD/$name/files
+f=$(file $input)
+zip7Type="7-zip archive"
+tarType="tar archive"
+bzip2Type="bzip2 compressed"
+gzipType="gzip compressed"
+zipType="Zip archive"
+rarType="RAR archive"
+zxType="XZ compressed data"
+
+if [[ "$f" == *"$zip7Type"* ]]; then
+	echo "7-zip"
+	echo "Trying: 7za e $input -o$PWD/files/"
+	7za e $input -o$PWD/$name/files
+fi
+
+if [[ "$f" == *"$tarType"* ]]
+then
+	echo "tar archive"
+	echo "Trying: tar xvf $input -C $PWD/files/"
+	tar -xvf $input -C $PWD/$name/files
+fi
+
+if [[ "$f" == *"$bzip2Type"* ]]
+then
+	echo "bzip2 compressed data"
+	echo "Trying: tar jxf $input -C $PWD/files/"
+	tar -jxf $input -C $PWD/$name/files
+fi
+
+if [[ "$f" == *"$gzipType"* ]]
+then
+	echo "gzip compressed data"
+	echo "Trying: tar xvzf $input -C $PWD/files/"
+	tar -xvzf $input -C $PWD/$name/files
+fi
+
+if [[ "$f" == *"$zipType"* ]]
+then
+	echo "Zip archive"
+	echo "Trying: unzip $input -d $PWD/files/"
+	unzip $input -d $PWD/$name/files > $PWD/unziplog.log
+fi
+
+if [[ "$f" == *"$rarType"* ]]
+then
+	echo "RAR archive"
+	echo "Trying: unrar e $input $PWD/files/"
+	unrar e $input $PWD/$name/files
+fi
+
+if [[ "$f" == *"$zxType"* ]]
+then
+	echo "xz compressed data"
+	echo "Trying: tar -xJf $input -C $PWD/files/"
+	tar xJf $input -C $PWD/$name/files
+fi
+find $PWD/$name/files -iname "1_*" -exec cat {} + > $PWD/$name/summ.txt
+find $PWD/$name/files -iname "3_*" -exec cat {} + > $PWD/$name/sequences.txt
+find $PWD/$name/files -iname "4_*" -exec cat {} + > $PWD/$name/gapped_aa.txt
+find $PWD/$name/files -iname "5_*" -exec cat {} + > $PWD/$name/aa.txt
+find $PWD/$name/files -iname "6_*" -exec cat {} + > $PWD/$name/junction.txt
+
+echo "summ.txt `cat $PWD/$name/summ.txt | wc -l`"
+echo "aa.txt `cat $PWD/$name/aa.txt | wc -l`"
+
+#python $dir/imgt_loader.py --summ $PWD/$name/summ.txt --aa $PWD/$name/aa.txt --junction $PWD/$name/junction.txt --output $output
+
+Rscript --verbose $dir/imgt_loader.r $PWD/$name/summ.txt $PWD/$name/sequences.txt $PWD/$name/aa.txt $PWD/$name/junction.txt $PWD/$name/gapped_aa.txt $output 2>&1