Mercurial > repos > kyost > atac_primer_tool
changeset 5:3cd53127a838 draft default tip
planemo upload commit 103cb51ec368438642504c3f98b76c363db478bb
| author | kyost |
|---|---|
| date | Fri, 04 May 2018 02:24:47 -0400 |
| parents | 72571a30f17b |
| children | |
| files | ATAC_Primer_Tool.ga concatenate_ordered.xml filter_bam.bash filter_bam.xml find_qPCR_regions.R find_qPCR_regions.xml make_window_bed.R make_window_bed.xml window_coverage.xml |
| diffstat | 2 files changed, 4 insertions(+), 8 deletions(-) [+] |
line wrap: on
line diff
--- a/filter_bam.bash Wed May 02 13:20:21 2018 -0400 +++ b/filter_bam.bash Fri May 04 02:24:47 2018 -0400 @@ -29,7 +29,7 @@ #create bed file of fragments for each bam file -samtools view -b -L $ext_peak_file $bam_file > $filtered_bam_file +samtools view -b -F 4 -L $ext_peak_file $bam_file > $filtered_bam_file rm name rm count
--- a/find_qPCR_regions.R Wed May 02 13:20:21 2018 -0400 +++ b/find_qPCR_regions.R Fri May 04 02:24:47 2018 -0400 @@ -1,8 +1,7 @@ ## Command to run tool: # Rscript --vanilla find_qPCR_regions.R o.coverage.bed f9.coverage.bed lib_sizes.txt cor_cutoff cov_cutoff output_file -# Set up R error handling to go to stderr -options(show.error.messages=F, error=function(){cat(geterrmessage(),file=stderr());q("no",1,F)}) +options(warn=-1) # Avoid crashing Galaxy with an UTF8 error on German LC settings loc <- Sys.setlocale("LC_MESSAGES", "en_US.UTF-8") @@ -58,9 +57,6 @@ a_split <- asplit[[i]] b_split <- bsplit[[i]] - if (sd(t(a_split[1,5:m])) == 0 || sd(t(b_split[-1,5:m])) == 0 ){ - next - } #calculate correlation between total peak reads and spanning reads in each window corb <- data.frame(cor(t(a_split[1,5:m]) @@ -83,7 +79,7 @@ #keep windows for which correlation with total peak reads is above threshold and #number of spanning reads is above cutoff - regions <- data[which(data$correlation > corr_cut & data3$averageReadDepth > cov_cut),] + regions <- data[which(data$correlation >= corr_cut & data3$averageReadDepth >= cov_cut),] regions <- regions[,2:4] #collapse windows to non-overlappin regions @@ -93,7 +89,7 @@ nregion <- 1 if (nrow(regions)>=2){ for (j in 2:nrow(regions)) { - if (as.numeric(regions[j,2]) <= (as.numeric(regions[j-1,3]))) { + if (as.numeric(regions[j,2]) < (as.numeric(regions[j-1,3]))) { collapsed_regions[nrow(collapsed_regions),3] <- regions[j,3] }else{ collapsed_regions <- rbind(collapsed_regions, regions[j,])