R语言实例操作分析GEO数据库甲基化芯片

x
用微信扫描二维码
分享至好友和朋友圈

  本文首发于“百味科研芝士”微信公众号,转载请注明:百味科研芝士,Focus科研人的百味需求。

  小伙伴们,上次为大家解读了一篇GEO甲基化芯片相关的SCI文献(Aberrantly methylated-diferentially expressed genes and pathways in colorectal cancer),今天,小编打算带领大家用R软件实例操作分析GEO甲基化芯片。作为目前最大的芯片数据库,GEO数据库提供给我们了海量的数据,但是,错综复杂的数据交织在一起,如何选择数据是摆在我们面前最重要的问题,读完今天这篇文章,我相信大家都能学会GEO甲基化芯片的分析。下面,就和大家一起跑一遍R,希望大家喜欢这篇文章!

  首先是GEO甲基化芯片的下载和预处理。进入GEO数据库主页(网址https://www.ncbi.nlm.nih.gov/geo/),可以直接打开链接进入,也可以通过NCBI官网的GEO DataSets或GEO Profiles进入(进入NCBI数据库后下来搜索框)。

  

  主页搜索框输入关键词“Methylation”,点击search,出来两行英语,一般选择第一行的数字,点进去。

  

  页面跳转到GEO DataSets,这和NCBI中直接进入是一样的,直接输入需要检索的肿瘤,或者如下图点击Advanced后,在高级检索中分别输入“Colorectal cancer”和“Methylation”,根据自己的需要选择合适的甲基化芯片。如果大家知道芯片的GSE号,也可以直接根据GSE号来检索。

  

  例如,本帖根据研究内容选择了“GSE29490”这张芯片。

  

  点击芯片标题,则可弹出这张芯片的全部注释信息

  

  下拉该页面,可以看到该芯片的探针平台信息,样本信息,以及矩阵文件(名字Series Matrix File(s)的TXT格式)和原始文件(TAR格式的文件)。一般选择下载矩阵文件,如果下载原始文件,需要我们自行整理文矩阵文件,还是比较麻烦的!

  

  点击Series Matrix File(s)后,选择文件的路径点击保存。

  

  对下载好的矩阵文件解压,使用EXCEL表格打开,如下图,其中感叹号开头的是注释文件,将其删除。

  

  将注释文件删除后,把EXCEL里面的矩阵粘贴在txt文档里,命名为M.txt。

  

  同时在EXCEL中建立两列以sample,group分组的表格,对样本进行分组,C为肿瘤组,T为病例组。将表格粘贴在txt文档里面,以“group.txt”命名。

  

  将M.txt和group.txt放在同一文件夹里,文件夹名就叫M吧,文件的准备已告一段落。

  接下来,就可以打开Rstdio了(和R软件运行一样,本质也是R软件,只是界面不同),做分析之前,需要安装甲基化芯片相关的包,这个过程一般比较慢,大概2h。而通过以下的代码就可以实现安装。对于这些包的安装和功能,我们可以参考Bioconductor网站。

  

  

  分析GEO甲基化芯片需要安装如下包:

  

  加载安装包,设置工作目录(注意R软件中工作目录需要用”\\”或“/”,不可用“\”来设置工作目录),工作目录直接设置M文件夹的路径,然后读取txt文件。

  

  

  对数据进行标准化处理并输出结果

  

   我们可以看出,标准化处理之前的箱线图的中位值未处于同一水平线上。

  

  而标准化处理后的箱线图中所有样本均处于同一水平线上,使各种实验条件下的测量可以相互比较,消除测量间的非实验差异。这样所有的样本就具有可比性了。

  

  对芯片进行质量控制(QC),这里需要输出两幅质量控制的图:DensityBean图,MDS图。

  

  DensityBean图:图中可以看出,峰值主要出现在0附近,说明这张芯片的甲基化水平较低。

  

  MDS图:选出样本中1000个变异最大的位点,观察肿瘤组和对照组样本分布情况,在我们实验中,对于这种分布明显偏离组内其他样本的样本要予以剔除,因为这样的样本对实验结果影响较大。

  

  甲基化位点的差异分析

  

  dmpFinder函数,其参数设置如下,对于这些函数的设置在R软件中输入“?dmpFinder”,运行代码即跳转到参数设置页面,非常方便。

  

  打开文件夹,即可看见dmpDiff命名的甲基化位点的表格,而我们挑选的差异甲基化位点的值则是根据q-val<0.05来判定的。

  

  甲基化差异位点制作热图:这里是根据M值(甲基化的率)来做差异甲基化位点的热图。

  

  

  甲基化差异区域分析

  

   运行代码后,得到名为dmrs的表格。

  

   我们在做甲基化差异区域注释时候需要五列数据:chr,start,end,Ref,Alt。其中Ref,Alt这两列数据缺失,需要手动用0补充。

  

   补充好了之后,就可以对甲基化芯片差异区域进行注释了。

  

   甲基化区域的注释,我们采用wANNOVAR网站(http://wannovar.wglab.org/),进入主页后需要输入邮箱(机构邮箱),这里小编也是在网上随便机构邮箱,这个网站有点不地道,不支持个人邮箱!当然,只是为了使用这个网站,不需要邮箱来接受消息。Sample identifer栏随便输入英文名称即可。将上述准备好的5列数据(chr,start,end,Ref,Alt)粘贴在第三个检索框里。下拉页面,在Input Fomat栏里选择ANNOVAR,然后点击Submit提交数据。

  

  

  大概几分钟的时间,注释结束,便可下载我们刚刚注释的甲基化区域了,其中第一行为外显子区域结果,第二行为整个基因上面的注释。

  

  打开文件,我们可以看到第6列是甲基化位于基因上的位置,第7列是甲基化区域所在的基因,第8行是对甲基化基因的注释。到这里,我们针对一张GEO甲基化分析已完成。

  

   当我们得到这些甲基化的基因后,可以对这些基因进行GO、KEGG富集分析,生存分析等,或者去联合基因表达的芯片进行分析。作为最大的芯片数据库,因GEO芯片来自用户的上传,我们需要对芯片的质量做质控,并且做标准化的处理,以尽可能消除对因实验条件不同造成的实验误差,这一点在TCGA数据库中则不需要。不需要花费经费,只需要大家动脑动手就能发SCI,哈哈,这样的GEO芯片来一打!

特别声明:本文为网易自媒体平台“网易号”作者上传并发布,仅代表该作者观点。网易仅提供信息发布平台。

跟贴 跟贴 5 参与 8
© 1997-2019 网易公司版权所有 About NetEase | 公司简介 | 联系方法 | 招聘信息 | 客户服务 | 隐私政策 | 广告服务 | 网站地图 | 意见反馈 | 不良信息举报

百味科研芝士

科研百味,玩转芝士(知识)

头像

百味科研芝士

科研百味,玩转芝士(知识)

7

篇文章

0

人关注

列表加载中...
请登录后再关注
x

用户登录

网易通行证/邮箱用户可以直接登录:
忘记密码