Is this procedure reasonable for SNP calling on RNA-seq data by GATK4?

AuraAura Member
edited October 2018 in Ask the GATK team

Hi, I have some RNA-seq samples from different individuals which belong to a species. I would like to combine these data as one bam file for SNP calling on species level by using GATK4. But I am not sure if it is reasonable?(The GATK Best Practices by GATK3 as reference)
Here is my steps:
1. STAR for alignment
2. Add read groups, sort, mark duplicates, and create index

#java -jar picard-tools-1.95/AddOrReplaceReadGroups.jar I=out.sam O=f.bam RGID=MN_S1B RGLB=library RGPL=platform RGPU=machine RGSM=S1B SORT_ORDER=coordinate
#java -jar picard-tools-1.95/MarkDuplicates.jar I=f.bam O=dedupped.bam  CREATE_INDEX=true VALIDATION_STRINGENCY=SILENT M=output.metrics

(I am not understand why need the output.metric because next steps don't apply it and this index also not apply next steps)
**
3. merge(How to set parameters for merge data from different individual?)**

#samtools merge [-nr] [-h] out.bam 1.bam 2.bam ...

**
4. create index****(I am not sure whether the input file is genome data or not)**

#samtools index genomic.fna

**
5. faidx(I am not sure whether the input file is genome data or not)**

#samtools faidx genomic.fa

**
6. 06_Split-Reassign**

#gatk SplitNCigarReads -R genomic.fna out.bam -O split.bam

**
7. Variant calling
8. Variant filtering**

In summary, I have the following question:
_
2. Is this procedure reasonable?
3. (Step2) I am not understand why need the output.metric because next steps don't apply it and this index also not apply next steps?
4. (Step 3) How to set parameters for merge data from different individual?
5. (Step4&5) I am not sure whether the input file is genome data or not?
Maybe bam file?
6. I don't do the Base Quality Score Recalibration (BQSR) for the lake of the known sites on the species. Does the lack of this step have an impact on the results? Can I replace it with known sites of similar species? _

Thanks

Post edited by shlee on
Tagged:

Comments

  • shleeshlee CambridgeMember, Broadie ✭✭✭✭✭

    Hi @Aura,

    I think perhaps you will find the RNA-Seq pipeline at https://github.com/gatk-workflows/gatk3-4-rnaseq-germline-snps-indels of interest. It is written in a workflow language called WDL.

    The MarkDuplicates metrics are for QC purposes, e.g. to gauge library complexity and to see the extent of duplicates. These metrics are of interest to sequencing centers, where they may guarantee a certain amount of non-duplicate library depth of coverage. For example, the Broad Genomics Platform guarantees 30X mean coverage for WGS and I believe this depth is that calculated after removing duplicates. This guarantees 30x independent observations towards conservative variant calling.

    I believe the GATK workflow does not support multi-sample calling. See https://software.broadinstitute.org/gatk/best-practices/workflow?id=11164. It states:

    This workflow is designed to operate on a set of samples one sample at a time; joint calling RNAseq is not supported.

    I think these two resources will help you most.

Sign In or Register to comment.