细胞周期分析
在单细胞测序分析中，植物（如苜蓿、拟南芥）的细胞周期分析与动物的主要区别在于**标志基因（Marker Genes）**的不同。由于 Seurat 内置的cc.genes是针对人类设计的，植物研究需要通过同源比对或查阅文献来构建个性化的基因集。
以下是根据你提供的代码总结的植物细胞周期分析思路与代码全集。
01. 目的与逻辑
目的
- 鉴定状态：识别每个细胞处于周期的哪个阶段（G1、S、G2M）。
- 消除噪音：细胞周期差异往往是异质性的主要来源。在关注细胞类型本身时，通过“回归（Regress out）”消除周期对聚类的干扰。
- 对比分析：比较不同处理（如 CK vs CD）下，各组织细胞的增殖活力差异。
逻辑步骤
- 同源映射：从拟南芥（经典模式植物）中提取已知的细胞周期标志基因，通过
Blastp映射到目标物种（如苜蓿）上。 - 基因集构建：将映射后的基因分为 S 期和 G2M 期列表。
- 评分赋值：利用
CellCycleScoring计算每个细胞的 S.Score 和 G2M.Score，并判定 Phase。 - 可视化验证：通过 PCA 观察周期基因是否能分开细胞，通过 UMAP 查看不同处理组的分布。
- 模块分析：利用
AddModuleScore细化观察 G1、S、G2、M 四个阶段的特征表达。 - 下游回归：在
ScaleData时移除周期影响，使下游分析聚焦于细胞类型。
02. 全部代码（整理版）
# 1. 环境准备与同源基因映射library(Seurat)library(tidyverse)library(patchwork)# 加载 Blastp 映射表并清理 IDblast_res<-read.table("ref_tar.m8.solar.cor.idAdd.rbh",header=FALSE)blast_res$V1<-gsub("\\.[0-9]+_ref$","",blast_res$V1)# 拟南芥 IDblast_res$V2<-gsub("\\.[^.]+$","",blast_res$V2)# 苜蓿 IDgene_map<-blast_res%>%select(At_ID=V1,Ms_ID=V2)%>%distinct(At_ID,.keep_all=TRUE)# 拟南芥阶段基因集 (拟南芥 -> 苜蓿)ath_genes_list<-list(S=c("AT2G28740","AT3G27060"),G2=c("AT1G76540","AT1G44110"),M=c("AT4G31840","AT3G25980"),G1=c("AT4G02060","AT2G29570","AT5G67260","AT5G54600","AT4G31700"))# 映射转换s_genes<-CaseMatch(search=gene_map$Ms_ID[match(ath_genes_list$S,gene_map$At_ID)],match=rownames(muxu))g2_genes<-CaseMatch(search=gene_map$Ms_ID[match(ath_genes_list$G2,gene_map$At_ID)],match=rownames(muxu))m_genes<-CaseMatch(search=gene_map$Ms_ID[match(ath_genes_list$M,gene_map$At_ID)],match=rownames(muxu))g1_genes<-CaseMatch(search=gene_map$Ms_ID[match(ath_genes_list$G1,gene_map$At_ID)],match=rownames(muxu))g2m_combined<-c(g2_genes,m_genes)# 2. 细胞周期评分与判定muxu<-CellCycleScoring(muxu,s.features=s_genes,g2m.features=g2m_combined)# 3. 可视化：CK 与 CD 对比 (Phase 分布)draw_phase_plot<-function(obj,ident,title_prefix){DimPlot(subset(obj,orig.ident==ident),group.by="Phase",order=TRUE,cols=c("G1"="#32CD32","S"="#00BFFF","G2M"="#FF4500","Undecided"="#80808044"))+labs(title=paste0(title_prefix,"_Phase"))+theme(plot.title=element_text(hjust=0.5))}p1<-draw_phase_plot(muxu,"CK","CK")p2<-draw_phase_plot(muxu,"CD","CD")(p1|p2)+plot_annotation(title="Cell Cycle Phase Comparison")# 4. 阶段模块评分 (AddModuleScore)phase_list<-list(S=s_genes,G2=g2_genes,M=m_genes,G1=g1_genes)for(p_nameinnames(phase_list)){feat_name<-paste0(p_name,"_mod")muxu<-AddModuleScore(muxu,features=list(phase_list[[p_name]]),name=feat_name)# 绘制 FeaturePlot 并保存 (代码略，参考前述 FeaturePlot 拼接逻辑)}# 5. 回归细胞周期影响muxu_regressed<-ScaleData(muxu,vars.to.regress=c("S.Score","G2M.Score"))# 回归后的 PCA 验证muxu_regressed<-RunPCA(muxu_regressed,features=c(s_genes,g2m_combined))DimPlot(muxu_regressed,reduction="pca",group.by="Phase")03. 代码解读
1. 核心映射逻辑 (convert_list&CaseMatch)
- 同源转换：通过
match函数将拟南芥的 AGI 编号转换为苜蓿的基因 ID。这是植物非模式生物分析的关键。 CaseMatch：Seurat 提供的函数，用于确保提取的基因确实存在于当前单细胞对象的矩阵中（处理大小写和拼写）。
2. 评分核心 (CellCycleScoring)
- 该函数计算两个得分：
S.Score和G2M.Score。 - 判定标准：如果两个得分均 < 0，细胞被判定为 G1 期；否则，哪个得分高就判定为哪个期。
- 注意：植物中 G1 基因虽然多，但 Seurat 默认只使用 S 和 G2M 基因来判定这三个阶段。
3. 细分阶段研究 (AddModuleScore)
CellCycleScoring只能粗略给出三个阶段。通过AddModuleScore并输入具体的 G2、M、G1 基因列表，可以更精细地观察。- 例如，在你的循环代码中，
p_name（如 “M_phase”）被计算为一个ModuleScore。这在观察细胞是否卡在有丝分裂期（M期）非常有用。
4. 颜色与可视化控制 (cols&patchwork)
#80808044：这是带 Alpha 通道的透明灰色。用于Undecided或背景细胞，使视觉重心集中在活跃分裂的细胞上。order = TRUE：在FeaturePlot和DimPlot中非常重要，它确保高表达（或关键相位）的细胞被画在顶层，不会被大量 G1 期的细胞遮盖。
5. 线性回归 (vars.to.regress)
- 在
ScaleData步骤中加入S.Score和G2M.Score。 - 物理意义：模型会计算每个基因表达量与周期得分的相关性，并扣除这部分贡献。回归后的 PCA 图中，细胞应按类型聚类，而不再按周期阶段（G1/S/G2M）呈圆环状分布。