snakemake使用小結

 首先在linux 裏配置conda

下載

wget https://mirrors.tuna.tsinghua.edu.cn/anaconda/archive/Anaconda3-5.3.1-Linux-x86_64.sh

chmod +x Anaconda3-5.3.1-Linux-x86_64.sh

bash Anaconda3-5.3.1-Linux-x86_64.sh

安裝完畢，若是忘記選擇yes，敲conda命令報錯「command not found" 加上source /root/anaconda3/etc/profile.d/conda.sh

conda env list   獲得 /root/anaconda3

export PATH=~/anaconda3/bin:$PATH

否則全局沒法使用conda命令，（可是重啓putty好像就無論用了,還不清楚緣由）

vs code能夠不安裝

 

安裝tree命令，yum install tree

tree -af 能夠查看樹形文件結構 

 

snakemake是純python的任務流程工具（基於python3），之前商業環境用過control-M

https://snakemake.readthedocs.io/en/stable/

 

首先作個變異檢測，就是和標準的序列作對比，有點相似於代碼的compare，用過Beyond Compare或svn和git的筒子們應該很熟悉了，

可是基因序列是個很是大的序列文件，在linux也沒有windows那樣簡便的圖形操做見面，並且，這種對比工做是大量重複的，須要腳本化。

cd snakemake-snakemake-tutorial-623791d7ec6d
conda env create --name snakemake-tutorial --file environment.yaml

--------------------------------------------------------------

export PATH=~/anaconda3/bin:$PATH
source activate snakemake-tutorial

--------------------------------------------------------------

 

# 退出當前環境
source deactivate

這裏使用到Samtools工具，具體使用方法能夠參考https://blog.csdn.net/g863402758/article/details/53081342

他是一個用於處理sam與bam格式的工具軟件，可以實現二進制查看、格式轉換、排序及合併等功能，

結合sam格式中的flag、tag等信息，還能夠完成比對結果的統計彙總。同時利用linux中的grep、awk等操做命令，

還能夠大大擴展samtools的使用範圍與功能。

conda install snakemake

 conda install samtools

 bowtie2和samtools都是對比工具，bowtie2暫時沒安裝，安裝方法先記錄下

sudo wget https://jaist.dl.sourceforge.net/project/bowtie-bio/bowtie2/2.3.4.1/bowtie2-2.3.4.1-linux-x86_64.zip

unzip bowtie2-2.3.4.1-linux-x86_64.zip

vi /etc/environment

添加 bin 目錄的路徑，並用 : 隔開

source /etc/enviroment 使配置生效

開始寫job腳本

期間一直出一個錯誤，說Command must be given as string after the shell keyword
檢查一下是否會出錯rule bwa_map:
    input:
        "data/genome.fa", "data/samples/A.fastq" output: "mapped_reads/A.bam" shell: """ bwa mem {input} | samtools view -Sb - > {output} """運行snakemake -np mapped_reads/A.bam

執行這個job，把-n去掉

能夠看到，生成了A.bam文件

將A改爲{sample}，在輸入命令的時候加上你的參數，自動匹配上了，（注意此時文件夾貌似只能有一個腳本文件），cp了一個好像報錯了rule bwa_map:
    input:
        "data/genome.fa", "data/samples/{sample}.fastq" output: "mapped_reads/{sample}.bam" shell: """ bwa mem {input} | samtools view -Sb - > {output} """

接下來，要作排序了，代碼最後一塊兒貼

可使用dag選項和dot命令對「規則的執行和依賴關係」進行可視化，

snakemake --dag sorted_reads/{A,B,C}.bam.bai | dot -Tpdf > dag.pdf  這個命令好像會報錯

snakemake --dag sorted_reads/{A,B,C}.bam.bai | dot -Tsvg > dag.svg

整合以前的BAM文件，作基因組變異識別
SAMPLES=["A","B","C"] rule bcftools_call: input: fa="data/genome.fa", bam=expand("sorted_reads/{sample}.bam", sample=SAMPLES), bai=expand("sorted_reads/{sample}.bam.bai", sample=SAMPLES) output: "calls/all.vcf" shell: "samtools mpileup -g -f {input.fa} {input.bam} | " "bcftools call -mv - > {output}"
其中expand是自動匹配變量求文件路徑的語法糖
檢查一下，snakemake -np calls/all.vcf

最後出report，以上都是在規則裏執行shell腳本，snakemake的一個優勢就是能夠在規則裏面寫Python腳本，只須要把shell改爲run，此外還不須要用到引號。

測試一下，snakemake -np report.html

畫出流程圖

snakemake --dag report.html | dot -Tsvg > final.svg

 

執行一下：snakemake -p report.html
能夠看到生成了報告文件

到此，還有

rule all：

log：

多線程thread：

-j 指定cpu核心

params：

加載configfile: "config.yaml"

這幾個功能沒有操做，留個之後有空再處理

最後，在新建一個snakemake項目時，都先用conda create -n 項目名 python=版本號建立一個全局環境，用於安裝一些經常使用的軟件，例如bwa、samtools、seqkit等。而後用以下命令將環境導出成yaml文件

conda env export -n 項目名 -f environment.yaml

之後再部署的時候，

只須要conda env create -f environment.yaml

這個過程相似於ghost系統，或者打包虛擬機相似

參考瞭如下網址，感謝！

https://www.jianshu.com/p/8e57fd2b81b2

http://pedagogix-tagc.univ-mrs.fr/courses/ABD/practical/snakemake/snake_intro.html