== SNP Calling ==
=== mpileup (SAMtools) ===
  * Create bcf file (binary version of vcf) from alignment file (bam)


{{{
  #Use Q to set base quality threshold
  samtools mpileup -Q 25 -ugf <reference.fasta> <file.bam> | bcftools view -bvcg - > accepted_hits.raw.bcf
  #mpileup BAQ computation might be too stringent, disabling BAQ may be needed
  samtools mpileup -B -Q 25 -ugf <reference.fasta> <file.bam> | bcftools view -bvcg - > accepted_hits.raw.bcf
}}}


  * Convert bcf to vcf (vairant call format) and filter, if needed

  
{{{
#filter using varFilter (from vcfutils)
  #use Q to set mapping quality
  #use d for minimum read depth
  bcftools view accepted_hits.raw.bcf | vcfutils.pl varFilter -d 5 -Q 20 >| accepted_hits.raw.vcf
}}}



== Filtering SNPs ==

  * Remove, or retrieve, known SNPs


{{{
  #To remove known SNPs, use intersectBed (bedtools).  Known SNPs can be downloaded from Ensembl or NCBI/UCSC (eg. dbSNP)
  intersectBed -a accepted_hits.raw.filtered.vcf -b Mus_musculus.NCBIM37.60.bed -wo > filteredSNPs.vcf
}}}

  
== Annotating SNPs ==
    * Determining the effect of the SNP

=== snpEff ===
 
{{{
 #Use appropriate i) organism  and ii) annotation (eg. UCSC, RefSeq, Ensembl)
  #Mouse: mm37 (UCSC/RefSeq), mm37.61 (Ensembl)
  #Human: hg37 (UCSC/RefSeq), hg37.61 (Ensembl)
  #Output is created in several files: an html summary file and text files with detailed information 
  snpEff -vcf4 mm37.61 accepted_hits.raw.filtered.vcf > accepted_hits.raw.filtered.snpEff

}}}

=== Resources from the Broad Institute ===

  * [http://www.broadinstitute.org/partnerships/education/broade/best-practices-variant-calling-gatk Best Practices For Variant Calling With The GATK] (workshop materials)