使用AWK处理bed文件：提取重组字段生成新bed文件

解决方案脚本 #

直接运行以下AWK命令即可完成需求：

BEGIN {
    FS = "\t"
    OFS = "\t"
}
{
    # 从第6列提取CAF属性值
    split($6, attrs, /;/)
    caf_str = ""
    for (i in attrs) {
        if (attrs[i] ~ /^CAF=/) {
            caf_str = substr(attrs[i], 5)
            break
        }
    }
    if (caf_str == "") next  # 无CAF属性的行直接跳过

    # 拆分碱基列表与CAF值列表为数组
    split($5, bases, /,/)
    split(caf_str, cafs, /,/)

    # 遍历输出每行对应单个碱基与CAF值
    for (j = 1; j <= length(bases); j++) {
        print $1, $2, $3, $4, bases[j], cafs[j]
    }
}

脚本说明 #

• BEGIN块：设置输入输出分隔符为制表符，严格匹配BED文件的格式要求
• CAF提取逻辑：将第6列按分号拆分为独立属性，遍历找到以CAF=开头的项，截取前缀后的字符串得到逗号分隔的CAF值列表
• 列拆分与匹配：分别把第5列的碱基、CAF值拆分为数组，利用数组索引一一对应碱基与CAF值
• 循环输出：遍历每个索引，输出原BED的前4列，加上单个碱基和对应的CAF值，生成重组后的行

示例演示 #

输入BED行 #

chr1    1000    1001    rs123    A,T,G    AF=0.1;CAF=0.8,0.1,0.1;DP=100

输出结果 #

chr1    1000    1001    rs123    A    0.8
chr1    1000    1001    rs123    T    0.1
chr1    1000    1001    rs123    G    0.1

注意事项 #

• 若原数据中碱基数量与CAF值数量不匹配，脚本会以较短数组的长度为输出上限，可自行添加长度校验逻辑
• 如需保留原第6列的其他属性，可在print语句中追加$6或指定属性
• 若不想跳过无CAF属性的行，可删除if (caf_str == "") next，并补充对应的默认值处理逻辑

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