这将是一个有点苛刻的问题,但到目前为止,由sankemake doku主演的《死亡》并没有产生预期的结果,我希望有更有经验的人能帮助我度过难关。
让我描述一下我所拥有的:
1:外壳上的三个程序,执行方式如下:
./program_1 $(cat file/from/dir_1/foo.fa) $(cat file/from/dir_2/bar.fa) > output.txt
./program_2 $(cat file/from/dir_1/foo.fa) $(cat file/from/dir_2/bar.fa) > output.txt
./program_3 $(cat file/from/dir_1/foo.fa) $(cat file/from/dir_2/bar.fa) > output.txt
是的,cat是必要的,它们使用文件中的字符串,而不是文件本身。
2:两个目录看起来像这样:
hsa-let-7f-5p.fa
hsa-let-7g-5p.fa
hsa-miR-100-3p.fa
hsa-miR-100-5p.fa
hsa-miR-101-2-5p.fa
(更多(
NM_000044.fa
NM_000059.fa
NM_000075.fa
NM_000088.fa
NM_000103.fa
(更多(
3:具有所有所需输出组合的.txt文件:
hsa-miR-29b-3p__NM_138473_programm_1.txt
hsa-miR-29b-3p__NM_138473_programm_2.txt
hsa-miR-29b-3p__NM_138473_programm_3.txt
hsa-miR-545-3p__NM_002332_programm_1.txt
hsa-miR-545-3p__NM_002332_programm_2.txt
hsa-miR-545-3p__NM_002332_programm_3.txt
(更多(
并不是两个目录中的所有文件都被使用,有些文件被多次使用,我不希望所有的组合,只希望指定一次。输出文件应该是独立的,并根据上面指定的.txt命名。最终的sankefile应该在集群上很好地并行化。
工作流程非常简单:
1.读取输出文件名
2.从两个目录中检索组合所需的输入
3.在所有3个程序上运行它们,并输出预定义的txt。通过将程序stdout管道传输到一个文件中。
完成
但是。。。如何告诉Snakemake?有人告诉我这很方便,但到目前为止运气不好。如果有问题,请提问。提前感谢(:
我根据我所知调整了下面的代码。我现在正在从文件中读出;程序";在顶部;dirs";在";规则一";拥有完整的路径。我还切换了shell命令中的两个输入,因为我很困惑,它们是这样的。是的,我知道文件名是次优的。
out = []
f = open("snakemake_output_small.txt", "r")
for line in f:
out.append(line.replace('n', ''))
# Get distinct filenames
hsa = set([x.split('__')[0] for x in out])
nm = [x.split('__')[1] for x in out]
nm = set([re.sub('_program_.*', '', x) for x in nm])
program = ['full/path/to/program_1', 'full/path/to/program_2', 'full/path/to/program_3']
# Force snakemake to use exactly these wildcard values
# i.e. do not match by regex
wildcard_constraints:
hsa= '|'.join([re.escape(x) for x in hsa]),
nm= '|'.join([re.escape(x) for x in nm]),
program= '|'.join([re.escape(x) for x in program]),
rule all:
input:
out,
rule one:
input:
hsa= 'full/path/to/dir1/{hsa}.fa',
nm= 'full/path/to/dir2/{nm}.fa',
output:
'{hsa}__{nm}_{program}.txt',
shell:
r"""
{wildcards.program} {input.nm} {input.hsa} > {output}
"""
现在它告诉我:
Building DAG of jobs...
MissingInputException in line 21 of Snakefile:
Missing input files for rule all:
hsa-miR-545-3p__NM_002332_program_3.txt
hsa-miR-29b-3p__NM_138473_program_1.txt
hsa-miR-29b-3p__NM_138473_program_1.txt
hsa-miR-545-3p__NM_002332_program_2.txt
hsa-miR-29b-3p__NM_138473_program_3.txt
hsa-miR-545-3p__NM_002332_program_2.txt
这是100%的用户错误,但我缺少什么。如果我不是完全失明,输入目录是正确的。
我也忘了提到,其中两个程序有自己的参数,必须包括在内。程序2的A-d。以及用于程序3的-P此处/a/parameter/文件。因此shell命令可能需要分开。Sry这一切都很混乱。
像以前一样,如果有问题,请提问。
我更喜欢将输出文件名以及输入和程序的组合放在csv文件中,用panda读取它,并使用数据帧来指导工作流程。在这里,我解析输出文件名以提取相关的输入,这在我看来是模糊的。
为了容纳传递给每个程序的选项,请使用字典,或者更好的是,使用yaml配置文件,并使用program
通配符查询该字典:
# Read this list from file:
out= ['hsa-miR-29b-3p__NM_138473_program_1.txt',
'hsa-miR-29b-3p__NM_138473_program_2.txt',
'hsa-miR-29b-3p__NM_138473_program_3.txt',
'hsa-miR-545-3p__NM_002332_program_1.txt',
'hsa-miR-545-3p__NM_002332_program_2.txt',
'hsa-miR-545-3p__NM_002332_program_3.txt']
hsa = set([x.split('__')[0] for x in out])
nm = [x.split('__')[1] for x in out]
nm = set([re.sub('_program_.*', '', x) for x in nm])
# Better to put this in a yaml file and read it with `--configfile progs.yaml`
programs = {'program_1':
{'path': '/path/to/program_1', 'opts': '-d foo -P bar'},
'program_2':
{'path': '/path/to/program_2', 'opts': '-P spam eggs'},
'program_3':
{'path': '/path/to/program_3', 'opts': ''}
}
wildcard_constraints:
hsa= '|'.join([re.escape(x) for x in hsa]),
nm= '|'.join([re.escape(x) for x in nm]),
program= '|'.join([re.escape(x) for x in programs.keys()]),
rule all:
input:
out,
rule one:
input:
hsa= 'dir1/{hsa}.fa',
nm= 'dir2/{nm}.fa',
output:
'{hsa}__{nm}_{program}.txt',
params:
path= lambda wc: programs[wc.program]['path'],
opts= lambda wc: programs[wc.program]['opts'],
shell:
r"""
{params.path} {params.opts} {input.hsa} {input.nm} > {output}
"""