如何使用grep提取scaffold序列行中非A/C/T/G/N的异常字符？

Got it！你需要从这种类FASTA格式的文件里，跳过scaffold标题行，只提取序列行里不属于A/C/T/G/N的异常字符，对吧？这里有几个实用的命令行方案，帮你快速搞定：

直接用grep组合命令，过滤标题行后提取异常字符：

grep -Pv "^>" your_sequence_file | grep -Po "[^ACTGN]+"

命令解释： #

• grep -Pv "^>" your_sequence_file：先筛掉所有以>开头的scaffold标题行，只保留序列内容行。-P启用Perl兼容正则（让匹配更灵活），-v表示反向匹配（排除符合条件的行）。
• grep -Po "[^ACTGN]+"：在序列行里，只提取连续的非ACTGN字符块。-o表示仅输出匹配到的部分，[^ACTGN]+匹配一个或多个不属于A/C/T/G/N的字符（如果你的序列包含小写的actgn，把正则改成[^ACTGNactgn]+即可兼容）。

针对你给出的示例，这个命令会输出：

YYYYYYYYYYYYYYYYYY
SSSSSSSSSS

如果想明确知道异常字符来自哪个scaffold，可以用awk把标题和异常字符关联起来：

awk '/^>/ {scaffold=$0; next} {while(match($0, /[^ACTGN]+/, arr)) {print scaffold ": " arr[0]; $0=substr($0, RSTART+RLENGTH)}}' your_sequence_file

命令解释： #

• 遇到以>开头的标题行时，先把该行内容存入scaffold变量，然后直接跳到下一行处理序列。
• 处理序列行时，循环查找所有连续的异常字符块，打印对应的scaffold名称和异常字符，同时把已经处理过的部分从当前行移除，确保找到所有异常块。

针对你的示例，这个命令会输出：

>scaffold3|size67203: YYYYYYYYYYYYYYYYYY
>scaffold4|size66423: SSSSSSSSSS

如果你的系统grep不支持-P选项（比如部分BSD系统），可以用POSIX正则实现：

grep -v "^>" your_sequence_file | grep -o "[^ACTGN]" | tr -d '\n'

这个命令会把所有异常字符连在一起输出（没有按块拆分），适合只需要收集所有异常字符的场景。

内容的提问来源于stack exchange，提问作者Reda