我有许多以 fastq.gz 结尾的 fastq 文件。
rep1.fastq.gz
rep2.fastq.gz
rep3.fastq.gz
rep4.fastq.gz
.....
我期望我的输出为
rep1.fastq.gz 23516782
rep2.fastq.gz 45126780
rep3.fastq.gz 67543908
rep4.fastq.gz 76425368
其中第 1 行显示我的每个输入文件,第 2 行显示每个文件中的序列数。
为了实现这一点,我编写了一个小的 bash 脚本来计算每个文件中的序列数,并将每个文件后面写入的数字作为输出
for sample in *.fastq.gz;do echo -en $sample "\t";(zcat $sample|wc -l)/4|bc ;done
我收到错误:-bash:意外标记“/4”附近出现语法错误
答案1
忽略每个序列 4 行的假设中可能出现的所有错误...上面显示的命令应使用以下格式
for file in *.fastq.gz; do echo -en $file "\t";echo "$(zcat $file| wc -l)"/4 |bc;done
答案2
你正在运行这个:
(zcat $sample|wc -l)/4|bc
/4
shell 无法理解其中唯一的内容,这就是它失败的原因。我怀疑你想要的是得到结果命令的值zcat $sample|wc -l
,然后打印该值并将/4
其传递给bc
.如果是这样,您$()
不仅需要而且()
还需要引用它:
echo "$(zcat $sample|wc -l)/4" | bc
所以这意味着:
for sample in *.fastq.gz; do
echo -en $sample "\t"; echo "$(zcat $sample|wc -l)/4" | bc
done
或者,更便携一点:
for sample in *.fastq.gz; do
printf '%s\t%s\n' "$sample" "$(echo "$(zcat "$sample" | wc -l)/4" | bc )"
done
或者,您可以在以下位置完成整个操作awk
:
for sample in *.fastq.gz; do
printf '%s\t' "$sample"
zcat "$sample" | awk '!(NR % 4){k++}END{print k}'
done
但是,请注意,fastq 格式的定义中没有任何内容表明文件每个序列只有 4 行。如果你很了解你的数据,你可以使用这种方法,但如果你需要处理任意 fastq 文件,你不能假设只有 4 行,最好使用专用工具。
您可能会发现这个问答很有趣:快速计算 fastq 文件中读取次数和碱基数量的方法?。
还有FASTQ文件格式规范澄清了您不能假设每个条目只有 4 行。也就是说,根据我过去 7 年在临床环境中处理人类 NGS 数据的丰富经验,我见过的每个文件每个样本只有 4 行。但我不处理长时间读取的数据,而且格式本身允许读取更多数据,因此这是需要考虑的事情。