Question: RNA-SEQ SNP calling
gravatar for babji2318
8 days ago by
babji23180 wrote:

I am comparing the effect of knocking down a gene in Aedes aegypti, I did a RNA-SEQ to understand expression profiles, I am also using this RNA-SEQ against the genome to check for variants I did a variant analysis using the GATK variant calling program, since GATK 4 doesnt support RNA-seq variant calling I used GATK 3. Used STAR for alignment I also did variant calling using samtools and bcftools pipeline

Is there a way to count the number of SNP's and INDEL's in my samples from the VCF files generated

I am looking for something that is going to say for example: Total SNP's detected: 3450, Total INDEL's detected: 460 Something along those lines so I can have a grasp of what is happening in the treatment samples

rna-seq snp variant indel • 89 views
ADD COMMENTlink modified 8 days ago by JC7.4k • written 8 days ago by babji23180

Please use the search function:

enter image description here

ADD REPLYlink written 8 days ago by ATpoint13k
gravatar for JC
8 days ago by
JC7.4k wrote:

A simple way is to parse and count the elements in the file, for example, to get how many variants are in your VCF, you can use grep command to count all lines which are not a header (marked with #):

grep -vc "#" file.vcf

Then you can filter all variant which looks like an SNP, I prefer to use Perl but you can use whatever you want:

perl -lane 'print if ($F[3] =~ m/^[ACGT]$/ and $F[4] =~ m/^[ACGT]$/)' | wc -l

In this case, we filter lines which have the fields #4 (reference call) and #5 (sample call) as a single DNA letter (A, C, G, or T), then we just count them with wc -l

Indels will be the subtraction of SNP from the total variants unless you are calling another type of variants, but GATK reports SNPs and short-Indels by default.

ADD COMMENTlink modified 8 days ago • written 8 days ago by JC7.4k
Please log in to add an answer.


Use of this site constitutes acceptance of our User Agreement and Privacy Policy.
Powered by Biostar version 2.3.0
Traffic: 1988 users visited in the last hour