linux - OSX 和 Linux 之间的 Grep 差异

Question

我编写了一个需要三个用户输入的脚本：

1. 文本文件的名称，其行对应于染色体和碱基对位置，例如extract.txt.
2. 要在文本文件中搜索并提取值的 VCF 文件的名称，即test.vcf.
3. 包含从中提取的行的新 VCF 文件的所需输出名称test.vcf

对于那些不熟悉 VCF 文件的人来说，有很多复杂的标题信息，但在这个问题上下文中唯一重要的标题行是以 . 开头的最终标题行#CHROM。这些文件之一的玩具示例如下所示（注意：在 VCF 文件中，这些列是制表符分隔的，但我不确定如何在此处格式化）：

#CHROM POS ID REF ALT QUAL FILTER INFO FORMAT
Chr06 5567134 . T C 999 PASS DP=13782;VDB=0.0302;AF1=0.7154;AC1=1312;DP4=1987,2142,4337,4552;MQ=39;FQ=999;PV4=0.49,1.5e-10,0,1 GT:PL:DP:SP:GQ
Chr06 5567140 . G A 999 PASS DP=13537;VDB=0.0304;AF1=0.7489;AC1=1374;DP4=1744,1858,4383,4621;MQ=39;FQ=999;PV4=0.8,0.068,0,1 GT:PL:DP:SP:GQ
Chr06 5567195 . G T 999 PASS DP=12016;VDB=0.0284;AF1=0.384;AC1=704;DP4=3311,4518,1537,1850;MQ=37;FQ=999;PV4=0.0026,0,0,1 GT:PL:DP:SP:GQ
Chr06 5567224 . TAAAAA TAGAAACAAAAA 999 PASS INDEL;DP=13190;VDB=0.0352;AF1=0.1229;G3=0.7854,0.1806,0.03405;HWE=1.71e-05;AC1=225;DP4=4930,5006,542,685;MQ=40;FQ=999;PV4=0.00035,1,0,1 GT:PL:DP:SP:GQ
Chr06 5567247 . T A 999 PASS DP=14383;VDB=0.03;AF1=0.1233;G3=0.7772,0.1986,0.02415;HWE=0.022;AC1=226;DP4=6484,6134,587,654;MQ=41;FQ=999;PV4=0.0062,9.3e-07,0,0.16 GT:PL:DP:SP:GQ
Chr06 5567444 . TAAAAAA TAAAAAAA 999 PASS INDEL;DP=10303;VDB=0.0235;AF1=0.6037;G3=0.1784,0.4375,0.3841;HWE=0.0162;AC1=1107;DP4=1996,2224,2158,2446;MQ=46;FQ=999;PV4=0.7,1,6.3e-30,1 GT:PL:DP:SP:GQ
Chr06 5567497 . AG A 999 PASS INDEL;DP=5243;VDB=0.028;AF1=0.07142;G3=0.873,0.08398,0.04299;HWE=0.000541;AC1=131;DP4=2297,2010,0,146;MQ=46;FQ=999;PV4=0,0,0.00097,0 GT:PL:DP:SP:GQ
Chr06 5567499 . TAAA TGGCCAAATAAGCCACTCAAAAGAAATACAGCCAAAAACATCTACAAA 999 PASS INDEL;DP=5243;VDB=0.0273;AF1=0.3508;G3=0.4662,0.2686,0.2652;HWE=1.95e-12;AC1=643;DP4=1993,1638,195,545;MQ=46;FQ=999;PV4=0,1,5.1e-18,0 GT:PL:DP:SP:GQ
Chr06 5567583 . G C 999 PASS DP=12372;VDB=0.0279;AF1=0.09276;AC1=170;DP4=4794,5928,512,569;MQ=42;FQ=999;PV4=0.095,1,5.5e-31,1 GT:PL:DP:SP:GQ
Chr06 5567628 . G T 999 PASS DP=12657;VDB=0.0244;AF1=0.1049;AC1=192;DP4=5230,6194,197,578;MQ=40;FQ=999;PV4=1.2e-29,9.8e-06,0,1 GT:PL:DP:SP:GQ

要提取的项目的 txt 文件的玩具示例如下所示：

Chr06\t5567140
Chr06\t5567583
Chr06\t5567224

我创建了这个 txt 文件以明确包含染色体和位置值之间的 VCF（制表符）分隔符。

当我在本地 OSX 计算机上测试它时，我有一个 bash 脚本。

这是我的 bash 脚本：

#!/bin/bash

# Reset getopts
OPTIND=1

usage="$(basename "$0") [-h] [-i] [-v] [-o] -- program to extract loci from VCF files

where:
    -h  show this help text
    -i  input text file containing loci to extract
      CHR and BP columns must be separated by appropriate delimiter
    -v  VCF file to extract loci from
    -o  output file name
"


# Read in variables
while getopts ":hi:v:o:" opt; do
  case "$opt" in
    h) echo "$usage"
   exit
   ;;
    i)
      input=$OPTARG
      ;;
    v)
      vcf=$OPTARG
      ;;
    o)
      output=$OPTARG
      ;;
    esac
done

echo
echo "File containing SNP positions: '$input'"
echo "VCF file to extract from:      '$vcf'"
echo "Output filename:               '$output'"


LC_ALL=C grep '#CHROM' $vcf > $output

while IFS= read -r line
do

LC_ALL=C grep "$line" $vcf >> $output

done < $input

当我将此文件上传到我的 Linux 服务器以大规模运行时，它会崩溃。while 循环中的 grep 命令 ieLC_ALL=C grep "$line" $vcf >> $output无法匹配较大 VCF 文件中的任何内容。

故障排除，我发现如果我更改文本文件以使条目不再Chr06\t5567140是，现在只是5567140. 这只是证明条目包含在我正在搜索的文件中，但是在扩展 CHR-tab-BP 模式以匹配时出现了一些问题。

显然，我可以更改所有文本文件并解决问题，但这会使我的工具不那么通用，我现在只是好奇为什么这个脚本在本地工作但在服务器上失败。

非常感谢任何帮助，或对如何改进我现有代码的评论。如果关于这个问题的任何内容不清楚，请告诉我。

谢谢！

编辑：

我得到了脚本工作。我只需要更改我的文本文件格式以使用实际的制表符而不是文字\t。using 在 OSX 上有效但在 Linux 上无效的事实\t表明 OSX grep 能够解释这个字符，而 GNU grep 则不能（正如@Sundeep 在下面的评论中指出的那样）。将此编辑放在这里，以防将来有人遇到类似问题。

标签： linuxmacosgrepvcf-variant-call-format