2017-08-28 61 views
0

正如你可以閱讀標題,我有興趣將存儲shell命令的結果並將其傳遞給另一個規則。存儲shell命令的結果

貝婁是我的規則:

SAMTOOLS = config["SAMTOOLS"] 
rule useDepth: 
    input: 
    depth = "{individual}_{chr}.fixmate.sort.rgmdup.bam.depth" 
    output: 
    tmpVCF = "{individual}_{chr}.vcf" 
    run: 
    depth = storage.fetch("chrDepth") 
    shell("echo {depth} | exit 1") 

rule calDepth: 
    input: 
    bam = "{individual}.fixmate.sort.rgmdup.bam" 
    output: 
    temp("{individual}_{chr}.fixmate.sort.rgmdup.bam.depth") 
    run: 
    import subprocess,shlex 
    depth=subprocess.check_output(shlex.split("{SAMTOOLS} depth -r {wildcards.chr} {input.bam} | awk '{{sum += $3}} END {{print sum/NR}}'"),shell=True) 
    storage.store("chrDepth", depth) 
    shell("echo \"Depth for {wildcards.chr} has been calculated\" > {output[0]}") 

我肯定在這裏接受,因爲1號出口的錯誤!但那只是爲了測試。

我試圖解決的錯誤是subprocess.check_output()中{SAMTOOLS}的值!

depth: 1: depth: {SAMTOOLS}: not found 
Error in job chrDepth while creating output file 
RuleException: 
Command '['{SAMTOOLS}', 'depth', '-r', '{wildcards.chr}', '{input.bam}', '|', 'awk', '{{sum += $3}} END {{print sum/NR}}']' 

要更多的信息,因爲不同勢用戶可能安裝samtools在不同的地方,我們做samtools的地址,通過CONFIGFILE配置。但是,在這裏我不能:

1)讀取{SAMTOOLS}的正確值!

2)使整個命令可以運行!

那麼,你能否告訴我是否有其他方式存儲/傳遞規則的輸出到另一個規則!?更特別的是,我如何增強snakemake來告訴shell {SAMTOOLS}可用。

謝謝!

回答

0

這是您設置訪問權限以用作Python變量。

SAMTOOLS = config["SAMTOOLS"] 

但你嘗試這裏訪問它,通過{} SAMTOOLS作爲Snakemake規則特定通配符:

depth=subprocess.check_output(shlex.split("{SAMTOOLS} depth -r {wildcards.chr} {input.bam} | awk '{{sum += $3}} END {{print sum/NR}}'"),shell=True) 

Snakemake通配符不被訪問的方式Python的變量相同。 此外,{SAMTOOLS}在這裏被作爲一個Snakemake通配符來訪問,但是你不會在規則的輸出中使用它作爲通配符。

假設{wildcards.chr}有效,並且{SAMTOOLS}調用是唯一未找到的通配符(不僅是第一個未知的通配符),我認爲您應該嘗試兩種方法之一。

沒有預分配:

depth=subprocess.check_output(shlex.split("config['SAMTOOLS'] depth -r {wildcards.chr} {input.bam} | awk '{{sum += $3}} END {{print sum/NR}}'"),shell=True) 

訪問它作爲Python變量作爲字符串(它是代表一個字符串的對象):

depth=subprocess.check_output(shlex.split(SAMTOOLS + " depth -r {wildcards.chr} {input.bam} | awk '{{sum += $3}} END {{print sum/NR}}'"),shell=True) 

最後,並且由於規則至少推薦它引入的規則耦合,有許多方法可以在Snakemake中通過規則傳遞變量,而且您已經在使用它了,但是,我不認爲這是需要的。正確的訪問和設計應該是足夠的。

Snakemake Tutorial FAQ: How to pass variables between rules

旁註

爲了消除跨規則的通過炭深度,同時也將其保存爲文件名的路徑,以及分離的規則,我強烈建議轉換chrDepth到命名通配符...

喜歡的東西...

rule useDepth: 
    input: 
    depth = "{individual}_{chr}_of_{chrDepth}.fixmate.sort.rgmdup.bam.depth" 
    output: 
    tmpVCF = "{individual}_{chr}_of{chrDepth}.vcf" 

但我不知道你是如何計算chrDepth。它讓我擔心你在所有這些規則之間傳遞它,而不僅僅依賴於良好的命名約定。它可能會不必要地耦合你的代碼,導致下游的問題和開銷。

+0

正如你所建議的使用SAMTOOLS以外的折扣是訣竅!但是,我想知道爲什麼我以前沒有收到這個錯誤。我正在使用{SAMTOOLS}或其他沒有通配符的程序,我的規則都沒有抱怨這種用法(即我在這些規則中使用shell)。 – khikho