在bioconductor官网链接是:
实际上跑一下cytofWorkflowbioconductor官网教程就足够了,我这里把他们的教程拓展一下,以一篇发表在nature medicine杂志的文章数据为例子,演示给大家全部的流程细节。
前面我们介绍了流式细胞术数据文件标准,FCS文件类型 ,来龙去脉。大家只需要了解最新版即可,是2020年9月3日Spidlen J等人在Cytometry A杂志上提出了的FCS 3.2,该修订版满足了一些新需求和建议,并结合了十年来整个细胞仪领域的进步。FCS 3.2规范完整版可以下载http://flowcyt.sf.net/fcs/fcs32.pdf
当然了,我们如果处理的是公共数据,通常是比较老的版本,不过,只要是FCS 3.0以上就ok。比如发表在nature medicine杂志的文章《Immune profiling of human tumors identifies CD73 as a combinatorial target in glioblastoma》:
它的 Data availability 部分就清晰地列出来了cytof数据和单细胞转录组数据存放的地方:
根据链接 http://flowrepository.org/id/FR-FCM-Z2B3 ,可下载到全部的FCS文件。我这里演示其中5个文件,如下:
3.6M Nov 6 10:20 WT-UnRx 3-G2_CD45pos.fcs
13M Nov 6 10:27 WT-UnRx-1_CD45pos.fcs
25M Nov 6 10:27 WT-UnRx-2_CD45pos.fcs
14M Nov 6 10:27 WT-UnRx-4_CD45pos.fcs
15M Nov 6 10:20 WT-UnRx-5_CD45pos.fcs
请自行下载哈!
只需要把全部的下载成功的FCS文件统一存放在一个文件夹,比如 paper_fcs文件夹,然后使用read.flowSet函数即可读取,示例代码如下:
rm(list = ls())
require(cytofWorkflow)
p1='paper_fcs/'
fs1=list.files(p1,'*fcs' )
fs1
samp <- read.flowSet(files = fs1,path = p1)
samp
使用read.flowSet函数全部的FCS文件后,会产生一个对象,这里面变量名是 samp 。如下所示:
> samp
A flowSet with 5 experiments.
column names:
As75Di Bi209Di Ce140Di Center Cs133Di Dy161Di Dy162Di Dy163Di Dy164Di Er166Di Er167Di Er168Di Er170Di Eu151Di Eu153Di Event_length Gd155Di Gd156Di Gd158Di Gd160Di Ho165Di Ir191Di Ir193Di Lu175Di Nd142Di Nd143Di Nd144Di Nd145Di Nd146Di Nd148Di Nd150Di Offset Os190Di Pb208Di Pd102Di Pd104Di Pd105Di Pd106Di Pd108Di Pd110Di Pr141Di Pt195Di Residual Sm147Di Sm149Di Sm152Di Sm154Di Sn120Di Tb159Di Tm169Di Xe131Di Y89Di Yb171Di Yb172Di Yb173Di Yb174Di Yb176Di beadDist Time
>
我演示的是读取了5个文件,如下:
3.6M Nov 6 10:20 WT-UnRx 3-G2_CD45pos.fcs
13M Nov 6 10:27 WT-UnRx-1_CD45pos.fcs
25M Nov 6 10:27 WT-UnRx-2_CD45pos.fcs
14M Nov 6 10:27 WT-UnRx-4_CD45pos.fcs
15M Nov 6 10:20 WT-UnRx-5_CD45pos.fcs
他们的抗体数量是一致的,所以这个read.flowSet函数没有报错。这个samp对象呢,有5个元素,每个元素也可以独立使用 exprs 函数获取其表达量信息,如下所示:
> exprs(samp[[1]])[1:6, 1:5]
As75Di Bi209Di Ce140Di Center Cs133Di
[1,] 0 0.00000000 0.0000000 639.176 0
[2,] 0 0.07410723 0.0000000 521.575 0
[3,] 0 0.00000000 0.1077218 716.634 0
[4,] 0 0.00000000 0.0000000 821.956 0
[5,] 0 0.00000000 0.0000000 684.997 0
[6,] 0 0.69866854 0.0000000 786.601 0
每个元素,都是通过索引获取:
> lapply(1:5, function(x) dim( exprs(samp[[x]]) ))
[[1]]
[1] 15841 59
[[2]]
[1] 58421 59
[[3]]
[1] 109203 59
[[4]]
[1] 62320 59
[[5]]
[1] 67497 59
cytofWorkflow流程后续分析,都是基于这个读取好的对象哦。可以看到,每个样品的细胞数量差异非常大!!!
其实R包 HDCytoData 就内置了一些cytof数据集哈, 不同数据集,需要不同的函数来下载,所以对网速要求比较高:
library(HDCytoData)
fs <- Bodenmiller_BCR_XL_flowSet()
# 如果网络不好,也可以自行下载
# 然后:loaded into R as a flowSet using read.flowSet() from the flowCore package
这个下载到的fs变量,就就等价于前面的使用read.flowSet函数全部的FCS文件后产生的对象:
> fs
A flowSet with 16 experiments.
column names:
Time Cell_length CD3(110:114)Dd CD45(In115)Dd BC1(La139)Dd BC2(Pr141)Dd pNFkB(Nd142)Dd pp38(Nd144)Dd CD4(Nd145)Dd BC3(Nd146)Dd CD20(Sm147)Dd CD33(Nd148)Dd pStat5(Nd150)Dd CD123(Eu151)Dd pAkt(Sm152)Dd pStat1(Eu153)Dd pSHP2(Sm154)Dd pZap70(Gd156)Dd pStat3(Gd158)Dd BC4(Tb159)Dd CD14(Gd160)Dd pSlp76(Dy164)Dd BC5(Ho165)Dd pBtk(Er166)Dd pPlcg2(Er167)Dd pErk(Er168)Dd BC6(Tm169)Dd pLat(Er170)Dd IgM(Yb171)Dd pS6(Yb172)Dd HLA-DR(Yb174)Dd BC7(Lu175)Dd CD7(Yb176)Dd DNA-1(Ir191)Dd DNA-2(Ir193)Dd group_id patient_id sample_id population_id
还有更多其它内置数据集,也可以使用如下所示的函数:
ehub <- ExperimentHub() # create ExperimentHub instance
ehub <- query(ehub, "HDCytoData") # find HDCytoData datasets
md <- as.data.frame(mcols(ehub)) # retrieve metadata table
帮助文档写的很清楚: