转载自yizhuo
基因、测序、分析
基因，生命的基本因素，是人类和其他生物的基础遗传物质。人有 23 对染色体，总共记录了大约 3Gb 个碱基（这里的 b 是 base，即碱基，可不是 bit，参考这里），每个位置上的碱基可能是 ATCG 中的一个。简单理解起来，就是有了这 3Gb 长的字符串，就能克隆一个你。基因测序，就是用化学和物理的方法，把你身体里这 3Gb 字符串检测出来。
当然，由于受当前测序技术的限制，我们并不能一次性测得一个完整的 3Gb 字符串，而是无数个 150bp 左右长度的小碎片。把这无数个小碎片重新组合还原成 3Gb 的长字符串的过程，叫全基因组组装。人类基因组计划干的就是这个组装拼图的事情，到了 2003 年，基本上算是拼完了。于是就有了一个标准的 3Gb 人类字符串，业界称其为『人类基因组参考序列』，也常被简称为『参考序列』。
基因测序分析，就是检测每个个体身上的这 3Gb 字符串，然后再跟标准的字符串（参考序列）做比对，来看是不是哪里发生了变异，是否有已知的遗传疾病风险。
现在随着基因测序技术的成熟，成本在飞速下降，所谓『旧时王谢堂前燕，飞入寻常百姓家』简直就是指日可待的事情。
于是越来越多的基因数据需要计算分析。
原来的分析计算流程
测序出来的数据是一大堆长度 150bp 的小碎片，但由于已经有了完整的人类参考序列的拼图，那么在这个拼图上寻找位置要比还原拼图的过程容易很多。多年来，从输入碎片数据比对到标准的 3Gb 参考序列，再到变异的检测，形成了如 BWA、Samtools/Picard、GATK 等为单机运行优化良好的业界公认软件。
将这些软件组合使用起来可以形成一个基因组数据分析的流程。华大基因基于阿里云 ECS 的 BGI Online 平台、七桥在 Amazon 和 Google 云上架构的 IGOR 分析平台，都有这样的业界金标准流程（下图为七桥，1.7G 左右的输入数据跑了 1 小时 48 分）。
但是上面的分析只是外显子测序分析（WES），所谓外显子就是那些会直接影响我们外在性状和健康与否的序列，它的长度只占了全基因组总长的 1%。因此需要测序的数据量比全基因组测序分析（WGS）少得多。使用 WES 的原因和意义，一方面是因为我们对这 3Gb 的人类基因组所蕴含的全面信息所知甚少，这 1% 的序列是我们目前能够很有效进行解读的位点，同时，也是由于它直接影响着我们的性状和健康，所以在很多疾病的研究中，通常只检测和分析这些位点和区域；另一方面也是因为全基因组测序目前成本相对较高，数据量巨大，单机计算时间基本要几天。
基于 MaxCompute 的方案
数据表达
来看基因数据的几个重要表达。
测序出来的原始输入数据 FastQ 格式（截取了两个 record），原本就是四行一个 record 的文本：
@ERR194147.1 HSQ1004:134:C0D8DACXX:1:1104:3874:86238/1
GGTTCCTACTTCAGGGTCATAAAGCCTAAATAGCCCACACGTTCCCCTTAAATAAGACATCACGATGGATCACAGGTCTATCACCCTATTAACCACTCACG
+
CC@FFFFFHHHHHJJJFHIIJJJJJJIHJIIJJJJJJJJIIGIJJIJJJIJJJIJIJJJJJJJJJJIJHHHHFFFDEEEEEEEEDDDCDDEEDDDDDDDDD
@ERR194147.2 HSQ1004:134:C0D8DACXX:2:2104:2852:75174/1
ACTTCAGGGTCATAAAGCCTAAATAGCCCACACGTTCCCCTTAAATAAGACATCACGATGGATCACAGGTCTATCACCCTATTAACCACTCACGGGAGCTC
+
@BBFDFFFHFFHHIHIJJJJFIHHFHFHJCIHFHIJJJJJJJIJIJIJJIIHIJJJJJJJBEGIGHIHGHHHEFCDFFEDEEDEEDDD?CCCDDDDDDDDC
AI 代码解读
测序产生的中间及最终结果的 VCF 格式（截取两个 record），也是一个多行表头，然后每行一个数据 record 的形式。
##fileformat=VCFv4.1
##...
##other headers
#CHROM POS ID REF ALT QUAL FILTER INFO
1 10019 rs376643643 TA T . . OTHERKG;R5;RS=376643643;RSPOS=10020;SAO=0;SSR=0;VC=DIV;VP=0x050000020001000002000200;WGT=1;dbSNPBuildID=138
1 10054 rs373328635 CAA C,CA . . NOC;OTHERKG;R5;RS=373328635;RSPOS=10055;SAO=0;SSR=0;VC=DIV;VP=0x050000020001000002000210;WGT=1;dbSNPBuildID=138
AI 代码解读
用表来表达描述完全没问题。当然，PairEnd 的 FastQ 总是成对出现，那么就编排成 8 列的表就好了。VCF 格式的多行表头也要稍微处理一下，不过这都不是事儿。
问题可分割
这块涉及到生物学的专业知识，要描述清楚比较困难。基本的想法，是将 3Gb 标准数据切成若干份，然后分布式执行每一份的计算。理论上，切割肯定会对计算的精度带来影响。但是影响多大要试过才知道。实践中，也有些办法来弥补这些影响：比如适当扩大 reference 的 interval，或者通过数据重叠来缓解。
简言之，跟单机计算结果对比，一致率 98.8%，这个数字可以接受。
MaxCompute Streaming 作业
如果数据如何在分布式系统里流动都想得很清楚的话，那么最大的障碍其实已经不在，剩下的就都是工程问题。
MaxCompute 支持 Streaming 作业，允许用户使用任何语言开发 MR 作业，开启隔离环境之后甚至可以用在公共云环境。用 Streaming 作业，可以非常方便的将 BWA、samtools、GATK 集成在一个 MR 作业当中：
这个作业，运行时间不到 3 个小时。
结语
或许很快，癌症患者可以当天拿到靶向药物是否有疗效的根据；新生病儿的早期诊断当天就可以排除是否遗传疾病；基因测序分析就像今天验血一样成为普通人司空见惯的体检手段。
或许将来，MaxCompute 可以更加发挥计算的能力，辅助挖掘 DNA 上目前 99% 还未知的领域。
我们期待这一天早些到来。
欢迎加入“数加·MaxCompute购买咨询”钉钉群（群号： 11782920）进行咨询，群二维码如下：
