ixxmu
commented
1 year ago 两个系统整理的转录组数据库：DEE2和recount3 by 生信菜鸟团
相信一些小伙伴在使用GEO上的转录组数据时常常会苦恼作者没有上传整理好的表达量矩阵，需要自己获得SRA或者fastq文件进行一些上游分析。这里我们将介绍两个经过系统整理的转录组数据库：DEE2和recount3，其中recount3在上一期推文中已经提到。
Digital expression explorer 2: a repository of uniformly processed RNA sequencing data
RNA测序是一种广泛用于研究基因表达的技术，但是对于非专业人士来说，处理和分析RNA测序数据可能会很困难。因此，本文介绍了Digital Expression Explorer 2（DEE2），这是一个在线工具，旨在帮助用户轻松地访问、搜索、下载和分析RNA测序数据。DEE2提供了来自9个重要模式生物的大量RNA测序数据，并且可以通过关键词或访问号搜索特定物种的RNA测序数据。此外，DEE2还提供了多种功能，例如基因和转录本水平的表达计数、质量控制信息和可视化工具等，以帮助用户更好地理解和分析RNA测序数据。因此，本文旨在为研究人员提供一个易于使用且功能强大的工具来处理和分析RNA测序数据。
DEE2是一个包含统一处理的RNA测序数据的存储库，提供基因水平和转录物水平的表达计数。它包含来自58万个RNA-seq数据集的5.3万亿个指定读数，包括大肠杆菌、酵母、拟南芥、蠕虫、果蝇、斑马鱼、大鼠、小鼠和人类。用户可以使用登录号和关键字搜索来快速识别感兴趣的数据集，还可以使用专门设计的R包进行编程访问数据,并且DEE2数据与edgeR或DESeq等统计软件包兼容。
此外，DEE2还提供了一个易于使用的Web界面，用户可以使用关键词或访问号搜索特定物种的RNA测序数据，并选择要下载的数据集。我们先来看看这个web界面https://dee2.io/
可以发现，DEE2是从NCBI SRA 获得公共数据并挖掘、统一处理的RNA-seq数据的存储库，并且提供了一个R包getDEE2
网页分为了5个部分:HOME ABOUT HELP BULK GENESETS
（1）HOME
在HOME页面中，作者提供了一个在线使用登录号和关键字搜索来快速识别感兴趣的数据集的功能
在使用登录号检索的输入栏，我们测试了网页提供的GSE63462,SRP070529,SRR401430 for A. thaliana .示例。虽然作者称“Search accession numbers Multiple accessions can be separated by commas.”，但最终只返回了SRP070529的序列结果。
单独通过GSE63462,SRR401430检索，发现也没有结果返回，通过人类GSE180123数据集登录号测试进行检索也没有结果返回
这可能说明，目前DEE2网页工具所涵盖的数据集并不是很多（作者有详细介绍如何进行质控过滤）。为了排除网页影响，我们使用作者提供的R包再次对作者提供的示例登录号GSE63462,SRP070529,SRR401430进行下载：
getDEE2 R package
进行R包的基本使用，并顺带对网页提供的GSE63462,SRP070529,SRR401430 for A. thaliana .示例再次进行检查
①下载R包
library(devtools)
# devtools::install_github("markziemann/dee2/getDEE2")
BiocManager::install("getDEE2")
library(getDEE2)
②寻找感兴趣的数据集
Try GSE63462,SRP070529,SRR401430 for A. thaliana.
获取athaliana的元数据
library(getDEE2)
mdat<-getDEE2Metadata("athaliana")
head(mdat)
看看GSE33569
mdat[which(mdat$GSE_accession %in% "GSE63462"),]
这么看来通过GSE63462确实找不到对应数据集，手动检查下载好的mdat数据集也发现GEO_seriers列全部都为NA，试一试其对应SRP登录号SRP050051
发现其实是有数据的！但是官网和R包直接通过GSE63462找不到，这说明DEE2并不支持多平台关联登录号检索
看看SRP070529
和网页结果一致
看看SRR401430
SRR401430也确实没有，检索其对应SRP登录号SRP010481
咄咄怪事，DEE2数据库里储存了该SRP数据却没有SRR401430通过NCBI SRA手动检查：GSE33808（SRR401413——SRR401430）
这可能与作者质控过滤或者数据管理有关
③Fetching DEE2 data using SRA run accession numbers
接下来的R包使用我们以SRP070529为例进行
mdat1<-mdat[which(mdat$SRP_accession %in% "SRP070529"),]
SRRlist<-as.vector(mdat1$SRR_accession)
SRRlist
前面说了这么多，大家最感兴趣的应该还是如何获取DEE2处理好SRA文件获得的表达矩阵, getDEE2函数就发挥了这样的功能
getDEE2(species,SRRvec,metadata,outfile="NULL",counts="GeneCounts")
该函数首先查询元数据以确保所请求的数据集存在，然后获取所请求的表达数据并构建一个SummarizedExperiment对象。
‘counts’参数控制提供的计数类型：
GeneCounts STAR基因水平计数（默认值）
TxCounts Kallisto转录本水平计数
Tx2Gene转录本计数聚合（求和）到基因水平。
如果定义了‘outfile’，则文件将被下载到指定路径。如果未定义，则文件将被下载到一个临时目录，并在使用后立即删除。
我们这里以STAR gene level counts (this is the default)为例
suppressPackageStartupMessages(library("SummarizedExperiment"))
x <- getDEE2("athaliana",SRRvec =SRRlist,metadata=mdat,counts="GeneCounts")
class(x)
这样得到了我们熟悉的eSets对象，接下来就可以使用exprs()来获取表达矩阵
当然也可以下载到本地持久化存储
x <- getDEE2("athaliana",SRRvec =SRRlist,metadata=mdat,counts="GeneCounts",outfile = 'DEE_count_data.zip')
# 读取之前下载好的DEE2 data 表达矩阵
myGeneCounts<-loadGeneCounts("DEE_count_data.zip")
head(myGeneCounts)
可以发现，读取本地文件直接为表达矩阵，而非eset对象
④接下来就是熟悉的转录组下游分析
DEE2 数据非常适合用于 edgeR、DESeq2 等基因表达和通路富集工具的下游分析
回到网页，在进行检索时发现，在使用某一方式检索时需要确保另一方式的输入框清空，比如我们这里，在使用了GEO登录号检索后在没有清除的情况下又使用了关键词检索，ERROR，这应该是一个bug，不应该这样，因为作者并不支持登录号和关键词同时检索
在网页中，以SRP070529登录号检索结果为例，
选择所需SRA runs, 点击 Get Counts
会得到一个压缩包，解压后得到这样的一个文件夹：
阅读README文件：
GeneCountMatrix.tsv和GeneInfo.tsv分别提供了由STAR（Gene-level mapping, Diagnose strandedness）产生的基因水平的表达计数矩阵和Ensembl gene与gene symbol的对应关系注释
网页会将所有选中的SRR文件整合：
TxCountMatrix.tsv和TxInfo.tsv文件分别提供了由Kallisto（Transcript-level mapping）产生的转录本水平的表达计数矩阵和Ensembl transcripts与gene symbol的对应关系注释
GeneCountMatrix.tsv文件提供了每个基因在每个样本中的表达计数，而TxCountMatrix.tsv文件提供了每个转录本在每个样本中的表达计数。此外，文件夹下还提供了meta元数据和相关QC数据。
Degust在线工具
DEE2官网提供了将检索的数据直接使用Degust进行在线分析的功能
这对自己不会转录组下游分析的同学很友好，我们也可以把自己的转录组数据传上去分析，属于黑盒操作
Degust (STAR gene counts)提供基因水平的表达计数，使用STAR软件进行比对和计数。 Degust (Kallisto tx counts)提供转录本水平的表达计数，使用Kallisto软件进行比对和计数。 Degust (Kallisto tx2gene counts)提供基因水平的表达计数，使用Kallisto软件将转录本水平的表达计数聚合到基因水平。
但具体使用效果却不理想
于是尝试在Degust上手动上传自己下载好的GeneCountMatrix.tsv
根据数据集样本信息，设置好分组和重复
initial select所选分组会在分析完成后作为默认contrast进行展示
（2） about
在"ABOUT"页面，作者提供了一些数据库相关信息
特别是在Data processing中，作者详细介绍了该数据库如何对获取的原始SRA数据进行处理，包括相关软件和参数以及软件使用的参考基因组数据来源
（3） HELP
在HELP页面作者回答了一些常见的问题（FAQ）以及放了一个油管学习视频
其中作者提到了关于如何处理目前不在DEE2数据库中的数据集
首先确认物种/登录号组合正确。然后检查序列数据是否可从SRA ftp站点获得。检查研究原始数据是否已发布且未受到封锁。如果数据集是最近添加的，那么您有几个选项。您可以在自己的服务器或云上运行docker或奇点映像。这里的github README.md https://github.com/markziemann/dee2 中有相关说明。或者联系我们，我们会添加它。
（4） BULK DATA
BULK DATA页面提供了各个物种汇总数据下载链接
根据HELP页面提出的一个常见问题，“The webserver is limited to 500 datasets but my project of interest contains more than that - can I still use DEE2?”
检查项目是否已经打包，为了方便起见，我们将通过DEE2流程完全处理过的runs(SRR)>200次的所有项目包含表达数据打包到zip文件中。如果项目已经被打包了，则可以在这里下载 https://dee2.io/huge/ 包含表达数据的ZIP文件。
如若项目没有被打包，则可以逐个提交 500 个运行 ID 下载。 也可以下载所有的bulk data，再根据需要进行过滤。
在下载地址 https://dee2.io/mx/中，
数据以长格式的表格形式提供，包括列：dataset，gene和count https://dee2.io/mx/README.txt
这种长格式的表格适合加载到数据库中，而不是宽格式矩阵。
STAR 计数（Gene-level mapping, Diagnose strandedness）文件以 'se.tsv.bz2' 为前缀，而 Kallisto 估计计数（Transcript-level mapping）文件以 'ke.tsv.bz2' 前缀。
质控指标也以长表格格式提供，包括列 'dataset'，'QC metric type' 和 'QC metric result'。
而对应的元数据（metadata）可以通过这里的链接 https://dee2.io/metadata/ 获得
（5） GENE SETS
在GENE SETS页面，作者提供了经过整理的基于特定疾病表型的基因SYMBOL信息，包括癫痫、糖尿病、心脏病和病毒感染（SARS / MERS / SARS-CoV-2），这些基因适用于富集分析，如GSEA和mitch，以帮助研究人员了解这些人类疾病。这个页面还处于“Under construction”状态。
作者还检查了这些基因集与MsigDB基因集是否存在不合理的高重叠，得分最高的为SARS-CoV-2与MsigDB中的一个数据集，Jaccard分数为0.39。
同时，该页面还提供了代码和方法 https://github.com/markziemann/dee2_gene_signatures 说明是如何生成这些基因集：
在GEO中寻找与该疾病相关的研究和数据集 研究中应设计实验组和对照组 研究的数据通过DEE2的流程处理后应该通过“QC”质控，QC为“FAIL”的数据集将被排除 实验被重复，也就是在通过QC质控后，每个条件下的样本数应仍大于2 在R软件中使用专门设计的流程分析对比，首先解析对照组和实验组数据，然后使用"getDEE2"包获取基因表达计数数据，排除QC总结为"FAIL"的数据，将数据从runs到实验进行聚合，排除每个样本平均每个基因表达数量小于10的基因，排除少于1000个表达基因的样本，如果在排除之后每组少于两个重复，则将对比分析废弃。 使用DESeq2版本1.28.1进行差异表达分析，FDR<0.05的基因分为上调和下调组，删除少于10个基因的集合，聚合集合并将其写成GMT格式。
recount3: summaries and queries for large-scale RNA-seq expression and splicing
recount3: uniformly processed RNA-seq https://rna.recount.bio/ 资源包含超过750,000个公开可用的人类和小鼠RNA测序样本。为了构建这个资源，作者们开发了一个新的分布式处理系统Monorail，并使用它来处理和汇总大量RNA-seq数据。具体来说，他们处理了来自GTEx v8、TCGA和Sequence Read Archive等多个来源的763,000多个人类和小鼠测序运行数据，并编写了R/Bioconductor软件包（如recount3和snapcount）来方便研究人员访问和分析这些数据。此外，作者还对数据进行了后续质控和分析，并探索了不同组织、不同条件下基因表达和剪接模式之间的关系。
recount3是ReCount项目的第三代，原始测序数据使用Monorail系统处理，生成了覆盖范围的bigWig文件和recount统一的文本文件。虽然这些原始输出文件可以通过AWS开放数据注册表以及IDIES SciServer获得，但为了方便统计分析，作者提供了recount3 R/Bioconductor程序包的接口，为每个研究构建基因、外显子和外显子-外显子连接数的RangedSummarizedExperiment R对象。
此外，snapcount R包使得可以基于查询访问recount3和recount2的数据。覆盖范围的bigWig文件可用于使用megadepth、derfinder R包和其他工具进行无注释表达分析。
可以通过shinyapps.io : https://jhubiostatistics.shinyapps.io/recount3-study-explorer/ 在线浏览recount3项目托管的数据
托管数据的一些特征上周的推文也有介绍：你还在花800块钱进行普通bulk转录组定量吗
这里我们尝试通过官网recount3 quick start guide: https://bioconductor.org/packages/release/bioc/vignettes/recount3/inst/doc/recount3-quickstart.html 进行recount3 R包使用，但是因为遇到了一些问题，并没有成功，在这里仅作记录供感兴趣的小伙伴学习讨论
一开始我使用的是R4.2.1 Bioconductor3.1.5 默认安装的为1.6.0版本
后升级了bioconductor，重新下载和官网recount3教程1.8.0版本一致但下面的问题一直存在：
于是访问了对应github项目查看issue，并检索了相关问题：
参考链接：recount3_url' is not a valid supported URL : https://github.com/LieberInstitute/recount3/issues/23
Quick recount3 : https://rna.recount.bio/docs/quick-access.html#quick-recount3
Recount3 error with default "recount3_url" option
https://duffel.rail.bio/recount3 https://www.idies.jhu.edu/recount3/data/human/homes_index https://sciserver.org/public-data/recount3/data https://registry.opendata.aws/recount/ https://support.bioconductor.org/p/9148071/
Unable to load data from duffel.rail : https://support.bioconductor.org/p/9148272/#9148276
[BUG] Is data server functional?
https://duffel.rail.bio/recount3/human/homes_index https://github.com/LieberInstitute/recount3/issues/29
发现此问题曾反复出现，并且就在两周前再次有人反馈相同问题
初步怀疑是data server那边的问题
总结：
DEE2和recount3都是比较新的数据库，作者们都通过自己搭建的流程对转录组数据进行了系统的过滤质控和整理。虽然DEE2在我们目前使用的过程中仍发现了一些问题，但其对于小白来说仍是非常友好的，如果大家没有在GEO直接获取到整理好的表达矩阵，不妨去DEE2试一试，recount3作为一个持续项目的最新成果仍值得期待，欢迎大家就这里记录的问题进行尝试和讨论。需要注意的是，这两个数据库目前都是正在维护建设状态，这篇推文以介绍为主，所记录的问题可能有时效性，在后面学习的过程中需要注意。
小编有话说：
经常看我们转录组专辑的同学可能发现这篇推文风格换了，谢谢曾老师的信任让我临时接手了这个专题，我作为一个小白也是刚刚学习生信，所以我的表述存在不规范、论述存在错误的地方欢迎大家批评指正，也希望以后和大家共同学习进步。
https://mp.weixin.qq.com/s/twrBE7bxNxI2kvpI4G9G8Q