vcf文件提取信息及排序

vcf文件示例

VCF文件简介：
VCF（Variant Call Format）文件是一种常用的存储基因组变异信息的文本文件格式。它由一系列的元数据行（以"#"开头）和数据行组成。下面是VCF文件的主要组成部分：

1. 元数据行（Metadata lines）：以"#"开头的行，用于描述VCF文件的相关信息，如文件版本、参考基因组等。其中，以"##"开头的行为注释行，用于描述或定义特定的元数据。例如，##fileformat用于指定文件的格式，##INFO用于定义变异的相关信息。
2. 列标题行（Header line）：以"#"开头的非注释行，用于描述数据列的名称和格式。其中，前九列（CHROM、POS、ID、REF、ALT、QUAL、FILTER、INFO、FORMAT）为必需的列，其他列用于存储样本相关的信息。
3. 数据行（Data lines）：由列标题行之后的每一行组成，用于存储具体的变异信息。每一列对应于一个样本，其中前九列是固定的，后续列用于存储每个样本的基因型信息。
4. INFO字段：位于每一行的INFO列中，用于存储与变异相关的附加信息。它以分号分隔的键值对形式表示，如"DP=20;AF=0.5"，其中DP表示测序深度，AF表示突变等位基因的频率。
5. FORMAT字段：位于列标题行的FORMAT列中，用于定义每个样本的基因型信息的格式。它以冒号分隔的关键字表示，如"GT:DP:GQ"，其中GT表示基因型，DP表示测序深度，GQ表示基因质量。
   VCF文件是广泛用于存储和共享基因组变异数据的标准格式，具有可读性和灵活性，方便进行基因组变异的注释和分析。

因为我的有不同的变异类型，所有用正则把它们区分出来

open IN,"<$ARGV[0]";
open OU,">$ARGV[1]";
while (){
if($_=~/^#/){next};
chomp;
my @list=split;
my $chr=@list[0],my $pos=@list[1];

my @r=@list[0]=~/(/d+)/g;
my $r=@r[0];


  # 使用正则表达式匹配冒号前的字符
my @match=@list[7] =~ /;(.*?)=/g;
my $match=@match[0];
# 替换字符VDB为SNP
$match =~ s/VDB/SNP/g;  
# 替换字符IS为INDEL
$match =~ s/IS/INDEL/g;

print OU "$chr/t$r/t$pos/t$match/n";
}
close IN;
close OU;

匹配的结果

这里我想看变异来自哪套染色体，所以要排序，来自A和B染色体的区分开，所以这里只排序第二和第三列的内容：

sort -k2,2n -k3,3n t1 >t2

排序后的结果

import os

# 设置目录路径
directory = './'

# 统计第1列含有字母A到I的数量的字典
letter_count = {'A': 0, 'B': 0, 'C': 0, 'D': 0, 'E': 0, 'F': 0, 'G': 0, 'H': 0, 'I': 0}

# 统计第一列相同时第四列中INDEL和SNP的数量的字典
indel_snp_count = {}

# 遍历目录中的文件
for filename in os.listdir(directory):
    filepath = os.path.join(directory, filename)
    
    # 只处理bed文件
    if filename.endswith('.bed'):
        # 统计第1列含有字母A到I的数量
        with open(filepath, 'r') as file:
            for line in file:
                fields = line.strip().split('/t')
                first_column = fields[0]
                fourth_column = fields[3]
                
                # 判断第1列是否含有字母A到I
                for letter in letter_count.keys():
                    if letter in first_column:
                        letter_count[letter] += 1
                        break
                    
                # 统计第一列相同时第四列中INDEL和SNP的数量
                if first_column in indel_snp_count:
                    if fourth_column == 'INDEL':
                        indel_snp_count[first_column]['INDEL'] += 1
                    elif fourth_column == 'SNP':
                        indel_snp_count[first_column]['SNP'] += 1
                else:
                    indel_snp_count[first_column] = {'INDEL': 0, 'SNP': 0}
                    if fourth_column == 'INDEL':
                        indel_snp_count[first_column]['INDEL'] = 1
                    elif fourth_column == 'SNP':
                        indel_snp_count[first_column]['SNP'] = 1

# 输出第1列含有字母A到I的数量统计结果
print("第1列含有字母A到I的数量统计：")
for letter, count in letter_count.items():
    print(f"{letter}: {count}")

# 输出第一列相同时第四列中INDEL和SNP的数量统计结果
print("第一列相同时第四列中INDEL和SNP的数量统计：")
for first_column, counts in indel_snp_count.items():
    print(f"第一列为{first_column}时：")
    print(f"INDEL数量：{counts['INDEL']}")
    print(f"SNP数量：{counts['SNP']}")
    print()