在提供的两篇研究文献中，非负矩阵分解（NMF）的应用主要与胰腺癌的分子亚型分析相关，具体计算逻辑结合文献内容及NMF基本原理总结如下：
一、文献中的NMF应用场景
在第二篇文献《Integrative analysis of spatial and single-cell transcriptome data from human pancreatic cancer》中，提到NMF用于分析胰腺癌的bulk转录组数据，以识别分子亚型。例如：
- 文中指出胰腺癌上皮亚簇 Ep_PIFO 的标记基因与先前研究中通过NMF鉴定的 Signature 10 高度相关，该Signature与纤毛功能相关。
- NMF在bulk数据中通过分解基因表达矩阵，将胰腺癌分为不同分子亚型（如经典型、基底样型等），为单细胞亚簇的功能注释提供参考。
二、NMF的核心计算原理（结合文献场景）
1. 数据输入与目标
- 输入数据：胰腺癌bulk转录组的基因表达矩阵 ( \mathbf{V} )（行=基因，列=样本，元素为非负表达值）。
- 分解目标：将 ( \mathbf{V} ) 近似分解为两个非负矩阵的乘积，即：
[
\mathbf{V} \approx \mathbf{W} \times \mathbf{H}
]
其中，( \mathbf{W} ) 为 基矩阵（行=基因，列=亚型，反映各亚型的特征基因表达模式），( \mathbf{H} ) 为 系数矩阵（行=亚型，列=样本，反映每个样本在各亚型上的权重）。
2. 非负性约束的意义
- 基因表达值本质为非负数，NMF的非负约束确保分解结果符合生物学意义（如亚型特征基因的表达量不会为负）。
- 相比PCA等方法，NMF更适合捕捉基因表达的“组合性”特征（如不同亚型的基因共表达模块）。
3. 迭代优化过程
- 初始化：随机生成非负矩阵 ( \mathbf{W} ) 和 ( \mathbf{H} )。
- 迭代更新：通过优化算法（如乘法更新规则、梯度下降）最小化重构误差，常用损失函数为：
[
\text{Loss} = |\mathbf{V} - \mathbf{W}\mathbf{H}|F^2 \quad (\text{ Frobenius范数})
]
每次迭代中，按以下规则更新矩阵元素（保持非负）：
[
W{ik} \leftarrow W_{ik} \times \frac{(\mathbf{V}\mathbf{H}T)_{ik}}{(\mathbf{W}\mathbf{H}\mathbf{H}T){ik}}, \quad H{kj} \leftarrow H_{kj} \times \frac{(\mathbf{W}T\mathbf{V})_{kj}}{(\mathbf{W}T\mathbf{W}\mathbf{H})_{kj}}
] - 终止条件：损失函数收敛或达到预设迭代次数。
4. 亚型识别与解释
- 分解后，( \mathbf{H} ) 中权重最高的亚型即代表样本的主要分子类型（如经典型、基底样型）。
- 通过分析 ( \mathbf{W} ) 中高权重基因的功能，可注释各亚型的生物学特征（如Ep_PIFO关联的纤毛功能相关基因）。
三、文献中NMF与单细胞分析的整合
- 验证单细胞亚簇的生物学意义：通过NMF在bulk数据中定义的亚型特征，可验证单细胞亚簇是否对应已知分子亚型。例如，Ep_PIFO簇的标记基因（如与纤毛相关的基因）在NMF的Signature 10中高度富集，提示其功能特异性。
- 桥接单细胞与bulk数据：NMF分解得到的亚型权重（( \mathbf{H} )）可与单细胞亚簇的比例进行关联分析，揭示肿瘤内异质性与整体亚型的关系。
四、NMF在癌症研究中的优势
- 鲁棒性：对高维、稀疏的基因表达数据具有较好的降维效果，且不受噪声影响。
- 可解释性：分解结果直接对应生物学上的“亚型”或“模块”，便于功能注释。
- 兼容性：可与单细胞分析结合，从bulk水平验证单细胞亚簇的临床意义（如预后相关性）。
总结
文献中NMF的计算核心是通过非负矩阵分解将bulk转录组数据降维为生物学可解释的分子亚型，其迭代优化过程确保分解结果符合基因表达的非负性和功能模块性。该方法在胰腺癌研究中用于定义分子亚型，并与单细胞亚簇的功能注释和预后分析相结合，为解析肿瘤异质性提供了关键工具。
在单细胞转录组和癌症研究中，非负矩阵分解（NMF）常用于分析基因表达数据以识别分子亚型或细胞状态。以下结合文献中的应用场景，提供Python实现NMF的完整流程，重点关注其在胰腺癌和食管癌研究中的应用逻辑。
一、Python实现NMF的核心库与数据准备
1. 关键库导入
import numpy as np
import pandas as pd
import matplotlib.pyplot as plt
from sklearn.decomposition import NMF
from sklearn.preprocessing import normalize
from sklearn.metrics import silhouette_score
import seaborn as sns
import scanpy as sc # 单细胞数据处理
2. 数据准备（以胰腺癌bulk转录组为例）
假设输入数据为基因表达矩阵（行=基因，列=样本），非负且已标准化：
# 加载bulk转录组数据（示例数据，实际需替换为真实数据）
# 数据格式：行为基因，列为样本，值为log2(TPM+1)
expression_matrix = pd.read_csv("pancreatic_cancer_expression.csv", index_col=0)
# 确保数据非负（NMF的基本要求）
if (expression_matrix < 0).any().any():
expression_matrix = expression_matrix - expression_matrix.min().min()
print("数据已调整为非负矩阵")
# 示例：若数据为单细胞矩阵，需先提取非免疫细胞或上皮细胞群体
# adata = sc.read_h5ad("single_cell_data.h5ad")
# non_immune_cells = adata[adata.obs['cell_type'] != 'immune_cell']
# expression_matrix = non_immune_cells.X.toarray()
二、NMF模型构建与参数优化
1. 模型训练与分解
def run_nmf(matrix, n_components, random_state=42):
"""执行NMF分解并返回基矩阵和系数矩阵"""
nmf = NMF(
n_components=n_components,
init='nndsvd', # 推荐初始化方法，提高收敛速度
solver='mu', # 乘法更新规则，适用于大规模数据
random_state=random_state,
max_iter=1000, # 最大迭代次数
tol=1e-4 # 收敛阈值
)
W = nmf.fit_transform(matrix) # 系数矩阵 (样本×亚型)
H = nmf.components_ # 基矩阵 (亚型×基因)
return W, H, nmf
# 示例：分解胰腺癌表达矩阵为k=2种亚型（经典型+基底样型）
W, H, nmf = run_nmf(expression_matrix, n_components=2)
# 计算重构误差（评估分解质量）
reconstruction_error = nmf.reconstruction_err_
print(f"NMF重构误差: {reconstruction_error}")
2. 确定最佳亚型数量（k值）
# 方法1：基于重构误差随k的变化
def find_optimal_k(matrix, max_k=10):
errors = []
for k in range(1, max_k+1):
_, _, nmf = run_nmf(matrix, k)
errors.append(nmf.reconstruction_err_)
plt.figure(figsize=(8, 6))
plt.plot(range(1, max_k+1), errors, 'o-')
plt.xlabel('亚型数量 (k)')
plt.ylabel('重构误差')
plt.title('NMF重构误差随k的变化')
plt.grid(True)
plt.show()
# 肘部法则：选择误差下降趋缓的k值
elbow_k = np.argmin(np.diff(errors) > 0) + 1
return elbow_k
optimal_k = find_optimal_k(expression_matrix)
print(f"最佳亚型数量: {optimal_k}")
# 方法2：结合生物学先验（如文献中已知胰腺癌分为2-3种亚型）
optimal_k = 2 # 假设根据文献预设k=2
三、NMF结果分析与生物学解释
1. 解析基矩阵（H）：识别亚型特征基因
def identify_subtype_markers(H, gene_names, n_top_genes=20):
"""从基矩阵中提取各亚型的 top 标记基因"""
subtype_markers = {}
for i in range(H.shape[0]):
# 按基因权重排序
gene_indices = np.argsort(H[i, :])[::-1]
top_genes = gene_names[gene_indices][:n_top_genes]
subtype_markers[f"Subtype_{i+1}"] = top_genes
return subtype_markers
# 示例：假设gene_names为表达矩阵的行索引
gene_names = expression_matrix.index.values
subtype_markers = identify_subtype_markers(H, gene_names)
# 可视化亚型特征基因热图
plt.figure(figsize=(12, 8))
top_genes = np.unique([gene for markers in subtype_markers.values() for gene in markers])
heatmap_data = expression_matrix.loc[top_genes].T
sns.heatmap(
heatmap_data,
cmap="virdis",
annot=False,
fmt=".1f",
cbar_kws={"label": "标准化表达量"}
)
plt.title("胰腺癌NMF亚型特征基因热图")
plt.show()
2. 解析系数矩阵（W）：样本亚型分类
# 将样本分配到主要亚型
sample_subtypes = np.argmax(W, axis=1)
subtype_proportions = pd.Series(sample_subtypes).value_counts(normalize=True)
print("各亚型样本比例:", subtype_proportions)
# 结合临床数据分析（如文献中的肿瘤分期、生存数据）
# 假设clinical_data为包含样本临床信息的DataFrame
# clinical_data['NMF_Subtype'] = sample_subtypes
# 后续可进行生存分析、亚型与临床特征的相关性分析
四、NMF与单细胞数据的整合分析
1. 单细胞亚簇与NMF亚型的关联
def correlate_sc_subclusters_with_nmf(
sc_expression, # 单细胞基因表达矩阵 (细胞×基因)
nmf_h, # NMF基矩阵 (亚型×基因)
sc_subclusters # 单细胞亚簇标签
):
"""计算单细胞亚簇与NMF亚型的相关性"""
# 计算单细胞亚簇的平均基因表达
cluster_avg_expr = {}
for cluster in np.unique(sc_subclusters):
cluster_cells = sc_expression[sc_subclusters == cluster]
cluster_avg_expr[cluster] = np.mean(cluster_cells, axis=0)
# 转换为DataFrame便于计算
cluster_avg_df = pd.DataFrame(cluster_avg_expr).T
cluster_avg_df.columns = sc_expression.columns
# 计算亚簇与NMF亚型的相关性
correlation_matrix = pd.DataFrame()
for i in range(nmf_h.shape[0]):
subtype_genes = np.argsort(nmf_h[i, :])[::-1][:50] # 取亚型前50个特征基因
subtype_expr = cluster_avg_df[subtype_genes].mean(axis=1)
correlation_matrix[f"Subtype_{i+1}"] = subtype_expr
# 计算相关性系数
correlation_coeffs = pd.DataFrame(
np.corrcoef(correlation_matrix.T),
index=correlation_matrix.index,
columns=correlation_matrix.index
)
return correlation_coeffs
# 示例：假设sc_data为单细胞表达矩阵，sc_clusters为细胞亚簇标签
# correlation = correlate_sc_subclusters_with_nmf(sc_data, H, sc_clusters)
# sns.heatmap(correlation, annot=True, cmap="coolwarm")
# plt.title("单细胞亚簇与NMF亚型的相关性")
2. 验证单细胞亚簇与NMF亚型的一致性
以文献中胰腺癌的Ep_PIFO亚簇为例，验证其与NMF的Signature 10相关性：
# 假设ep_pifo_genes为Ep_PIFO亚簇的标记基因
ep_pifo_genes = ["PIFO", "KRT6B", "TRIM54", ...] # 替换为实际标记基因
# 计算Ep_PIFO基因在NMF亚型中的富集程度
pifo_gene_indices = [np.where(gene_names == gene)[0][0] for gene in ep_pifo_genes if gene in gene_names]
pifo_enrichment = np.mean(H[:, pifo_gene_indices], axis=1)
# 可视化Ep_PIFO与各亚型的关联
plt.bar(range(1, len(pifo_enrichment)+1), pifo_enrichment)
plt.xlabel("NMF亚型")
plt.ylabel("Ep_PIFO基因平均权重")
plt.title("Ep_PIFO亚簇标记基因在NMF亚型中的富集")
plt.show()
五、进阶应用：结合空间转录组数据
1. NMF亚型的空间分布分析
def analyze_nmf_spatial_distribution(
spatial_expression, # 空间转录组表达矩阵 (spot×基因)
nmf_h, # NMF基矩阵
spot_coordinates # spot的空间坐标 (x, y)
):
"""分析NMF亚型在空间上的分布"""
# 计算每个spot的亚型权重
spot_subtype_weights = NMF(n_components=nmf_h.shape[0]).fit_transform(spatial_expression)
# 转换为DataFrame并添加坐标
spot_subtype_df = pd.DataFrame(spot_subtype_weights, columns=[f"Subtype_{i+1}" for i in range(nmf_h.shape[0])])
spot_subtype_df["x"] = spot_coordinates[:, 0]
spot_subtype_df["y"] = spot_coordinates[:, 1]
# 可视化亚型空间分布（以Subtype_1为例）
plt.figure(figsize=(10, 8))
plt.scatter(
spot_subtype_df["x"],
spot_subtype_df["y"],
c=spot_subtype_df["Subtype_1"],
cmap="viridis",
s=50,
alpha=0.8
)
plt.colorbar(label="Subtype_1权重")
plt.title("NMF亚型在空间转录组中的分布")
plt.show()
六、代码应用场景与文献关联
-
胰腺癌分子亚型分析
在第二篇文献中，NMF用于将胰腺癌分为经典型和基底样型，对应单细胞亚簇Ep_TRIM54（经典型）和Ep_KRT6A（基底样型）。通过上述代码可复现：- 分解bulk转录组数据，识别亚型特征基因（如经典型的GATA6，基底样型的KRT6B）。
- 验证单细胞亚簇的标记基因是否在NMF亚型中富集（如Ep_PIFO与纤毛相关Signature 10的关联）。
-
食管癌免疫微环境分析
第一篇文献中虽未明确提及NMF，但可类比应用：- 对食管癌bulk转录组进行NMF分解，识别与免疫抑制相关的亚型（如富含Treg、TEx细胞的亚型）。
- 结合单细胞数据，分析NMF亚型与免疫细胞亚群（如tDC、Treg）的空间共定位关系。
七、注意事项
-
数据预处理：
- 确保输入数据非负，若为原始计数数据，建议使用log(CPM+1)或标准化处理。
- 对于高维数据，可先通过PCA降维（保留80%方差）再进行NMF，提高计算效率。
-
生物学验证：
- NMF亚型需结合标记基因的功能注释（如GO/KEGG富集分析），避免纯数学分解导致的生物学误解。
- 参考文献中的已知亚型（如胰腺癌的经典-基底样二分法），调整k值以匹配生物学先验。
-
单细胞-空间整合：
- 若需复现文献中Ep_VGLL1亚簇的过渡特性，可通过NMF亚型权重与单细胞亚簇比例的相关性分析，验证其“中间状态”。
通过上述代码框架，可实现从NMF模型构建到生物学解释的全流程分析，紧密结合两篇文献中胰腺癌和食管癌的研究场景，为肿瘤微环境的分子亚型解析提供量化工具。
扫一扫
被折叠的 条评论
为什么被折叠?
到【灌水乐园】发言
