Skip to content

Commit

Permalink
Error checking in R scripts
Browse files Browse the repository at this point in the history
  • Loading branch information
richardmleggett committed Feb 3, 2017
1 parent 2041e3a commit f131185
Show file tree
Hide file tree
Showing 6 changed files with 335 additions and 284 deletions.
4 changes: 2 additions & 2 deletions Dockerfile
Original file line number Diff line number Diff line change
Expand Up @@ -19,8 +19,8 @@ RUN cd /usr ; git clone https://github.com/TGAC/NanoOK
ENV NANOOK_DIR="/usr/NanoOK"
RUN echo "export PATH=/usr/NanoOK/bin:${PATH}" >> ~/.bashrc
RUN Rscript -e "install.packages('ggplot2', repos='https://cran.ma.imperial.ac.uk/')"
RUN Rscript -e "install.packages('ggmap, repos='https://cran.ma.imperial.ac.uk/')"
RUN Rscript -e "install.packages('plyr, repos='https://cran.ma.imperial.ac.uk/')"
RUN Rscript -e "install.packages('ggmap', repos='https://cran.ma.imperial.ac.uk/')"
RUN Rscript -e "install.packages('plyr', repos='https://cran.ma.imperial.ac.uk/')"
RUN Rscript -e "install.packages('scales', repos='https://cran.ma.imperial.ac.uk/')"
RUN Rscript -e "install.packages('gridExtra', repos='https://cran.ma.imperial.ac.uk/')"
RUN Rscript -e "install.packages('reshape', repos='https://cran.ma.imperial.ac.uk/')"
60 changes: 38 additions & 22 deletions bin/nanook_plot_lengths.R
Original file line number Diff line number Diff line change
Expand Up @@ -35,37 +35,53 @@ for (t in 1:3) {

# Count vs length
filename_lengths <- paste(analysisdir, "/", "all_",type,"_lengths.txt", sep="");
filename_kmers <- paste(analysisdir, "/", "all_",type,"_kmers.txt", sep="");

if (file.exists(filename_lengths)) {
data_lengths = read.table(filename_lengths, col.name=c("name", "length"))

if (format=="png") {
lengths_png <- paste(graphsdir, "/", "all_",type,"_lengths.png", sep="");
png(lengths_png, width=1200, height=800)
print(ggplot(data_lengths, aes(x=data_lengths$length), xlab="Length") + geom_histogram(binwidth=1000, fill=colourcode) + xlab("Length") +ylab("Count") + scale_x_continuous(limits=c(0, 35000)) + ggtitle(type) + theme(text = element_text(size=textsize)) + theme(plot.margin = unit(c(0.02,0.02,0.04,0.02), "npc")) + theme(axis.title.x=element_text(vjust=-xvjust)) + theme(axis.title.y=element_text(vjust=yvjust)))
if (nrow(data_lengths) > 1) {
if (format=="png") {
lengths_png <- paste(graphsdir, "/", "all_",type,"_lengths.png", sep="");
png(lengths_png, width=1200, height=800)
print(ggplot(data_lengths, aes(x=data_lengths$length), xlab="Length") + geom_histogram(binwidth=1000, fill=colourcode) + xlab("Length") +ylab("Count") + scale_x_continuous(limits=c(0, 35000)) + ggtitle(type) + theme(text = element_text(size=textsize)) + theme(plot.margin = unit(c(0.02,0.02,0.04,0.02), "npc")) + theme(axis.title.x=element_text(vjust=-xvjust)) + theme(axis.title.y=element_text(vjust=yvjust)))
} else {
lengths_pdf <- paste(graphsdir, "/", "all_",type,"_lengths.pdf", sep="");
pdf(lengths_pdf, width=6, height=4)
print(ggplot(data_lengths, aes(x=data_lengths$length), xlab="Length") + geom_histogram(binwidth=1000, fill=colourcode) + xlab("Length") +ylab("Count") + scale_x_continuous(limits=c(0, 35000)) + ggtitle(type) + theme(text = element_text(size=textsize)) + theme(plot.margin = unit(c(0.02,0.02,0.04,0.02), "npc")) + theme(axis.title.x=element_text(vjust=-xvjust)) + theme(axis.title.y=element_text(vjust=yvjust)))
}
garbage <- dev.off()
} else {
lengths_pdf <- paste(graphsdir, "/", "all_",type,"_lengths.pdf", sep="");
pdf(lengths_pdf, width=6, height=4)
print(ggplot(data_lengths, aes(x=data_lengths$length), xlab="Length") + geom_histogram(binwidth=1000, fill=colourcode) + xlab("Length") +ylab("Count") + scale_x_continuous(limits=c(0, 35000)) + ggtitle(type) + theme(text = element_text(size=textsize)) + theme(plot.margin = unit(c(0.02,0.02,0.04,0.02), "npc")) + theme(axis.title.x=element_text(vjust=-xvjust)) + theme(axis.title.y=element_text(vjust=yvjust)))
cat("WARNING: No data in ", filename_lengths, "\n");
}
garbage <- dev.off()

# Number of perfect 21mers verses length scatter
filename_kmers <- paste(analysisdir, "/", "all_",type,"_kmers.txt", sep="");
data_alignments = read.table(filename_kmers, header=TRUE)
} else {
cat("WARNING: Couldn't find", filename_lengths, "\n");
}

# Number of perfect 21mers verses length scatter
if (file.exists(filename_kmers)) {
data_kmers = try(read.table(filename_kmers, header=TRUE), silent=TRUE)

if (format=="png") {
kmers_png <- paste(graphsdir, "/", "all_",type,"_21mers.png", sep="");
png(kmers_png, width=1200, height=800)
print(ggplot(data_alignments, aes(x=data_alignments$Length, y=data_alignments$nk21), xlab="Read length") + geom_point(shape=pointshape, size=pointsize, alpha=pointalpha, color=colourcode) + xlab("Read length") +ylab("Number of perfect 21mers") + ggtitle(type) + theme(text = element_text(size=textsize)) + scale_x_continuous(breaks=seq(0, 40000, 4000)) + scale_y_continuous(breaks=seq(0, 400, 20)) + theme(plot.margin = unit(c(0.02,0.02,0.04,0.02), "npc")) + theme(axis.title.x=element_text(vjust=-xvjust)) + theme(axis.title.y=element_text(vjust=yvjust)))
#grid.edit("geom_point.points", grep = TRUE, gp = gpar(lwd = pointwidth))
if (inherits(data_kmers, "try-error")) {
cat("WARNING: Couldn't read", filename_kmers,"\n");
} else {
kmers_pdf <- paste(graphsdir, "/", "all_",type,"_21mers.pdf", sep="");
pdf(kmers_pdf, width=6, height=4)
print(ggplot(data_alignments, aes(x=data_alignments$Length, y=data_alignments$nk21), xlab="Read length") + geom_point(shape=pointshape, alpha=pointalpha, color=colourcode) + xlab("Read length") +ylab("Number of perfect 21mers") + ggtitle(type) + theme(text = element_text(size=textsize)) + scale_x_continuous(breaks=seq(0, 40000, 4000)) + scale_y_continuous(breaks=seq(0, 400, 20)) + theme(plot.margin = unit(c(0.02,0.02,0.04,0.02), "npc")) + theme(axis.title.x=element_text(vjust=-xvjust)) + theme(axis.title.y=element_text(vjust=yvjust)))
if (nrow(data_kmers) > 1) {
if (format=="png") {
kmers_png <- paste(graphsdir, "/", "all_",type,"_21mers.png", sep="");
png(kmers_png, width=1200, height=800)
print(ggplot(data_kmers, aes(x=data_kmers$Length, y=data_kmers$nk21), xlab="Read length") + geom_point(shape=pointshape, size=pointsize, alpha=pointalpha, color=colourcode) + xlab("Read length") +ylab("Number of perfect 21mers") + ggtitle(type) + theme(text = element_text(size=textsize)) + scale_x_continuous(breaks=seq(0, 40000, 4000)) + scale_y_continuous(breaks=seq(0, 400, 20)) + theme(plot.margin = unit(c(0.02,0.02,0.04,0.02), "npc")) + theme(axis.title.x=element_text(vjust=-xvjust)) + theme(axis.title.y=element_text(vjust=yvjust)))
#grid.edit("geom_point.points", grep = TRUE, gp = gpar(lwd = pointwidth))
} else {
kmers_pdf <- paste(graphsdir, "/", "all_",type,"_21mers.pdf", sep="");
pdf(kmers_pdf, width=6, height=4)
print(ggplot(data_kmers, aes(x=data_kmers$Length, y=data_kmers$nk21), xlab="Read length") + geom_point(shape=pointshape, alpha=pointalpha, color=colourcode) + xlab("Read length") +ylab("Number of perfect 21mers") + ggtitle(type) + theme(text = element_text(size=textsize)) + scale_x_continuous(breaks=seq(0, 40000, 4000)) + scale_y_continuous(breaks=seq(0, 400, 20)) + theme(plot.margin = unit(c(0.02,0.02,0.04,0.02), "npc")) + theme(axis.title.x=element_text(vjust=-xvjust)) + theme(axis.title.y=element_text(vjust=yvjust)))
}
garbage <- dev.off()
} else {
cat("WARNING: No data in ", filename_kmers, "\n");
}
}
garbage <- dev.off()
} else {
cat("WARNING: Couldn't find", filename_lengths, "\n");
cat("WARNING: Couldn't find", filename_kmers, "\n");
}
}
Loading

0 comments on commit f131185

Please sign in to comment.