求助：如何用R的ChemmineR处理Excel格式CHEMBL分子数据集做聚类

前提准备 #

确保你的Excel数据集包含SMILES字符串列（CHEMBL数据自带该字段，用于描述分子结构），同时保留蛋白质ID、IC50等你需要的关联字段。

步骤1：读取Excel数据 #

安装并加载专门处理Excel的readxl包：

install.packages("readxl")
library(readxl)
# 替换为你的Excel文件路径
df <- read_excel("your_molecules.xlsx")
# 检查SMILES列是否正常读取
head(df$smiles)

步骤2：将SMILES转成ChemmineR兼容的分子对象 #

用rcdk解析SMILES，再转换为ChemmineR核心格式SDFset：

library(rcdk)
library(ChemmineR)

# 提取SMILES列并转为字符向量
smiles_vec <- as.character(df$smiles)
# 解析SMILES，过滤无效分子结构
mols <- parse.smiles(smiles_vec)
valid_idx <- sapply(mols, function(x) !is.null(x)) # 标记有效分子
mols <- mols[valid_idx]
df <- df[valid_idx, ] # 同步过滤原数据框

# 转换为SDFset，并给分子设置唯一ID（用蛋白质ID或行号）
sdfset <- as(mols, "SDFset")
cid(sdfset) <- df$protein_id # 假设Excel中蛋白质ID列名为"protein_id"

步骤3：计算分子指纹（聚类核心依据） #

采用常用的Morgan指纹进行结构编码：

# 计算1024位的Morgan指纹
fps <- calcFingerprints(sdfset, type = "Morgan", fpType = "bit", size = 1024)
# 转换为矩阵格式，方便后续聚类计算
fp_matrix <- fp2bit(fps)

步骤4：执行聚类分析 #

用层次聚类实现直观的结构分组，也可替换为k-means等方法：

# 计算指纹间的Jaccard距离（适配二进制指纹）
dist_mat <- dist(fp_matrix, method = "binary")
# 用Ward法进行层次聚类
hc <- hclust(dist_mat, method = "ward.D2")

# 绘制树状图查看聚类结果
plot(hc, labels = cid(sdfset), cex = 0.7, main = "Molecular Clustering Tree")

# 切割聚类：比如分成5个组，可根据需求调整k值
cluster_groups <- cutree(hc, k = 5)
# 将聚类结果合并回原数据框
df$cluster <- cluster_groups

步骤5：提取每个聚类的常见结构 #

用rcdk识别聚类内的共同子结构：

# 按聚类分组拆分SDFset
cluster_sdfs <- split(sdfset, cluster_groups)

# 遍历每个聚类，提取共同子结构（至少2个分子才执行匹配）
common_substructures <- lapply(cluster_sdfs, function(subset_sdf) {
  if(length(subset_sdf) >= 2) {
    findCommonSubstructure(subset_sdf[[1]], subset_sdf[[2]], match.params = list(type = "smarts"))
  } else {
    subset_sdf[[1]] # 单个分子直接返回自身结构
  }
})

# 将共同子结构转为SMILES格式输出，方便查看
lapply(common_substructures, function(x) get.smiles(x))

关键注意事项 #

• 提前清理Excel中的缺失SMILES值，否则会导致解析失败
• 聚类参数可灵活调整：比如换用MACCS指纹、修改指纹长度、切换Tanimoto距离、更改聚类合并方法（如complete）
• 可结合IC50数据做后续分析：统计各聚类的IC50均值、分布，关联结构特征与活性表现

内容的提问来源于stack exchange，提问作者Haifa Hassanie