泛癌分析方法与流程¶

# === 泛癌分析：TCGA 数据获取与突变全景图 ===

# 1. 使用 TCGAbiolinks 下载数据（通过 R）
# 或使用 GDC API / cBioPortal

import pandas as pd                    # 数据处理
import matplotlib.pyplot as plt        # 绑图
import seaborn as sns                  # 统计绑图

# 2. 读取泛癌突变数据（MAF格式合并）
# 假设已下载并合并了多癌种MAF
maf = pd.read_csv("pancancer_mutations.maf",
                   sep="\t",           # 制表符分隔
                   usecols=["Hugo_Symbol", "Variant_Classification",
                            "Tumor_Sample_Barcode", "cancer_type"],
                   low_memory=False)

# 3. 计算每个基因在每种癌症中的突变频率
def calc_mutation_freq(maf_df):
    """计算基因-癌种突变频率矩阵"""
    # 每个癌种的总样本数
    samples_per_cancer = maf_df.groupby("cancer_type")[
        "Tumor_Sample_Barcode"].nunique()             # 各癌种样本数

    # 每个基因在每个癌种中突变的样本数
    gene_cancer = maf_df.groupby(
        ["Hugo_Symbol", "cancer_type"])[
        "Tumor_Sample_Barcode"].nunique().reset_index()  # 突变样本数
    gene_cancer.columns = ["gene", "cancer", "mut_samples"]

    # 计算频率
    gene_cancer["freq"] = gene_cancer.apply(
        lambda r: r["mut_samples"] / samples_per_cancer[r["cancer"]],
        axis=1)                        # 突变频率 = 突变样本数 / 总样本数
    return gene_cancer

freq_df = calc_mutation_freq(maf)

# 4. 找出跨癌种高频突变基因（在>=5个癌种中频率>5%）
gene_cancer_count = freq_df[freq_df["freq"] > 0.05].groupby("gene")[
    "cancer"].nunique()                # 每个基因突变频率>5%的癌种数
pan_cancer_genes = gene_cancer_count[
    gene_cancer_count >= 5             # 至少5个癌种
].sort_values(ascending=False)
print("=== 泛癌高频突变基因 ===")
print(pan_cancer_genes.head(20))       # 打印前20个

# 5. 绘制 oncoplot 风格的热图
top_genes = pan_cancer_genes.head(15).index.tolist()  # 前15基因
pivot = freq_df[freq_df["gene"].isin(top_genes)].pivot(
    index="gene", columns="cancer", values="freq").fillna(0)

plt.figure(figsize=(16, 8))            # 设置大小
sns.heatmap(pivot,                     # 热图
            cmap="YlOrRd",             # 黄-红配色
            annot=True,                # 显示数值
            fmt=".0%",                 # 百分比格式
            linewidths=0.5)            # 网格线
plt.title("泛癌基因突变频率热图")       # 标题
plt.tight_layout()                     # 紧凑布局
plt.savefig("pancancer_oncoplot.png", dpi=150)  # 保存

# === 使用 TCGAbiolinks 下载 TCGA 数据 ===
library(TCGAbiolinks)                  # TCGA数据下载包

# 下载某癌种的突变数据
query <- GDCquery(                     # 构建查询
  project = "TCGA-BRCA",              # 乳腺癌项目
  data.category = "Simple Nucleotide Variation",  # 突变数据
  data.type = "Masked Somatic Mutation",          # 体细胞突变
  workflow.type = "Aliquot Ensemble Somatic Variant Merging and Masking"
)
GDCdownload(query)                     # 下载数据
maf <- GDCprepare(query)              # 整理为MAF格式

# 使用 maftools 分析
library(maftools)                      # MAF分析包
laml <- read.maf(maf = maf)           # 读取MAF
oncoplot(maf = laml, top = 20)        # 绘制oncoplot（前20基因）