GWAS提取显著位点信息

提取位点，要确定显著SNP上下游多长，来查找基因，这就需要计算LD衰减距离：LD衰减图绘制，然后根据上下游去和gff文件合并，把区间内的基因找到，这就找到目标基因了。

下面，用一个例子，来介绍一下操作的方法：

需要基因注释信息，SNP位点信息

对于SV，要先计算SV覆盖的范围：
计算vcf文件中SV的长度 - 简书 (jianshu.com)

###第一列为染色体，第三列为类型，第四五列为起止位置，如果第三列类型为gene，则打印染色体，起止位置
use strict;
open IN,"<$ARGV[0]";


open OU,">$ARGV[1]";
while(){
    chomp;
    my @vcf=split//t/;
    if ($vcf[2] eq 'gene') {
        print OU join("/t", $vcf[0], $vcf[3], $vcf[4]), $vcf[8]"/n";
    }
}

close IN;
close OU;

D_GW

nip.ge

使用bedtools提取信息

bedtools intersect -a D_GW -b nip.ge -loj >dgw

以下是一个Python脚本，用于从命令行中传入文件名，并提取读取文件中第七列以"ID=gene："开头的内容，并截取第一个":"之前的字符，并去重：

```python
import re
import sys

# 从命令行获取文件名
filename = sys.argv[1]

# 读取文件内容
with open(filename, 'r') as file:
    lines = file.readlines()

# 提取第七列ID=gene：后面及第一个；前面的字符
result = set()
for line in lines:
    columns = line.split('/t')
    if len(columns) >= 7:
        match = re.search(r'ID=gene:(.*?);', columns[6])
        if match:
            result.add(match.group(1))

# 输出结果
for item in result:
    print(item)

匹配的结果信息，-1表示没有匹配