Question: RNA-SEQ SNP calling
gravatar for babji2318
4 months 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 • 408 views
ADD COMMENTlink modified 4 months ago by JC7.9k • written 4 months ago by babji23180

Please use the search function:

enter image description here

ADD REPLYlink written 4 months ago by ATpoint17k
gravatar for JC
4 months ago by
JC7.9k 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 4 months ago • written 4 months ago by JC7.9k
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: 1610 users visited in the last hour