By liyupeng, 31 March, 2026
Forums

基因重命名脚本

第四版(convert4.sh)

  1. 基因ID默认前缀改为TGA,前缀而已,可以直接拆开脚本进行调整
  2. 添加5端和3端的注释内容,就是区别分开,用于基于jbrower的可视化框架使用
  3. 不再对cds进行计数
  4. 在gene行上面插入空行,让注释行变得好更整洁一些
  5. 示例命令行,i指定输入注释文件,单个染色体的注释,c指定染色体序列编号,t指定需要添加的固定前缀的后半部本,之类是因为雌树雄树什么对的,所以添加F代表雌树,o指定输出文件名称
  6. 这个脚本是针对evm的整合结果和pasa更新后的注释文件设计的,所以还是具有一定通用性。
bash convert4.sh -i Chr01_sorted.gff -c 01G -t F -o Chr01_rename.gff3
#!/bin/bash

# 默认值
input_file=""
chromosome=""
suffix=""
output_file=""

# 解析命令行参数
while getopts "i:c:t:o:" opt; do
    case ${opt} in
        i) input_file=${OPTARG} ;;
        c) chromosome=${OPTARG} ;;
        t) suffix=${OPTARG} ;;
        o) output_file=${OPTARG} ;;
        \?) 
            echo "Usage: $0 [-i input_file] [-c chromosome] [-t suffix] [-o output_file]"
            exit 1
            ;;
    esac
done

# 检查输入文件是否指定
if [ -z "$input_file" ]; then
    echo "Error: Input file not specified. Use -i option."
    exit 1
fi

# 如果没有手动指定 chromosome,通过输入文件名解析染色体编号
if [ -z "$chromosome" ]; then
    base_filename=$(basename "$input_file")
    if [[ $base_filename =~ ^Chr([0-9]+)\. ]]; then
        chrom_num=$(printf "%02d" "${BASH_REMATCH[1]}")
        chromosome="${chrom_num}G"
    else
        echo "Error: Could not extract chromosome number from file name."
        exit 1
    fi
fi

# 如果没有指定输出文件名,则默认使用输入文件名加上 "_convert.gff" 后缀
if [ -z "$output_file" ]; then
    output_file="${input_file%.gff}_convert.gff"
fi

# 清除输入文件中的注释行并处理 GFF 文件
awk '!/^#/' "$input_file" | awk -v chromosome="$chromosome" -v suffix="$suffix" '
BEGIN { gene_id = 1 }
$3 == "gene" {
    id = sprintf("TGA%s%s%04d0", suffix, chromosome, gene_id);
    sub(/ID=[^;]+/, "ID=" id);
    sub(/Name=[^;]+/, "Name=" id);
    print;
    gene_map[gene_id] = id;  # 记录基因 ID
    gene_id++;
}
$3 == "mRNA" {
    parent_id = gene_map[gene_id-1];
    id = sprintf("%s-mRNA%d", parent_id, mRNA_count[gene_id-1]+1);
    sub(/ID=[^;]+/, "ID=" id);
    sub(/Parent=[^;]+/, "Parent=" parent_id);
    sub(/Name=[^;]+/, "Name=" id);
    print;
    mRNA_map[parent_id] = id;  # 记录 mRNA ID
    mRNA_count[gene_id-1]++;
}
$3 != "gene" && $3 != "mRNA" {
    parent_id = mRNA_map[gene_map[gene_id-1]];
    sub(/ID=[^;]+/, "ID=" parent_id);
    sub(/Parent=[^;]+/, "Parent=" parent_id);
    print;
}
' | awk -F'\t' '
BEGIN { OFS="\t"; }
{
    if ($3 == "exon" || $3 == "CDS") {
        split($9, id_arr, ";");
        split(id_arr[1], id_parts_arr, "=");
        id = id_parts_arr[2];

        if (!seen[id]) {
            exon_count[id] = 0;
            cds_count[id] = 0;
            seen[id] = 1;
        }

        if ($3 == "exon") {
            exon_count[id]++;
            new_id = id ":exon" exon_count[id];
        } else if ($3 == "CDS") {
            cds_count[id]++;
            new_id = id ":CDS" cds_count[id];
        }

        sub(/ID=[^;]+/, "ID=" new_id, $9);
        print $1, $2, $3, $4, $5, $6, $7, $8, $9;
    } else {
        print $0;
    }
}
' > "$output_file"

echo "Processing complete. Output saved to $output_file"