change keep markdup reads in bam file to don't markdup

Author: tuan-tt

README.md

@@ -2,16 +2,19 @@
 
 Developed by BioTuring (www.bioturing.com), <i><b>mdup</b></i> is a tool that preprocess cloud-read data (read has barcode). <i><b>mdup</b></i> will do: 
 
-- Markdup PCR duplication, remove not primary reads, duplicate reads, secondary alignment, unmapped reads.
+- Remove duplicate reads, remove not primary reads, secondary alignment, unmapped reads.
 - Detect molecule by clustering reads have same barcode into group.
 - Get stats about sequencing and GEM performance.
 
+Two reads are consider duplicate if they share same mapped position, mapped target, cigar,
+mate info (if paired-end).
+
 ## Install
 
 ```shell
 git clone https://github.com/bioturing/mdup.git
 cd mdup
-sh build.sh
+bash build.sh
 ```
 
 ## Usage
@@ -32,7 +35,7 @@ Optional arguments:
   -t INT                number of threads [default: 1]
   -o DIR                output directory [default: "./mdup_out/"]
   -g FILE               reference file that generated bam file (for better stats)
-  -k                    keep all record from BAM file, turn on duplicate bit flag instead.
+  -k                    don't mark duplicate.
 ```
 
 # Contacts

argument.c

@@ -16,7 +16,7 @@ void print_usage()
 	__VERBOSE("  -t INT                number of threads [4]\n");
 	__VERBOSE("  -o DIR                output directory [\"./mdup_out/\"]\n");
 	__VERBOSE("  -g FILE               reference file that generated bam file (better statictis)\n");
-	__VERBOSE("  -k                    keep all record from BAM file, turn on duplicate bit flag instead.\n");
+	__VERBOSE("  -k                    don't mark duplicate\n");
 	__VERBOSE("\n");
 	__VERBOSE("This tool will generate some file in output directory:\n");
 	__VERBOSE("  output.bam            BAM file after processed\n");
@@ -33,7 +33,7 @@ void get_args(int argc, char *argv[])
 	args.is_remove = true;
 	args.reference = NULL;
 
-	if (argc < 3) {
+	if (argc < 2) {
 		print_usage();
 		exit(0);
 	}

argument.h

@@ -3,7 +3,7 @@
 
 #include "attr.h"
 
-#define VERSION			"1.2"
+#define VERSION			"1.3"
 
 struct prog_args {
 	int n_thread;

bam.c

@@ -37,8 +37,6 @@ void read_bam_unmapped(struct bam_inf_t *bam_inf, struct stats_t *stats)
 			__ERROR("Unmapped read doesn't have not primary flag or suplementary flag!");
 		}
 		get_basic_stats(b, stats);
-		if (!args.is_remove)
-			sam_write1(out_bam_f, bam_inf->b_hdr, b);
 	}
 
 	bam_destroy1(b);
@@ -71,8 +69,6 @@ void read_bam_target(struct bam_inf_t *bam_inf, int id, struct stats_t *stats)
 		if (b->core.flag & (FLAG_NOT_PRI | FLAG_SUPPLEMENT)) {
 			duplicate_try_process(b->core.pos, stats,
 					      out_bam_f, bam_inf->b_hdr);
-			if (!args.is_remove)
-				sam_write1(out_bam_f, bam_inf->b_hdr, b);
 			continue;
 		}
 
@@ -82,8 +78,6 @@ void read_bam_target(struct bam_inf_t *bam_inf, int id, struct stats_t *stats)
 		if (b->core.n_cigar == 0) {
 			duplicate_try_process(b->core.pos, stats,
 					      out_bam_f, bam_inf->b_hdr);
-			if (!args.is_remove)
-				sam_write1(out_bam_f, bam_inf->b_hdr, b);
 			continue;
 		}
 

build.sh

@@ -1,13 +1,24 @@
 #!/bin/bash
 
-git submodule init && git submodule update
+RED="\033[0;31m\033[1m"
+NC="\033[0m" # No Color
 
-dir="$(pwd)/htslib/"
+if [ "$BASH_VERSION" = '' ]; then
+    printf "${RED}Please run by command: bash build.sh${NC}\n"
+    exit 1
+fi
 
-cd htslib && autoheader && autoconf
-
-./configure --prefix=$dir --disable-lzma --disable-libcurl --disable-bz2
+DIR="$( cd "$( dirname "${BASH_SOURCE[0]}" )" && pwd )"
+cd ${DIR}
+git submodule update --init --recursive || exit 1
 
-make -j 4 && make install && cd ../
+# build htslib
+DIR="$(pwd)/htslib/"
+cd htslib && autoheader && autoconf
+./configure --prefix=${DIR} --disable-lzma --disable-libcurl --disable-bz2 || exit 1
+make || exit 1
+make install || exit 1
+cd ../
 
-make
+# build mdup
+make || exit 1

duplicate.c

@@ -78,9 +78,11 @@ void duplicate_process(struct stats_t *stats, samFile *out_bam_f, bam_hdr_t *b_h
 		for (i = u; i < v; ++i) {
 			if (i == pos)
 				continue;
-			align[i].b->core.flag |= FLAG_DUPLICATE;
-			if (!args.is_remove)
+			if (!args.is_remove) {
 				sam_write1(out_bam_f, b_hdr, align[i].b);
+				coverage_add(align[i].b, stats);
+				mlc_insert(align[i].bx_id, align[i].b, stats);
+			}
 		}
 
 		sam_write1(out_bam_f, b_hdr, align[pos].b);