写在前面
最近在用的包经常涉及到SummarizedExperiment格式的文件，不知道大家有没有遇到过。🤒
一开始觉得这种格式真麻烦，后面搞懂了之后发现真是香啊，爱不释手！~😜
什么是SummarizedExperiment
这种class主要包括了以下几个内容：👇
assay(), 包含一个或多个矩阵, 如行为基因名,列为样本;colData(), 对列的注释，格式为DataFrame, 如样本信息的描述;rowData()和（或）rowRanges(), 对列的注释，如rowRanges()描述基因坐标，rowData()描述每个基因差异分析后的logFC和pvalue;metadata()， 描述整个object的list;
用到的包
rm(list = ls())
library(SummarizedExperiment)
library(tidyverse)
library(RColorBrewer)
library(airway)
library(edgeR)
处理SummarizedExperiment对象
4.1 示例数据
这里我们用到airway包内的示例数据，讲解一下如何操作。🧐
这个data的基本研究设计是，用地塞米松处理人气道平滑肌细胞后进行RNA-seq。🤠
data(airway, package="airway")
se <- airway
se
4.2 查看colData
colData包含了样本或表型信息，返回的格式为DataFrame。🥰
colData(se)
4.3 提取colData的指定列
se$cell
4.4 查看列名和行名
我们看一下行名和列名。😉
colnames(se)
head(rownames(se))
4.5 查看表达矩阵名
assayNames(se)
4.6 查看表达矩阵
一个SummarizedExperiment格式的object是可以包含多个assay的。
assays(se)
4.7 查看指定assay
head(assay(se, "counts"))
4.8 rowRanges或granges
接下里是重中之重了，SummarizedExperiment允许代表不同特征的rowRanges (或granges)数据。🤩
length(rowRanges(se))
dim(se)
这里我们可以看到行特征对应了很多注释信息，这样我们在操作的时候就更加方便调取了。🤩
rowRanges(se)
4.9 获取start信息
start(rowRanges(se))
对于这种IRanges 对象， 你也可以直接使用start()函数获取，其他常见的函数还有end，width。🤩
start(se)
4.10 提取制定对象
如果我们只想获取制定条件下的SummarizedExperiment对象，可以用subsetByOverlaps()函数，或者直接使用GRanges[List]。🤓
gr <- GRanges(seqnames = "1", ranges = IRanges(start = 1, end = 10^7))
subsetByOverlaps(airway, gr)
手动创建SummarizedExperiment
5.1 读入数据
这里我准备了样本数据和counts矩阵两个文件，大家跟着我一起试一下吧。
pdata <- read.csv("./SummarizedExperiment/airway-sample-sheet.csv")
counts <- read.csv("./SummarizedExperiment/airway-read-counts.csv")
5.2 整理数据并创建SummarizedExperiment
pdata <- column_to_rownames(pdata, "Run")
counts <- column_to_rownames(counts, "Run")
se_juan <- SummarizedExperiment(t(counts), colData = pdata, rowRanges = )
se_juan
5.3 准备rowData
我们再试着把rowData加进之前的SummarizedExperiment里。🤠
这里我们用一下EnsDb.Hsapiens.v86包来获取基因的各种信息，如染色体位置、起止位点、类型、id等等，这个包以后我们再具体讲怎么用。
输出的文件为Granges，完美匹配。😁
library(EnsDb.Hsapiens.v86)
edb <- EnsDb.Hsapiens.v86
filter <- rownames(se_juan)
genes <- genes(edb)
genes <- genes[genes$gene_id %in% filter]
head(genes)
5.4 添加rowData
这里需要说一下，有的基因没有具体的位点信息等，可能和版本有关系，以后我们再讲怎么处理。🥰
rowData(se_juan) <- genes
se_juan
rowData(se_juan)
小练习
我们做个小练习，试试画个基因平均表达的boxplot吧， 还要取一下log哦。😏
assay(se_juan) %>%
log() %>%
boxplot(col = colorRampPalette(brewer.pal(8, "Set2"))(8))
<center>最后祝大家早日不卷!~</center>
需要示例数据的小伙伴，在公众号回复
SummarizedExperiment获取吧！点个在看吧各位~ ✐.ɴɪᴄᴇ ᴅᴀʏ 〰
<center> <b>📍 往期精彩 <b> </center>
📍 <font size=1>🤩 WGCNA | 值得你深入学习的生信分析方法！~</font>
📍 <font size=1>🤩 ComplexHeatmap | 颜狗写的高颜值热图代码！</font>
📍 <font size=1>🤥 ComplexHeatmap | 你的热图注释还挤在一起看不清吗！？</font>
📍 <font size=1>🤨 Google | 谷歌翻译崩了我们怎么办！？（附完美解决方案）</font>
📍 <font size=1>🤩 scRNA-seq | 吐血整理的单细胞入门教程</font>
📍 <font size=1>🤣 NetworkD3 | 让我们一起画个动态的桑基图吧~</font>
📍 <font size=1>🤩 RColorBrewer | 再多的配色也能轻松搞定！~</font>
📍 <font size=1>🧐 rms | 批量完成你的线性回归</font>
📍 <font size=1>🤩 CMplot | 完美复刻Nature上的曼哈顿图</font>
📍 <font size=1>🤠 Network | 高颜值动态网络可视化工具</font>
📍 <font size=1>🤗 boxjitter | 完美复刻Nature上的高颜值统计图</font>
📍 <font size=1>🤫 linkET | 完美解决ggcor安装失败方案（附教程）</font>
📍 <font size=1>......</font>
本文由mdnice多平台发布
