tcga-bulk-data-preprocessing-with-omicverse

使用 OmicVerse 预处理 TCGA 批量数据

概述

使用此技能从 GDC 下载内容加载 TCGA 数据、构建标准化表达矩阵、附加临床元数据，并通过 ov.bulk.pyTCGA 运行生存分析。

说明

1. 收集所需下载内容

确认用户已从 GDC Data Portal 下载以下三个项目：

• gdc_sample_sheet.<date>.tsv — 样本表导出文件
• 已解压的 gdc_download_xxxxx/ 目录，其中包含表达数据压缩包
• clinical.cart.<date>/ 目录，其中包含临床 XML/JSON 文件

2. 初始化 TCGA 辅助工具

import omicverse as ov
import scanpy as sc
ov.plot_set()

aml_tcga = ov.bulk.pyTCGA(sample_sheet_path, download_dir, clinical_dir)
aml_tcga.adata_init()  # Builds AnnData with raw counts, FPKM, and TPM layers

3. 保存并重新加载

aml_tcga.adata.write_h5ad('data/ov_tcga_raw.h5ad', compression='gzip')

# To reload later:
new_tcga = ov.bulk.pyTCGA(sample_sheet_path, download_dir, clinical_dir)
new_tcga.adata_read('data/ov_tcga_raw.h5ad')

4. 初始化元数据和生存信息

aml_tcga.adata_meta_init()   # Gene ID → symbol mapping, patient info
aml_tcga.survial_init()      # NOTE: "survial" spelling — see Critical API Reference below

5. 运行生存分析

# Single gene
aml_tcga.survival_analysis('MYC', layer='deseq_normalize', plot=True)

# All genes (can take minutes for large gene sets)
aml_tcga.survial_analysis_all()  # NOTE: "survial" spelling

6. 导出结果

aml_tcga.adata.write_h5ad('data/ov_tcga_survival.h5ad', compression='gzip')

关键 API 参考

重要：方法名称拼写不一致

pyTCGA API 中有意保留了拼写不一致。其中两个方法使用了“survial”（漏掉了“v”），另一个方法使用了正确的“survival”：

方法	拼写	用途
survial_init()	survial（没有“v”）	初始化生存元数据列
survival_analysis(gene, layer, plot)	survival（正确）	单基因 Kaplan-Meier 曲线
survial_analysis_all()	survial（没有“v”）	对所有基因进行生存显著性检验

# CORRECT — use the exact method names as documented
aml_tcga.survial_init()                    # "survial" — no 'v'
aml_tcga.survival_analysis('MYC', layer='deseq_normalize', plot=True)  # "survival" — correct
aml_tcga.survial_analysis_all()            # "survial" — no 'v'

# WRONG — these will raise AttributeError
# aml_tcga.survival_init()                 # AttributeError! Use survial_init()
# aml_tcga.survival_analysis_all()         # AttributeError! Use survial_analysis_all()

生存分析方法

survival_analysis() 执行 Kaplan-Meier 分析：

1. 使用中位数作为阈值，将患者分为高表达组和低表达组
2. 计算 log-rank 检验 p 值，以评估显著性
3. 如果 plot=True，则绘制带置信区间的生存曲线

图层选择很重要：建议使用 layer='deseq_normalize'，因为 DESeq2 标准化会校正文库大小和组成偏差，使不同样本间的表达量具有可比性。另一种选择是使用经过 TPM 标准化的值：layer='tpm'。

防御性验证模式

import os

# Before pyTCGA init: verify all paths exist
for name, path in [('sample_sheet', sample_sheet_path),
                    ('downloads', download_dir),
                    ('clinical', clinical_dir)]:
    if not os.path.exists(path):
        raise FileNotFoundError(f"TCGA {name} path not found: {path}")

# After adata_init(): verify expected layers were created
expected_layers = ['counts', 'fpkm', 'tpm']
for layer in expected_layers:
    if layer not in aml_tcga.adata.layers:
        print(f"WARNING: Missing layer '{layer}' — check if TCGA archives are fully extracted")

# Before survival analysis: verify metadata is initialized
if 'survial_init' not in dir(aml_tcga) or aml_tcga.adata.obs.shape[1] < 5:
    print("WARNING: Run adata_meta_init() and survial_init() before survival analysis")

故障排除

• AttributeError: 'pyTCGA' object has no attribute 'survival_init'：使用拼写错误的名称 survial_init()（漏掉“v”）。survial_analysis_all() 也是如此。请参阅上面的关键 API 参考。
• adata_meta_init() 期间出现 KeyError：表达矩阵中的基因 ID 与预期格式不匹配。TCGA 使用 ENSG ID；该方法会在内部将其映射为基因符号。请确保压缩包来自同一次 GDC 下载。
• 生存图为空或 p 值为 NaN：临床 XML 文件缺少日期字段（days_to_death、days_to_last_follow_up）。请检查 clinical.cart.* 目录中是否包含完整的 XML 文件，而不只是元数据 JSON 文件。
• survial_analysis_all() 运行非常缓慢：该方法会逐一检验每个基因。对于包含约 20,000 个基因的基因组，预计需要 5–15 分钟。可以先筛选出感兴趣的基因。
• 样本表列不匹配：确认 TSV 使用制表符分隔，且表头行符合 GDC 的预期格式。如果列名不同，请从 GDC 重新下载。
• 缺少 deseq_normalize 图层：该图层会在 adata_meta_init() 期间创建。如果缺失，请重新运行元数据初始化步骤。

示例

• “读取我的 TCGA OV 下载内容，初始化元数据，并使用 DESeq 标准化计数绘制 MYC 生存曲线。”
• “重新加载已保存的 AnnData 文件，附加生存注释，并导出更新后的 .h5ad。”
• “对所有基因运行生存分析，并保存富集后的数据集。”

参考资料

• 教程笔记本：t_tcga.ipynb
• 可快速复制粘贴的命令：reference.md