最近在实习，第一个工作是这个，现在这里记录一些随笔，之后会系统总结
不同方法的其他流程部分怎么分析可以看20年综述怎么处理这部分。
答：
对于 Seurat 2、Harmony、MNN Correct、fastMNN 和 limma，使用 Seurat 2 包执行标准化、缩放和高可变基因 (HVG) 选择的数据预处理步骤。 对于 Seurat 3 批量校正，使用了包中的相应功能。 对于其他方法，采用了它们各自推荐的预处理管道。
Q：Seurat 选择哪个，现在Seurat4已经出来了，这些方法还用Seurat2的进行预处理吗？Seurat2、3和4有什么区别？
一、单细胞测序整体流程与去批次效应在流程中哪部分
1. 下机后进行预处理得到表达矩阵
2. 降维，聚类，可视化
3. 聚类和差异分析
- step1: 创建对象
- step2: 质量控制
- step3: 表达量的标准化和归一化
- step4: 去除干扰因素(多个样本整合)
- step5: 判断重要的基因
- step6: 多种降维算法
- step7: 可视化降维结果
- step8: 多种聚类算法
- step9: 聚类后找每个细胞亚群的标志基因
- step10: 继续分类
https://github.com/theislab/single-cell-tutorial
python 从第二部分开始
1. load data
# Set up data loading
#Data files
sample_strings = ['Duo_M1', 'Duo_M2', 'Jej_M1', 'Jej_M2', 'Il_M1', 'Il_M2']
sample_id_strings = ['3', '4', '5', '6', '7', '8']
file_base = '../data/Haber-et-al_mouse-intestinal-epithelium/GSE92332_RAW/GSM283657'
exp_string = '_Regional_'
data_file_end = '_matrix.mtx'
barcode_file_end = '_barcodes.tsv'
gene_file_end = '_genes.tsv'
cc_genes_file = '../Macosko_cell_cycle_genes.txt'
# First data set load & annotation
#Parse Filenames
sample = sample_strings.pop(0)
sample_id = sample_id_strings.pop(0)
data_file = file_base+sample_id+exp_string+sample+data_file_end
barcode_file = file_base+sample_id+exp_string+sample+barcode_file_end
gene_file = file_base+sample_id+exp_string+sample+gene_file_end
#Load data
adata = sc.read(data_file, cache=True)
adata = adata.transpose()
adata.X = adata.X.toarray()
barcodes = pd.read_csv(barcode_file, header=None, sep='\t')
genes = pd.read_csv(gene_file, header=None, sep='\t')
#Annotate data
barcodes.rename(columns={0:'barcode'}, inplace=True)
barcodes.set_index('barcode', inplace=True)
adata.obs = barcodes
adata.obs['sample'] = [sample]*adata.n_obs
adata.obs['region'] = [sample.split("_")[0]]*adata.n_obs
adata.obs['donor'] = [sample.split("_")[1]]*adata.n_obs
genes.rename(columns={0:'gene_id', 1:'gene_symbol'}, inplace=True)
genes.set_index('gene_symbol', inplace=True)
adata.var = genes
# Loop to load rest of data sets
for i in range(len(sample_strings)):
#Parse Filenames
sample = sample_strings[i]
sample_id = sample_id_strings[i]
data_file = file_base+sample_id+exp_string+sample+data_file_end
barcode_file = file_base+sample_id+exp_string+sample+barcode_file_end
gene_file = file_base+sample_id+exp_string+sample+gene_file_end
#Load data
adata_tmp = sc.read(data_file, cache=True)
adata_tmp = adata_tmp.transpose()
adata_tmp.X = adata_tmp.X.toarray()
barcodes_tmp = pd.read_csv(barcode_file, header=None, sep='\t')
genes_tmp = pd.read_csv(gene_file, header=None, sep='\t')
#Annotate data
barcodes_tmp.rename(columns={0:'barcode'}, inplace=True)
barcodes_tmp.set_index('barcode', inplace=True)
adata_tmp.obs = barcodes_tmp
adata_tmp.obs['sample'] = [sample]*adata_tmp.n_obs
adata_tmp.obs['region'] = [sample.split("_")[0]]*adata_tmp.n_obs
adata_tmp.obs['donor'] = [sample.split("_")[1]]*adata_tmp.n_obs
genes_tmp.rename(columns={0:'gene_id', 1:'gene_symbol'}, inplace=True)
genes_tmp.set_index('gene_symbol', inplace=True)
adata_tmp.var = genes_tmp
adata_tmp.var_names_make_unique()
# Concatenate to main adata object
adata = adata.concatenate(adata_tmp, batch_key='sample_id')
if 'gene_id-1' in adata.var.columns:
adata.var['gene_id'] = adata.var['gene_id-1']
adata.var.drop(columns = ['gene_id-1', 'gene_id-0'], inplace=True)
adata.obs.drop(columns=['sample_id'], inplace=True)
adata.obs_names = [c.split("-")[0] for c in adata.obs_names]
adata.obs_names_make_unique(join='_')
#Assign variable names and gene id columns
adata.var_names = [g.split("_")[1] for g in adata.var_names]
adata.var['gene_id'] = [g.split("_")[1] for g in adata.var['gene_id']]
# Annotate the data sets
print(adata.obs['region'].value_counts())
print('')
print(adata.obs['donor'].value_counts())
print('')
print(adata.obs['sample'].value_counts())
2. 预处理
2.1 QC
数据质量控制可分为细胞质量控制和基因质量控制。评估细胞质量的典型质量指标包括分子计数（UMI）、表达基因数和线粒体计数分数。由于细胞中核mRNA的比例很低，所以高比例的线粒体读数可以指示细胞压力。应当注意的是，高线粒体RNA部分也可以是指示呼吸升高的生物信号。
通过可视化选择合适的筛选阈值，一般建议开始选择我不非常严格的阈值，观察分析结果后视情况更新阈值。代码见jupyter。
2.2 标准化
教程中用scran包，几年前最好的方法。
到目前为止，数据仅可用作计数矩阵。计数代表在scRNA-seq实验中捕获的分子。由于并非细胞中的所有mRNA分子都被捕获，因此细胞间检测到的总计数存在变异性，这是由细胞中的分子数量和采样造成的。由于我们不能假设所有细胞都含有相同数量的分子（细胞大小可能有很大差异），因此我们必须估计细胞中最初存在的分子数量。事实上，我们并没有估计分子的确切数量，而是估计与分子的真实数量成比例的细胞特异性因子。这些被称为尺寸因子。标准化表达式值通过将测量的计数除以单元格的大小因子来计算。
基于差异测试（BeatheVieth的个人通信）和批量校正的标准化方法的比较[Buettner等人，2019]，scran包中实现的标准化法表现最佳。该方法需要粗聚类输入，以提高大小因子估计性能。因此，我们使用简单的预处理方法，以低分辨率对数据进行聚类，以获得尺寸因子估计的输入。基本的预处理包括假设所有大小因子都相等（库大小归一化为每百万CPM计数），并对计数数据进行对数变换。
2.3 批次更正
执行批次校正以调整加载的6个样品的批次效应。由于来自样品和上皮区域的批次效应是重叠的，因此对该批次效应的校正也将部分地回归出区域之间的差异。我们考虑到这一点，以便对数据进行最佳聚类。这种方法也有助于找到分化轨迹，但我们返回到非批量校正数据，用于差异测试和计算标记基因。
本文使用ComBat，一种古老的处理方式，是用于微阵列，效果不好
# ComBat batch correction
sc.pp.combat(adata, key='sample')
2.4 高度可变基因提取
我们提取高度可变基因（HVG）以进一步降低数据集的维数，并仅包括信息量最大的基因。在整个数据集中差异很大的基因是数据中潜在生物变异的信息。由于我们只想捕捉这些基因中的生物变异，我们在标准化和批量校正后选择了高度可变的基因。HVG用于聚类、轨迹推断和降维/可视化，而完整数据集用于计算标记基因、差异测试、细胞周期评分和可视化数据上的表达值。
在这里，我们使用标准技术从10X基因组预处理软件CellRanger中提取高度可变的基因。通常选择1000至5000个基因。在这里，我们提取了前4000个最可变的基因用于进一步处理。如果已知特定的重要基因，人们可以评估需要多少高度可变的基因才能包括所有或大部分这些基因。
2.5 可视化
略
2.6 细胞周期打分
已经调查并纠正了数据中已知的技术变化来源（例如批次、计数深度）。可以解释这些数据的已知生物变异来源是细胞周期。这里，Macosko等人，Cell 161（2015）的基因列表用于对数据中的细胞周期效应进行评分，并按细胞周期阶段对细胞进行分类。该文件可以在单单元教程github存储库中找到，也可以从论文的补充材料中获取。
请注意，由于基因列表是为人类HeLa细胞生成的，所以基因名称用小写字母开头，以映射到相应的小鼠基因。在将此脚本修改为您自己的数据时，必须考虑到来自除老鼠以外的物种的数据。
我们对adata中的整批校正数据集进行细胞周期评分。
3 下游分析
3.1 Clustering.
3.2 marker 基因和 cluster 注释
3.3 subclustering
3.4 成分分析
略
标签：数据分析,tmp,单细胞,测序,sample,file,gene,adata,id From： https://www.cnblogs.com/agis/p/16833545.html