引言:在前面我们了解了如何使用tcgabiolinks检索并获取tcga数据库的公开数据。今天小编就用前面涉及到的代码,下载今天数据准备需要用到的tcga样本数据。
一、数据下载阶段
第一步:gdcquery()筛选我们需要的数据,tcgabiolinks包下载tcga数据进行表达差异分析-肝癌案例
library("tcgabiolinks")
query <- gdcquery(project = "tcga-lihc",
data.category = "transcriptome profiling",
data.type = "gene expression quantification",
workflow.type = "htseq - counts")
上图为通过tcga gdc链接中根据筛选条件查看的符合要求结果。下图为通过gdcquery()函数中传入对应的参数得到的结果。两者对比,我们可以发现,两者是一模一样的。说明代码执行正确。前面一期中,我们有详细谈及 gdcquery,可做参考。
samplesdown <- getresults(query,cols=c("cases"))
#getresults(query, rows, cols)根据指定行名或列名从query中获取结果,此处用来获得样本的barcode
# 此处共检索出424个barcodes
getresults()中用到的参数:
参数用法query
来自gdcquery的结果rows用于指定特定的行cols用于指定特定的列
# 从samplesdown中筛选出tp(实体肿瘤)样本的barcodes
# tcgaquery_sampletypes(barcode, typesample)
# tp代表primary solid tumor;nt-代表solid tissue normal(其他组织样本可参考学习文档)
##此处共检索出371个tp样本barcodes
datasmtp <- tcgaquery_sampletypes(barcode = samplesdown,
typesample = "tp")
# 从samplesdown中筛选出nt(正常组织)样本的barcode
#此处共检索出50个nt样本barcodes
datasmnt <- tcgaquery_sampletypes(barcode = samplesdown,
typesample = "nt")
tcgaquery_sampletypes中的参数详解:
参数用法barcodetcga中的barcodes列表typesample用于指定筛选哪种类型的组织样本,如肿瘤组织“tp”,正常组织“nt”
补充tcga中的组织样本类型:
tpprimary solid tumortmmetastatictrrecurrent solid tumortamadditional metastatictbprimary blood derived cancer-peripheral bloodthochuman tumor original cellstrbmrecurrent blood derived cancer-bone marrowtbm primary blood derived cancer-bone marrowtapadditional-new primarynb blood derived normal ntsolid tissue normalnbcbuccal cell normal???nebvebv immortalized normalnbmbone marrow normal
###设置barcodes参数,筛选符合要求的371个肿瘤样本数据和50正常组织数据
querydown <- gdcquery(project = "tcga-lihc",
data.category = "transcriptome profiling",
data.type = "gene expression quantification",
workflow.type = "htseq - counts",
barcode = c(datasmtp, datasmnt))
#barcode参数:根据传入barcodes进行数据过滤
上图为 querydown<-gdcquery()的结果,仅选择了选择371个正常组织和50个肿瘤组织样本。
第二步:gdcdownload()下载gdcquery()得到的结果
# 下载数据,默认存放位置为当前工作目录下的gdcdata文件夹中。
gdcdownload(querydown,method = "api", directory = "gdcdata",
files.per.chunk = 10)
#method ;"api"或者"client"。"api"速度更快,但是容易下载中断。
#directory:下载文件的保存地址。default: gdcdata。
#files.per.chunk = null:使用api下载大文件的时候,可以把文件分成几个小文件来下载,可以解决下载容易中断的问题。
gdcdownload(query = querydown)
说明:由于小编前面已经下载过该tcga数据,所以这里显示的是421个文件已存在。如果还没有下载的话,可能需要根据自己的网速等待一些时间。
显示这样的结果,就算下载成功啦!文件默认保存在 rstudio默认路径下的gdcdata中。前面就是我们利用第一期知识进行数据下载环节,权当温习功课吧——接下来我们就开始此期的数据处理~~
二、数据处理
第三步:gdcprepare()将前面gdcquery()的结果准备成r语言可处理的se(summarizedexperiment)文件。
#读取下载的数据并将其准备到r对象中,在工作目录生成(save=true)lihc_case.rda文件
# gdcprepare():prepare gdc data,准备gdc数据,使其可用于r语言中进行分析
dataprep1 <- gdcprepare(query = querydown, save = true, save.filename =
"lihc_case.rda")
gdcprepare()中的参数:
参数用法query来自gdcquery的结果save是否将结果保存为rdata object,默认为truesave.filename文件名,如果没有设置,系统将默认设置directory文件数据的文件夹,默认为“gdcdata”summarizedexperiment是否生成summarizedexperiment对象,默认true
第四步:tcgaanalyze_preprocessing()对数据进行预处理:使用spearman相关系数去除数据中的异常值
# 去除dataprep1中的异常值,dataprep1数据中含有肿瘤组织和正常组织的数据
# tcgaanalyze_preprocessing(object, cor.cut = 0, filename = null,
width = 1000, height = 1000, datatype = names(assays(object))[1])
# 函数功能描述:array array intensity correlation (aaic) and correlation boxplot to define outlier
dataprep2 <- tcgaanalyze_preprocessing(object = dataprep1,
cor.cut = 0.6,
datatype = "htseq - counts")
#将预处理后的数据dataprep2,写入新文件“lihc_dataprep.csv”
write.csv(dataprep2,file = "lihc_dataprep.csv",quote = false)
这里将生成一个array-array intensity correlation(aaic)相关性热图,如下:
tcgaanalyze_preprocessing()中的参数:
参数用法object来自tcgaprepare的结果cor.cut设置阈值,根据样本中各个样本之间的spearman相关系数进行过滤。默认为0filename设置生成图片文件的名称,默认为preprocessingoutput.pngwidth生成图片的宽度?? height生成图片的高度datatype描述rangedsummarizedexperiment 数据类型的字符串
第五步:tcgatumor_purity()筛选肿瘤纯度大于60%的肿瘤barcodes
# tcgatumor_purity(barcodes, estimate, absolute, lump, ihc, cpe),使用来自5种方法的5个估计值作为阈值对tcga样本进行过滤,这5个值是estimate, absolute, lump, ihc, cpe,这里设置cpe=0.6(cpe是派生的共识度量,是将所有方法的标准含量归一化后的均值纯度水平,以使它们具有相等的均值和标准差)
#筛选肿瘤纯度大于等于60%的样本数据
puritydata <- tcgatumor_purity(colnames(dataprep1), 0, 0, 0, 0, 0.6)
# filtered 为被过滤的数据, pure_barcodes是我们要的肿瘤数据
purity.lihc<-puritydata$pure_barcodes
normal.lihc<-puritydata$filtered
filtered 为被过滤的数据(为正常组织的数据barcodes), pure_barcodes是我们要的肿瘤样本barcodes。
第六步:将肿瘤表达矩阵与正常组织表达矩阵合并,进行基因注释
#获取肿瘤纯度大于60%的340个肿瘤组织样本+50个正常组织样本,共计390个样本
puried_data <-dataprep2[,c(purity.lihc,normal.lihc)]
第七步:进行表达矩阵基因注释
#基因注释,需要加载“summarizedexperiment”包,“summarizedexperiment container”每个由数字或其他模式的类似矩阵的对象表示。行通常表示感兴趣的基因组范围和列代表样品。
#if (!requirenamespace("biocmanager", quietly = true))
install.packages("biocmanager")
#biocmanager::install("summarizedexperiment") #没有的需要执行下载代码
library("summarizedexperiment")
rowdata(dataprep1) #传入数据dataprep1必须为summarizedexperiment对象
# dataframe with 56512 rows and 3 columns
# ensembl_gene_id external_gene_name original_ensembl_gene_id
# <character> <character> <character>
# ensg00000000003 ensg00000000003 tspan6 ensg00000000003.13
# ensg00000000005 ensg00000000005 tnmd ensg00000000005.5
# ensg00000000419 ensg00000000419 dpm1 ensg00000000419.11
# ensg00000000457 ensg00000000457 scyl3 ensg00000000457.12
#将结果写入文件“puried.lihc.cancer.csv”
rownames(puried_data)<-rowdata(dataprep1)$external_gene_name
write.csv(puried_data,file = "puried.lihc.csv",quote = false)
第八步:进行表达矩阵标准化和过滤,得到用于差异分析的表达矩阵
`tcgaanalyze_normalization()`使用edaseq软件包标准化mrna转录本和mirna。
#tcgaanalyze_normalization()执行edaseq包中的如下功能:
1. edaseq::newseqexpressionset
2. edaseq::withinlanenormalization
3. edaseq::betweenlanenormalization
4. edaseq::counts
datanorm <- tcgaanalyze_normalization(tabdf = puried_data,
geneinfo = geneinfo,
method = "gccontent")
tcgaanalyze_normalization中的参数:
参数用法tabdfrnaseq表达矩阵,行代表基因,列代表样本geneinfo关于genelength和gccontent的20531个基因的矩阵,“geneinfoht”和“geneinfo”可选。method选择标准化的方法,基于’gccontent’ 或 ’genelength’的标准化方法可选
#将标准化后的数据再过滤,去除掉表达量较低(count较低)的基因,得到最终的数据
datafilt <- tcgaanalyze_filtering(tabdf = datanorm,
method = "quantile",
qnt.cut = 0.25)
str(datafilt)
#num [1:13083, 1:340] 274 2432 60347 1012 1947 ...
#- attr(*, "dimnames")=list of 2
# ..$ : chr [1:13083] "a1bg" "a1cf" "a2m" "a4galt" ...
# ..$ : chr [1:390] "tcga-dd-aad5-01a-11r-a41c-07" "tcga-dd-a4no-01a-11r-a28v-07" "tcga-ep-a2ka-01a-11r-a180-07" "tcga-dd-aacp-01a-11r-a41c-07" ...
tcgaanalyze_filtering()中的参数:
参数用法tabdf数据框或者矩阵,行代表基因,列代表来自tcga的样本method用于过滤较低count数的基因的方法,有’quantile’, ’varfilter’, ’filter1’, ’filter2’qnt.cut选择均值作为过滤的阈值
最后将过滤后的数据写入文件“tcga_lihc_final.csv”,就得到我们用于后续差异分析的表达文件:
write.csv(datafilt,file = "tcga_lihc_final.csv",quote = false)
#保留的是390个样本(前340肿瘤,后50正常组织)
今天的数据预处理就讲到这里,接下来我们将分享:数据分析(差异表达分析、富集分析和聚类分析等)。如果你喜欢的话,就加入我们一起挖数据吧~~
物联网患者防走失系统的产品简介与系统构成
石墨烯在平面微型电容器中的应用进展与展望
赛梅斯凯参与奇瑞雄狮生态联盟 助推5G智能驾驶
Hifn容量优化卡成功支持Bull公司StoreWay Vi
Analog Devices ADMV4530 Ka频段上变频器在贸泽开售
如何使用TCGAbiolinks进行数据预处理?
测土配肥设备有哪些,性能是怎样的
74ls290构成31进制计数器电路图文详解
高效率、降压型μModule稳压器打造真正的大功率密度解决方案
2017年旺宏以在NOR型快闪存储器的市占率约30%成为全球霸主
磁珠是电感的一种吗 电感与磁珠的区别有哪些 磁珠和电感怎么区分
易华录高分通过“NAST-PT23001”软件性能效率测试能力验证
Starlink星座通信建模仿真分析
百度云区块链解决方案的推出,将为更多的企业提供区块链技术服务
纳芯微推出全新带弥勒钳位功能的高可靠性隔离单管驱动NSi6601M
重新审视基于FPGA的原型设计
Google服务大面积瘫痪,智能家居设备也无法幸免于难
雷达物位计的5个优点
能源危机下辅助驾驶算力技术路线展望
自动检测设备项目的设计与实现 深圳麦逊电子有限公司