公共数据一般会以h5ad形式进行保存，格式跟R会有些出入，本文旨在对数据进行处理，转换常用的格式，包括基因转换，矩阵处理等操作

以allen 小鼠的h5ad数据为例

import matplotlib.pyplot as plt import pandas as pd import scanpy as sc import numpy as np import anndata as ad import seaborn as sns import harmonypy as hm import leidenalg import os sc.settings.set_figure_params(dpi=100, facecolor="white") STR_mous1 = sc.read_h5ad("/data/input/Files/scRNA/mouse_allen/WMB-10Xv3/WMB-10Xv3-PAL-raw.h5ad") STR_mous2 = sc.read_h5ad("/data/input/Files/scRNA/mouse_allen/WMB-10Xv3/WMB-10Xv3-STR-raw.h5ad")

目前矩阵所使用的数据是用Ensembl 作为基因名称，一般我们gene- symbol作为基因名称，所以要对其进行转换，转换存在名称重复的问题，所以写了一段函数进行重复值处理及基因转换

import numpy as np import pandas as pd import anndata as ad from scipy import sparse def fast_collapse(adata, new_index_column, mode='sum'): """ 极速版：使用稀疏矩阵乘法进行基因合并。 比 Pandas GroupBy 快 10 to 100 倍。 """ print(f"🚀 正在使用矩阵乘法加速合并 (模式: {mode})...") # 1. 准备新名字 # 填充 NaN，防止报错 new_gene_names = adata.var[new_index_column].astype(str).fillna("Unknown") # 将新名字转换为 Categorical 类型 # codes 将作为新矩阵的列索引 (0, 1, 2, ...) # categories 将作为新的 var_names new_categories = pd.Categorical(new_gene_names) new_var_names = new_categories.categories row_indices = np.arange(adata.n_vars) col_indices = new_categories.codes # 2. 构建“映射矩阵” (Aggregator Matrix) # 形状: (n_old_vars, n_new_unique_vars) # data: 全是 1 # 这是一个极其稀疏的矩阵 aggregator = sparse.coo_matrix( (np.ones(adata.n_vars), (row_indices, col_indices)), shape=(adata.n_vars, len(new_var_names)) ).tocsr() # 转为 CSR 方便乘法 # 3. 执行矩阵乘法 (核心加速步骤) # (n_cells, n_old) @ (n_old, n_new) -> (n_cells, n_new) X_new = adata.X @ aggregator # 4. 如果是 mean 模式，需要除以每个新基因合并了多少个旧基因 if mode == 'mean': # 计算每个新基因包含了多少个旧基因 # 相当于对 aggregator 列求和 counts = np.array(aggregator.sum(axis=0)).flatten() # 避免除以 0 (理论上不应该发生) counts[counts == 0] = 1 # 利用广播机制除以数量 # 稀疏矩阵乘法只能做乘法，所以我们乘 1/counts # 创建对角矩阵进行缩放 diag_weights = sparse.diags(1 / counts) X_new = X_new @ diag_weights # 5. 重建 var 表 # 这种方法没法像之前的代码那样分 "unique" 和 "dup" 处理 # 所以我们统一用 groupby 取第一条作为元数据 print(" 正在重建元数据...") old_var = adata.var.copy() old_var['old_index'] = old_var.index # Groupby 取第一条 (这是这一步里唯一慢的地方，但也比全量数据快得多) new_var = old_var.groupby(new_index_column).first() # 确保顺序与 categories 一致 new_var = new_var.loc[new_var_names] # 6. 构建新对象 new_adata = ad.AnnData( X=X_new, obs=adata.obs.copy(), var=new_var ) print(f"✅ 完成！从 {adata.n_vars} 个基因合并为 {new_adata.n_vars} 个。") return new_adata # === 使用方法 === # adata_allen_final = fast_collapse(adata_allen_v2, 'feature_name', mode='sum')

转换成功之后对两个数据进行合并

######第一种方法########## import scanpy as sc import pandas as pd # 假设你已经合并了 combined_adata = sc.concat( [hpf_mous1_final, hpf_mous2_final], join="outer", # 取并集 fill_value=0 ) ######第二种方法########## obs1 = data_public_clas.obs.copy() obs2 = adata_claustrum.obs.copy() combined_obs = pd.concat([obs1, obs2], axis=0) combined_adata = sc.concat( [data_public_clas, adata_claustrum], join="inner", # 基因取交集 # label="dataset", # 给不同来源打标签（可选） fill_value=0 # 这个只影响表达矩阵 ) combined_obs = combined_obs.reindex(combined_adata.obs_names) combined_adata.obs = combined_obs combined_adata.X

最后保存即可

combined_adata.write("/data/work/human_data/05.Cross_species_comparson/mouse.publci.STR.h5ad",compression= "lzf")

然后进行h5ad的转换

library(sceasy) library(reticulate) library(Seurat) sceasy::convertFormat(obj = "/data/work/human_data/05.Cross_species_comparson/marmoset.publci.STR.h5ad", from="anndata",to="seurat", outFile = '/data/work/human_data/05.Cross_species_comparson/marmoset.publci.STR.rds')

转化之后存在一个问题，小鼠基因名称是小写，所以还需要对小鼠基因名称进行大写转换

# 1. 提取原始计数矩阵 (Sparse Matrix) # 我们只提取 counts，因为 data 和 scale.data 可以重新算，这样最干净 raw_counts <- GetAssayData(mouse.STR, assay = "RNA", slot = "counts") # 2. 修改矩阵的行名 (变成全大写) # 这一步是对矩阵操作，是合法的 rownames(raw_counts) <- toupper(rownames(raw_counts)) # 3. 使用改名后的矩阵，重新创建一个 Seurat 对象 # 注意：我们要把原来的 meta.data (分组信息等) 完美继承过来 mouse.STR.new <- CreateSeuratObject( counts = raw_counts, meta.data = mouse.STR@meta.data, project = "mouse.STR_Upper" # 可选：给个新名字 )

tips：

公共数据有些注释只存在少许细胞，可对其进行初步过滤

FilterRareCellTypes <- function(seurat_obj, group.by, threshold = 100) { # 1. 安全检查：确认列名存在 if (!group.by %in% colnames(seurat_obj@meta.data)) { stop(glue::glue("❌ 错误：在 meta.data 中找不到列名 '{group.by}'")) } # 2. 提取分类标签向量 # 使用 [[ ]] 提取比 $ 更灵活，支持变量名 labels <- seurat_obj@meta.data[[group.by]] # 3. 统计每个类型的细胞数量 # table() 会自动忽略 NA，这是一个好特性 counts <- table(labels) # 4. 筛选出保留和剔除的类型 types_to_keep <- names(counts)[counts >= threshold] types_to_remove <- names(counts)[counts < threshold] # 5. 打印顾问报告 (让用户知道发生了什么) if (length(types_to_remove) > 0) { message(paste0("⚠️ 正在剔除 {length(types_to_remove)} 个稀有类型 (阈值 < {threshold}):")) print(paste0(types_to_remove, " (n=", counts[types_to_remove], ")")) } else { message("✅ 没有类型被剔除，所有类型细胞数均达标。") return(seurat_obj) } # 6. 【核心步骤】物理锁定细胞 ID # 只有标签属于“保留名单”的细胞才会被选中 cells_to_keep <- rownames(seurat_obj@meta.data)[labels %in% types_to_keep] # 7. 执行提取 (使用最稳的 cells 参数) new_obj <- subset(seurat_obj, cells = cells_to_keep) # 8. 返回结果 return(new_obj) }

BGI.STR.filtered <- FilterRareCellTypes( seurat_obj = BGI.STR, group.by = "celltype", # 变量形式传入列名 threshold = 100 # 变量形式传入阈值 )

又学习到一个好用的代码，🐮🐮🐮

import scanpy as sc import glob import os import scanpy as sc import scipy.io import pandas as pd #读取文件所在路径 glob.glob("/data/input/Files/luoting1/human_data/download/*h5ad") import pandas as pd import scipy.io #主要函数export_files def export_files(adata, key): scipy.io.mmwrite(f"/data/work/organized_samples/public_marmoset/fengguoping/amygdala_{key}_matrix.mtx", adata.X.T) features_df = pd.DataFrame({ 'gene_id': adata.var_names, 'gene_name': adata.var['gene_name'] }) features_df.to_csv( f"/data/work/organized_samples/public_marmoset/fengguoping/amygdala_{key}_features.tsv", sep="\t", index=False, header=False ) adata.obs[['donor_id', 'structure', 'region', 'cell_type']].to_csv( f"/data/work/organized_samples/public_marmoset/fengguoping/amygdala_{key}_metadata.csv" ) micro = sc.read_h5ad("/data/input/Files/luoting1/human_data/download/MIC.h5ad") micro_sub = micro[micro.obs.region.isin(['amyg'])] export_files(micro_sub, "micro")

太简洁了

然后是转换成rds

library(Seurat) library(Matrix) library(data.table) library(dplyr) library(patchwork)

# microglia counts <- Matrix::readMM("./data/download/micro_matrix.mtx") genes <- read.table("./data/download/micro_features.tsv", stringsAsFactors = FALSE) rownames(counts) <- genes$V1 metadata <- read.csv("./data/download/micro_metadata.csv", row.names = 1) colnames(counts) <- rownames(metadata) micro <- CreateSeuratObject(counts = counts, meta.data = metadata) micro