如何在Bash中对Fasta文件指定序列执行反向互补？

针对你需要仅对头部末尾为"C"的Fasta序列做反向互补的需求，下面是一个可直接通过管道使用的awk脚本，支持处理多行序列、兼容大小写碱基：

awk '
BEGIN {
    comp["A"]="T"; comp["T"]="A"; comp["C"]="G"; comp["G"]="C";
    comp["a"]="t"; comp["t"]="a"; comp["c"]="g"; comp["g"]="c";
}
/^>/ {
    if (prev_seq != "") {
        if (need_revcomp) {
            len = length(prev_seq);
            revcomp = "";
            for (i=len; i>=1; i--) {
                revcomp = revcomp comp[substr(prev_seq, i, 1)];
            }
            print revcomp;
        } else {
            print prev_seq;
        }
        prev_seq = "";
        need_revcomp = 0;
    }
    need_revcomp = (substr($0, length($0)) == "C") ? 1 : 0;
    print $0;
    next;
}
{ prev_seq = prev_seq $0 }
END {
    if (prev_seq != "") {
        if (need_revcomp) {
            len = length(prev_seq);
            revcomp = "";
            for (i=len; i>=1; i--) {
                revcomp = revcomp comp[substr(prev_seq, i, 1)];
            }
            print revcomp;
        } else {
            print prev_seq;
        }
    }
}
'

工作逻辑 #

1. BEGIN块：预先定义所有大小写碱基的互补映射，确保兼容不同格式的序列。
2. 头部行处理：
   • 遇到>开头的头部时，先输出上一条积累完成的序列（根据标记决定是否反向互补）。
   • 检查当前头部的最后一个字符是否为"C"，设置need_revcomp标记。
   • 直接输出当前头部。
3. 序列行处理：将多行序列拼接成单个字符串，避免拆分处理的问题。
4. END块：处理最后一条序列（循环结束后没有后续头部触发输出，需要单独处理）。

使用方式 #

直接将该脚本接在输出Fasta的管道后即可，比如：

# 处理本地Fasta文件
cat your_sequences.fasta | awk '...' > aligned_ready.fasta

# 处理其他命令的标准输出
your_fasta_generator_command | awk '...'

示例验证 #

输入 #

seq1_C
ATCGGATC
seq2_T
GGATCCGA
seq3_C
TTAAGCGC

输出 #

seq1_C
GATCCGAT
seq2_T
GGATCCGA
seq3_C
GCGCTTAA

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