problem regarding mitochondria analysis
REQUIRED for all errors and issues: No read was passed after PrintReads.
Hello, I'm currently trying to perform mtDNA analysis. I sequenced 2 amplicons covered the whole mitochondrial genome using Miseq.
After trimming, I perform mapping to hg38 using bwa.
bwa index ref.fa
bwa mem ref.fa read1.fq read2.fq > aln.pe.sam
samtools view -bS aln.pe.sam > aln.pe.bam
samtools sort aln.pe.bem -o aln.pe.sorted.bam
samtools index aln.pe.sorted.bam
Then, I used gatk PrintReads as mentioned in best practice.
But no read was passed.
samtools view -c chrM.bam
0
Using samtools, reads mapped to chrM were successfully output.
samtools view -b aln.pe.sorted.bam chrM > chrM.samtools.bam
Is it due to my input file being mitochondria genome only, rather than WGS?
If so, which step should I start in best practice of mitochondrial short variant discovery?
Best regards,
a) GATK version used:gatk 4.5.0.0
b) Exact command used:
gatk PrintReads -I aln.pe.sorted.bam -L chrM -O chrM.bam
c) Entire program log:
Using GATK jar /lustre7/home/user/gatk-4.5.0.0/gatk-package-4.5.0.0-local.jar
Running:
java -Dsamjdk.use_async_io_read_samtools=false -Dsamjdk.use_async_io_write_samtools=true -Dsamjdk.use_async_io_write_tribble=false -Dsamjdk.compression_level=2 -jar /lustre7/home/user/gatk-4.5.0.0/gatk-package-4.5.0.0-local.jar PrintReads -I aln.pe.sorted.bam -L chrM -O chrM2.bam
02:01:06.459 INFO NativeLibraryLoader - Loading libgkl_compression.so from jar:file:/lustre7/home/user/gatk-4.5.0.0/gatk-package-4.5.0.0-local.jar!/com/intel/gkl/native/libgkl_compression.so
02:01:06.742 INFO PrintReads - ------------------------------------------------------------
02:01:06.745 INFO PrintReads - The Genome Analysis Toolkit (GATK) v4.5.0.0
02:01:06.745 INFO PrintReads - For support and documentation go to https://software.broadinstitute.org/gatk/
02:01:06.745 INFO PrintReads - Executing as user@at138 on Linux v5.15.0-87-generic amd64
02:01:06.745 INFO PrintReads - Java runtime: Java HotSpot(TM) 64-Bit Server VM v21.0.5+9-LTS-239
02:01:06.746 INFO PrintReads - Start Date/Time: November 19, 2024, 2:01:06 AM JST
02:01:06.746 INFO PrintReads - ------------------------------------------------------------
02:01:06.746 INFO PrintReads - ------------------------------------------------------------
02:01:06.746 INFO PrintReads - HTSJDK Version: 4.1.0
02:01:06.747 INFO PrintReads - Picard Version: 3.1.1
02:01:06.747 INFO PrintReads - Built for Spark Version: 3.5.0
02:01:06.747 INFO PrintReads - HTSJDK Defaults.COMPRESSION_LEVEL : 2
02:01:06.747 INFO PrintReads - HTSJDK Defaults.USE_ASYNC_IO_READ_FOR_SAMTOOLS : false
02:01:06.747 INFO PrintReads - HTSJDK Defaults.USE_ASYNC_IO_WRITE_FOR_SAMTOOLS : true
02:01:06.747 INFO PrintReads - HTSJDK Defaults.USE_ASYNC_IO_WRITE_FOR_TRIBBLE : false
02:01:06.747 INFO PrintReads - Deflater: IntelDeflater
02:01:06.748 INFO PrintReads - Inflater: IntelInflater
02:01:06.748 INFO PrintReads - GCS max retries/reopens: 20
02:01:06.748 INFO PrintReads - Requester pays: disabled
02:01:06.748 INFO PrintReads - Initializing engine
02:01:06.983 INFO IntervalArgumentCollection - Processing 16569 bp from intervals
02:01:06.984 INFO PrintReads - Done initializing engine
02:01:07.038 INFO ProgressMeter - Starting traversal
02:01:07.038 INFO ProgressMeter - Current Locus Elapsed Minutes Reads Processed Reads/Minute
02:01:08.060 INFO PrintReads - 167439 read(s) filtered by: WellformedReadFilter
02:01:08.061 INFO ProgressMeter - unmapped 0.0 0 0.0
02:01:08.061 INFO ProgressMeter - Traversal complete. Processed 0 total reads in 0.0 minutes.
02:01:08.371 INFO PrintReads - Shutting down engine
[November 19, 2024, 2:01:08 AM JST] org.broadinstitute.hellbender.tools.PrintReads done. Elapsed time: 0.03 minutes.
Runtime.totalMemory()=1191182336
-
I tried to the next step, nut error was occured.
samtools view -b aln.pe.sorted.bam chrM > samtools.chrM.bam
gatk RevertSamSpark -I samtools.chrM.bam -O reverted.samtools.chrM.bam
java -jar picard.jar MergeBamAlignment -ALIGNED samtools.chrM.bam -UNMAPPED reverted.samtools.chrM.bam -O merge_alignments.bam -R hg38_v0_chrM_Homo_sapiens_assembly38.chrM.shifted_by_8000_bases.fasta16:02:21.871 INFO NativeLibraryLoader - Loading libgkl_compression.so from jar:file:picard.jar!/com/intel/gkl/native/libgkl_compression.so
[Tue Nov 19 16:02:21 JST 2024] MergeBamAlignment --UNMAPPED_BAM reverted.samtools.chrM.bam --ALIGNED_BAM samtools.chrM.bam --OUTPUT merge_alignments.bam --REFERENCE_SEQUENCE hg38_v0_chrM_Homo_sapiens_assembly38.chrM.shifted_by_8000_bases.fasta --ADD_PG_TAG_TO_READS true --PAIRED_RUN true --CLIP_ADAPTERS true --IS_BISULFITE_SEQUENCE false --ALIGNED_READS_ONLY false --MAX_INSERTIONS_OR_DELETIONS 1 --ATTRIBUTES_TO_REVERSE OQ --ATTRIBUTES_TO_REVERSE U2 --ATTRIBUTES_TO_REVERSE_COMPLEMENT E2 --ATTRIBUTES_TO_REVERSE_COMPLEMENT SQ --READ1_TRIM 0 --READ2_TRIM 0 --ALIGNER_PROPER_PAIR_FLAGS false --SORT_ORDER coordinate --PRIMARY_ALIGNMENT_STRATEGY BestMapq --CLIP_OVERLAPPING_READS true --HARD_CLIP_OVERLAPPING_READS false --INCLUDE_SECONDARY_ALIGNMENTS true --ADD_MATE_CIGAR true --UNMAP_CONTAMINANT_READS false --MIN_UNCLIPPED_BASES 32 --MATCHING_DICTIONARY_TAGS M5 --MATCHING_DICTIONARY_TAGS LN --UNMAPPED_READ_STRATEGY DO_NOT_CHANGE --VERBOSITY INFO --QUIET false --VALIDATION_STRINGENCY STRICT --COMPRESSION_LEVEL 5 --MAX_RECORDS_IN_RAM 500000 --CREATE_INDEX false --CREATE_MD5_FILE false --help false --version false --showHidden false --USE_JDK_DEFLATER false --USE_JDK_INFLATER false
[Tue Nov 19 16:02:22 JST 2024] Executing as user@at138 on Linux 5.15.0-87-generic amd64; Java HotSpot(TM) 64-Bit Server VM 21.0.5+9-LTS-239; Deflater: Intel; Inflater: Intel; Provider GCS is available; Picard version: Version:3.3.0-SNAPSHOT
INFO 2024-11-19 16:02:22 SamAlignmentMerger Processing SAM file(s): [samtools.chrM.bam]
[Tue Nov 19 16:02:22 JST 2024] picard.sam.MergeBamAlignment done. Elapsed time: 0.00 minutes.
Runtime.totalMemory()=2164260864
To get help, see http://broadinstitute.github.io/picard/index.html#GettingHelp
Exception in thread "main" htsjdk.samtools.SAMFormatException: SAM validation error: ERROR::INVALID_FLAG_SUPPLEMENTARY_ALIGNMENT:Record 1, Read name M02699:70:000000000-LR2P7:1:1101:8983:4984, Supplementary alignment flag should not be set for unmapped read.
at htsjdk.samtools.SAMUtils.processValidationErrors(SAMUtils.java:460)
at htsjdk.samtools.BAMFileReader$BAMFileIterator.advance(BAMFileReader.java:865)
at htsjdk.samtools.BAMFileReader$BAMFileIterator.<init>(BAMFileReader.java:836)
at htsjdk.samtools.BAMFileReader$BAMFileIterator.<init>(BAMFileReader.java:824)
at htsjdk.samtools.BAMFileReader.getIterator(BAMFileReader.java:527)
at htsjdk.samtools.SamReader$PrimitiveSamReaderToSamReaderAdapter.iterator(SamReader.java:495)
at picard.sam.AbstractAlignmentMerger.mergeAlignment(AbstractAlignmentMerger.java:374)
at picard.sam.SamAlignmentMerger.mergeAlignment(SamAlignmentMerger.java:181)
at picard.sam.MergeBamAlignment.doWork(MergeBamAlignment.java:375)
at picard.cmdline.CommandLineProgram.instanceMain(CommandLineProgram.java:281)
at picard.cmdline.PicardCommandLine.instanceMain(PicardCommandLine.java:105)
at picard.cmdline.PicardCommandLine.main(PicardCommandLine.java:115)Is this due to supplementary alignments?
Thanks,
-
Hi man Sac
We recommend mapping mitochondria reads to mitochondria genome only rather than the whole genome fasta file. It may be better if you can try that way. Our recommended mitochondria workflows only map reads to mitochondria genome.
Regards.
Please sign in to leave a comment.
2 comments