跨物种对比的人类特异性基因挖掘
在理解了细胞类型的演化关系后,下一步通常是挖掘那些在人类中表达显著上调、而在其他物种(包括非人灵长类和啮齿类)中表达较低的基因。
1. 核心逻辑:全物种动态比对策略
为了避免筛选出的基因在小鼠中也高表达(即避免“灵长类特异性下调”基因),我们采用全物种严谨比对逻辑:
- 动态比对:自动识别数据中除 Human 以外的所有物种。
- 绝对占优:要求 对所有 个参考物种同时成立。
- 稳健排序:计算人类与各物种差异的最小值,以此作为特异性评分进行排序。
2. 人类绝对特异性基因挖掘函数
library(Seurat)library(ComplexHeatmap)library(circlize)#' 挖掘人类绝对特异性基因 (全物种动态比对版)#' @param obj 跨物种整合的 Seurat 对象#' @param target_cell_type 目标细胞类型名称#' @param celltype_col meta.data 中细胞类型所在的列名#' @param species_col meta.data 中物种信息所在的列名#' @param lfc_threshold LogFC 差异阈值,默认 0.5#' @param min_expression 人类中的最小表达量阈值FindHumanSpecificGenes_Strict<-function(obj,target_cell_type,celltype_col="cell_type",species_col="species",lfc_threshold=0.5,min_expression=0.1){message(paste0(">>> 正在启动全物种比对,挖掘 [",target_cell_type,"] 的人类特异基因..."))# 1. 动态提取特定细胞类型 (使用 [[ ]] 访问软编码参数)cells_to_keep<-colnames(obj)[obj@meta.data[[celltype_col]]==target_cell_type]if(length(cells_to_keep)==0)stop("未找到指定细胞类型。")sub_obj<-subset(obj,cells=cells_to_keep)# 2. 计算伪大块均值 (Pseudo-bulk)avg_data<-AverageExpression(sub_obj,group.by=species_col,assay="RNA",slot="data")$RNAif(!"Human"%in%colnames(avg_data))stop("数据中未检测到 'Human' 分组。")# 3. 动态计算 Human 与所有其他物种的差异other_species<-setdiff(colnames(avg_data),"Human")message(paste(">>> 正在比对 Human 与:",paste(other_species,collapse=", ")))# 对每一个物种,计算 Human - Speciesdiff_list<-lapply(other_species,function(sp){return(avg_data[,"Human"]-avg_data[,sp])})# 4. 严谨筛选条件:Human 必须大于【所有】其他物种# 使用 Reduce 确保基因在所有对比中都满足阈值is_human_spec<-Reduce(`&`,lapply(diff_list,function(d)d>lfc_threshold))expressed_in_human<-avg_data[,"Human"]>min_expression target_genes<-rownames(avg_data)[is_human_spec&expressed_in_human]# 5. 稳健排序:取 Human 与其他物种间“最小差异值”# 这样排在前面的基因在所有对比中都是绝对高表达min_diffs<-apply(as.matrix(do.call(cbind,diff_list)),1,min)target_genes<-target_genes[order(min_diffs[target_genes],decreasing=TRUE)]message(paste(">>> 筛选完成,共找到",length(target_genes),"个绝对特异基因。"))return(target_genes)}3. 配套可视化:演化轨迹热图
该函数支持按照预设的演化轴顺序排列物种,并自动进行 Z-score 标准化。
#' 绘制人类特异性进化热图 (全参数版)PlotEvolutionHeatmap<-function(obj,target_cell_type,features,celltype_col="cell_type",species_col="species"){# 1. 提取数据并计算均值cells_to_keep<-colnames(obj)[obj@meta.data[[celltype_col]]==target_cell_type]sub_obj<-subset(obj,cells=cells_to_keep)avg_mat<-AverageExpression(sub_obj,group.by=species_col,assay="RNA",slot="data")$RNA# 2. 准备绘图矩阵:自动适配进化顺序ideal_order<-c("Mouse","Marmoset","Macaque","Human")plot_species<-intersect(ideal_order,colnames(avg_mat))# 截取有效基因valid_features<-intersect(features,rownames(avg_mat))plot_mat<-avg_mat[valid_features,plot_species,drop=FALSE]# 3. Z-score 标准化 (Row-scale)plot_mat_scaled<-t(scale(t(plot_mat)))# 4. 绘制热图Heatmap(plot_mat_scaled,name="Z-score",cluster_columns=FALSE,# 保持演化轴顺序cluster_rows=FALSE,# 保持特异性分值排序show_row_names=TRUE,col=colorRamp2(c(-1.5,0,1.5),c("#4575b4","white","#d73027")),column_title=paste0("Absolute Human Specificity: ",target_cell_type),row_names_gp=gpar(fontsize=8),column_names_rot=45)}4. 实战示例
# 1. 在指定列名下挖掘 MLI1 细胞的人类特异基因hs_genes<-FindHumanSpecificGenes_Strict(seurat_obj,target_cell_type="MLI1",celltype_col="celltype.V1",species_col="species_std")# 2. 展示前 50 个基因的跨物种演化模式PlotEvolutionHeatmap(seurat_obj,target_cell_type="MLI1",features=head(hs_genes,50),celltype_col="celltype.V1",species_col="species_std")