单细胞学习

学习原因：闲着无聊，自己瞎看看

数据读取，合并 （合并的原因，这个不用说了吧，多个样本要一起看，当然了，一个样本的话，合并鸡毛啊）

samples_file <- list.file("./path") # 输入的文件的路径，每个文件夹下的文件只有  barcodes.tsv、genes.tsv、matrix.mtx 
seurat_list <- lapply(samples,function(sample_name){
  folder=file.path('./raw_data',sample_name,sep='/') # 别说这个不能用，这个就是那三个文件所在的文件夹
  CreateSeuratObject(counts = Read10X(folder),# 本次读取的就是这个10X格式的，其他格式的还没学，没空
  project = sample_name, #添加一个标签，这样那个数据是哪个样本的，一目了然
  min.cells = 3,  #这个参数是要过滤掉至少在3个细胞中表达了这个基因，为什么是3 ,官网推荐的, 个人理解就是，一个基因没有几个细胞表达，可能就是系统错误，或者是因为表达太少没什么意义
  min.features = 200 )# 表示一个细胞至少表达200个基因，理由和上面差不多，表达太少，细胞有问题，细胞碎片啥的
})
# 简简单单的读取数据 	，这个有啥问题的，就看看自己的路径啥的，文件的格式，剩下的就是你抄错了

合并数据

# 一个数据的，合并个屁啊
sec.merge <- seurat_list[[1]]
# 多个数据的合并
y <- c()
for (i in 2:length(seurat_list)){
   y <- append(y,seurat_list[[i]])
}
sec.merge <- merge(seurat_list[[1]],y,project="merge")
# ok ,已经合并了啊，那部分不理解的，应该没有了吧，就是将剩下的数据合并到第一个数据中

合并之后的数据就要开始质控了，因为数据不是同一个批次的数据，可能出现什么批次效应啊，什么线粒体、红细胞啥的都要过滤掉

# 为啥去除线粒体、红细胞？怎么去？
#线粒体是因为线粒体会干扰分群,线粒体占比太高，细胞可能凋亡或坏死。线粒体是独立遗传的,不是染色体上基因控制的。
#去除红细胞是因为红细胞无细胞核,基因少,与后续研究无联系
#基因数太少(barcode 是细胞碎片),很高的话就是多个细胞,也需要去掉
#网上查的 [线粒体 人：MT 小鼠：mt]、[核糖体 人：RP[SL] 小鼠：rp[sl]]、[血红 人："^HB[^(P)]"小鼠："^hb[(^p)]"]
object[["pMT"]] <- PercentageFeatureSet(object,pattern = "^MT-") # 一般情况下人的线粒体(^MT-)，小鼠的(^mt)，计算线粒体基因百分比
# 如果不知道是人的样本还是小鼠的，这个代码执行之后会有提示，All cells have the same value of pMT 
# 红细胞这个怎么算的，有问题，需要好好看注释！！！错了活该，因为我这个红细胞的是我抄的
HB.genes <- c("HBA1","HBA2","HBB","HBD","HBE1","HBG1","HBG2","HBM","HBQ1","HBZ")
HB_m <- match(HB.genes, rownames(scRNA@assays$RNA)) 
HB.genes <- rownames(scRNA@assays$RNA)[HB_m] 
HB.genes <- HB.genes[!is.na(HB.genes)] 
scRNA[["percent.HB"]]<-PercentageFeatureSet(scRNA, features=HB.genes) 
# 我是这样过滤的，这个个人感觉好像也没什么问题
object[["percent.HB"]]<-PercentageFeatureSet(object, pattern = "^HB[^(P)]")

先定义点颜色，画图用的，问为啥的，看啥啊

my36colors <- c('#E5D2DD', '#53A85F', '#F1BB72', '#F3B1A0', '#D6E7A3', '#57C3F3', '#476D87', 
                '#E95C59', '#E59CC4', '#AB3282', '#23452F', '#BD956A', '#8C549C', '#585658',
                '#9FA3A8', '#E0D4CA', '#5F3D69', '#C5DEBA', '#58A4C3', '#E4C755', '#F7F398',
                '#AA9A59', '#E63863', '#E39A35', '#C1E6F3', '#6778AE', '#91D0BE', '#B53E2B',
                '#712820', '#DCC1DD', '#CCE0F5', '#CCC9E6', '#625D9E', '#68A180', '#3A6963','#968175')

my27colors <- c('#CD1E24','#262C68','#1E843F','#8499CA','#84278B','#EF7A2A','#F4DC00',
                '#684797','#8CCEDD','#0B6E78','#EA7678','#CEA063','#16FF32','#1CFECE',
                '#7B1315','#393939','#85C25D','#38B4A2','#947EB5','#CC4826','#3283FE',
                '#B968A5','#DCBAD3','#F8AA18','#BDD240','#A1C6DB',"#66B0FF")

my63colors <- c(my27colors,my36colors)
# 就是单纯的调色啊，画图用的，没别的事儿

画个图，一个小提琴图、一个DimPlot 图

plotDim <- function(object, reduction, group.by, color, prefix){
  if(group.by == 'orig.ident'){
    pdf.file <- paste(prefix,reduction,'samples.pdf',sep='_')
    png.file <- paste(prefix,reduction,'samples.png',sep='_')
  }else{
    pdf.file <- paste(prefix,reduction,'cluster.pdf',sep='_')
    png.file <- paste(prefix,reduction,'cluster.png',sep='_')
  }
  plot1<-DimPlot(object, reduction = reduction, group.by=group.by, label = F, pt.size = 1.2, cols=color)
  pdf(pdf.file,w=10,h=8)
  print(plot1)
  dev.off()
  png(png.file,w=1000,h=800)
  print(plot1)
  dev.off()
}

plotVln <- function(object ,features, color, prefix){
  pdf.file <- paste(prefix,features,'vln.pdf',sep='_')
  png.file <- paste(prefix,features,'vln.png',sep='_')
  p <- VlnPlot(object , features = features,pt.size = 0)+ NoLegend()+scale_fill_manual(values = color) 
  pdf(pdf.file,width=12,height=8) ## 注意画图输出
  print(p)
  dev.off()
  p1 <- VlnPlot(object , features = features,pt.size = 0)+ NoLegend()+scale_fill_manual(values = color)
  png(png.file,width=1200,height=800)
  print(p1)
  dev.off()
}
# 两个图都要保存，格式是保存 pdf 和 png
# 这两个图没啥要点，直接用就好了，详细的细节等着我学会了之后再说吧~

#我自己的数据哈，反正我画出来就这样，怎么过滤，过滤条件等我问问大神
#线粒体表达量：绝大多数小于20% <20 |
#血红蛋白表达量：常规值 小于5% <5
#Feature_rna和count_rna的关系，Feature_rna是值每个barcode当中基因检测到的基因数，count_rna统计的是每个barcode中识别到的基因数，每个基因可以检出好几条RNA
#实际理解就是把尖儿上的那些线以上的部分去除掉就行
# nfeature_RNA这个图就是反应样本中的每个细胞表达的基因数量，表达过高就可能是(可能个屁)，新手理解成就是得了，双细胞或者多细胞，表达过低就是空的液滴或者环境的RNA，去掉就完了
#ncount_rna

**

PCA降维

归一化 [基因水平上的]

**
消除不同样本之间的技术差异和批次效应，以及使得不同细胞之间的表达值在数量级上可比较。

#为什么归一化，归一化什么，怎么归一化，不同的方式有啥结果差异，归一化的方式有哪些
#去除测序深度给每个细胞中的表达量带来的影响 ，如果两个细胞的表达量相同，测序越深，表达量受测序深度影响而变化
# 归一化的之后：一个基因的表达量就显而易见了，容易找到高表达量基因
#归一化，归的是技术误差，同一基因在不同样本之间有可比性；降低数据的离散程度，将离散的数据集中。
#常用的归一化方式是:log
#1.对表达量的范围有要求  2.表达量稳定，无极值  3.数据不符合正态
"""这个大概就这么个意思，其他的方式还有啥的，我得找大神问问，肯定还有""" 
object <- NormalizeData(object,normalization.method = "LogNormalize") # 归一化就这么操作，缩放因子

---

获得高变基因

#高变基因[细胞间的基因表达趋于相似，表达量的多少差异不大，这一步为了查找表达差异比较大的基因]
这是个稀疏矩阵，大多数基因都是无用基因，寻找高变基因
频率为0，方差为1 ，基因变化大小权重保持一致

object <- FindVariableFeatures(object, selection.method = "vst", nfeatures = 2000) # 一般默认2000个基因
# selection.method的方式：vst:基于方差稳定化的方式, dispersion:基于离散度的方法，"mean.var.plot":使用平均方差图

object <- ScaleData(object)

#--------------------------------------------------------------------------------

以上内容为数据的基本处理，剩下的就是降维 ，降维的目的是减少特征属性的个数，确保特征之间是相互独立的；简化计算，便于发现规律

#在我的认知里面就是基础分析就到降维，然后就是画图，看见一簇一簇的细胞，具体怎么分，我也不知道，继续看吧~

#PCA  样本少，细胞少，数据简单，有效的可视化，样本的异质性高，就完犊子了
#t-SNE，效果比PCA好，有效的展示细胞间的异质性，当细胞数据多，运行慢
# umap,比之前的两种更好，更好的展示细胞的异质性，更好的展示细胞的亚型

#主成分分析，数据简化，降噪、去冗余
object <- RunPCA(object,npcs=30) # 保存的主成分的数量，npcs值越大，获得基因越多，当然基因越多，噪音越多			

PCA线性降维
筛选标准
主成分累计贡献大于90% ，

去除批次效应，harmony，不同样本的细胞均匀分布在不同的细胞中 ，抑制性比较弱的聚在一起，抑制性强的不做改变

非线性降维 umap 和TSNE
使相似的细胞聚在一起，不相似的分开