samtools depth for multiple files
1
0
Entering edit mode
4.8 years ago
evelyn ▴ 230

Can samtools depth be used for finding out sequence depths for individual sorted bam files? I am not trying to get one value for all the files altogether.

 samtools depth ls sorted.bam |  awk '{sum+=$3} END { print "Average = ",sum/NR}' > log.txt

I want output (log.txt) as:

File name         Average 
1sorted.bam        --
2sorted.bam        --
alignment • 2.2k views
ADD COMMENT
0
Entering edit mode
4.8 years ago
 samtools depth  -H src/jvarkit-git/src/test/resources/S*.bam | awk -F '\t' '(NR==1) {split($0,header);N=0.0;next;} {N++;for(i=3;i<=NF;i++) a[i]+=int($i);} END { for(x in a) print header[x], a[x]/N;}'

src/jvarkit-git/src/test/resources/S2.bam 7.55906
src/jvarkit-git/src/test/resources/S3.bam 7.55906
src/jvarkit-git/src/test/resources/S4.bam 7.55717
src/jvarkit-git/src/test/resources/S5.bam 7.55722
src/jvarkit-git/src/test/resources/S1.bam 7.55722
ADD COMMENT
0
Entering edit mode

Thanks @Pierre! However, option -H does not work for me. When I tried running without -H, it gives the output as

0 5.4
1 4.4

It does not show file names instead it shows 0 and 1.

ADD REPLY
0
Entering edit mode

you version of samtools is too old. http://www.htslib.org/doc/samtools.html

$ samtools  depth

Usage: samtools depth [options] in1.bam [in2.bam [...]] 
Options: (...)
   -H                  print a file header (...) 

$ samtools    Program: samtools (Tools for alignments in the SAM format) Version:
1.9-74-gf69e678 (using htslib 1.9-160-g10e6de3)
ADD REPLY
0
Entering edit mode

I am using samtools/1.9

ADD REPLY
0
Entering edit mode

which version of samtools has -H for depth? I am still getting depth: invalid option -- 'H'

ADD REPLY
0
Entering edit mode

Thanks, I am new to bioinformatics and I am having some trouble installing this code. I am using mac.

[htslib]$ ./configure
checking for gcc... gcc
checking whether the C compiler works... yes
checking for C compiler default output file name... a.out
checking for suffix of executables... 
checking whether we are cross compiling... no
checking for suffix of object files... o
checking whether we are using the GNU C compiler... yes
checking whether gcc accepts -g... yes
checking for gcc option to accept ISO C89... none needed
checking for ranlib... ranlib
checking for grep that handles long lines and -e... /usr/bin/grep
checking for C compiler warning flags... -Wall
checking for pkg-config... /usr/bin/pkg-config
checking pkg-config is at least version 0.9.0... yes
checking for special C compiler options needed for large files... no
checking for _FILE_OFFSET_BITS value needed for large files... no
checking for _LARGEFILE_SOURCE value needed for large files... no
checking shared library type for unknown-Linux... plain .so
checking how to run the C preprocessor... gcc -E
checking for egrep... /usr/bin/grep -E
checking for ANSI C header files... yes
checking for sys/types.h... yes
checking for sys/stat.h... yes
checking for stdlib.h... yes
checking for string.h... yes
checking for memory.h... yes
checking for strings.h... yes
checking for inttypes.h... yes
checking for stdint.h... yes
checking for unistd.h... yes
checking for stdlib.h... (cached) yes
checking for unistd.h... (cached) yes
checking for sys/param.h... yes
checking for getpagesize... yes
checking for working mmap... yes
checking for gmtime_r... yes
checking for fsync... yes
checking for drand48... yes
checking whether fdatasync is declared... yes
checking for fdatasync... yes
checking for library containing log... -lm
checking for zlib.h... yes
checking for inflate in -lz... yes
checking for library containing recv... none required
checking for bzlib.h... yes
checking for BZ2_bzBuffToBuffCompress in -lbz2... yes
checking for lzma.h... no
checking for lzma_easy_buffer_encode in -llzma... no
configure: error: liblzma development files not found

The CRAM format may use LZMA2 compression, which is implemented in HTSlib
by using compression routines from liblzma <http://tukaani.org/xz/>.

Building HTSlib requires liblzma development files to be installed on the
build machine; you may need to ensure a package such as liblzma-dev (on Debian
or Ubuntu Linux), xz-devel (on RPM-based Linux distributions or Cygwin), or
xz (via Homebrew on macOS) is installed; or build XZ Utils from source.

Either configure with --disable-lzma (which will make some CRAM files
produced elsewhere unreadable) or resolve this error to build HTSlib.
ADD REPLY
0
Entering edit mode
ADD REPLY

Login before adding your answer.

Traffic: 1064 users visited in the last hour
Help About
FAQ
Access RSS
API
Stats

Use of this site constitutes acceptance of our User Agreement and Privacy Policy.

Powered by the version 2.3.6