refseqID转换为geneID（2018-06-05）

简单点lili已关注
0.2492018.06.05 13:39:02字数 191阅读 3,260

1 refseqID转换成geneID

下载gene2refseq.gz 

wget -c ftp://ftp.ncbi.nih.gov/gene/DATA/gene2refseq.gz

2 提取相似物种的基因

gene2refseq.csv文件较大，有4G多，所以将blast比对过的物种中的提取出来，在进行ID转换。

  • 首先查找物种的分类号
    Glycine max 3847
    Cicer arietinum 3827
    Medicago truncatula 3880
    vitis japonicus 29760
    Homo sapiens 9606
    Citrus sinensis 2711
    Theobroma cacao 3641
    Zea mays 4577

  • grep命令(例) 

    grep ^9606 gene2refseq1 | grep ^38 gene2refseq1 | …

  • cut命令合并（应该不会写脚本） 

    cut  gene2refseq* > gene.csv

(得到的两个文件：gene.csv和db_blast3.csv)

3 R语言-merge函数

merge函数的声明:

merge( x, y, by = intersect(names(x), names(y)),
by.x = by, by.y = by, all = FALSE, all.x = all, all.y = all,
sort = TRUE, suffixes = c(".x",".y"),

incomparables = NULL, ...） 

#merge过程中两个文件匹配行一致，不然会报错。

DF1<-read.csv("lyr/rna-seq/blast/geneID/gene.csv") #读取gene.csv
DF2<-read.csv("lyr/rna-seq/blast/geneID/db_blast3.csv") #读取db_blast3.csv
dim(DF1) #看一下表格维度
dim(DF2)
merge(DF1,DF2,by ="RNA_nucleotide_accession.version", all.y = TRUE)
data1<-(merge(DF1,DF2,by ="RNA_nucleotide_accession.version", all.y = TRUE)) #将merge结果写道data中
write.csv(data1, file = "result_data.csv", quote = FALSE, row.names = FALSE) #输出文件为result_data.csv







	接着看result_data.csv中基因的功能。