我前面有一篇介绍TCGA数据库中miRNA数据的下载与整理的文章，在这篇文章中，我们下载的数据是
miRNAExpressionQuantification数据。
我整理得到的数据结果是这样的。
这里我们可以发现，miRNA的前体可能对应多个成熟的miRNA，比如hsa-let-7a-1，有两个对应的成熟体，MIMAT(hsa-let-7a-5p)和MIMAT(hsa-let-7a-3p)。这里的值是对所有成熟体miRNA求和的结果。
如果我们想要得到成熟体miRNA的数据，我们需要下载IsoformExpressionQuantification的数据。
下载后的数据处理方式和前面文章差不多。下载后的数据，我们打开单个文件看一下，比miRNAExpressionQuantification数据多了2列。在miRNA_region列中就有是否成熟的信息。我们可以参考前面的文章进行转换为我们常用的名称，比如：hsa-let-7a-5p。
下面我们来看一下，这些数据的差异。
加载一些需要用到的R包
rm(list=ls())options(stringsAsFactors=F)#作者：DoubleHelix，loudlibrary(rjson)#没有安装事先安装library(dplyr)#没有安装事先安装library(tidyr)library(limma)#没有安装事先安装library(stringr)#没有安装事先安装library(miRBaseVersions.db)
处理一下json文件。
###---1.处理json文件----------------------jsonFile-fromJSON(file="./metadata.cart.-03-10.json")filesNameToBarcode-data.frame(filesName=c(),TCGA_Barcode=c())for(iin1:length(jsonFile)){TCGA_Barcode-jsonFile[][["associated_entities"]][[1]][["entity_submitter_id"]]file_name-jsonFile[][["file_name"]]filesNameToBarcode-rbind(filesNameToBarcode,data.frame(filesName=file_name,TCGA_Barcode=TCGA_Barcode))}rownames(filesNameToBarcode)-filesNameToBarcode[,1]
这个和前文中一样的。
这里读入其中一个文件：
filepath-dir(path="./raw",pattern=".isoforms.quantification.txt$",full.names=TRUE,recursive=T)oneSampExp-read.table(filepath[1],header=T,check.names=F)
选择我们需要的列，并重新命名。
oneSampExp-oneSampExp[grep("mature",oneSampExp$miRNA_region),]oneSampExp-separate(oneSampExp,into=c("State","ACCESSION"),col="miRNA_region",sep=",")oneSampExp-oneSampExp[c("miRNA_ID","ACCESSION","read_count","reads_per_million_miRNA_mapped")]colnames(oneSampExp)-c("miRNA_ID","ACCESSION","count","RPM")
在众多前体下，这里获得个成熟体。据我自己所知，目前已知的miRNA成熟体应该有多，这里一个样本只检测到多，是因为不是所有miRNA在一个样本中都能检测到，你把所有的样本加起来统计大概也就多个。
acc-oneSampExp$ACCESSIONacc-acc[!duplicated(acc)]length(acc)
我们可以通过miRBaseVersions.db包进行转换。
ano-select(miRBaseVersions.db,keys=acc,keytype="VW-MIMAT-22.0",columns=c("ACCESSION","NAME"))head(ano)
与前面的数据融合。
data-merge(ano,oneSampExp,by="ACCESSION")head(data)
我们只提取hsa-let-7a-1的RPM数据。求和是。
temp-data[data$miRNA_ID=="hsa-let-7a-1",]sum(temp$RPM)
这里只是一个病人的数据，该病人的ID是：TCGA-AA--01A-01T--13
tempPath-unlist(strsplit(filepath[1],"/"))filename-tempPath[length(tempPath)]filesNameToBarcode[filename,2]
利用这个ID，我们加载前面通过文章得到的数据。提取该病人的RPM数据。当然，个别miRNA可能会有点出入，值没多大区别。
load("miRNAdata.Rdata")rpm-miRNAdata[["miRNA_RPM_Expdata"]][["Exp"]]rpm["hsa-let-7a-1",filesNameToBarcode[filename,2]]%%as.numeric()
这里的值就是一样的。也就是说，我们通过miRNAExpressionQuantification数据得到的是miRNA前体的数据。一个miRNA前体可能有多个成熟体，一个成熟体可能来自多个前体。
当然，还有其他处理方式，比如利用TCGAbiolinks包，关于TCGAbiolinks包我在很久前的一篇文章中有介绍，同时，也有差异分析的教程，但是现在这个包下载数据是真的感人，完全看GDC的心情，而且下载很慢，所以我是自己下载数据处理，我也发了自己下载提取数据进行差异分析的文章，当然还有GDCRNATools包。这个包下载数据的时候有时候也会挂，开VPN有时候可以解决。但我们可以自己下载数据后用该包的函数提取数据。如果你想偷懒，那么你可以在一些数据库中直接下载，比如：GDACFirehose数据库，以及UCSC数据库，但从UCSC下载的数据和我们前面处理的一样，是miRNAExpressionQuantification的数据，而不是IsoformExpressionQuantification。当然还是TCGAbiolinks包好用。
最后说明，上面我只是读入一个样本的数据来进行讲解，你可以写一个循环，把所有文件遍历读入，处理合并就行，和前文介绍的差不多，那样我们会得到如下这样的数据。
代码在文章，以前付费过的，你回复文章中的关键词，重新获取。
经典栏目
预览时标签不可点收录于话题#个上一篇下一篇