查看: 138|回复: 2

[动植物重测序] 学员来稿|全基因组关联分析(GWAS)学习笔记分享(二)

[复制链接]
  • TA的每日心情

    前天 14:14
  • 签到天数: 74 天

    连续签到: 1 天

    [LV.6]常住居民II

    管理员

    Rank: 15Rank: 15Rank: 15Rank: 15Rank: 15

    主题
    164
    奥币
    329
    积分
    2590
    注册时间
    2018.4.19
    在线时间
    482 小时

    推广达人宣传达人


    发表于 2019.1.8 09:34:43 | 显示全部楼层 |阅读模式
    表型处理
    上一节(戳这里)给大家讲解了以vcftools为例如何安装GWAS所需要的软件,今天给大家讲解一下,我们的表型处理。可能你会想,表型需要处理什么,不就直接拿来做就可以了?打住,你这种思想是很危险的,一个准确的表型是进行GWAS分析的前提!!!

    我们在刚进入实验室的时候,可能干的第一件事情就是帮师兄师姐输数据,那我们如何检验我们是否由于手抖把数据输入错误了呢?最简单的办法就是把我们的数据绘制成一张频率分布直方图,看看有没有明显的偏离。 今天给大家带来一个用R快速绘制GWAS表型数据的方法,废话不多说直接上代码。 我们的表型数据是长这个样子滴:



    ## 首先第一步是设置工作路径,你的表型文件应该放在这个路径下
    setwd("E:/LX")
    ##第二部将表型文件读入,这里我们将表型文件命名为lecture06_trait
    trait <- read.table("lecture06_trait.txt", header = T,check.names = F, sep = "\t")
    ## subset data (extract 2009_OH only)第三步提取数据,这里将表型数据中的重复1,年份2009,地点为OH的表型提取出来赋值为trait2009_OH
    trait2009_OH <- subset(trait, Rep == 1 & Year == 2009 &Loc == "OH")
    ## histogram plot and return data results to hist_data object绘制频率分布直方图,trait2009_OH$Brix为绘制Brix这个表型, col ="red"为直方图的颜色 xlab = "Brix"横轴名字,main ="Histogram of Brix"横坐标名字。hist_data <- hist(trait2009_OH$Brix, col = "red", xlab= "Brix", main = "Histogram of Brix")

    下图就是我们绘制出的图片效果了,当然你喜欢黑色,你可以把col改为black,你喜欢蓝色,可以改为blue,这个都随您了。



    很幸运,我们这组表型值,没有偏离值。那偏离之后会有什么严重的影响呢,给大家举个栗子,下图是我们在实际分析中遇到的情况,表型绘制完频率分布直方图后如下图,那下图这种表型分布会对我们GWAS产生什么影响呢,答案是产生极强的假阳性。





    所以,我们在做GWAS之前,一定要再三确定我们的表型准确性。
    那我们如何快速的去判断表型是否有问题,如何移除并替换异常值并且对表型数据进行一定的统计分析?再或者您有一份多年多点的数据,如何消除环境效应对GWAS结果的影响,是多年分开分析?还是合在一起分析? 更详细的GWAS学习课程,尽在基迪奥的Omicshare课堂,链接http://www.omicshare.com/class/

    基因型(vcf文件)处理

    上一节给大家讲了如何对我们考察的表型进行处理,我们都知道做GWAS最主要的就是表型文件、基因型文件,那我们已经讲了表型文件的处理,这一期为大家讲解如何用vcftools对vcf文件进行处理。 一个vcf文件刚出来的时候是很青涩的,要经过各种各样的磨练,才可以成长为我们最终可堪大用的最终版本,那有哪些磨练呢,首先第一步,就是要为文件添加id,先给大家看一下不成熟的vcf文件和成熟的vcf文件是长什么样子的。






    那我们为什么要给vcf文件添加id呢,在做完GWAS之后你会得到下图这个表:


    如果你不添加id 你得到的都是一堆...那你又如何知道你找到的显著位点标记是什么呢?有问题就有解决方式,这里我们写了一个perl脚本给大家。 命令行格式如下:
    Perl 脚本 你的vcf文件 新的vcf文件名字perl lecture06_07_add_id.pllecture06_genotype.vcf lecture06_genotype.id.vcf



    运行之后你ls -t你就会发现生成一个新的vcf文件:



    我们less -S lecture06_genotype.id.vcf这个文件,你会发现id已经加好啦~


    那这个文件就可以进行过滤啦。这个时候vcf文件已经成长为完全体了,那如何过度到我们可以做GWAS分析的最终形态呢,需要进行几步过滤,首先按照质量值进行过滤,因为vcf文件会有一些低质量的SNP,会对GWAS结果产生一定的影响,然后按照miss和maf进行过滤,在GWAS过程中,稀有等位基因的突变也会对结果产生一定的影响,第三,如果您想只保留SNP或者只保留indel,也可以用vcftools进行过滤。 你可以到vcftools官网上去寻找您需要的相应命令。网址:https://vcftools.github.io/man_latest.html



    当然,对基因型的处理还有更多的方法,如果您想对GWAS有更深入的了解以及学习,也可以观看我们基迪奥推出的课程,从原理到技术所有的内容尽在其中。教程链接:http://www.omicshare.com/class/。 当然,大家也可在手机端通过基迪奥公众号底部菜单栏 视频教程 观看。



    本文作者:可乐加冰啊

                   

    本帖子中包含更多资源

    您需要 登录 才可以下载或查看,没有帐号?立即注册

    x
    回复

    使用道具 举报

  • TA的每日心情
    吃饭
    12 小时前
  • 签到天数: 368 天

    连续签到: 19 天

    [LV.9]以坛为家II

    版主

    Rank: 10Rank: 10Rank: 10

    主题
    61
    奥币
    7162
    积分
    2835
    注册时间
    2017.9.21
    在线时间
    254 小时

    突出贡献优秀版主论坛元老


    发表于 2019.1.8 14:57:26 | 显示全部楼层
    优秀!
    回复

    使用道具 举报

    该用户从未签到

    草履虫

    Rank: 2

    主题
    0
    奥币
    23
    积分
    4
    注册时间
    2019.1.9
    在线时间
    0 小时

    发表于 2019.1.9 09:17:17 | 显示全部楼层
    求GWAS的原始格式的基因型数据啊!就是对照组案例组的基因型编码为0,1,2的数据格式,或者未编码的数据都行啊!计算机专业要写这方向的论文,没有原始数据头大啊!自己模拟总不是个办法
    回复 支持 反对

    使用道具 举报

    您需要登录后才可以回帖 登录 | 立即注册

    本版积分规则

    快速回复 返回顶部 返回列表