【直播】我的基因组49:Y染色体的SNV不能用常规流程来找?

在上一次直播中,我们说到了一个不符合我们的认知的问题。就是我的全基因组测序数据里找到的SNV的纯合杂合比例失衡,这着实让我非常纠结。在朋友圈大量求助中,肿瘤所的朋友非常热心的帮我检查了她手头的几百个外显子测序样品,给了我下面这个表格,我简单的截取一部分。从这张表格中可以看到,女性样本X染色体的纯合杂合比例符合我们的认知。

不过,我更好奇女性样本的Y染色体SNV(虽然理论上女性是不可能有Y染色体的)。而且我真正想看的是男性样本的性染色体,在朋友电脑里面只有sort好的bam文件,没有vcf直接统计。所以我就借了朋友的电脑亲自上阵来统计这些指标,把所有她已有的外显子测序文件批量统计一下:

为了节省时间,我就用了bcftools来做SNP-calling,批量统计的代码如下:

  1. ls /media/cancer_path/*bam |while read id
  2. do
  3. file=$(basename $id )
  4. sample=${file%%.*}
  5. echo $sample
  6. samtools mpileup -r X -ugf /media/software/bwa/human_g1k_v37.fasta $id | bcftools call -vmO z -o $sample.chrX.vcf.gz
  7. samtools mpileup -r Y -ugf /media/software/bwa/human_g1k_v37.fasta $id | bcftools call -vmO z -o $sample.chrY.vcf.gz
  8. echo "chrX"
  9. zcat $sample.chrX.vcf.gz |perl -alne 'next if /^#/;/DP=(\d+);/;print if $1>20'|grep -v "^#" |cut -f 10|cut -d":" -f 1|sort |uniq -c
  10. echo "chrY"
  11. zcat $sample.chrY.vcf.gz |perl -alne 'next if /^#/;/DP=(\d+);/;print if $1>20'|grep -v "^#" |cut -f 10|cut -d":" -f 1|sort |uniq -c
  12. done

得到的统计表格我稍微进行了整理了(左边是男性,右边是女性):

假设朋友给我提供的性别与样本对应表格是准确无误的!

那么她提供的样本中:女性的X染色体的杂合数量远高于X的纯合。这合情合理,而且女性没有Y染色体,但是X,Y有同源区域,所以女性样本仍然会有Y染色体的SNV,也符合情理,毕竟比例很小嘛。

而她提供的男性样本数据里面出现我现在全基因组数据结果相同的困惑,明明男性只有一条X和一条Y染色体,那么上面的SNV应该是纯合的,但是这里面都是杂合的多于纯合的。跟我面临的情况一模一样!

对此,我提出了几个假设:

1.就是人类的X,Y染色体同源区域太多了,即使是PE150的建库测序策略也无法保证reads正确的匹配到参考基因组应有的位置。

2.参考基因组在这两条染色体本来就是模糊不清。

3.我们常规的SNV calling流程在,X,Y染色体上面,准确率很有限!

既然我已经在大样本里面验证了这个现象,那么可以暂时排除是公司把我的样本弄错了那个假设啦!

接下来,我就需要详细解释我自己提出的3个假设咯!

同时在这里向朋友圈给我提出各种建议的朋友表示衷心的感谢!

下面是大家的建议列表的部分摘抄:

不是说男性的就一定都是纯合的,只是男女比例不一样。这在之前的gwas中也可以观察到。甚至有可以导致性别完全误判的个例基因组型。

xy是绝大部分是同源的,这个现象正常。再有看下突变比例分布,0/1什么的说明不了太多问题。选uniq的方法是什么?最差的结果是,你的数据是混lane测的,污染了。

可能1: 女性样品污染;

可能2: 搜索gametologys evolution;

可能3: chr x link gene tends to be duplicate more .

男性中X,Y上出现0/1的情况主要是同源区域导致,这个可以从这些0/1突变所在区域发现,这些突变强烈富集,主要集中在几个同源区域。但是X,Y上1/1的突变就分布均匀很多了

对了,有朋友反映用我的samtools和bcftools代码报错,我看了一下,只是因为他们的samtools和bcftools没有升级到最新版,所以给大家提醒一下:

  1. ## Download and install samtools
  2. ## http://samtools.sourceforge.net/
  3. ## http://www.htslib.org/doc/samtools.html
  4. cd ~/biosoft
  5. mkdir samtools && cd samtools
  6. wget https://github.com/samtools/samtools/releases/download/1.3.1/samtools-1.3.1.tar.bz2
  7. tar xvfj samtools-1.3.1.tar.bz2
  8. cd samtools-1.3.1
  9. ./configure --prefix=/home/jianmingzeng/biosoft/myBin
  10. make
  11. make install
  12. ~/biosoft/myBin/bin/samtools --help
  13. ~/biosoft/myBin/bin/plot-bamstats --help
  14. cd htslib-1.3.1
  15. ./configure --prefix=/home/jianmingzeng/biosoft/myBin
  16. make
  17. make install
  18. ~/biosoft/myBin/bin/tabix
  19. ## Download and install bcftools
  20. ## http://www.htslib.org/download/
  21. ## http://www.htslib.org/doc/bcftools-1.0.html
  22. cd ~/biosoft
  23. mkdir bcftools && cd bcftools
  24. wget https://github.com/samtools/bcftools/releases/download/1.3.1/bcftools-1.3.1.tar.bz2
  25. tar xvfj bcftools-1.3.1.tar.bz2
  26. cd bcftools-1.3.1
  27. make
  28. cp bcftools /home/jianmingzeng/biosoft/myBin
  29. ~/biosoft/myBin/bin/bcftools --helpniq -c
  30. done

文:Jimmy

图文编辑:吃瓜群众

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

原文发表时间:2017-01-22

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

发表于

我来说两句

0 条评论
登录 后参与评论

相关文章

来自专栏Y大宽

Cytoscape插件3:Enrichment Map(1)

早期的基因列表解释依赖于选择一系列高得分的基因,然后建立相当主观奇怪的关系。富集分析是一个自动的,基于严格的统计学的方法来分析和解释很大的基因列表,使用的是先验...

712
来自专栏生信技能树

基因表达调控的顺式作用因子 (CREs) 了解一下

基因表达调控 基因调控是现代分子生物学研究的中心课题之一。因为要了解动植物生长发育规律、形态结构特征及生物学功能,就必须搞清楚基因表达在时间和空间上的调控机制,...

2744
来自专栏镁客网

Broad研究所公布两大研究成果,共同提高CRISPR编辑的精准度 | 黑科技

1190
来自专栏生信技能树

【直播】我的基因组56:探索遗传起源

首先,节日快乐!在这个众人狂欢的节日里,我冷静了冷静,听说知识量储备差太多的人做不了朋友,于是默默的搬起了板凳专心学习。 ? 昨天我们看了看千人基因组计划的公共...

2767
来自专栏生信技能树

比较不同单细胞转录组数据寻找features方法

挑选到的跟feature相关的基因集,有点类似于在某些组间差异表达的基因集,都需要后续功能注释。 背景介绍 单细胞转录组测序的确可以一次性对所有细胞都检测到上千...

4299
来自专栏生信技能树

RNA-seq老司机领读转录组结题报告

其实大家更关心的是数据处理问题,为此我们在前期已经推送过两篇相关内容,如果还没看过的朋友可别落下。 WGS,WES,RNA-seq组与ChIP-seq之间的异同...

3585
来自专栏生信宝典

TCGA数据库在线使用

最近做培训时整理的一部分TCGA相关数据库的使用总结。在线数据库更新改版都比较快,使用时需要参照最新的线上数据教程。不过癌症相关的数据库操作起来也都比较类似,输...

2815
来自专栏生信技能树

【直播】我的基因组50:从测序深度和位点间距来看SNV分布情况

今天的我们,还是继续探究那一个困扰我这么久的问题。为什么我作为堂堂正正的男性,明明X,Y染色体都只有一条,可是却测到了那么多的杂合突变的问题。 在之前,我们在Q...

4758
来自专栏企鹅号快讯

推荐 38个miRNA数据库

-1.psRNATarget http://plantgrn.noble.org/psRNATarget/ 一个在线的植物小RNA的靶基因预测工具,发表在核酸研...

2117
来自专栏SAP最佳业务实践

SAP最佳业务实践:SD–客户寄售(119)-1业务概览

用途 寄售货物是指属于贵公司但存储在客户地点的货物。客户无需支付这些货物的货款,直到货物从寄售库存中移除。另外,客户可以时常退回不需要的寄售货物。 优点 ...

3436

扫描关注云+社区