Я пытаюсь настроить конвейер snakemake, который снова выравнивает и читает эталонный геном. Правило, которое индексирует геном, имеет параметр, которому требуется максимальная длина считываний, которые будут сопоставлены с геномом. Сейчас я вручную установил это в файле config.yaml, но было бы неплохо установить это в зависимости от того, какие входные данные передаются конвейеру. Например, запустив небольшой фрагмент кода awk, который находит максимальную длину во всех файлах fastq.
Конвейер может принимать любое заданное количество образцов для сопоставления с геномом, но следует использовать только самую длинную длину чтения в образце для построения индекса (который я бы предпочел сделать только один раз)
Есть ли способ сделать это?
Мое текущее правило выглядит так:
rule index_reference_genome:
"""
Index reference genome.
"""
input:
ref_fasta=config['ref_fasta'],
ref_gtf=config['ref_gtf']
output:
directory("{outdir}/references/star")
threads: 40
params:
len=int(config['max_read_len'])-1,
#star=config['star']
log: "{outdir}/logs/star/index_reference_genome.log"
benchmark: "{outdir}/benchmarks/star/index_reference_genome.log"
conda: config['conda_env']
shell:
"""
if [ ! -d {output} ] ; then mkdir {output}; fi;
STAR \
--runThreadN {threads} \
--runMode genomeGenerate \
--genomeDir {output} \
--genomeFastaFiles {input.ref_fasta} \
--sjdbGTFfile {input.ref_gtf} \
--sjdbOverhang {params.len} &> {log}; #--outTmpDir {output}/STAR_tmp
"""