单细胞day3
原创
1 加载数据
还是放在原来的文件夹下没有动
load("obj.Rdata") #载入数据
> p1 = DimPlot(seu.obj,reduction = "umap",label = T)
> p1
> rm(list = ls()) #清空列表区域省的乱
> library(Seurat)#加载这个Seurat这个包
> load("obj.Rdata") #载入数据
> p1 = DimPlot(seu.obj,reduction = "umap",label = T) #昨天画图的口令，有一个label = true
> p1
这张图的每一个点都是一个细胞，同一个颜色的点被认为时一类细胞，那末到底是什么细胞呢，可以通过marker基因进行分析。
2marker 基因
可以把marker基因理解为特征基因，特征基因在某一类细胞中表达量比较高，在其它类细胞中表达量低。
理论上，上调和下调都算marker基因，但我们会设置 only.pos = TRUE，只关注表达量变高的。
FindAllMarkers，又是一个限速步骤，随着细胞数量的增大会更慢，它会计算出所有细胞簇的marker基因，比如计算0这一簇时，就是0 vs 剩下的1-10，1-10全部算一个组。我理解的就是找某一类和其他剩余类别的不同，工作量不小。
还有一个参数min.pct，min.pct = 0.25的意思是，以0簇为例，就是说要么在0这一簇里25%以上的细胞中表达，要么在剩下的全部簇里25%以上的细胞中表达。这样可以去掉一些表达比例太少（就是表达占比不高）的基因，比例太少的话不能代表这个簇的整体情况，其实很好理解就是你这个在这个群体里占比比较高才可以说你是marker基因。
Seurat开发了一个新的方法加快FindAllMarkers的计算速度，在初次运行时它会给出提示让你安装这个包，我的给了提示证明我并没有安装过，可以安装一下来提提速
(没有提示就说明已经安装过啦)。
2.1 R语言知识补充github的R包安装
有一些包你可能特别想用，但是不在各种管理网站上，而是放在开发者的github仓库里。这样的R包也是可以安装和使用的。
两种安装方法
1 在线安装,注意安装的时presto这个包
if(!require(devtools))install.packages("devtools")
if(!require(presto))devtools::install_github("immunogenomics/presto",upgrade = F,dependencies = T)2 本地下载一下
在github网站下载压缩包放在工作目录下,然后运行这个口令，这样安装可能会有一些依赖包安不上，那还得在进行安装，注意看报错
if(!require(presto))devtools::install_github("presto-master.zip",upgrade = F,dependencies = T)2.2 取每个细胞簇logFC排名前n的marker基因。logFC用于衡量两个样本或条件下基因表达水平的相对变化
> mks = markers %>% group_by(cluster) %>% top_n(n = 2, wt = avg_log2FC)
> g = unique(mks$gene)n=2就是每个簇取前2,这个可以自己设置，看情况吧，然后为什么要g = unique(mks$gene)这个呢，因为可能有些基因在不同簇中都被记为marker基因，有重复会报错，因此需要去除一下。
3可视化
3.1 热图
> library(ggplot2)
> DoHeatmap(seu.obj, features = g) + NoLegend()+
+ scale_fill_gradientn(colors = c("#2fa1dd", "green","#f87669"))
Scale for fill is already present.
Adding another scale for fill, which will replace the existing scale.
3.2 气泡图
DotPlot(seu.obj, features = g,cols = "RdYlBu") +
RotatedAxis()
3.3 小提琴图，每一个基因都会有一个小提琴，所以少画点
VlnPlot(seu.obj, features = g[1:6]) #提取了前六个
3.4 Featureplot，本质上还是umap图，只是换了一种着色方式，一张图对应一个基因，展示了该基因在所有细胞里的表达情况。
3.5 峰峦图,峰峦图，或者是山脊图，其实就是密度图
RidgePlot(seu.obj, features = g[1:2])
4 注释亚群
手动注释很复杂，手动注释需要自行查找每种可能的细胞类型已知的marker基因，可以从文献里查，可以从数据库里查，真的很复杂得自己看，还是自动注释把，自动注释也有问题：没有把组织、疾病类型等背景知识考虑在内 > 现成的参考数据只有人类和小鼠的 > 切换不同的参考数据就会给出不同的注释结果
4.2自动注释
常用的自动注释R包是singleR，使用的参考数据是celldex里面的
花花老师的有7个，但是我有15个很奇怪。解密了是更新啦
library(celldex)
library(SingleR)
ls("package:celldex")
z = "ref_Human_all.RData" #花花老师用的是f = "ref_BlueprintEncode.RData"
if(!file.exists(z)){
ref <- celldex::ref_Human_allData() #上面用啥这块改啥保留RData前面的内容，放在Data前面
save(ref,file = z)
}
ref <- get(load(z))
library(BiocParallel)
library(SingleR) #报错显示没找到，得加载一下
scRNA = seu.obj
test = scRNA@assays$RNA$data
pred.scRNA <- SingleR(test = test,
ref = ref,
labels = ref$label.main,
clusters = scRNA@active.ident)
pred.scRNA$pruned.labels #看一下分类
new.cluster.ids <- pred.scRNA$pruned.labels
names(new.cluster.ids) <- levels(scRNA)
library(Seurat) #"RenameIdents"是Seurat里面的
scRNA <- RenameIdents(scRNA,new.cluster.ids)
p3 <- DimPlot(scRNA, reduction = "umap",label = T,pt.size = 0.5) + NoLegend()
p3 #没操作手动注释的所以没有p2
跟花花老师的不一样是因为我用了不同的参考数据集。yeah.
原创声明：本文系作者授权腾讯云开发者社区发表，未经许可，不得转载。
如有侵权，请联系 cloudcommunity@tencent.com 删除。
原创声明：本文系作者授权腾讯云开发者社区发表，未经许可，不得转载。
如有侵权，请联系 cloudcommunity@tencent.com 删除。
腾讯云开发者
