Mercurial > repos > davidvanzessen > plotting_merged
changeset 0:edcce37b28b8 draft
Uploaded
| author | davidvanzessen |
|---|---|
| date | Wed, 25 Sep 2013 09:03:29 -0400 |
| parents | |
| children | b1706c96df48 |
| files | RScript.r plotting_merged.xml r_wrapper.sh tool_dependencies.xml |
| diffstat | 4 files changed, 150 insertions(+), 0 deletions(-) [+] |
line wrap: on
line diff
--- /dev/null Thu Jan 01 00:00:00 1970 +0000 +++ b/RScript.r Wed Sep 25 09:03:29 2013 -0400 @@ -0,0 +1,103 @@ +options( show.error.messages=F, error = function () { cat( geterrmessage(), file=stderr() ); q( "no", 1, F ) } ) + +args <- commandArgs(trailingOnly = TRUE) + +inFile = args[1] +outFile = args[2] + +if (!require("gridExtra")) { +install.packages("gridExtra", repos="http://cran.xl-mirror.nl/") +} +library (gridExtra) +if (!require("ggplot2")) { +install.packages("ggplot2", repos="http://cran.xl-mirror.nl/") +} +require(ggplot2) +if (!require("plyr")) { +install.packages("plyr", repos="http://cran.xl-mirror.nl/") +} +require(plyr) + +test = read.table(inFile, sep="\t", header=TRUE) + +test$Top.V.Gene = gsub("[*]([0-9]+)", "", test$Top.V.Gene) +test$Top.D.Gene = gsub("[*]([0-9]+)", "", test$Top.D.Gene) +test$Top.J.Gene = gsub("[*]([0-9]+)", "", test$Top.J.Gene) + +test$VDJCDR3 = do.call(paste, c(test[c("Top.V.Gene", "Top.D.Gene", "Top.J.Gene","CDR3.Seq.DNA")], sep = ":")) + +PROD = test[test$VDJ.Frame != "In-frame with stop codon" & test$VDJ.Frame != "Out-of-frame" & test$CDR3.Found.How != "NOT_FOUND" , ] + +NONPROD = test[test$VDJ.Frame == "In-frame with stop codon" | test$VDJ.Frame == "Out-of-frame" | test$CDR3.Found.How == "NOT_FOUND" , ] + +PRODF = PROD[ -1] + +#unique(PRODF[duplicated(PRODF),]) +#length(row.names(PRODF[duplicated(PRODF),])) + +#length(row.names(PRODF)) +PRODF = unique(PRODF) +#length(row.names(PRODF)) + +PRODFV = ddply(PRODF, c("Sample", "Top.V.Gene"), function(x) summary(x$VDJCDR3)) +PRODFV$Length = as.numeric(PRODFV$Length) +Total = 0 +Total = ddply(PRODFV, .(Sample), function(x) data.frame(Total = sum(x$Length))) +PRODFV = merge(PRODFV, Total, by.x='Sample', by.y='Sample', all.x=TRUE) +PRODFV = ddply(PRODFV, c("Sample", "Top.V.Gene"), summarise, relFreq= (Length*100 / Total)) + +PRODFD = ddply(PRODF, c("Sample", "Top.D.Gene"), function(x) summary(x$VDJCDR3)) +PRODFD$Length = as.numeric(PRODFD$Length) +Total = 0 +Total = ddply(PRODFD, .(Sample), function(x) data.frame(Total = sum(x$Length))) +PRODFD = merge(PRODFD, Total, by.x='Sample', by.y='Sample', all.x=TRUE) +PRODFD = ddply(PRODFD, c("Sample", "Top.D.Gene"), summarise, relFreq= (Length*100 / Total)) + +PRODFJ = ddply(PRODF, c("Sample", "Top.J.Gene"), function(x) summary(x$VDJCDR3)) +PRODFJ$Length = as.numeric(PRODFJ$Length) +Total = 0 +Total = ddply(PRODFJ, .(Sample), function(x) data.frame(Total = sum(x$Length))) +PRODFJ = merge(PRODFJ, Total, by.x='Sample', by.y='Sample', all.x=TRUE) +PRODFJ = ddply(PRODFJ, c("Sample", "Top.J.Gene"), summarise, relFreq= (Length*100 / Total)) + +V = c("v.name\tchr.orderV\nIGHV7-81\t1\nIGHV3-74\t2\nIGHV3-73\t3\nIGHV3-72\t4\nIGHV3-71\t5\nIGHV2-70\t6\nIGHV1-69\t7\nIGHV3-66\t8\nIGHV3-64\t9\nIGHV4-61\t10\nIGHV4-59\t11\nIGHV1-58\t12\nIGHV3-53\t13\nIGHV3-52\t14\nIGHV5-a\t15\nIGHV5-51\t16\nIGHV3-49\t17\nIGHV3-48\t18\nIGHV3-47\t19\nIGHV1-46\t20\nIGHV1-45\t21\nIGHV3-43\t22\nIGHV4-39\t23\nIGHV3-35\t24\nIGHV4-34\t25\nIGHV3-33\t26\nIGHV4-31\t27\nIGHV4-30-4\t28\nIGHV4-30-2\t29\nIGHV3-30-3\t30\nIGHV3-30\t31\nIGHV4-28\t32\nIGHV2-26\t33\nIGHV1-24\t34\nIGHV3-23\t35\nIGHV3-22\t36\nIGHV3-21\t37\nIGHV3-20\t38\nIGHV3-19\t39\nIGHV1-18\t40\nIGHV3-15\t41\nIGHV3-13\t42\nIGHV3-11\t43\nIGHV3-9\t44\nIGHV1-8\t45\nIGHV3-7\t46\nIGHV2-5\t47\nIGHV7-4-1\t48\nIGHV4-4\t49\nIGHV4-b\t50\nIGHV1-3\t51\nIGHV1-2\t52\nIGHV6-1\t53") +tcV = textConnection(V) +Vchain = read.table(tcV, sep="\t", header=TRUE) +PRODFV = merge(PRODFV, Vchain, by.x='Top.V.Gene', by.y='v.name', all.x=TRUE) +close(tcV) + +D = c("v.name\tchr.orderD\nIGHD1-1\t1\nIGHD2-2\t2\nIGHD3-3\t3\nIGHD6-6\t4\nIGHD1-7\t5\nIGHD2-8\t6\nIGHD3-9\t7\nIGHD3-10\t8\nIGHD4-11\t9\nIGHD5-12\t10\nIGHD6-13\t11\nIGHD1-14\t12\nIGHD2-15\t13\nIGHD3-16\t14\nIGHD4-17\t15\nIGHD5-18\t16\nIGHD6-19\t17\nIGHD1-20\t18\nIGHD2-21\t19\nIGHD3-22\t20\nIGHD4-23\t21\nIGHD5-24\t22\nIGHD6-25\t23\nIGHD1-26\t24\nIGHD7-27\t25") +tcD = textConnection(D) +Dchain = read.table(tcD, sep="\t", header=TRUE) +PRODFD = merge(PRODFD, Dchain, by.x='Top.D.Gene', by.y='v.name', all.x=TRUE) +close(tcD) + + +J = c("v.name\tchr.orderJ\nIGHJ1\t1\nIGHJ2\t2\nIGHJ3\t3\nIGHJ4\t4\nIGHJ5\t5\nIGHJ6\t6") +tcJ = textConnection(J) +Jchain = read.table(tcJ, sep="\t", header=TRUE) +PRODFJ = merge(PRODFJ, Jchain, by.x='Top.J.Gene', by.y='v.name', all.x=TRUE) +close(tcJ) + +pV = ggplot(PRODFV) +pV = pV + geom_bar( aes( x=factor(reorder(Top.V.Gene, chr.orderV)), y=relFreq, fill=Sample), stat='identity', position="dodge") + theme(axis.text.x = element_text(angle = 90, hjust = 1)) +#pV + +pD = ggplot(PRODFD) +pD = pD + geom_bar( aes( x=factor(reorder(Top.D.Gene, chr.orderD)), y=relFreq, fill=Sample), stat='identity', position="dodge") + theme(axis.text.x = element_text(angle = 90, hjust = 1)) +#pD + +pJ = ggplot(PRODFJ) +pJ = pJ + geom_bar( aes( x=factor(reorder(Top.J.Gene, chr.orderJ)), y=relFreq, fill=Sample), stat='identity', position="dodge") + theme(axis.text.x = element_text(angle = 90, hjust = 1)) +#pJ + +#plotall = grid.arrange(pV, arrangeGrob(pD, pJ, ncol=2), ncol=1, widths=c(1,1.2)) +#plotall +#ggsave(outFile, dpi=125) + + +png(outFile,width = 1920, height = 1200) +print(plotall = grid.arrange(pV, arrangeGrob(pD, pJ, ncol=2), ncol=1, widths=c(1,1.2))) +dev.off() + +
--- /dev/null Thu Jan 01 00:00:00 1970 +0000 +++ b/plotting_merged.xml Wed Sep 25 09:03:29 2013 -0400 @@ -0,0 +1,15 @@ +<tool id="plotting_merged" name="Plotting the merged data" version="1.0"> + <description>Creating a graph of the merged reports.</description> + <command interpreter="bash"> + r_wrapper.sh $in_file $out_file + </command> + <inputs> + <param name="in_file" format="tabular" type="data" label="Data to Process" /> + </inputs> + <outputs> + <data format="png" name="out_file" /> + </outputs> + <help> + Creating a graph of the merged and parsed igblast reports. + </help> +</tool>
--- /dev/null Thu Jan 01 00:00:00 1970 +0000 +++ b/r_wrapper.sh Wed Sep 25 09:03:29 2013 -0400 @@ -0,0 +1,4 @@ +#!/bin/bash +dir="$(cd "$(dirname "$0")" && pwd)" +$R_LOCATION/bin/Rscript --verbose $dir/RScript.r $1 "$2.png" 2>&1 +mv -f "$2.png" $2
--- /dev/null Thu Jan 01 00:00:00 1970 +0000 +++ b/tool_dependencies.xml Wed Sep 25 09:03:29 2013 -0400 @@ -0,0 +1,28 @@ +<?xml version="1.0"?> +<tool_dependency> + <set_environment version="1.0"> + </set_environment> + + <package name="igBlastn" version="1.0.0"> + <install version="1.0"> + <actions> + <action type="download_by_url" target_filename="R-3.0.1.tar.gz"> + http://cran.rstudio.com/src/base/R-3/R-3.0.1.tar.gz + </action> + <action type="shell_command">./configure --with-x=no</action> + <action type="shell_command">make</action> + + <action type="move_directory_files"> + <source_directory>R-3.0.1</source_directory> + <destination_directory>$INSTALL_DIR/R-3.0.1</destination_directory> + </action> + + <action type="set_environment"> + <environment_variable name="R_LOCATION" action="set_to">$INSTALL_DIR/R-3.0.1</environment_variable> + </action> + </install> + <readme> + Download and configure R. + </readme> + </package> +</tool_dependency>