GATK PrintRead Producing Empty File

Question

[BAGeno]

Hi,

I am running gatk command on human chr21 which I downloaded from UCSC. I have performed alignment with bwa. Then sorted, marked duplicates and added read groups with Picard tool. After I perform base quality score recalibration with gatk. But when I performed printreads step it produced file only which header, there were no reads in it. Please tell me on which step am I making mistake?

Here are BQSR command

/home/hpcc/tools/gatk-4.0.8.1/./gatk --java-options -Xmx8g  BaseRecalibrator -R hg19.fa -I  chr21.sorted.dedup.RG.bam --known-sites 1000G_phase1.indels.hg19.vcf -O chr21_BQRC.grp  --read-validation-stringency SILENT

and here is printreads command

/home/hpcc/tools/gatk-4.0.8.1/./gatk --java-options -Xmx8g PrintReads -R hg19.fa -I chr21.sorted.dedup.RG.bam   -O chr21_BQRC.bam  --read-validation-stringency SILENT

This is output of samtools flagstat

9 + 0 in total (QC-passed reads + QC-failed reads)
8 + 0 secondary
0 + 0 supplementary
0 + 0 duplicates
9 + 0 mapped (100.00% : N/A)
0 + 0 paired in sequencing
0 + 0 read1
0 + 0 read2
0 + 0 properly paired (N/A : N/A)
0 + 0 with itself and mate mapped
0 + 0 singletons (N/A : N/A)
0 + 0 with mate mapped to a different chr
0 + 0 with mate mapped to a different chr (mapQ>=5)

This is the message of GATK.

09:58:05.908 INFO  PrintReads - ------------------------------------------------------------
09:58:05.909 INFO  PrintReads - The Genome Analysis Toolkit (GATK) v4.0.8.1
09:58:05.909 INFO  PrintReads - For support and documentation go to https://software.broadinstitute.org/gatk/
09:58:05.910 INFO  PrintReads - Executing as hpcc@hpcc on Linux v4.4.0-142-generic amd64
09:58:05.911 INFO  PrintReads - Java runtime: OpenJDK 64-Bit Server VM v1.8.0_191-8u191-b12-2ubuntu0.16.04.1-b12
09:58:05.911 INFO  PrintReads - Start Date/Time: March 19, 2019 9:58:05 AM PKT
09:58:05.912 INFO  PrintReads - ------------------------------------------------------------
09:58:05.912 INFO  PrintReads - ------------------------------------------------------------
09:58:05.913 INFO  PrintReads - HTSJDK Version: 2.16.0
09:58:05.914 INFO  PrintReads - Picard Version: 2.18.7
09:58:05.914 INFO  PrintReads - HTSJDK Defaults.COMPRESSION_LEVEL : 2
09:58:05.914 INFO  PrintReads - HTSJDK Defaults.USE_ASYNC_IO_READ_FOR_SAMTOOLS : false
09:58:05.914 INFO  PrintReads - HTSJDK Defaults.USE_ASYNC_IO_WRITE_FOR_SAMTOOLS : true
09:58:05.914 INFO  PrintReads - HTSJDK Defaults.USE_ASYNC_IO_WRITE_FOR_TRIBBLE : false
09:58:05.914 INFO  PrintReads - Deflater: IntelDeflater
09:58:05.914 INFO  PrintReads - Inflater: IntelInflater
09:58:05.914 INFO  PrintReads - GCS max retries/reopens: 20
09:58:05.915 INFO  PrintReads - Using google-cloud-java fork https://github.com/broadinstitute/google-cloud-java/releases/tag/0.20.5-alpha-GCS-RETRY-FIX
09:58:05.915 INFO  PrintReads - Initializing engine
09:58:06.690 INFO  PrintReads - Done initializing engine
09:58:06.716 INFO  ProgressMeter - Starting traversal
09:58:06.718 INFO  ProgressMeter -        Current Locus  Elapsed Minutes       Reads Processed     Reads/Minute
09:58:07.691 INFO  PrintReads - 9 read(s) filtered by: WellformedReadFilter

09:58:07.693 INFO  ProgressMeter -             unmapped              0.0                     0              0.0

[Pierre Lindenbaum]

what is the output of samtools flagstat for chr21.sorted.dedup.RG.bam , what is are the log messages of gatk ?

Answer 1

[Pierre Lindenbaum]

the output of samtools flagstat says that you almost have no read ( 1 read + 8 supplementary) in your bam.

Furthermore, its seems that there is a problem with those reads as they've been fitered-out by : WellformedReadFilter