Kraken报告格式详解
说实话，我第一次拿到Kraken报告的时候，盯着那几列数字看了半天，心里直犯嘀咕——这到底是什么意思？后来做项目遇到过一次数据解析错误，才发现读懂报告格式是正确使用Pavian的前提。这一章我就带你把这6列数据彻底弄明白。
为什么理解Kraken报告是关键一步
你想想看，Pavian本身不产生分类结果，它只是一个可视化工具。喂给它什么数据，它就展示什么图表。如果你连输入格式都搞不清楚，后面导出表格时很容易误解数据。
举个例子，我遇到过一位同事，他把"该层级的读段数"当成了"总读段数"，结果物种相对丰度算出来全是错的。嗯，这种坑能避就避吧。
六列数据的来龙去脉
Kraken报告的每一行都有6列，用制表符隔开。没有列名，没有额外注释。直接上结构：
| 列号 | 字段 | 类型 | 含义 |
|---|---|---|---|
| 1 | 百分比 | float | 该分类的读段占全部读段的百分比 |
| 2 | 总读段数 | int | 该分类及其所有子分类的读段总数 |
| 3 | 该层级读段数 | int | 刚好分类到该层级的读段数（不含子分类） |
| 4 | 分类层级 | string | 分类学级别：D、P、C、O、F、G、S、U等 |
| 5 | NCBI分类ID | int | NCBI Taxonomy数据库中的唯一标识 |
| 6 | 物种名称 | string | 科学名称，包含层级前缀（如"Escherichia coli"） |
下面这个SVG图把6列结构和它们之间的关系画清楚了，建议你保存下来当参考。
百分比列 — 看起来简单，暗藏玄机
第一列的百分比，是该分类的读段占所有读段（包括未分类）的比例。注意，这里的分母是总读段数，不是已分类读段数。这一点经常被人忽略。
我遇到过一个案例：某样本的未分类读段占80%，剩下的20%里大肠杆菌占了5%。在报告里大肠杆菌的百分比显示的是5%吗？不对，显示的是1%（因为分母是总读段数）。你想想看，如果没搞清这个，后面的丰度分析全偏了。
两列读段数 — 最容易搞混的地方
第2列和第3列，乍一看都是数字，但含义完全不同。
- 第2列（总读段数）：从该分类往下，所有子分类的读段加起来。说白了，就是这棵分类子树上所有读段的总和。
- 第3列（该层级读段数）：刚好停在这个分类上的读段。它没有被细分到更具体的子分类。
举个例子：
1.20 8000 2000 G 561 Escherichia
0.90 6000 6000 S 562 Escherichia coli这里，"Escherichia"属的总读段数是8000，但其中只有2000个读段刚好停在了"属"这个层级，剩下的6000个读段被进一步分到了"Escherichia coli"这个种。说白了，第3列是"停在这里的"，第2列是"从这里往下全部的"。
分类层级与ID — 生物信息学的核心
第4列是层级代码，我用过最常见的就这几个：
| 代码 | 全称 | 示例 |
|---|---|---|
| D | Domain（域） | Bacteria |
| P | Phylum（门） | Proteobacteria |
| C | Class（纲） | Gammaproteobacteria |
| O | Order（目） | Enterobacterales |
| F | Family（科） | Enterobacteriaceae |
| G | Genus（属） | Escherichia |
| S | Species（种） | Escherichia coli |
| U | Unclassified（未分类） | — |
| -1 | 无法分类 | — |
第5列是NCBI分类ID，一个数字ID，对应NCBI Taxonomy数据库里的唯一分类节点。我个人习惯在导出表格时保留这一列，方便后续跟其他数据库做关联查询。
物种名称列 — 不是你想的那样简单
最后一列是物种的科学名称。但有个细节值得注意：Kraken报告里的名称是带层级前缀的。
比如，你看到的是这样：
0.50 3000 3000 S 9606 Homo sapiens名称就是"Homo sapiens"，没有额外前缀。但有些老版本的Kraken报告会在名称前面加一些特殊字符或空格，用来表示缩进层级。如果你用文本编辑器打开报告，可能会看到类似这样的格式：
0.12 5000 1200 S 12345 Escherichia coli
0.05 2000 2000 G 561 Escherichia注意"Escherichia"前面多了两个空格。这不是错误，是Kraken用来表示缩进关系的方式。Pavian在解析时会自动处理这些缩进，但如果你自己写脚本解析，记得先 strip 一下。
line.strip().split('\t') 来分割每一行。这样可以去除首尾空白，避免被缩进空格坑到。
实战：快速检验一个Kraken报告
拿到一个Kraken报告后，我通常先做三件事来验证格式是否正确：
- 统计行数 — 用 wc -l 看有多少行。一般几百到几千行是正常的。
- 检查列数 — 用 awk 检查每行是否都有6列。
- 看第一行 — 确认第一列是数字（百分比），第六列是字符串。
嗯，这里可以用一个简单命令来验证：
# 查看报告的前3行，检查格式
head -3 sample.kreport
# 输出示例：
0.12 5000 1200 S 12345 Escherichia coli
0.08 3200 3200 S 9606 Homo sapiens
0.05 2000 1800 G 561 Escherichia
# 检查每行的列数是否都是6
awk '{print NF}' sample.kreport | sort | uniq -c
# 输出应该是：
# 100 6
# 表示所有100行都是6列如果发现有行不是6列，千万别直接塞进Pavian。我曾经遇到过文件名里带空格导致解析错位的情况，查了好久才发现是数据预处理的问题。
dos2unix sample.kreport 转换一下。
好了，这章的内容不少，但都是硬知识。把Kraken报告的6列结构吃透，后面用Pavian导出表格时你就能清楚地知道每个数据是从哪里来的，不会出现"这个数字怎么对不上"的尴尬。