Hi everybody, I'm trying to carry out some VC analysis, but I'm not very advanced at it, so I am a bit struggle... For my RNA-seq samples, I've run the following command lines:
#1 Mark duplicates
ls *bam | xargs -I {} -n 1 java -jar picard.jar MarkDuplicates I={} O={}.dedup M={}.metrix.txt
#2 Index file
ls *bam | xargs -I {} -n 1 java -jar picard.jar BuildBamIndex I={} O={}.bam.bai
#3 Split'N'Trim and reassign mapping qualities
ls *bam | xargs -I {} -n 1 java -jar GenomeAnalysisTK.jar -T SplitNCigarReads -R GRCh38.p10.genome.fa -I {} -o {}.split.bam -rf ReassignOneMappingQuality -RMQF 255 -RMQT 60 -U ALLOW_N_CIGAR_READS
But then, I get the following error:
##### ERROR MESSAGE: SAM/BAM/CRAM dedup.bam is malformed. Please see https://software.broadinstitute.org/gatk/documentation/article?id=1317for more information. Error details: SAM file doesn't have any read groups defined in the header. The GATK no longer supports SAM files without read groups
##### ERROR ------------------------------------------------------------------------------------------
Any ideas about how to solve it? Thanks!
1 answer
Error details: SAM file doesn't have any read groups defined in the header. The GATK no longer supports SAM files without read groups
Your sam file doesn't have read groups. This is a common issue with GATK, as the program now wants the reads to be grouped for analysis purposes.
You can assign them manually, or you can use Picard's tool AddOrReplaceReadGroups. Third possibility, if this tool of Picard doesn't work, is to make your own script to perform this task, as I had to do (this is mine: https://github.com/MatteoSchiavinato/Utilities/blob/master/add-read-groups.py ).
Make stuff up, the IDs and such are arbitrary.
Both Devon and Genomax are right, because - at least my personal experience - it depends on your analysis.
Are you running an analysis composed of different replicates/accessions/genomes/species, or are the data coming from one source and one replicate?
- Many sources: add read groups that make sense (in that case my script might not work as you desire)
- One source: make stuff up, like using the sample name and that's it (my script will likely do the job - not 100% guaranteed, cause it worked with my data but haven't tested it extensively)
Log in to answer this question.
-
How to solve the error 'ERROR::MATE_NOT_FOUND:Found xxx unpaired mates' when run SamToFastq
written by Lila M
-
SNP Calling using Genome Analysis Toolkit
written by Fid_o •
-
GATK multiple samples
written by evelyn
-
Getting only Indels and no SNPs after variant calling variant calling using HaplotypeCaller
written by Susmita Mandal
-
Using SSrGE for SNV calling
written by landscape95
-
Base Recabiliration Filters out 88% of the reads
written by Sharon
-
Regarding GATK for snp analysis
written by DL •
-
Malformed BAM file produced by BWA
written by joreamayarom
-
How to improve the pipeline of call snp in RNA-seq?
written by liuhui •
-
Badly formed genome loc: Contig '20' does not match any contig in the GATK sequence dictionary deri…
written by rubic
Have a look at one of the files, do you see anything obviously wrong with it?
Not at all with
but nothing comes up with:
does it means that my BAM file doesn't have read group and sample information?
yes
PAA* because they make it harder to read your post. In this case, it's probably "variant calling" what you are talking about. While VC may be a CUA** for you, this is not not necessarily the case for the rest of the users here.
* Please Avoid Abbreviations
** Commonly Used Abbreviation