From 460f8534c5d4e6410bee946734a34e0d71b9c1e5 Mon Sep 17 00:00:00 2001
From: rlyu <rlyu@svi.edu.au>
Date: Mon, 21 Jun 2021 15:48:13 +1000
Subject: [PATCH] add cell barcode option
---
src/sscocaller.nim | 17 +++++++++++------
1 file changed, 11 insertions(+), 6 deletions(-)
diff --git a/src/sscocaller.nim b/src/sscocaller.nim
index 9471c6b..83d7b2a 100755
--- a/src/sscocaller.nim
+++ b/src/sscocaller.nim
@@ -59,13 +59,13 @@ proc getTrans(pos1:int64,pos2:int64,cmPmb=0.1): float =
proc countAllele(rec:Variant, ibam:Bam, maxTotalReads:int,
minTotalReads:int,
chrom:string, mapq: int,
- barcodeTable:OrderedTableRef,minbsq:int,bulkBam:bool): Table[string,allele_expr] =
+ barcodeTable:OrderedTableRef,minbsq:int,bulkBam:bool,barcodeTag:string): Table[string,allele_expr] =
var alleleCountTable = initTable[string,allele_expr]()
var rec_alt:char
var total_reads=0
rec_alt = rec.ALT[0][0]
for aln in ibam.query(chrom = chrom,start = rec.POS.cint-1, stop = rec.POS.cint):
- var cbt = tag[string](aln, "CB")
+ var cbt = tag[string](aln, barcodeTag)
var currentCB: string
if cbt.isNone:
# echo "no cb "
@@ -127,7 +127,7 @@ proc countAllele(rec:Variant, ibam:Bam, maxTotalReads:int,
proc sscocaller(threads:int, vcff:string, barcodeFile:string,
bamfile:string, out_dir:string, mapq:int,
minbsq:int, mintotal:int, maxtotal:int, mindp:int, maxdp:int,
- thetaREF:float, thetaALT:float, cmPmb:float,s_Chrs:seq): int =
+ thetaREF:float, thetaALT:float, cmPmb:float,s_Chrs:seq,barcodeTag:string): int =
var
ibam:Bam
@@ -202,7 +202,8 @@ proc sscocaller(threads:int, vcff:string, barcodeFile:string,
minbsq=minbsq,
maxTotalReads = maxtotal,
minTotalReads = mintotal,
- bulkBam = bulkBam)
+ bulkBam = bulkBam,
+ barcodeTag = barcodeTag)
if alleleCountTable.len==0: continue
var rec_alt:char
@@ -404,7 +405,7 @@ proc sscocaller(threads:int, vcff:string, barcodeFile:string,
return 0
when(isMainModule):
- let version = "0.2.0"
+ let version = "0.2.1"
var doc = format("""
$version
@@ -424,6 +425,7 @@ Arguments:
Options:
-t --threads <threads> number of BAM decompression threads [default: 4]
+ -cb --cellbarcode <cellbarcode> the cell barcode tag, by default it is CB
-MQ --minMAPQ <mapq> Minimum MAPQ for read filtering [default: 20]
-BQ --baseq <baseq> base quality threshold for a base to be used for counting [default: 13]
-CHR --chrom <chrom> the selected chromsome (whole genome if not supplied,separate by comma if multiple chroms)
@@ -459,9 +461,12 @@ Options:
thetaREF:float
thetaALT:float
cmPmb:float
+ barcodeTag="CB"
if($args["--threads"] != "nil"):
threads = parse_int($args["--threads"])
+ if($args["--cellbarcode"] != "nil"):
+ barcodeTag = $args["--barcodeTag"]
else: threads = 4
if($args["--minDP"] != "nil"):
mindp = parse_int($args["--minDP"])
@@ -530,5 +535,5 @@ Options:
#var args = commandLineParams()
discard sscocaller(threads, vcff, barcodeFile, bamfile,
out_dir, mapq, minbsq, mintotal,
- maxtotal, mindp, maxdp, thetaREF, thetaALT, cmPmb,s_Chrs)
+ maxtotal, mindp, maxdp, thetaREF, thetaALT, cmPmb,s_Chrs,barcodeTag)
--
GitLab