SNPfilter 0.0.5
A handy little tool for filtering SNPs.
Requires Python: >=3.5

# SNPfilter
A handy little tool for filtering SNPs

## Installation
```
conda install -c bioconda samtools bwa
pip install snpfilter
```

## Usage

**prepare**: Use BWA and SAMTools to map the sequencing data, generate BAM files and BCF files, etc.

```
usage: SNPfilter prepare [-h] [-t THREADS] [-q MIN_MQ] [-Q MIN_BQ] sample_id reference R1 R2

Prepare work environment

positional arguments:
  sample_id             sample id
  reference             reference genome in fasta format, should be indexed by bwa
  R1                    fastq file 1
  R2                    fastq file 2

optional arguments:
  -h, --help            show this help message and exit
  -t THREADS, --threads THREADS
                        num of threads
  -q MIN_MQ, --min-MQ MIN_MQ
                        skip alignments with mapQ smaller than INT
  -Q MIN_BQ, --min-BQ MIN_BQ
                        skip bases with baseQ/BAQ smaller than INT
```

**qcfilter**: filtering SNPs from BCF file with min depth and min variant frequency

```
usage: SNPfilter qcfilter [-h] [-d MIN_DEPTH] [-v MIN_VARIANT_FREQUENCY] [-b BACKGROUND_SITE] sample_id input_bcf

filtering snp from bcf file with min depth and min variant frequency

positional arguments:
  sample_id             sample id
  input_bcf             input bcf file

optional arguments:
  -h, --help            show this help message and exit
  -d MIN_DEPTH, --min_depth MIN_DEPTH
                        min depth
  -D MAX_DEPTH, --max_depth MAX_DEPTH
                        max depth                        
  -v MIN_VARIANT_FREQUENCY, --min_variant_frequency MIN_VARIANT_FREQUENCY
                        min variant frequency
  -b BACKGROUND_SITE, --background_site BACKGROUND_SITE
                        A comma-separated list of bcf files, loci that appear in these bcf files will be filtered out
```

**codefilter**: filtering SNPs based on whether they cause changes in coding amino acids

```
usage: SNPfilter codefilter [-h] sample_id input_bcf reference gff_file

filtering SNPs based on whether they cause changes in coding amino acids

positional arguments:
  sample_id   sample id
  input_bcf   input bcf file
  reference   reference genome in fasta format
  gff_file    gff file for reference genome

optional arguments:
  -h, --help  show this help message and exit
```

## Example
1. Indexing of the reference genome using BWA
```
bwa index reference.fasta
```

2. Trim the sequencing data 
```
zcat sample1_R1.fastq.gz | bioawk -cfastx '{print "@"$name"\n"substr($seq, 16)"\n+"$name"\n"substr($qual, 16)}' | gzip > sample1_R1.trimmed.fastq.gz
zcat sample1_R2.fastq.gz | bioawk -cfastx '{print "@"$name"\n"substr($seq, 16)"\n+"$name"\n"substr($qual, 16)}' | gzip > sample1_R2.trimmed.fastq.gz
zcat sample2_R1.fastq.gz | bioawk -cfastx '{print "@"$name"\n"substr($seq, 16)"\n+"$name"\n"substr($qual, 16)}' | gzip > sample2_R1.trimmed.fastq.gz
zcat sample2_R2.fastq.gz | bioawk -cfastx '{print "@"$name"\n"substr($seq, 16)"\n+"$name"\n"substr($qual, 16)}' | gzip > sample2_R2.trimmed.fastq.gz
zcat sample3_R1.fastq.gz | bioawk -cfastx '{print "@"$name"\n"substr($seq, 16)"\n+"$name"\n"substr($qual, 16)}' | gzip > sample3_R1.trimmed.fastq.gz
zcat sample3_R2.fastq.gz | bioawk -cfastx '{print "@"$name"\n"substr($seq, 16)"\n+"$name"\n"substr($qual, 16)}' | gzip > sample3_R2.trimmed.fastq.gz
```

3. Use **prepare** to generate bam files and BCF files
```
SNPfilter prepare -t 10 sample1 reference.fasta sample1_R1.trimmed.fastq.gz sample1_R2.trimmed.fastq.gz
SNPfilter prepare -t 10 sample2 reference.fasta sample2_R1.trimmed.fastq.gz sample2_R2.trimmed.fastq.gz
SNPfilter prepare -t 10 sample3 reference.fasta sample3_R1.trimmed.fastq.gz sample3_R2.trimmed.fastq.gz
```

4. Use **qcfilter** to filter SNPs with a relaxed threshold（d=2, D=40, v=0.3）
```
SNPfilter qcfilter -d 2 -D 40 -v 0.3 sample1 sample1.bcf
SNPfilter qcfilter -d 2 -D 40 -v 0.3 sample2 sample2.bcf
SNPfilter qcfilter -d 2 -D 40 -v 0.3 sample3 sample3.bcf
```

5. Use **qcfilter** to filter SNPs with a strict threshold（d=5, D=40, v=0.9, other sample BCF file as background)
```
SNPfilter qcfilter -d 5 -D 40 -v 0.9 -b sample2.d2.v0.30.bcf,sample3.d2.D40.v0.30.bcf sample1 sample1.d2.D40.v0.30.bcf
SNPfilter qcfilter -d 5 -D 40 -v 0.9 -b sample1.d2.v0.30.bcf,sample3.d2.D40.v0.30.bcf sample2 sample2.d2.D40.v0.30.bcf
SNPfilter qcfilter -d 5 -D 40 -v 0.9 -b sample1.d2.v0.30.bcf,sample2.d2.D40.v0.30.bcf sample3 sample3.d2.D40.v0.30.bcf
```

6. Use **codefilter** to filter SNPs that do not cause changes in coding amino acids
```
SNPfilter codefilter sample1 sample1.d5.v0.90.bcf reference.fasta reference.gff
SNPfilter codefilter sample2 sample2.d5.v0.90.bcf reference.fasta reference.gff
SNPfilter codefilter sample3 sample3.d5.v0.90.bcf reference.fasta reference.gff
```