| Previous changeset 2:5c4bca9dd4a2 (2019-01-14) Next changeset 4:3ffca09599ca (2020-10-16) |
|
Commit message:
"planemo upload for repository https://github.com/galaxyproject/tools-iuc/tree/master/tools/scpipe commit 8007f71281553ddfa45e6f8e1172952d956bb000" |
|
modified:
scpipe.R |
| b |
| diff -r 5c4bca9dd4a2 -r 7397e6badc11 scpipe.R --- a/scpipe.R Mon Jan 14 08:06:47 2019 -0500 +++ b/scpipe.R Thu Jun 11 07:18:37 2020 -0400 |
| [ |
| b'@@ -1,4 +1,8 @@\n-options( show.error.messages=F, error = function () { cat( geterrmessage(), file=stderr() ); q( "no", 1, F ) } )\n+err_foo <- function() {\n+ cat(geterrmessage(), file = stderr());\n+ q("no", 1, F)\n+}\n+options(show.error.messages = F, error = err_foo)\n \n # we need that to not crash galaxy with an UTF8 error on German LC settings.\n loc <- Sys.setlocale("LC_MESSAGES", "en_US.UTF-8")\n@@ -18,170 +22,171 @@\n })\n \n option_list <- list(\n- make_option(c("-bam","--bam"), type="character", help="BAM file"),\n- make_option(c("-fasta","--fasta"), type="character", help="Genome fasta file"),\n- make_option(c("-exons","--exons"), type="character", help="Exon annotation gff3 file"),\n- make_option(c("-organism","--organism"), type="character", help="Organism e.g. hsapiens_gene_ensembl"),\n- make_option(c("-barcodes","--barcodes"), type="character", help="Cell barcodes csv file"),\n- make_option(c("-read1","--read1"), type="character", help="Read 1 fastq.gz"),\n- make_option(c("-read2","--read2"), type="character", help="Read 2 fastq.gz"),\n- make_option(c("-samplename","--samplename"), type="character", help="Name to use for sample"),\n- make_option(c("-bs1","--bs1"), type="integer", help="Barcode start in Read 1"),\n- make_option(c("-bl1","--bl1"), type="integer", help="Barcode length in Read 1"),\n- make_option(c("-bs2","--bs2"), type="integer", help="Barcode start in Read 2"),\n- make_option(c("-bl2","--bl2"), type="integer", help="Barcode length in Read 2"),\n- make_option(c("-us","--us"), type="integer", help="UMI start in Read 2"),\n- make_option(c("-ul","--ul"), type="integer", help="UMI length in Read 2"),\n- make_option(c("-rmlow","--rmlow"), type="logical", help="Remove reads with N in barcode or UMI"),\n- make_option(c("-rmN","--rmN"), type="logical", help="Remove reads with low quality"),\n- make_option(c("-minq","--minq"), type="integer", help="Minimum read quality"),\n- make_option(c("-numbq","--numbq"), type="integer", help="Maximum number of bases below minq"),\n- make_option(c("-stnd","--stnd"), type="logical", help="Perform strand-specific mapping"),\n- make_option(c("-max_mis","--max_mis"), type="integer", help="Maximum mismatch allowed in barcode"),\n- make_option(c("-UMI_cor","--UMI_cor"), type="integer", help="Correct UMI sequence error"),\n- make_option(c("-gene_fl","--gene_fl"), type="logical", help="Remove low abundant genes"),\n- make_option(c("-max_reads","--max_reads"), type="integer", help="Maximum reads processed"),\n- make_option(c("-min_count","--min_count"), type="integer", help="Minimum count to keep"),\n- make_option(c("-metrics_matrix","--metrics_matrix"), type="logical", help="QC metrics matrix"),\n- make_option(c("-keep_outliers","--keep_outliers"), type="logical", help="Keep outlier cells"),\n- make_option(c("-report","--report"), type="logical", help="HTML report of plots"),\n- make_option(c("-rdata","--rdata"), type="logical", help="Output RData file"),\n- make_option(c("-nthreads","--nthreads"), type="integer", help="Number of threads")\n- )\n-\n-parser <- OptionParser(usage = "%prog [options] file", option_list=option_list)\n-args = parse_args(parser)\n-\n-bam = args$bam\n-fa_fn = args$fasta\n-anno_fn = args$exons\n-fq_R1 = args$read1\n-fq_R2 = args$read2\n-read_structure = list(\n- bs1 = args$bs1, # barcode start position in fq_R1, -1 indicates no barcode\n- bl1 = args$bl1, # barcode length in fq_R1, 0 since no barcode present\n- bs2 = args$bs2, # barcode start position in fq_R2\n- bl2 = args$bl2, # barcode length in fq_R2\n- us = args$us, # UMI start position in fq_R2\n- ul = args$ul # UMI length\n+ make_option(c("-bam", "--bam"), type = "character", help = "BAM file"),\n+ make_option(c("-fasta", "--fasta"), type = "character", help = "Genome fasta file"),\n+ make_option(c("-exons", "--exons"), type = "character", help = "Exon annotation gff3 file"),\n+ make_option(c("-organism", "--organism"'..b' bc_len = read_structure$bl2,\n+ max_reads = args$max_reads,\n+ min_count = args$min_count,\n+ max_mismatch = args$max_mis\n )\n }\n } else {\n- aligned_bam = file.path(out_dir, bam)\n- barcode_anno=args$barcodes\n+ aligned_bam <- file.path(out_dir, bam)\n+ barcode_anno <- args$barcodes\n }\n \n print("Assigning reads to exons")\n-sc_exon_mapping(aligned_bam, mapped_bam, anno_fn, bc_len=read_structure$bl2, UMI_len=read_structure$ul, stnd=args$stnd)\n+sc_exon_mapping(aligned_bam, mapped_bam, anno_fn, bc_len = read_structure$bl2, UMI_len = read_structure$ul, stnd = args$stnd)\n \n print("De-multiplexing data")\n-sc_demultiplex(mapped_bam, out_dir, barcode_anno, has_UMI=has_umi)\n+sc_demultiplex(mapped_bam, out_dir, barcode_anno, has_UMI = has_umi)\n \n print("Counting genes")\n-sc_gene_counting(out_dir, barcode_anno, UMI_cor=args$UMI_cor, gene_fl=args$gene_fl)\n+sc_gene_counting(out_dir, barcode_anno, UMI_cor = args$UMI_cor, gene_fl = args$gene_fl)\n \n print("Performing QC")\n sce <- create_sce_by_dir(out_dir)\n pdf("plots.pdf")\n plot_demultiplex(sce)\n if (has_umi) {\n- p = plot_UMI_dup(sce)\n+ p <- plot_UMI_dup(sce)\n print(p)\n }\n-sce = calculate_QC_metrics(sce)\n-sce = detect_outlier(sce)\n-p = plot_mapping(sce, percentage=TRUE, dataname=args$samplename)\n+sce <- calculate_QC_metrics(sce)\n+sce <- detect_outlier(sce)\n+p <- plot_mapping(sce, percentage = TRUE, dataname = args$samplename)\n print(p)\n-p = plot_QC_pairs(sce)\n+p <- plot_QC_pairs(sce)\n print(p)\n dev.off()\n \n print("Removing outliers")\n if (is.null(args$keep_outliers)) {\n- sce = remove_outliers(sce)\n+ sce <- remove_outliers(sce)\n gene_counts <- counts(sce)\n- write.table(data.frame("gene_id"=rownames(gene_counts), gene_counts), file="gene_count.tsv", sep="\\t", quote=FALSE, row.names=FALSE)\n+ write.table(data.frame("gene_id" = rownames(gene_counts), gene_counts), file = "gene_count.tsv", sep = "\\t", quote = FALSE, row.names = FALSE)\n }\n \n if (!is.null(args$metrics_matrix)) {\n- metrics <- colData(sce, internal=TRUE)\n- write.table(data.frame("cell_id"=rownames(metrics), metrics), file="metrics_matrix.tsv", sep="\\t", quote=FALSE, row.names=FALSE)\n+ metrics <- colData(sce, internal = TRUE)\n+ write.table(data.frame("cell_id" = rownames(metrics), metrics), file = "metrics_matrix.tsv", sep = "\\t", quote = FALSE, row.names = FALSE)\n }\n \n if (!is.null(args$report) & (!is.null(fa_fn))) {\n print("Creating report")\n- create_report(sample_name=args$samplename,\n- outdir=out_dir,\n- r1=fq_R1,\n- r2=fq_R2,\n- outfq=combined_fastq,\n- read_structure=read_structure,\n- filter_settings=filter_settings,\n- align_bam=aligned_bam,\n- genome_index=fasta_index,\n- map_bam=mapped_bam,\n- exon_anno=anno_fn,\n- stnd=args$stnd,\n- fix_chr=FALSE,\n- barcode_anno=barcode_anno,\n- max_mis=args$max_mis,\n- UMI_cor=args$UMI_cor,\n- gene_fl=args$gene_fl,\n- organism=args$organism,\n- gene_id_type="ensembl_gene_id")\n+ create_report(sample_name = args$samplename,\n+ outdir = out_dir,\n+ r1 = fq_r1,\n+ r2 = fq_r2,\n+ outfq = combined_fastq,\n+ read_structure = read_structure,\n+ filter_settings = filter_settings,\n+ align_bam = aligned_bam,\n+ genome_index = fasta_index,\n+ map_bam = mapped_bam,\n+ exon_anno = anno_fn,\n+ stnd = args$stnd,\n+ fix_chr = FALSE,\n+ barcode_anno = barcode_anno,\n+ max_mis = args$max_mis,\n+ UMI_cor = args$UMI_cor,\n+ gene_fl = args$gene_fl,\n+ organism = args$organism,\n+ gene_id_type = "ensembl_gene_id")\n }\n \n-if (!is.null(args$rdata) ) {\n- save(sce, file = file.path(out_dir,"scPipe_analysis.RData"))\n+if (!is.null(args$rdata)) {\n+ save(sce, file = file.path(out_dir, "scPipe_analysis.RData"))\n }\n \n sessionInfo()\n' |