Правило snakemake вызывает скрипт оболочки, но завершается после первой команды - PullRequest
0 голосов
/ 11 ноября 2018

У меня есть сценарий оболочки, который хорошо работает, если я просто запускаю его из командной строки. Когда я вызываю это из правила внутри змеиного мастера, оно терпит неудачу.

Скрипт запускает цикл for для файла идентификаторов и использует их для извлечения последовательностей из файла fastq, за которым следует многократное выравнивание последовательностей, и достигается консенсус.

Вот сценарий. Я поместил туда несколько эхо-операторов, и по какой-то причине он не вызывает команды. Это останавливается на заявлении grep.

Я попытался добавить set + o pipefail; в правиле, но это тоже не работает.

 #!/bin/bash


function Usage(){
    echo -e "\
    Usage: $(basename $0) -r|--read2 -l|--umi-list -f|--outfile \n\ 
    where: ... \n\
    " >&2
        exit 1
}

# Check argument count
[[ "$#" -lt 2 ]] && Usage
# parse arguments
while [[ "$#" -gt 1 ]];do
    case "$1" in
        -r|--read2)
            READ2="$2"
            shift
            ;;
        -l|--umi-list)
            UMI="$2"
            shift
            ;;
        -f|--outfile)
            OUTFILE="$2"
            shift
            ;;
        *)
            Usage
            ;;
    esac
    shift
done

# Set defaults


# Check arguments
[[ -f "${READ2}" ]] || (echo "Cannot find input file ${READ2}, exiting..." >&2; exit 1)
[[ -f "${UMI}" ]] || (echo "Cannot find input file ${UMI}, exiting..." >&2; exit 1)

#Create output directory
OUTDIR=$(dirname "${OUTFILE}")
[[ -d "${OUTDIR}" ]] || (set -x; mkdir -p "${OUTDIR}")


# Make temporary directories
TEMP_DIR="${OUTDIR}/temp"
[[ -d "${TEMP_DIR}" ]] || (set -x; mkdir -p "${TEMP_DIR}")


#RUN consensus script

for f in $( more "${UMI}" | cut -f1);do
NAME=$(echo $f)
grep  "${NAME}" "${READ2}" | cut -f1 -d ' ' | sed 's/@M/M/' > "${TEMP_DIR}/${NAME}.name"
echo subsetting reads
seqtk  subseq "${READ2}" "${TEMP_DIR}/${NAME}.name" | seqtk  seq  -A > "${TEMP_DIR}/${NAME}.fasta"

~/software/muscle3.8.31_i86linux64  -msf -in "${TEMP_DIR}/${NAME}.fasta" -out "${TEMP_DIR}/${NAME}.muscle.fasta"

echo make consensus
~/software/EMBOSS-6.6.0/emboss/cons -sequence "${TEMP_DIR}/${NAME}.muscle.fasta" -outseq "${TEMP_DIR}/${NAME}.cons.fasta"

  sed -i 's/n//g' "${TEMP_DIR}/${NAME}.cons.fasta"
  sed -i "s/EMBOSS_001/${NAME}.cons/" "${TEMP_DIR}/${NAME}.cons.fasta"
done
cat "${TEMP_DIR}/*.cons.fasta" > "${OUTFILE}"

Правило змеиного мейкера:

rule make_consensus:
    input:
        r2=get_extracted,
        lst="{prefix}/{sample}/reads/cell_barcode_umi.count" 
    output:
        fasta="{prefix}/{sample}/reads/fasta/{sample}.R2.consensus.fa"
    shell:
        "sh ./scripts/make_consensus.sh -r {input.r2} -l {input.lst} -f {output.fasta}"

Редактировать сообщения об ошибках Snakemake Я изменил некоторые пути на нейтральный путь к файлу

RuleException:
CalledProcessError in line 29 of ~/user/scripts/consensus.smk:
Command ' set -euo pipefail;  sh ./scripts/make_consensus.sh -r  ~/user/file.extracted.fastq -l ~/user/cell_barcode_umi
.count -f ~/user/file.consensus.fa ' returned non-zero exit status 1.
  File "~/user/scripts/consensus.smk", line 29, in __rule
_make_consensus
  File "~/user/miniconda3/lib/python3.6/concurrent/futures/thread.py", line 56, in run
Shutting down, this might take some time.
Exiting because a job execution failed. Look above for error message

Если есть лучшие способы сделать это, чем использовать оболочку для цикла, пожалуйста, дайте мне знать!

спасибо!

Редактировать

Скрипт работал автономно: первый grep

grep  AGGCCGTTCT_TGTGGATG R_extracted/wgs_5_OL_debug.R2.extracted.fastq | cut -f1 -d ' ' | sed 's/@M/M/' > ./fasta/temp/AGGCCGTTCT_TGTGGATG.name 

Скрипт проходил через змеиную фабрику: первые два оператора grep

grep  :::::::::::::: R_extracted/wgs_5_OL_debug.R2.extracted.fastq | cut -f1 -d ' ' | sed 's/@M/M/' > ./fasta/temp/::::::::::::::.name

Я сейчас пытаюсь выяснить, откуда взялись те :::: в змейке. Все идеи приветствуются

1 Ответ

0 голосов
/ 12 ноября 2018

Останавливается на операторе grep

Я предполагаю, что команда grep в make_consensus.sh ничего не захватывает.grep возвращает код выхода 1 в таких случаях, и ненулевой статус выхода распространяется на создание змеи.(см. также Обработка ошибки SIGPIPE в snakemake )

Слабая связь ... Существует несоответствие между shebang make_consensus.sh, который говорит, что сценарий должен быть выполнен с bash (#!/bin/bash) и фактическое выполнение с использованием sh (sh ./scripts/make_consensus.sh).(На практике это не должно иметь никакого значения, поскольку в любом случае sh, вероятно, перенаправляется на bash)

...