Parameters for BS-Snper

[yaaminiv]

TL;DR: Is there a good resource for understanding parameters used in SNP calling? BS-Snper has minimal documentation.

Context

Currently running BS-Snper in this Jupiter notebook with minimal arguments specified! On the off-chance I get code to execute properly the first time, I want to understand some of the other BS-Snper options I can use when calling SNPs.

The full list of options is in their README. Below are parameters I'm unfamiliar with + their explanations from the README:

--minhetfreq: Threshold of frequency for calling heterozygous SNP
--minhomfreq: Threshold of frequency for calling homozygous SNP
--minquali: Threshold of base quality
--minread2: Minimum mutation reads number
--errorate: Minimum mutation rate
--mapvalue: Minimum read mapping value

Are there certain thresholds that are standard practice for calling SNPs (i.e., a minimal mutation reads or mutation rate that's used in the literature)? I can do my own search but figured I'd ask first.

[yaaminiv]

@laurahspencer Do you have any insight into parameters since you're working on something similar?

[laurahspencer]

I don't have any insight into these settings - I used the tool HaplotypeCaller from gatk to call variants, and didn't mess with their settings. You could check that tool's page to see what their default settings are?

When it comes to filtering SNPs, should it be helpful, I found this site useful when hard-filtering SNPs.

[yaaminiv]

Related question! I ran this code:

!perl /Users/Shared/Apps/BS-Snper-master/BS-Snper.pl \
--input /Volumes/web/spartina/project-oyster-oa/Haws/bismark-2/*sorted.bam \
--fa /Volumes/web/spartina/project-oyster-oa/data/Cg-roslin/cgigas_uk_roslin_v1_genomic-mito.fa \
--output SNP-candidates.txt \
--methcg CpG-meth-info.tab \
--methchg CHG-meth-info.tab \
--methchh CHH-meth-info.tab \
--mincover 5 \
> SNP-results.vcf 2> ERR.log

I thought the wildcard would allow for multiple samples to be processed, but looks like BS-Snper only processes one file at a time. Is that normal for SNP calling? If I have multiple samples I want to consider when calling SNPs, do I call SNPs for each sample separately, then concatenate the data?

[kubu4]

The documentation indicates it only expects a single BAM file. You can concatenate all of your BAMs into a single file with samtools merge.

Assuming each individual BAM file has the appropriate @RG header, and assuming that BS-Snper understands how to process BAM files properly, then BS-Snper should be able to figure everything out from the merged BAM.

A lot of assumptions, I know, but since BS-Snper is specifically written to identify SNPs and it only accepts a single BAM file as input, I would I expect it to properly process a merged BAM file.

EDITED: Changed formatting to remove accidental, automatic GitHub user-tagging of @RG.

[mgavery]

I ran samples one at a time and merged the C/T SNP files afterward (as a bed file). As a note, I was only using the BS-Snper data to filter potential C/T SNPs from my methylation data (if a locus was called a SNP in any sample, it was removed from further analysis)

[kubu4]

@mgavery's point is a good one.

What are you planning on doing with the data after running it through BS-Snper? Ultimately, that will be the best guide on how to use it.

[yaaminiv]

I want to do something similar to what Mac described: filter potential C/T SNPs from methylation data. In that case, I guess it does make sense to identify SNPs in individual samples, since a SNP in one sample shouldn't be considered.

@mgavery Do you have a script/code for BS-Snper anywhere? I'd like to compare parameters!

[mgavery]

This is just an example script. I messed around with the --minqual parameter (here it's set to 15). I ended up using default parameters though

EpiMet_Sp_BSsnper_minqual15.txt

[ksil91]

One of my fave snp filtering papers is this one:
https://onlinelibrary.wiley.com/doi/full/10.1111/mec.14792

@yaaminiv if you are just trying to identify C/T snps to remove them, you probably don't need to worry too much about making sure you have high-confidence snps (that is, not due to genotyping error). But you are right, you probably only want to remove SNPs found in at least 2 individuals. You can use BS-snper to make a vcf file for each individual separately, then combine the vcfs with vcfmerge in vcftools . You can then use vcftools to filter your combined vcf to snps with a Minor Allele Count > 3 (--mac 3).

[nclowell]

I'm very late to the party but just wanted to add I also like the paper @ksil91 shared and that the pipelines I've used don't have these parameters except for minimum mapping quality per base, for which I've used 20. Hope you're feeling better about it a couple of weeks out @yaaminiv!