【直播】我的基因组63:wegene芯片跟二代测序的简单比较

首先要说的是,wegene公司效率很赞。

从我送样品(2017-2-20)到拿到测序数据(2017-3-10),总共耗时不到20天~~~

前面我们已经介绍过我的全基因组测序结果如何模拟成芯片数据,里面详细说明了芯片数据和高通量测序的基因型区别

【直播】我的基因组59:把我的数据伪装成23andme或wegene的芯片数据

今天我终于拿到了自己的两个数据,可以比较一下这两个平台的基因型的一致性了~

首先看看wegene芯片数据咯:

还是需要简单的探索一下他们提供的芯片覆盖的染色体情况,还有给出的基因型情况。

对此,我专门录制了一个视频来说明wegene芯片基因型跟二代测序的基因型比较,我视频里面有秀探索的方法,还有代码。

我这里就说结果吧,他们的芯片包括了22条常染色体,还有X、Y、MT信息。每条染色体以及对应的探针个数如下:

然后基因型除了ATCG的组合,还有DD,DI,II这样的表明是插入或者缺失,在比较的时候,我把它们剔除掉了。

然后看看我自己的vcf:

好吧,很明显,两个文件都是有dbSNP的,所以需要用到一个公共数据如下: ~/annotation/variation/human/dbSNP/dbsnp.pos

(希望这么一点点的数据不会暴露我的隐私,唉,为了这个直播,我也是贡献了不少了)

简单的写一个脚本,就好啦~

ln ~/data/project/myGenome/fastq/variation/autochr.highQuali.dbsnp.vcf jmzeng.vcf
cat  jmzeng_wegene.txt |grep -v "^#" |grep -v "\-\-" |grep -v "DD" |grep -v "II" |grep -v "DI" |grep -v -w 'X'|grep -v -w "Y" |grep -v -w "MT"  >jmzeng_wegene.dbsnp.txt
cat jmzeng_wegene.dbsnp.txt ~/annotation/variation/human/dbSNP/dbsnp.pos |perl -lane '{if(/^rs/){$h{$F[0]}=$F[3]}else{ if(exists $h{$F[2]}) {if($h{$F[2]} eq "$F[3]$F[3]"){$gt='wild'}elsif($h{$F[2]} eq "$F[4]$F[4]"){$gt='hom'}else{$gt='het'}print "$F[2]\t$gt" } }}' >jmzeng_wegene.genotype
cat jmzeng_wegene.genotype  jmzeng.vcf |perl -alne '{next if/^#/;$h{$F[0]}=$F[1] if/^rs/; if(exists $h{$F[2]}){@tmp=split/:/,$F[9];$h{$F[2]}.="\t$tmp[0]"}}END{foreach(sort keys %h){print "$_\t$h{$_}"}}' >merge.genotype
cut -f 2-3 merge.genotype |sort |uniq -c

看得懂代码的朋友应该知道,我把X,Y,MT染色体探针都去掉了,不想比较这个,因为我本身自己二代测序的全基因组vcf文件的准确性不敢保证。还有就是我只比较了wild,het,hom,其实并没有考虑,如果同一个位点,在芯片和测序结果是突变成了不同的碱基的情况,如果要考虑,会更复杂一点,下次再讲吧。

结果如下:

中间文件:

首先过滤掉wegene数据里面的头文件还有那些没有被dbSNP数据库收录的位点,还有没有被检测到的位点,还有X,Y,MT染色体探针。这样过滤了3万多位点。

然后根据dbSNP数据库文件把wegene的芯片基因型转换成wild,het,hom,因为我的vcf文件里面没有记录的就是wild,记录0/1的就是het,记录1/1的就是hom的mutation

很明显可以看到,我的全基因组测序数据跟自己芯片数据一致性非常好!

含有2的那些不需要管,是一个位点有多种突变形式,本来就不准确。

在总共的534239位点里面,有266728个位点,在芯片测序结果和WGS的VCF里面都是野生型。有107091都是纯合突变,有151228是杂合突变!

也就是说,基因型一致性达到98.27942%

既然一致性这么好,我非常肯定我的测序没有问题啦~~~

接下来我们就可以毫无顾忌的进行各种各样的解读啦!!!

友情提示:该视频请在wifi下观看,土豪请随意!

可能腾讯视频播放比较模糊,可以考虑下载到本地观看:

链接:https://pan.baidu.com/s/1c1obJg 密码:zs2w

最后,点击阅读原文可以直达wegene官网查看他家产品详情,还有,预告一下,下周我们会放出一定数量的wegene优惠券,希望大家跟我一样来探索自己的基因组数据,感兴趣的朋友可以在下面留言。

文:Jimmy

图文编辑:吃瓜群众

原文发布于微信公众号 - 生信技能树(biotrainee)

原文发表时间:2017-03-12

本文参与腾讯云自媒体分享计划,欢迎正在阅读的你也加入,一起分享。

发表于

我来说两句

0 条评论
登录 后参与评论

相关文章

来自专栏大数据

看看那些不在gnomAD数据库出现的常见人群变异位点是什么

前面我们说到了对3784343个的SNP位点来说,3353921个因为人群频率大于了0.05会被过滤掉,还剩下430304值得我好好研究一下。 那么,现在就开始...

20010
来自专栏SDNLAB

SDN实战团分享(三十):解读DC中的overlay与underlay

企业在上云的时候,一般不会抛弃现有的物理服务器与物理网络设备,而选择完全的虚拟化环境。其原因有如下几点:1. 保护存量投资,进行增量部署;2. 一些特殊类型的工...

4956
来自专栏生信宝典

生信宝典之傻瓜式 (二) 如何快速查找指定基因的调控网络

我是谁?我在哪儿?我在查什么? 在信息爆炸的时代,相信很多小伙伴在查文章时会因信息量太大而抓狂。今天带来一款设计简洁、功能全面的基因功能查询工具,助你事半功倍,...

2086
来自专栏腾讯大讲堂的专栏

昨夜,你可是朋友圈人生赢家?

听说昨晚,全人类都在微信朋友圈做同一件事——找广告,由此引发了如下惨绝(dou)人寰(bi)的一幕幕: ? ? ? ? ? ? ? ? ? ? ? ? 文章转载...

20510
来自专栏轮子工厂

Linux系统的前世今生

上世纪六十年代,人们还在用批处理计算机,也就是一次性给一批任务到计算机,然后等待结果,中途不能和计算机进行交互,而且准备作业需要耗费大量时间。于是1965年,贝...

792
来自专栏ionic3+

入驻简书一周年

不知不觉一年了,还记得简书写下第一篇文章的契机,是作为某技术群管理员,想帮助新人,又不想重复回答相同问题,同时也为了描述清楚和便于查看。

1073
来自专栏Python爬虫与算法进阶

大佬,我代码哪错了?

问题无处不在 我: “大佬,帮我看看这个问题错在哪了?” 大佬: “你的代码呢、你的错误提示呢?” 我: “好的,我发给你” 大佬: “。。。 再见” 留下一脸...

3468
来自专栏生信宝典

UCSC XENA - 集大成者(TCGA, ICGC)

TCGA有自己的一批工具,ICGC也有自己的网站,但好的资源都是要整合起来,整合越多越好(虽然事实不一定如此,但有这个想法的人不少),用着才更方便。这就靠今天介...

6143
来自专栏Golang语言社区

【Golang语言社区】--无线视频监控传输原理

监视前端一体摄像机安装在圆形防护罩内,内带全方位云台。摄制的图像转换成视频信号传输到微波发射机的调制端,微波发射机将其加载到载波上,经微波天线定向辐射到监控中心...

4559
来自专栏生信技能树

用gnomDB数据库对个人vcf变异文件进行过滤

直播我的基因组前面的上游分析到此为止了,这里是一个分界线,经过孜孜不倦的探索挖掘我已经拿到了我个人基因组跟hg19参考基因组的全部差异位点,而且可以肯定方法学上...

4337

扫码关注云+社区