如何用Bash从聚类文件提取含指定字符串且有多序列的簇编号

需求说明 #

从基于100%序列一致性的聚类文件中，筛选出同时满足以下两个条件的簇编号：

1. 簇内存在序列ID包含指定搜索字符串（例如crgi）的条目；
2. 簇内除了含该搜索字符串的序列外，还有其他不含该字符串的序列。

输入文件格式示例 #

>Cluster 94
0       416aa, >crgi_XP_019918436.1... *
>Cluster 95
0       415aa, >TRINITY_DN19027_c1_g1_i3_77... *
>Cluster 96
0       415aa, >TRINITY_DN31043_c0_g1_i1_68... *
>Cluster 97
0       414aa, >TRINITY_DN4649_c0_g1_i1_24... *
1       414aa, >TRINITY_DN4649_c0_g1_i2_29... at 100.00%
2       414aa, >TRINITY_DN4649_c0_g1_i3_23... at 100.00%
3       414aa, >TRINITY_DN4649_c0_g1_i4_23... at 100.00%
4       414aa, >TRINITY_DN4649_c0_g1_i5_27... at 100.00%
>Cluster 107
0       410aa, >TRINITY_DN9528_c0_g1_i1_30... *
1       410aa, >TRINITY_DN9528_c0_g1_i2_36... at 100.00%
2       404aa, >crgi_XP_011414097.1... at 100.00%
>Cluster 109
0       410aa, >crgi_XP_011450995.2... *

预期输出（搜索字符串为crgi时） #

Clusters in file1.ids containing 'crgi':
107

原脚本问题分析 #

原脚本存在两个核心错误：

1. 簇编号提取错误：原awk命令$0 ~ search {print $1}提取的是序列行的索引数字（如0、2），而非簇编号；
2. 簇内序列判断错误：grep -c ">$search_string"统计的是整个文件中匹配字符串的行数，而非当前簇内的内容，无法准确判断簇内是否存在其他序列。

正确实现脚本 #

使用awk一次性完成簇的遍历、条件判断和结果输出，无需嵌套循环，效率更高：

#!/bin/bash

# 定义要搜索的字符串
search_string="crgi"

# 遍历所有.ids文件
for file in *.ids; do
    # 用awk处理文件，提取符合条件的簇
    result=$(awk -v search="$search_string" '
        BEGIN {
            current_cluster = ""
            has_target = 0
            has_other = 0
        }
        # 匹配簇开头行，更新当前簇编号，并重置标记
        /^>Cluster/ {
            # 如果上一个簇符合条件，输出编号
            if (current_cluster != "" && has_target == 1 && has_other == 1) {
                print current_cluster
            }
            current_cluster = $2
            has_target = 0
            has_other = 0
            next
        }
        # 处理序列行
        {
            # 判断当前序列是否包含搜索字符串
            if ($0 ~ search) {
                has_target = 1
            } else {
                has_other = 1
            }
        }
        # 文件末尾检查最后一个簇
        END {
            if (current_cluster != "" && has_target == 1 && has_other == 1) {
                print current_cluster
            }
        }
    ' "$file")

    # 如果有结果，输出格式内容
    if [ -n "$result" ]; then
        echo "Clusters in $file containing '$search_string':"
        echo "$result"
        echo
    fi
done

脚本逻辑说明 #

1. 初始化变量：在BEGIN块中定义当前簇编号、是否含目标字符串、是否含其他序列的标记；
2. 处理簇开头行：每次遇到>Cluster行时，先检查上一个簇是否符合条件，符合则输出簇编号，然后重置当前簇和标记；
3. 处理序列行：对每一行序列，判断是否包含搜索字符串，更新对应标记；
4. 文件末尾处理：检查最后一个簇是否符合条件，避免遗漏；
5. 输出结果：对每个文件，若有符合条件的簇，按要求格式输出。

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