前面我们在介绍TCGA数据库数据挖掘的时候,课程中使用了人了所有miRNA的ID号。
TCGA数据库介绍及数据挖掘
课程网址:
https://ke.qq.com/course/package/37633
课程代码中加载了mirbase.rds这个文件,里面保存了人的所有miRNA的成熟体ID和miRNA名字。很多粉丝问这个文件是如何得到的,另外如果miRBase数据库中更新了人的miRNA数据,个数变多了,怎么才能保证这个文件是最新的。
代码语言:javascript复制#加载mirbase.rds文件,里面保存了人的所有miRNA的成熟体ID和miRNA名字
load("mirbase.rds")
其实,前面小编就用视频给大家介绍过,如何使用Excel来提取人的所有的miRNA的ID号,可能大家觉得比较麻烦。能不能把这一部分也整合到R代码中。
接下来小编就给大家讲讲如何使用R来从miRBase数据库中下载人的最新的miRNA注释信息,然后使用R来出来提取所有的miRNA的ID号。对miRBase这个数据库还不了解的小伙伴,请猛戳下面链接。
☞miRBase数据库介绍及miRNA数据下载
代码语言:javascript复制#miRBase数据库中人的miRNA的注释文件
link="https://www.mirbase.org/ftp/CURRENT/genomes/hsa.gff3"
#保存到本地的文件名
file="hsa.gff3"
#下载注释文件并保存到本地的hsa.gff3中
download.file(link,file)
#读取hsa.gff3的内容,跳过#开始的行
mir=read.table("hsa.gff3",comment.char = "#",sep="t",stringsAsFactors = F)
#第三列为miRNA的行包含成熟体信息,具体在第九列
mature=mir[mir$V3=="miRNA",9]
#根据;Alias=,;Name=,;Derives_from=来拆分第九列的内容
#提取拆分开的向量中的第二和三个元素,MIMAT0027618 hsa-miR-6859-5p
#转置之后,强制转换成数据框,去除重复
human_mirs=data.frame(unique(t(sapply(strsplit(mature,";.*?=",fixed=F),"[",2:3))))
#将miRNA的ID号和名字保存到mirbase.rds中
saveRDS(human_mirs,file="mirbase.rds")
#读取mirbase.rds中的内容,可以赋给任意变量名
mirbase=readRDS("mirbase.rds")
#查看前几行
head(mirbase)
这段代码中用了saveRDS和readRDS这样一对函数来保存和读取数据。前面小编还给大家介绍过R中另外一对函数save和load,同样可以用来保存和读取数据。
☞R的save,load函数和 .rda文件
使用这段代码能够保证,大家分析时用到的人的miRNA肯定是最新最全的。其实也没有必要每次分析之间都去下载hsa.gff3这个文件,处理一遍。只需要先去看下这个文件的表头,看看miRBase的版本和时间。你会发现其实目前最新版本的数据还是2018年的,已经有四年没有更新了。