關於 mates 文件與單曲文件的“xargs”錯誤
我目前正在使用 HISAT2 並嘗試使用 xargs 在輸入多個樣本時讓我的生活更輕鬆。
所以我有一個文本文件“samples.txt”,其中每個樣本名稱用空格分隔 -
ERR199044 ERR188104 ERR188234 ERR188245 ...
我目前的命令行輸入看起來像這樣 -
> cat ./samples.txt| xargs -I {} sh -c "./hisat2-2.1.0/hisat2 -p 8 --dta -x ./indexes/chrX_tran -1 ./samples/{}_chrX_1.fastq.gz -2 ./samples/{}_chrX_2.fastq.gz -S ./map/{}_chrX.sam"
對於每個輸入樣本名稱,我要使用的命令行輸出應該像這樣格式化:
./hisat2-2.1.0/hisat2 -p 8 --dta -x ./indexes/chrX_tran -1 ./samples/ERR199044_chrX_1.fastq.gz -2 ./samples/ERR199044_chrX_2.fastq.gz -S ./map/ERR199044_chrX.sam
我收到一條錯誤消息,上面寫著 -
Warning: Same mate file "./chrX_data/samples/ERR199044" appears as argument to both -1 and -2 Extra parameter(s) specified: "ERR188104_chrX_1.fastq.gz", "ERR188104_chrX_2.fastq.gz", "ERR188104_chrX.sam" Note that if <mates> files are specified using -1/-2, a <singles> file cannot also be specified. Please run bowtie separately for mates and singles. Error: Encountered internal HISAT2 exception (#1) Command: /mnt/c/Alex/Lab_Files/RNAseq/tools/hisat2-2.1.0/hisat2-align-s --wrapper basic-0 -p 8 --dta -x ./chrX_data/indexes/chrX_tran -S ./map/ERR199044 -1 ./chrX_data/samples/ERR199044 -2 ./chrX_data/samples/ERR199044 ERR188104_chrX_1.fastq.gz ERR188104_chrX_2.fastq.gz ERR188104_chrX.sam (ERR): hisat2-align exited with value 1
所以問題看起來是忽略“{}”之後的所有內容,使兩個不同的文件看起來相同,並且 HISAT2 停止工作。
我不清楚的是“mates”和“singles”之間的區別是什麼,如果有任何方法可以解決這個問題,所以我可以輸入相同的範例名稱並讓 Unix 理解它指定了多個具有該名稱的不同文件在裡面?
謝謝!
xargs 版本:
為此,您需要在空白處分割輸入文件的每一行(這是’
xargs
預設行為),然後使用. 或者,您可以讓 shell 腳本遍歷“$@”。xargs``sh``-n 1
您不能
-I {}
在此處使用,因為這會導致xargs
一次讀取輸入文件一行。即使您將分隔符設置為帶有 的空格-d ' '
,當它開始讀取下一個輸入行時,您也會在每個輸入行的末尾收到錯誤。幸運的是,您根本不需要使用
-I {}
- 您已經只是將輸入的 echo 字sh
作為單獨的參數附加到命令行的末尾,這是xargs
沒有 -I 的預設行為。在 shell 腳本中,您通過其位置參數(即 as
$1
)引用參數,就像在任何 shell 腳本中一樣。您可以$1
在 shell 腳本中隨意使用。您還需要為命令提供參數 0(sh 程序的任意名稱 - 字元串“sh”很方便 - 但如果您希望在 中容易找到它,您可以使用任何名稱
ps
)sh -c '...script...'
。它必須是參數之後的第一個'...script...'
參數,並且將$0
在 shell 腳本中。所以,你需要做這樣的事情(沒有 for 循環):
xargs -n 1 sh -c \ './hisat2-2.1.0/hisat2 -p 8 --dta -x ./indexes/chrX_tran \ -1 "./samples/$1_chrX_1.fastq.gz" -2 "./samples/$1_chrX_2.fastq.gz" \ -S "./map/$1_chrX.sam"' sh < ./samples.txt
或(使用 for 循環):
xargs sh -c ' for f in "$@"; do \ ./hisat2-2.1.0/hisat2 -p 8 --dta -x ./indexes/chrX_tran \ -1 "./samples/${f}_chrX_1.fastq.gz" -2 "./samples/${f}_chrX_2.fastq.gz" \ -S "./map/${f}_chrX.sam" done' sh < ./samples.txt
for 循環版本會更快,因為它不需要執行
sh
多次(每個“單詞”一次)。它執行
sh
的次數越少越好,限制為 shell 的最大命令行長度(在現代系統上約為 2MB)。除非samples.txt
非常大(超過 200,000 個條目),這意味著它只會執行sh
一次。shell while-read循環:
像這樣的工作不需要 xargs 。以下將在 bash 中工作,可能還有其他支持數組
-a
和read
.while read -a words; do for f in "${words[@]}"; do ./hisat2-2.1.0/hisat2 -p 8 --dta -x ./indexes/chrX_tran \ -1 "./samples/${f}_chrX_1.fastq.gz" \ -2 "./samples/${f}_chrX_2.fastq.gz" \ -S "./map/${f}_chrX.sam" done done < samples.txt
這會將每個輸入行的每個單詞讀入一個 bash 數組,然後(使用
for
循環遍歷數組中的每個單詞)執行帶有適當參數的 hisat2 程序。但是,請參閱:為什麼使用 shell 循環處理文本被認為是不好的做法?
awk 版本:
awk '{ for (i=1;i<=NF;i++) { printf "./hisat2-2.1.0/hisat2 -p 8 --dta -x ./indexes/chrX_tran -1 \"./samples/%s_chrX_1.fastq.gz\" -2 \"./samples/%s_chrX_2.fastq.gz\" -S \"./map/%s_chrX.sam\"\n", $i, $i,$i; } }' ./samples.txt | sh
請注意,這會將 awk 的輸出通過管道傳輸到 sh 以執行。如果沒有進入 sh 的管道,輸出如下所示:
./hisat2-2.1.0/hisat2 -p 8 --dta -x ./indexes/chrX_tran -1 "./samples/ERR199044_chrX_1.fastq.gz" -2 "./samples/ERR199044_chrX_2.fastq.gz" -S "./map/ERR199044_chrX.sam" ./hisat2-2.1.0/hisat2 -p 8 --dta -x ./indexes/chrX_tran -1 "./samples/ERR188104_chrX_1.fastq.gz" -2 "./samples/ERR188104_chrX_2.fastq.gz" -S "./map/ERR188104_chrX.sam" ./hisat2-2.1.0/hisat2 -p 8 --dta -x ./indexes/chrX_tran -1 "./samples/ERR188234_chrX_1.fastq.gz" -2 "./samples/ERR188234_chrX_2.fastq.gz" -S "./map/ERR188234_chrX.sam" ./hisat2-2.1.0/hisat2 -p 8 --dta -x ./indexes/chrX_tran -1 "./samples/ERR188245_chrX_1.fastq.gz" -2 "./samples/ERR188245_chrX_2.fastq.gz" -S "./map/ERR188245_chrX.sam"
這也將快速執行,因為 sh 腳本在 awk 腳本列印時執行每一行。
或者,您可以使用
sprintf
將命令行放入變數,而不是printf
標準輸出。然後你可以使用awk
’system()
函式直接執行它,類似於下面的 perl 範例:perl 版本:
perl -lane ' foreach $f (@F) { system(qw(./hisat2-2.1.0/hisat2 -p 8 --dta -x ./indexes/chrX_tran), -1, "./samples/${f}_chrX_1.fastq.gz", -2, "./samples/${f}_chrX_2.fastq.gz", "-S", "./map/${f}_chrX.sam"); };' ./samples.txt
這使用
perl
’system()
函式,因此它直接執行命令,而不需要通過管道傳遞到sh
.對於測試執行,請在引號運算符
echo
之後立即添加單詞和空格。qw(
順便說一句,添加程式碼以檢查文件是否存在,和/或測試每次執行是否
./hisat2-2.1.0/hisat2
成功或在此 perl 版本中對輸出進行後處理比在 shell 或 awk 版本中更容易 - 特別是如果它被寫為一個獨立的腳本,而不是一個單行。例如:#!/usr/bin/perl -w use strict; while(<>) { foreach my $f (split) { my $f1 = "./samples/${f}_chrX_1.fastq.gz"; my $f2 = "./samples/${f}_chrX_2.fastq.gz"; my $sam = "./map/${f}_chrX.sam"; if (!(-r $f1 && -r $f2 && -r $sam)) { warn "Missing or unreadable file for $f\n"; next }; my $rc = system( qw(echo ./hisat2-2.1.0/hisat2 -p 8 --dta -x ./indexes/chrX_tran), -1, $f1, -2, $f2, '-S', $sam ); if ($rc) { warn "hisat2 returned non-zero exit code for $f: $rc\n"; }; } }