152_微生物组源追踪¶

一句话概述¶

微生物组源追踪(Source Tracking)利用贝叶斯模型或矩阵分解方法，从已知来源环境(sources)的微生物群落组成出发，定量估计未知样本(sink)中各来源的贡献比例，广泛应用于污染溯源、传播追踪和生态学研究。

核心知识点总览¶

知识点	说明
源追踪概念	估计sink样本中来自各source的微生物贡献比例
SourceTracker2	基于Gibbs采样的贝叶斯源追踪方法
FEAST	基于EM算法的快速源追踪(比SourceTracker快100倍)
Source-Sink模型	Source为已知来源，Sink为待分解的目标样本
Unknown source	未被已知source解释的微生物比例
OTU/ASV表	输入数据：微生物分类单元的丰度矩阵
应用场景	水源污染、医院感染溯源、人体微生物来源
方法对比	速度、精度、对稀有taxa的敏感度差异

各步骤详解¶

第一步：源追踪的基本概念¶

白话解释： 想象你在一条河里发现了一些"可疑"的微生物。上游有一个养殖场、一个城市污水口、一片自然湿地。源追踪就是通过比较河水样本和三个上游来源的微生物"指纹"，估计河水中有多少比例来自养殖场、多少来自污水、多少来自自然环境。

技术原理： 假设sink样本的微生物组成是各source的混合： $$P(sink) = \sum_{k=1}^{K} \pi_k \cdot P(source_k) + \pi_{unknown} \cdot P(unknown)$$

π_k：source k的贡献比例 (Σπ = 1)
P(source_k)：source k的微生物群落分布
P(unknown)：未知来源的均匀分布

代码概念：

import numpy as np

# 概念演示
# 3个source + 1个unknown
source_profiles = np.array([
    [0.4, 0.3, 0.2, 0.1, 0.0],  # 土壤
    [0.1, 0.1, 0.3, 0.3, 0.2],  # 污水
    [0.2, 0.4, 0.1, 0.2, 0.1],  # 肠道
])
unknown_profile = np.array([0.2, 0.2, 0.2, 0.2, 0.2])  # 均匀

# 真实混合比例
true_mix = np.array([0.3, 0.4, 0.2, 0.1])  # 土壤30%, 污水40%, 肠道20%, unknown10%

# 生成sink样本
all_profiles = np.vstack([source_profiles, [unknown_profile]])
sink_expected = true_mix @ all_profiles
print("Sink profile:", sink_expected)

第二步：SourceTracker2分析¶

白话解释： SourceTracker2使用贝叶斯方法，通过Gibbs采样来估计每条序列最可能来自哪个source。就像一个侦探逐条检查"证据"(每条序列)，根据各source的"档案"(微生物组成)判断它最可能来自哪里。

技术原理： - 使用Dirichlet-Multinomial模型 - Gibbs采样逐个分配每条序列到source - 输出：每个sink样本中各source的比例 + 未知比例 - 支持LOO(leave-one-out)交叉验证评估准确性

代码示例：

# === 安装 ===
pip install sourcetracker

# === 输入文件准备 ===
# 1. OTU/ASV表 (BIOM格式或TSV)
# 行: 样本, 列: OTU/ASV
# 需要包含source和sink样本

# 2. Metadata文件 (mapping file)
# SampleID  SourceSink  Env
# Source1   source      Soil
# Source2   source      Soil
# Source3   source      Sewage
# Sink1     sink        River
# Sink2     sink        River

# === Python API ===
import sourcetracker
import pandas as pd
import numpy as np
from biom import load_table

# 加载数据
otu_table = pd.read_csv("otu_table.tsv", sep="\t", index_col=0)
metadata = pd.read_csv("metadata.tsv", sep="\t", index_col=0)

# 分离source和sink
source_mask = metadata["SourceSink"] == "source"
sink_mask = metadata["SourceSink"] == "sink"

source_data = otu_table.loc[source_mask]
sink_data = otu_table.loc[sink_mask]
source_envs = metadata.loc[source_mask, "Env"]

# 运行SourceTracker2
# 命令行方式(推荐)：
# sourcetracker2 gibbs \
#   -i otu_table.biom \
#   -m mapping.txt \
#   --source_sink_column SourceSink \
#   --source_column_value source \
#   --sink_column_value sink \
#   --source_category_column Env \
#   -o results/ \
#   --jobs 4

# === 命令行运行 ===
sourcetracker2 gibbs \
  -i otu_table.biom \
  -m mapping.txt \
  --source_sink_column SourceSink \
  --source_column_value source \
  --sink_column_value sink \
  --source_category_column Env \
  --loo \
  -o results/ \
  --draws_per_restart 10 \
  --burnin 100 \
  --restarts 5 \
  --jobs 8

# 输出文件:
# results/
#   mixing_proportions.txt    : 各sink中source比例
#   mixing_proportions_stds.txt : 标准差
#   source_predictions/       : LOO验证结果

第三步：FEAST快速源追踪¶

白话解释： FEAST(Fast Expectation-mAximization for Source Tracking)比SourceTracker2快100倍以上，特别适合大数据集。它用EM算法替代了Gibbs采样，数学上等价但计算效率高得多。

技术原理： - 基于期望最大化(EM)算法 - 将source组成建模为Dirichlet-Multinomial - 通过交替更新mixing比例(E步)和source参数(M步) - 支持多个source环境的嵌套结构

代码示例：

# === 安装FEAST ===
# devtools::install_github("cozygene/FEAST")
library(FEAST)

# === 准备输入数据 ===
# OTU表: 行=样本, 列=OTU/ASV (count数据)
otu_table <- read.csv("otu_counts.csv", row.names = 1)

# Metadata
metadata <- data.frame(
  Env = c(rep("Soil", 5), rep("Sewage", 5), rep("Gut", 3), 
          rep("River", 4)),  # 环境类型
  SourceSink = c(rep("Source", 13), rep("Sink", 4)),  # source/sink标记
  id = 1:17,
  row.names = rownames(otu_table)
)

# === 运行FEAST ===
# FEAST逐个sink运行
results_list <- list()

for (i in which(metadata$SourceSink == "Sink")) {
  # 当前sink
  sink_sample <- otu_table[i, , drop = FALSE]

  # 对应的source
  source_indices <- which(metadata$SourceSink == "Source")
  source_samples <- otu_table[source_indices, ]
  source_envs <- metadata$Env[source_indices]

  # 组装FEAST输入
  combined <- rbind(sink_sample, source_samples)
  meta_feast <- data.frame(
    Env = c("Sink", source_envs),
    SourceSink = c("Sink", rep("Source", length(source_indices))),
    id = c(i, source_indices),
    row.names = rownames(combined)
  )

  # 运行
  feast_result <- FEAST(
    C = as.matrix(combined),
    metadata = meta_feast,
    different_sources_flag = 1,  # 1=不同source有不同Env
    dir_path = "feast_output/",
    outfile = paste0("sink_", i)
  )

  results_list[[rownames(otu_table)[i]]] <- feast_result
}

# === 结果整合与可视化 ===
library(ggplot2)
library(tidyr)

# 读取FEAST结果
results_df <- do.call(rbind, lapply(names(results_list), function(sink) {
  res <- results_list[[sink]]
  data.frame(Sink = sink, Source = names(res), Proportion = as.numeric(res))
}))

# 堆叠柱状图
ggplot(results_df, aes(x = Sink, y = Proportion, fill = Source)) +
  geom_bar(stat = "identity") +
  theme_minimal() +
  scale_fill_brewer(palette = "Set3") +
  labs(title = "FEAST Source Tracking Results",
       y = "Proportion", x = "Sink Sample") +
  theme(axis.text.x = element_text(angle = 45, hjust = 1))

第四步：方法对比与选择指南¶

对比表：

特征	SourceTracker2	FEAST
算法	Gibbs采样(贝叶斯)	EM算法
速度	慢(分钟到小时)	快(秒到分钟)
精度	高	高(与ST2相当)
Unknown估计	内置	内置
大数据集	较慢	适合
不确定性估计	提供后验分布/SD	仅点估计
语言	Python	R
稀有taxa敏感度	较好	中等
适用场景	小数据集/需要置信区间	大数据集/快速筛查

选择建议：

choose_source_tracker <- function(n_samples, need_CI, compute_time_limit) {
  if (n_samples > 500 | compute_time_limit == "minutes") {
    return("FEAST - 大数据集或需快速结果")
  }
  if (need_CI) {
    return("SourceTracker2 - 需要不确定性估计")
  }
  return("两者皆可, FEAST更快")
}

第五步：模拟数据验证与结果评估¶

白话解释： 源追踪结果靠不靠谱？用模拟数据验证——已知真实混合比例的情况下，看工具能不能正确还原。

代码示例：

# === 模拟混合数据 ===
simulate_mixed_community <- function(source_profiles, mixing_props, n_reads = 10000) {
  # source_profiles: K个source的相对丰度矩阵(K × P)
  # mixing_props: 混合比例向量(长度K+1, 最后一个是unknown)

  K <- nrow(source_profiles)
  P <- ncol(source_profiles)

  # 各source贡献的reads数
  n_from_source <- round(mixing_props[1:K] * n_reads)
  n_unknown <- n_reads - sum(n_from_source)

  # 从各source抽样
  sink_counts <- rep(0, P)
  for (k in 1:K) {
    if (n_from_source[k] > 0) {
      sink_counts <- sink_counts + rmultinom(1, n_from_source[k], source_profiles[k,])
    }
  }
  # Unknown (均匀分布)
  if (n_unknown > 0) {
    sink_counts <- sink_counts + rmultinom(1, n_unknown, rep(1/P, P))
  }

  return(as.vector(sink_counts))
}

# 创建模拟数据
set.seed(42)
n_taxa <- 200
n_sources <- 3
source_profiles <- matrix(rgamma(n_sources * n_taxa, shape = 0.5), 
                          nrow = n_sources, ncol = n_taxa)
source_profiles <- source_profiles / rowSums(source_profiles)

# 模拟多个sink
true_props <- matrix(c(
  0.5, 0.3, 0.1, 0.1,
  0.1, 0.6, 0.2, 0.1,
  0.3, 0.3, 0.3, 0.1
), nrow = 3, byrow = TRUE)

simulated_sinks <- t(apply(true_props, 1, function(p) {
  simulate_mixed_community(source_profiles, p, n_reads = 50000)
}))

# 运行源追踪并与真实比例比较
# ... (运行FEAST/SourceTracker2)
# 评估指标: MAE, RMSE, Pearson相关

第六步：实际案例——水环境污染溯源¶

代码示例：

library(FEAST)
library(ggplot2)
library(ggalluvial)

# === 水环境源追踪案例 ===
# Source类别: 
#   - 人粪便(Human fecal)
#   - 动物粪便(Animal fecal: 牛/鸡/猪)
#   - 土壤(Soil)
#   - 污水处理厂出水(WWTP effluent)
# Sink: 下游河水样本

# 结果可视化 - 冲积图(Alluvial plot)
results_long <- data.frame(
  Sink = rep(c("Site_A", "Site_B", "Site_C", "Site_D"), each = 5),
  Source = rep(c("Human", "Cattle", "Poultry", "Soil", "Unknown"), 4),
  Proportion = c(
    0.35, 0.15, 0.05, 0.25, 0.20,  # Site A
    0.45, 0.10, 0.08, 0.15, 0.22,  # Site B
    0.20, 0.30, 0.15, 0.20, 0.15,  # Site C
    0.10, 0.05, 0.03, 0.45, 0.37   # Site D
  )
)

# 堆叠柱状图
ggplot(results_long, aes(x = Sink, y = Proportion, fill = Source)) +
  geom_bar(stat = "identity", width = 0.7) +
  scale_fill_manual(values = c("Human" = "#E41A1C", "Cattle" = "#377EB8",
                                "Poultry" = "#4DAF4A", "Soil" = "#984EA3",
                                "Unknown" = "grey70")) +
  theme_minimal() +
  labs(title = "Microbial Source Tracking - River Sites",
       y = "Source Proportion", x = "Sampling Site") +
  theme(legend.position = "right")

实战命令¶

# === 环境安装 ===
# SourceTracker2
pip install sourcetracker

# FEAST (R)
Rscript -e 'devtools::install_github("cozygene/FEAST")'

# === SourceTracker2 完整运行 ===
sourcetracker2 gibbs \
  -i feature_table.biom \
  -m sample_metadata.txt \
  --source_sink_column SourceSink \
  --source_column_value source \
  --sink_column_value sink \
  --source_category_column Env \
  --loo \
  --draws_per_restart 10 \
  --burnin 100 \
  --restarts 5 \
  --jobs 16 \
  -o sourcetracker_results/

# === QIIME2集成 (如果在QIIME2生态中) ===
qiime sourcetracker2 gibbs \
  --i-feature-table table.qza \
  --m-sample-metadata-file metadata.tsv \
  --p-source-sink-column SourceSink \
  --p-source-column-value source \
  --p-sink-column-value sink \
  --p-source-category-column Env \
  --o-mixing-proportions proportions.qza

面试常问点¶

Q1：微生物组源追踪的基本原理是什么？¶

A：源追踪假设sink样本(未知来源)的微生物群落是多个已知source群落的混合。通过概率模型(如Dirichlet-Multinomial)估计各source对sink的贡献比例。核心是将sink中每个taxa的丰度分解为各source贡献的加权和，加上一个"unknown"成分代表未被已知source解释的部分。

Q2：SourceTracker2和FEAST的算法区别是什么？¶

A： SourceTracker2使用Gibbs采样(MCMC方法)逐个将序列分配到source，可以获得后验分布和不确定性估计，但计算慢。FEAST使用EM算法迭代优化mixing比例和source参数，收敛快(通常几十次迭代)，速度比SourceTracker2快约100倍，但只输出点估计。两者精度相当。

Q3："Unknown"比例该如何解读？¶

A： Unknown代表sink中无法被任何已知source解释的微生物比例。高Unknown可能说明：(1) 存在未采样的source环境；(2) source采样不够充分(样本量少)；(3) sink中有独特的环境特异性微生物。实践中，Unknown > 30%可能提示需要补充更多source类型或增加source样本量。

Q4：源追踪对输入数据有什么要求？¶

A： (1) Source样本数量：每个source类别建议至少5-10个代表性样本(越多越稳定)；(2) 测序深度：建议rarefaction到合理深度(如10,000 reads/sample)，确保source和sink深度一致；(3) 分类分辨率：通常使用genus或ASV级别；(4) Source覆盖度：应尽量覆盖sink中可能的所有来源。

Q5：源追踪在公共卫生中有哪些应用？¶

A： (1) 水体粪便污染溯源(人粪vs动物粪vs污水处理厂)；(2) 医院感染传播链追踪(患者→环境→患者)；(3) 食品安全中的污染来源鉴定；(4) 新生儿肠道菌群来源分析(母体阴道/皮肤/口腔/环境)；(5) 室内微生物组来源(人体/室外/建筑材料)。

易错点¶

1. Source样本代表性不足¶

错误： 每个source类别只有1-2个样本。 正确做法： 每个source至少5-10个样本以捕获组内变异性。Source样本量不足会导致source profile估计不准确，下游所有比例都不可靠。

2. 未做rarefaction或深度归一化¶

错误： Source测序深度50000, Sink深度5000，直接分析。 正确做法： 先rarefaction到相同深度，或使用相对丰度。深度差异会系统性偏倚结果。

3. 遗漏重要source¶

错误： 水体源追踪只考虑人粪和动物粪，忽略了土壤和植物来源。 正确做法： 尽可能覆盖所有潜在source。遗漏的source的贡献会被错误分配给已有source或inflating unknown。

4. 对同一环境的不同批次未考虑时间变异¶

错误： 用夏季采集的source数据去追踪冬季的sink。 正确做法： 微生物群落有季节性变化。Source采样应尽可能与sink时间匹配，或多季节采样覆盖时间变异。

5. 将比例当作绝对丰度¶

错误： "Site A有35%来自人粪便"解读为"Site A有大量人粪便微生物"。 正确做法： 比例反映的是群落组成的来源，而非绝对生物量。高比例可能只是因为其他source贡献极低。需要结合指示微生物(如HF183)的qPCR来验证绝对量。

补充知识¶

其他源追踪工具¶

mSourceTracker: 基于机器学习的源追踪(SVM/RF)
BayesianSourceTracking: R包，贝叶斯非参数方法
STENSL: 基于稀疏tensor的源追踪
指示微生物qPCR: 非测序方法，如HF183(人特异性)、Rum2Bac(反刍动物)

实验设计建议¶

Source样本：每个环境类别10+个样本
Sink样本：根据研究问题设计采样点
测序深度：至少10,000-20,000 reads/sample
分类数据库：SILVA/Greengenes2, 保持一致性
阴性对照：含空白提取和测序对照，排除污染

与其他微生物组分析的整合¶

源追踪 → 差异丰度分析：鉴定各source特异性指示taxa
源追踪 → 网络分析：不同source微生物的互作模式
源追踪 → 功能预测：各source贡献的功能潜力差异

152_微生物组源追踪¶

一句话概述¶

核心知识点总览¶

各步骤详解¶

第一步：源追踪的基本概念¶

第二步：SourceTracker2分析¶

第三步：FEAST快速源追踪¶

第四步：方法对比与选择指南¶

第五步：模拟数据验证与结果评估¶

第六步：实际案例——水环境污染溯源¶

实战命令¶

面试常问点¶

Q1：微生物组源追踪的基本原理是什么？¶

Q2：SourceTracker2和FEAST的算法区别是什么？¶

Q3："Unknown"比例该如何解读？¶

Q4：源追踪对输入数据有什么要求？¶

Q5：源追踪在公共卫生中有哪些应用？¶

易错点¶

1. Source样本代表性不足¶

2. 未做rarefaction或深度归一化¶

3. 遗漏重要source¶

4. 对同一环境的不同批次未考虑时间变异¶

5. 将比例当作绝对丰度¶

补充知识¶

其他源追踪工具¶

实验设计建议¶

与其他微生物组分析的整合¶

📚 相关文章推荐