python 处理snp的vcf文件，统计snp在基因的intron、exon还是上游、下游还是不在基因及基因附近

最新推荐文章于 2024-05-16 21:04:33 发布
最新推荐文章于 2024-05-16 21:04:33 发布
阅读量4.5k 收藏 15
点赞数 3
分类专栏： 基因组 SNP Gff3 文章标签： 生物信息学 python
版权声明：本文为博主原创文章，遵循 CC 4.0 BY-SA 版权协议，转载请附上原文出处链接和本声明。
本文链接：https://blog.csdn.net/qq_26012913/article/details/110390087
版权
基因组 同时被 3 个专栏收录
2 篇文章 0 订阅
订阅专栏
SNP
1 篇文章 0 订阅
订阅专栏
Gff3
1 篇文章 0 订阅
订阅专栏

1、c处理vcf文件，初步统计snp位置
位置信息有5种

down_stream是gene下游2k以内，up_stream 是gene上游2k以内，gene_exon是snp在外显子内，gene_intron是gene在内含子内，uninfluncial是gene不在上诉所有区域

SnpType.py
内容：

#这个脚本是用来看snp哪些是在基因中（分为exon和intron），哪些是在基因的上游、下游
#最终输出结果为：scaffold position 位置 geneid
#输入文件1：snpeff结果文件（/public/home/wangwen_lab/zhangjiexiong/Ageing/SNP/snpeff/snpeffC_F/02.runC/mor_imp_strictestfilter_Snp_AncRef_selected_anno4.vcf）
#输入文件2：gff3文件（/public/home/wangwen_lab/zhangjiexiong/Ageing/SNP/genome.gff3）

import sys,re
file1 = sys.argv[1]
file2 = sys.argv[2]
file3 = sys.argv[3]
f1 = open(file1,'r')
f2 = open(file2,'r')
f3 = open(file3,'w')
flag = 'fuck'
f1dick = {}
arr1 = []
for i in f1:
        if re.match("#",i):
                continue
        a = i.split("\t")
        if flag != a[0] and flag is not 'fuck':
                f1dick[flag] = arr1
                flag = a[0]
                arr1 = []
        #       print(a[0])
                arr1.append(a[1])
        elif flag is 'fuck':
                flag = a[0]
        else:
                arr1.append(a[1])
f1dick[flag] = arr1
#for i in f1dick:
#        print(i,f1dick[i])

sc = "fuck"
f2dick1 = {}
f2dick2 = {}
f2arr1 = []
for i in f2:
        if re.match("M",i) == None:
                f2dick1[geneid] = f2arr1
                continue
        a = i.split("\t")
        if a[2] == 'gene':
                f2arr1 = []
                f2arr1.append(a[3:5])
                if sc != a[0]:
                        f2dick2[sc] = f2dick1
                        f2dick1 = {}
                        sc = a[0]
                geneid = a[8].split(";")[0].split("=")[1]
        if a[2] == 'exon':
                f2arr1.append(a[3:5])
f2dick1[geneid] = f2arr1
f2dick2[sc] = f2dick1
del f2dick2['fuck']
f2.close()
for i in f1dick.keys():
        for j in f1dick[i]:#取到scaffold内的所有snp位点
                try:
                        taq = 0
                        for k in f2dick2[i].keys():#对gff3文件中所有的scaffold为i的gene进行遍历
                                a = f2dick2[i][k]
                                amax = max(a[0])
                                amin = min(a[0])
                                flag = 0
                                if int(amin) <= int(j) <= int(amax):#判断在gene中后判断在exon还是在intron
                                        for u in a[1:]:
                                                if  int(min(u))<int(j)<int(max(u)):
                                                        flag = "gene_exon"
                                        if flag != "gene_exon":
                                                flag = "gene_intron"
                                elif int(min(a[0]))-2000 <= int(j) <=int(min(a[0])):#判断是不是在gene下游
                                        flag = 'down_stream'
                                elif int(max(a[0]))+2000 >= int(j) >= int(max(a[0])):#判断是不是gene上游
                                        flag = "up_stream"
                                else:
                                        pass
                                if flag != 0:
                                        f3.write(i+"\t"+j+"\t"+k+"\t"+flag+"\t"+max(a[0])+"\t"+min(a[0])+"\n")
                                        taq = 1
                        if taq == 0:
                                f3.write(i+"\t"+j+"\t"+"\t"+"Uninfluncial"+"\n")

                except KeyError:
                        print("KeyError",sys.exc_info())
f3.close()

输入文件为vcf和gff3
vcf内容如下：

##contig=<ID=MI_M04M24_Scaffold001_2344856_2346921_target_region_2345149_2346921,length=2066>
##contig=<ID=MI_M04M24_Scaffold018_52187_53324_target_region_52187_53324,length=1138>
##contig=<ID=MI_M04M24_Scaffold037_441019_444340_target_region_441019_444340,length=3322>
##contig=<ID=MI_M04M24_Scaffold041_338612_340340_target_region_338612_340099,length=1729>
##contig=<ID=MI_M04M24_Scaffold041_340322_342654_target_region_340322_342654,length=2333>
##contig=<ID=MI_M04M24_Scaffold114_32039_35224_target_region_32039_35224,length=3186>
##reference=file:///public/home/wangwen_lab/zhangjiexiong/Ageing/SNP/genome.fasta
##source=SelectVariants
#CHROM  POS     ID      REF     ALT     QUAL    FILTER  INFO    FORMAT  Mor_importuna
Morimp01GS001   1387    .       C       T       969.80  PASS    AC=2;AF=1.00;AN=2;BaseQRankSum=1.56;ClippingRankSum=-4.910e-01;DP=28;ExcessHet=
Morimp01GS001   1418    .       C       T       1427.78 PASS    AC=2;AF=1.00;AN=2;BaseQRankSum=-1.073e+00;ClippingRankSum=-2.400e-02;DP=43;Exce
Morimp01GS001   8939    .       A       C       282.77  PASS    AC=1;AF=0.500;AN=2;BaseQRankSum=1.67;ClippingRankSum=-5.330e-01;DP=13;ExcessHet
Morimp01GS001   9251    .       C       T       1538.77 PASS    AC=1;AF=0.500;AN=2;BaseQRankSum=0.460;ClippingRankSum=0.631;DP=88;ExcessHet=3.0
Morimp01GS001   12787   .       T       G       2514.77 PASS    AC=1;AF=0.500;AN=2;BaseQRankSum=2.25;ClippingRankSum=0.383;DP=109;ExcessHet=3.0
Morimp01GS001   13324   .       C       G       2177.77 PASS    AC=1;AF=0.500;AN=2;BaseQRankSum=2.17;ClippingRankSum=1.50;DP=88;ExcessHet=3.010
Morimp01GS001   13446   .       G       C       1158.77 PASS    AC=1;AF=0.500;AN=2;BaseQRankSum=0.994;ClippingRankSum=-1.990e-01;DP=58;ExcessHe
Morimp01GS001   13866   .       A       C       2370.77 PASS    AC=1;AF=0.500;AN=2;BaseQRankSum=-1.094e+00;ClippingRankSum=-6.080e-01;DP=100;Ex
Morimp01GS001   14181   .       T       C       1455.77 PASS    AC=1;AF=0.500;AN=2;BaseQRankSum=0.710;ClippingRankSum=-6.500e-01;DP=71;ExcessHe
Morimp01GS001   14320   .       A       G       1473.77 PASS    AC=1;AF=0.500;AN=2;BaseQRankSum=1.42;ClippingRankSum=-8.810e-01;DP=64;ExcessHet
Morimp01GS001   14344   .       C       T       1471.77 PASS    AC=1;AF=0.500;AN=2;BaseQRankSum=1.48;ClippingRankSum=-1.880e-01;DP=69;ExcessHet
Morimp01GS001   14606   .       T       G       2329.77 PASS    AC=1;AF=0.500;AN=2;BaseQRankSum=-1.280e-01;ClippingRankSum=-1.583e+00;DP=94;Exc
Morimp01GS001   14869   .       A       G       2209.77 PASS    AC=1;AF=0.500;AN=2;BaseQRankSum=-4.340e-01;ClippingRankSum=-2.854e+00;DP=99;Exc
Morimp01GS001   14897   .       C       T       2048.77 PASS    AC=1;AF=0.500;AN=2;BaseQRankSum=0.603;ClippingRankSum=1.03;DP=88;ExcessHet=3.01
Morimp01GS001   15056   .       A       T       1930.77 PASS    AC=1;AF=0.500;AN=2;BaseQRankSum=0.870;ClippingRankSum=-3.810e-01;DP=88;ExcessHe
Morimp01GS001   15384   .       A       G       1690.77 PASS    AC=1;AF=0.500;AN=2;BaseQRankSum=-1.241e+00;ClippingRankSum=0.431;DP=86;ExcessHe
Morimp01GS001   15564   .       C       T       2023.77 PASS    AC=1;AF=0.500;AN=2;BaseQRankSum=0.198;ClippingRankSum=0.014;DP=85;ExcessHet=3.0
Morimp01GS001   15729   .       C       T       2458.77 PASS    AC=1;AF=0.500;AN=2;BaseQRankSum=

gff3如下：

Morimp01GS001   AUGUSTUS        gene    16797   18718   0.81    +       .       ID=MIM04M26Gene00001;Name=MIM04M26Gene00001
Morimp01GS001   AUGUSTUS        mRNA    16797   18718   0.81    +       .       ID=MIM04M26Gene00001.t1;Parent=MIM04M26Gene00001
Morimp01GS001   AUGUSTUS        exon    16797   16949   .       +       .       ID=MIM04M26Gene00001.t1.exon1;Parent=MIM04M26Gene00001.t1
Morimp01GS001   AUGUSTUS        CDS     16797   16949   0.82    +       0       ID=MIM04M26Gene00001.t1.CDS1;Parent=MIM04M26Gene00001.t1
Morimp01GS001   AUGUSTUS        exon    17024   17219   .       +       .       ID=MIM04M26Gene00001.t1.exon2;Parent=MIM04M26Gene00001.t1
Morimp01GS001   AUGUSTUS        CDS     17024   17219   0.99    +       0       ID=MIM04M26Gene00001.t1.CDS2;Parent=MIM04M26Gene00001.t1
Morimp01GS001   AUGUSTUS        exon    17298   18718   .       +       .       ID=MIM04M26Gene00001.t1.exon3;Parent=MIM04M26Gene00001.t1
Morimp01GS001   AUGUSTUS        CDS     17298   18718   0.98    +       2       ID=MIM04M26Gene00001.t1.CDS3;Parent=MIM04M26Gene00001.t1

Morimp01GS001   AUGUSTUS        gene    45970   47158   0.55    -       .       ID=MIM04M26Gene00002;Name=MIM04M26Gene00002
Morimp01GS001   AUGUSTUS        mRNA    45970   47158   0.55    -       .       ID=MIM04M26Gene00002.t1;Parent=MIM04M26Gene00002
Morimp01GS001   AUGUSTUS        exon    47073   47158   .       -       .       ID=MIM04M26Gene00002.t1.exon1;Parent=MIM04M26Gene00002.t1
Morimp01GS001   AUGUSTUS        CDS     47073   47158   0.6     -       0       ID=MIM04M26Gene00002.t1.CDS1;Parent=MIM04M26Gene00002.t1
Morimp01GS001   AUGUSTUS        exon    46887   46999   .       -       .       ID=MIM04M26Gene00002.t1.exon2;Parent=MIM04M26Gene00002.t1
Morimp01GS001   AUGUSTUS        CDS     46887   46999   0.92    -       1       ID=MIM04M26Gene00002.t1.CDS2;Parent=MIM04M26Gene00002.t1
Morimp01GS001   AUGUSTUS        exon    46644   46784   .       -       .       ID=MIM04M26Gene00002.t1.exon3;Parent=MIM04M26Gene00002.t1
Morimp01GS001   AUGUSTUS        CDS     46644   46784   0.91    -       2       ID=MIM04M26Gene00002.t1.CDS3;Parent=MIM04M26Gene00002.t1
Morimp01GS001   AUGUSTUS        exon    45970   46559   .       -       .       ID=MIM04M26Gene00002.t1.exon4;Parent=MIM04M26Gene00002.t1
Morimp01GS001   AUGUSTUS        CDS     45970   46559   0.92    -       2       ID=MIM04M26Gene00002.t1.CDS4;Parent=MIM04M26Gene00002.t1

输入结果为c.txt
对结果再进行排序

sort -k1 -k2 -V c.txt > snp_type_F.list

snp_type_F.list结果如下:

Morimp01GS001   1418            Uninfluncial
Morimp01GS001   1476            Uninfluncial
Morimp01GS001   8939            Uninfluncial
Morimp01GS001   9251            Uninfluncial
Morimp01GS001   12787           Uninfluncial
Morimp01GS001   13324           Uninfluncial
Morimp01GS001   13446           Uninfluncial
Morimp01GS001   13866           Uninfluncial
Morimp01GS001   14181           Uninfluncial
Morimp01GS001   14320           Uninfluncial
Morimp01GS001   14344           Uninfluncial
Morimp01GS001   14606           Uninfluncial
Morimp01GS001   14869   MIM04M26Gene00001       down_stream     18718   16797
Morimp01GS001   14897   MIM04M26Gene00001       down_stream     18718   16797
Morimp01GS001   15056   MIM04M26Gene00001       down_stream     18718   16797
Morimp01GS001   15384   MIM04M26Gene00001       down_stream     18718   16797
Morimp01GS001   15564   MIM04M26Gene00001       down_stream     18718   16797
Morimp01GS001   15729   MIM04M26Gene00001       down_stream     18718   16797
Morimp01GS001   15747   MIM04M26Gene00001       down_stream     18718   16797
Morimp01GS001   15943   MIM04M26Gene00001       down_stream     18718   16797
Morimp01GS001   15981   MIM04M26Gene00001       down_stream     18718   16797
Morimp01GS001   16351   MIM04M26Gene00001       down_stream     18718   16797
Morimp01GS001   16374   MIM04M26Gene00001       down_stream     18718   16797
Morimp01GS001   16376   MIM04M26Gene00001       down_stream     18718   16797
Morimp01GS001   16472   MIM04M26Gene00001       down_stream     18718   16797
Morimp01GS001   16496   MIM04M26Gene00001       down_stream     18718   16797
Morimp01GS001   16640   MIM04M26Gene00001       down_stream     18718   16797
Morimp01GS001   16761   MIM04M26Gene00001       down_stream     18718   16797
Morimp01GS001   17077   MIM04M26Gene00001       gene_exon       18718   16797
Morimp01GS001   17340   MIM04M26Gene00001       gene_exon       18718   16797
Morimp01GS001   17470   MIM04M26Gene00001       gene_exon       18718   16797
Morimp01GS001   17478   MIM04M26Gene00001       gene_exon       18718   16797
Morimp01GS001   17963   MIM04M26Gene00001       gene_exon       18718   16797
Morimp01GS001   18077   MIM04M26Gene00001       gene_exon       18718   16797
Morimp01GS001   18153   MIM04M26Gene00001       gene_exon       18718   16797
Morimp01GS001   18448   MIM04M26Gene00001       gene_exon       18718   16797
Morimp01GS001   18449   MIM04M26Gene00001       gene_exon       18718   16797
Morimp01GS001   18798   MIM04M26Gene00001       up_stream       18718   16797
Morimp01GS001   18802   MIM04M26Gene00001       up_stream       18718   16797
Morimp01GS001   18970   MIM04M26Gene00001       up_stream       18718   16797
Morimp01GS001   19137   MIM04M26Gene00001       up_stream       18718   16797
Morimp01GS001   19292   MIM04M26Gene00001       up_stream       18718   16797
Morimp01GS001   19367   MIM04M26Gene00001       up_stream       18718   16797
Morimp01GS001   19394   MIM04M26Gene00001       up_stream       18718   16797
Morimp01GS001   19401   MIM04M26Gene00001       up_stream       18718   16797
Morimp01GS001   19446   MIM04M26Gene00001       up_stream       18718   16797
Morimp01GS001   19785   MIM04M26Gene00001       up_stream       18718   16797
Morimp01GS001   19832   MIM04M26Gene00001       up_stream       18718   16797
Morimp01GS001   19885   MIM04M26Gene00001       up_stream       18718   16797

次文件position存在重复（eg: 既在A gene内，又在B gene的上游）
处理一下
display_SnpType.py
内容如下：

#这个是用来处理SnpType.py
import sys,re,math
file1 = sys.argv[1]
file2 = sys.argv[2]
f1 = open(file1,'r')
f2 = open(file2,'w')
la = [0,0]
a = []
fuck  = []
arr =[]
arrline = []
for i in f1:
        if re.search("Uninfluncial",i):
                f2.write(i)
        else:
                line = i.strip()
                a = line.split("\t")
                if a[1] == la[1]:
                        arr.append(la[4])
                        arr.append(la[5])
                        arrline.append(la)
                else:
                        if la[1] == 0:
                                la = a
                                continue
                        arrline.append(la)
                        arr.append(la[4])
                        arr.append(la[5])
                        for j in arr:
                                c = abs(int(j)-int(la[1]))
                                fuck.append(c)
                        t = min(fuck)
                        zjx = fuck.index(t)
                        x = zjx//2
                        fuck = []
                        for zz in arrline[x]:
                                f2.write(zz+"\t")
                        f2.write("\n")
                        arrline = []
                        arr = []
                la = a
arrline.append(la)
arr.append(la[4])
arr.append(la[5])
for j in arr:
        c = abs(int(j)-int(la[1]))
        fuck.append(c)
t = min(fuck)
zjx = fuck.index(t)
x = zjx//2
fuck = []
for zz in arrline[x]:
        f2.write(zz+"\t")
f2.write("\n")
arrline = []
arr = []
f1.close()
f2.close()

python display_SnpType.py snp_type_F.list snp_type_F_clear.list

结果和处理之前差不多，就是根据距离筛选了与SNP的position关联最强的gene，一个position对应一个gene

离子回旋
关注
  • 3
    点赞
  • 踩
  • 15
    收藏
    觉得还不错? 一键收藏
  • 8
    评论
专栏目录
vcf_stats:独立 HTML5 中 VCF 文件的统计信息
06-13
vcf_stats 独立 HTML5 中 VCF 文件的统计信息 例子 跑步 输入 file1.vcf file2.vcf filen.vcf pyhton vcf_stats.py file1.vcf file2.vcf filen.vcf 输出 file1.vcf.json file2.vcf.json filen.vcf.json python vcf_stats_report.py file1.vcf.json file2.vcf.json filen.vcf.json 输出 vcf_stats_report.py_report.html open vcf_stats_report.py_report.html in your browser and select the desired file to see the report.
SNP密度分布模式
weixin_30780649的博客
03-08 1700
1. window=100k，step=2k 统计每个window的snp密度，然后用mixtools的normalmixEM（两个组分的混合模型）统计snp的分布模式。 R command： library(mixtools) SNPdensity=read.table("snp.density.file") mixmdl=normalmixEM(SNPdensity) plot(mi...
生信软件19 - vcftools高级用法技巧合辑
最新发布
LittleComputerRobot的博客
05-16 1165
vcftools高级用法
根据VCF文件探测差异SNP并设计引物序列
青笋的博客
08-04 633
本文将分享一种基于R语言从VCF文件快速获得引物设计序列的方法，用于从变异位点的vcf文件中寻找两样品的差异位点，并寻找参考基因组指定位置上下游区间，设计Marker引物序列。：如果有两个样品的测序数据，并通过GATK等上游分析得到了变异位点信息，现在我想要找任意两个样品的差异SNP，并提取该位置上下游50bp序列，用于设计引物，应该怎么做？上述R脚本提供了一种方法来识别两个材料序列之间的差异位点，并设计marker引物序列，可用于生物学研究，有助于解放双手，早点下班哈哈哈哈。8. 迭代获取参考序列信息。
VCF变异结果文件详解
qq_36608036的博客
01-16 9500
看懂变异记录结果文件（VCF） VCF 做过DNA重测序，群体遗传进化，BSA，GWAS等项目的人都会遇到VCF文件，这个文件记录了全基因组的变异信息，如果不懂VCF文件就无法进行后续分析。 VCF文件介绍: 做过DNA重测序，群体遗传进化，BSA，GWAS等项目的人都会遇到VCF文件，这个文件记录了所有样品基因组中所有位置变异（主要包括SNP和InDel）信息。后续几乎所有的分析内容都是基于此文...
从SNP_VCF文件提取SNV
weixin_47707171的博客
03-15 2136
文章目录前期了解文件分割(样本少可以不做)数据注释文件整合(样本少可以不做)注释的过滤 前期了解 VCF是测序文件的一种格式，详细记录测序样本的SNP信息，有很多介绍vcf文件内容的，在这里就不详细描述了。以下的工作主要对人类样本的SNP过滤出SNV的操作进行描述。 我收到的样本是这样子的： 共有24个，分别是24条染色体的.vcf.gz文件，每一条染色体的vcf文件都包含多个人的样本。 文件分割(样本少可以不做) VCF测序文件可能包含有多个样本，所以需要执行bcftools脚本进行分割，首先把每一条染
如何使用vcftoos呢？CHATGPT告诉你
Jasmine_X_C的博客
03-03 1739
这个命令会从 input.vcf 中统计 SNP、INDEL、MNPs 和其他变异类型的数量，并将结果保存到 output.frq.count。这个命令会从 input.vcf 中统计每个样本的杂合度和缺失率，并将结果保存到 output.het 和 output.imiss。这个命令会从 input.vcf 中添加每个位点的全局和样本群体的频率信息，并将结果保存到 output.frq。这个命令会从 input.vcf 中保留所有杂合位点，并将结果保存到 output.vcf。如何使用vcftoos呢？
基因型文件转换为vcf文件.rar_vcf中./. 基因型_基因型文件转换为vcf文件
09-15
总之，基因型文件到VCF的转换是生物信息学工作流程中的关键步骤，它允许研究人员在统一的标准格式下处理和分析遗传数据。熟悉这一过程，掌握相关的工具和技巧，对于提升数据分析效率和准确性至关重要。
python通用读取vcf文件的类（复制粘贴即可用）
09-17
Python中的VCF文件处理主要涉及基因组数据的分析，通常用于存储变异呼叫格式的数据。在生物信息学领域，VCF（Variant Call Format）是一种标准格式，用于存储基因组序列中的变异信息，如单核苷酸多态性（SNP）、插入...
vcf2DIYABC.py:将vcf文件转换为基因型数据矩阵以输入DIYABC程序的脚本
05-22
vcf2DIYABC.py目的转换DIYABC输入文件中的VCF文件以获取SNP数据的工具要求python 2.x 分支python3.x上提供了与python 3.x兼容的版本。 Windows用户：如果您的计算机上未安装python，则此脚本无法运行。 我强烈建议您...
vcf2phylip:将VCF格式的SNP转换为PHYLIP，NEXUS，二进制NEXUS或FASTA比对以进行系统发育分析
05-05
vcf2phylip 将VCF格式的SNP转换为PHYLIP，NEXUS，二进制NEXUS或FASTA比对以进行系统发育分析 简要描述;简介 该脚本以VCF文件作为输入，并将使用SNP基因型创建PHYLIP（松弛版本），FASTA，NEXUS或二进制NEXUS格式的系统发育分析矩阵。 对于杂合SNP，已达成共识，并将IUPAC核苷酸歧义码写入最终矩阵（ces），允许任何倍性水平并自动检测。 该代码针对大型VCF矩阵（数百个样本和数百万个基因型）进行了优化，例如，在我们的测试中，它在约27分钟内处理了20GB VCF（约300万个SNP x 650个个体）。 该脚本的初始版本仅产生了PHYLIP矩阵，但现在我们添加了其他流行的格式，包括使用BEAST中的SNAPP插件运行SNP分析的二进制NEXUS文件（仅适用于二倍体基因型）。 此外，您可以为每个SNP选择最少数量的样本，以控制丢失数
论文研究 - 全基因组适应中单核苷酸多态性（SNP）生物学的数学建模
05-25
我们的双盲方法的统计能力在标记具有特定环境参数的多个人群中基于SNP信息的电位的数学函数相关性时得到证明。 此外，我们的方法为常见变异生物学的新发现提供了见识。 本文演示了种群多样性的生物物理模型的功能，...
怎么对一组SNP 数据进行统计（频率、哈温平衡检验）
邓飞----育种数据分析之放飞自我
02-28 1万+
怎么对一组SNP 数据进行统计（频率、哈温平衡检验）library(SNPassoc)Loading required package: haplo.stats Loading required package: survival Loading required package: mvtnorm Loading required package: parallel data(SNPs)查看SNPs
GWAS数据分析流程—SNP、Indel注释
weixin_44482816的博客
09-28 4287
SNP和Indel注释的两种方法：annovar、snpeff 推荐使用annovar
生物基因数据文件——vcf格式详解
热门推荐
vickyleexy's blog
04-24 6万+
转载于：1.什么是vcf文件VCF是用于描述SNP，INDEL和SV结果的文本文件。在GATK软件中得到最好的支持，当然SAMtools得到的结果也是VCF格式，和GATK的CVF格式有点差别。2.VCF的主体结构##fileformat=VCFv4.0 ##FILTER= ##FORMAT= ##FORMAT= ##FORMAT= ##FORMAT= ##FORMAT= ##INFO= ##IN
2021.06.08|提取、比较各样品vcf文件中snp突变频率
06-08 5769
目录摘要环境与方法使用代码分析结果总结 摘要 接到一个wgs项目，要帮助客户统计vcf文件中snp突变频率，比较两个样品的突变位点。这个工作在上一个项目中是手动处理的，当时参考序列短，突变位点少。这次经过比对后，发现了有个样品有上万个snp位点，肯定不能用手动处理的方式。因此，写了一个脚本来统计各个样品的突变频率。 需要统计的信息 包括染色体，突变位置，参考位点，各样品突变位点，突变率（AD杂合位点覆盖度/DP总覆盖度） 环境与方法 python 3.7 R version 3.6.1 使用代码 统
VCF格式详解
VincentLuo91的博客
05-23 1万+
1. 什么是VCF？ VCF是用于描述SNP，INDEL和SV结果的文本文件。在GATK软件中得到最好的支持，当然SAMtools得到的结果也是VCF格式，和GATK的VCF格式有点差别。 2. VCF的主体结构 先给出一个VCF文件的范例： ##fileformat=VCFv4.0 ##FILTER= ##FORMAT= ##FORMAT= ##FORMAT= ##FORMA
基因组分析：VCF文件中位点提取以及R软件中的分析
sunfishtan的博客
09-23 1万+
VCF文件时基因组分析中最常见的文件类型，有时需要从中提取部分信息进行后续分析，在《vcf、plink文件格式互转》中我们已经提及了SNPs的提取方法： #在file.txt中， snp名字作为一列，无header，输出格式为vcf vcftools --gzvcf my.vcf --snps snps.txt --recode --recode-INFO-all --out filter.snp 在实际应用中，还有许多变化，这里列出两个例子： 1.根据染色体（CHR）和物理位置（PS）筛选SNP 有些情

“相关推荐”对你有帮助么？

  • 非常没帮助
  • 没帮助
  • 一般
  • 有帮助
  • 非常有帮助
提交
评论 8
添加红包

请填写红包祝福语或标题

个

红包个数最小为10个

元

红包金额最低5元

当前余额3.43元 前往充值 >
需支付：10.00元
成就一亿技术人!
领取后你会自动成为博主和红包主的粉丝 规则
hope_wisdom
发出的红包
实付元
使用余额支付
点击重新获取
扫码支付
钱包余额 0

抵扣说明：

1.余额是钱包充值的虚拟货币，按照1:1的比例进行支付金额的抵扣。
2.余额无法直接购买下载，可以购买VIP、付费专栏及课程。

余额充值