标签：
主页：github: PacificBiosciences/FALCON
Falcon是一组通过快速比对长reads，从而来consensus和组装的工具。
Falcon工具包是一组简单的代码集合，我使用它们来研究单倍体和二倍体基因组的高效组装算法。
为了提高计算速度，它有一些后台代码是使用C来实现的，为了方便一些简单的前端是用Python编写的。
Falcon不是一个傻瓜的组装工具（除了很小的基因组），为了得到最好的结果，你可能需要了解各种分布式计算系统和一些基本的基因组组装理论。FAQ页面有一些常见的问题及其回答。
主要有如下几个步骤：
简单翻译：
每个步骤使用不同的命令行工具，实现不同的算法来完成这项工作。同样，每一步的计算需求也大不一样。假定用户已经拥有合理数量的计算资源。
例如，在合理的时间内装配100M大小的基因组，可能需要至少32核cpu 和 128 Gb RAM。代码是按照集群计算环境来编写的。需要一个作业队列来进行长时间的脚本工作，以及cpu-rich的计算工作队列。
一个包含了一组sub-reads的文件，fc_run.py脚本可以驱动工作流来管理和检查数据依赖和提交每一步的工作，从给定的数据中生成一个组装的草图。
fc_run.py是工作流驱动脚本，需要运行在一台计算机上，允许持续长时间地工作一段时间来完成整个装配过程。它需要一个配置文件fc_local.cfg作为单个输入。原始序列数据的输入文件包含在配置文件中。
配置文件可用于控制计算资源的使用和重要算法参数，根据输入数据集达到最有组装的效果。然而，没有万能的方式来猜测各种具体项目的最合适的选项和参数。参数调优的最好方法是了解一些组装理论，一点点实现方法，这样才能正确的了解更改某些选项对组装结果的影响。快速浏览reads长度分布、整体覆盖和调整相应的选项也是非常重要。
这里我们会详细介绍基因组分层组装策略 以及 fc_run.py 各个参数的含义
这个版本的 Falcon 的 overlapping 是用 Gene Myers‘ Daligner 的修改版本来做的（参考 github 的 forked version）。主要的修改是adapting some I/O for the downstream bioinformatics process。可以通过git diff来查看有哪些修改。
git diff #查看版本之间差异input_fofn 选项指向了包含所有输入数据的文件。来自Daligner的fasta2DB 在 fc_run.py 会被调用（I/O集中器，会在你执行 fc_run.py的电脑节点出运行，如果这在你的服务器集群中是个问题，你必须修改代码来打包相关的操作到一个脚本中，以便提交到你的任务管理系统中）。
input_fofn #指向了包含所有输入数据的文件这个版本的 fc_run.py 支持从错误校正的reads中运行组装程序（可以跳过错误校正环节，直接开始组装）
input_type = preads #输入为 all error-corrected reads，会跳过错误校正，直接开始最终的重叠群组装步骤
input_type = raw #原始数据你需要决定cutoff的长度，一般，选择一定的阈值，使得你能得到15x 到 20x的覆盖度是比较好的。 然而，如果计算资源充足，校正reads你有其他用途，你可以设置lower length cutoff，来获得更多的error corrected reads。
最终组装的时候，更多的reads不意味着更好的组装结果，因为会有噪音。
有时，尝试不同的length_cutoff_pr是有意义的，较初次overlapping step for error correction而言，后期计算更cheap。
我们的策略是选择一个短的length_cutoff 进行一次计算。之后， 我们使用不同的length_cutoff_pr来获得更好的组装结果。
length_cutoff #控制错误校正过程中的cutoff
length_cutoff_pr #控制之后组装重叠群步骤的cutoffpa_concurrent_jobs选项控制由fc_run.py 提交的 number of concurrent jobs。sge_option_da 和 sge_option_la控制作业队列和daligner jobs 的 number of slots。daligner 使用的线程数默认是4。然而，依据集群设置、计算节点的内存数量，你可以使用超过4的slots。为了选择最好的数量，你最好咨询你们实验室的HPC专家，做一些前期的测试。
pa_concurrent_jobs #fc_run.py提交的当前作业数
sge_option_da #作业队列
sge_option_la #daligner作业的slots数最终运行作业的数量取决于你怎么"splits"序列数据库，仔细阅读Gene Myers‘s blog ( http://dazzlerblog.wordpress.com ) ，了解如何设置 pa_DBsplit_option 和 pa_HPCdaligner_option选项。一般，对大型基因组而言，在pa_DBsplit_option中，你应该使用-s400 (400Mb sequence per block)。这样作业的数量会变少，但是持续时间会延长。如果，作业调度系统限制了一个作业可以运行的时间，你应该设置a smaller number for the -s option。
另一个影响作业总数的参数是pa_HPCdaligner_option中的 -dal option。它决定how many blocks are compared to each in single jobs。大的 -dal 给出了大的作业，但是数量少；反之。小的，作业小，但是提交的数量多。
这个工作流中，daligner生成的trace point没有被用到（为了高效，应该使用trace point，但是首先你的知道如何正确的pull out它们）。pa_HPCdaligner_option 中的 -s1000 使得trace point变稀疏，以此节约磁盘空间。我们使用-l1000忽视短于1kb的reads。
pa_DBsplit_option #-s400 for large genomes smaller number for the -s
pa_HPCdaligner_option #
-dal # determines how many blocks are compared to each in single jobs;
-s1000 #
-l1000 #daligner的输出时一组.las文件，它包含了reads间相互比对的信息。Such information is dumped as sequences for error correction by a binary executable LA4Falcon to fc_consensus.py.（什么意思？） fc_consensus.py生成了consensus（这部分是用C写的）。
fc_consensus.py有许多选项，可以使用falcon_sense_option来控制它。大部分情况下，--min_cov 和 --max_n_read是最重要的选项。--min_cov控制了什么时候a seed read gets trimmed or broken due to low coverage（什么意思？），--max_n_read puts a cap on the number of reads used for error correction（什么意思？）。
在高度重复的基因组汇总，设置小的 --max_n_read，来确保 consensus code does not waste time aligning repeats。最长的合适的overlaps 被用来进行correction ，来减少the probability of collapsed repeats。
你可以使用cns_concurrent_jobs来设置提交给任务管理系统的当前作业的最大数量。
falcon_sense_option #
--min_cov #when a seed read gets trimmed or broken due to low coverage???
--max_n_read #
cns_concurrent_jobs #这部分与最初的overlapping 十分相似，只有部分修改，daligner只以本地的raw reads作为默认值。
fc_run.py generates a fasta file of error-corrected reads where the fasta header is parse-able by daligner.（什么意思？）
下面的参数控制这个步骤的计算过程：
sge_option_pda = -pe smp 8 -q jobqueue
sge_option_pla = -pe smp 2 -q jobqueue
ovlp_concurrent_jobs = 32
ovlp_DBsplit_option = -x500 -s50
ovlp_HPCdaligner_option = -v -dal4 -t32 -h60 -e.96 -l500 -s1000这个设置与第一步的overlapping 相平行，主要的区别是ovlp_HPCdaligner_option的-e 选项。错误率小了很多，因此我们期待 much higher correlation between the p-reads。
不是所有的overlaps 都是独立的，所以利用某些filtering step来减少计算和组装图的复杂性是可能的。例如：如果一个read被全部包含在另一个read里，那么它们的overlap不会提供额外的构建基因组的信息。由于重叠关系的传递属性（transitive property），许多overlap 信息可以被简单的推测出来。
事实上，构建contigs 的第一步就是删除 transitive reducible edges（？？？），这意味着，我们只需要best n overlaps 在5‘ or 3‘ 端。overlap_filtering_setting option的--bestn parameter参数可以被用来控制每个read报告的maximum overlap。
overlap_filtering_setting
--bestn #control the maximum overlap reported for each read另一个有用的启发式，仅保留有平均 5‘ and 3‘ 覆盖的reads。这是因为，如果一个read以重复结束，它会有更高的重复。这种reads不提供价值，不能解决相关的重复问题。我们可以过滤掉它们，希望能有reads（跨越重复，在两端有正常的独特的anchors）。同样，如果一个reads一端的覆盖度太低，那它可能有太多的测序错误。这种 reads产生了 "spurs"在组装图中，通常都会将它们过滤掉。--max_cov 和 --min_cov就是用来过滤那些有太高或者太低的reads。
--max_cov #过滤太高的reads
--min_cov #过滤太低的reads过滤的脚本也允许过滤一些 "split" reads。如果一个read在5‘ 和 3‘端有非常不相等的覆盖度，这也是某一段是重复的信号。 --max_diff可以用来过滤掉一端比另一端的覆盖哦高很多的reads。
--max_diff #过滤掉一端比另一端的覆盖哦高很多的reads使用这些参数的正确的数量是多少？真的很难设置正确，如果 error corrected reads的总覆盖度比知道的length cut off 并且合理的高 (e.g. greater than 20x)，那么将平均的覆盖度min_cov设置为5，max_cov设置为3倍将会很安全， max_diff设置为2倍。然而，在低覆盖度的情况下，将min_cov设为1或2将会更好。一个名为fc_ovlp_stats.py的帮助脚本可以帮助dump（丢弃）。 A helper script called fc_ovlp_stats.py can help to dump the number of the 3‘ and 5‘ overlap of a given length cutoff, you can plot the distribution of the number of overlaps to make a better decision.（？？？）
你也可以设置max_diff和min_cov为很高的值，来避免过滤，在某些特殊的场合。
这个过滤的过程会显著过滤掉高度拷贝重复的信息。也就是，这些重复会被过滤掉，并且不会在最终的组装结果中显示。如果你对这些重复感兴趣，总之，单独处理这些重复会更加有效和实用。如何解决重复序列这个问题具有非常重要的意义。
鉴于overlapping 数据，string graph 由 fc_ovlp_to_graph.py来创建，使用默认的参数。 fc_ovlp_to_graph.py生成一些文件代表组装后的最终的 string graph。最终的ctg_path包含了每一个contig的图的信息。一个 contig是一个simple paths 和 compound paths的线性的path。Compound paths are those subgraph that is not simple but have unique inlet and outlet after graph edge reduction.
They can be induced by genome polymorphism or sequence errors. By explicitly encoding such information in the graph output, we can examine the sequences again to classify them later.
create draft contigs的最后一步就是find a single path of each contig graph，以及generate sequences accordingly。在a contig graph is not a simple path的情况下，我们寻找 end-to-end path that has the most overlapped bases。这被称为 primary contig。图中的每一个path，如果它与 primary one不同，我们会construct the associated contig。在the associated contigs are induced by sequencing error的情况下，the identity of the alternative contig and the primary contig will be high ( > 99% identity most of time)。在there are true structure polymorphism的情况下，there are typically bigger difference between the associated contigs and the primary contigs。
fc_graph_to_contig.py脚本，takes the sequence data and graph output to construct contigs。它会生成所有的 associated contigs at this moment。将来，一些后处理程序，来删除errors中的some of the associated contigs将会被开发出来。（正在开发）
daligner由pa_HPCdaligner_option 和ovlp_HPCdaligner_option控制。
为了限制内存，我们可以使用-M选项。对于人类基因组组装，我们测试 -M 32，每个daligner使用32G RAM。其他可能性正在调查中。
更多的daligner的细节，参见： dazzlerblog。
代码设计工作在单一目录。工作目录的典型布局是这样的：
$ ls -l
total 56
drwxr-xr-x 84 jchin Domain Users 8192 Nov 30 12:30 0-rawreads
drwxr-xr-x 18 jchin Domain Users 4096 Nov 30 12:33 1-preads_ovl
drwxr-xr-x 2 jchin Domain Users 4096 Nov 30 12:44 2-asm-falcon
-rwxr-xr-x 1 jchin Domain Users 1041 Nov 30 12:13 fc_run.cfg
-rw-r--r-- 1 jchin Domain Users 378 Nov 29 23:20 input.fofn
drwxr-xr-x 2 jchin Domain Users 4096 Nov 30 12:13 scripts
drwxr-xr-x 2 jchin Domain Users 24576 Nov 30 12:33 sge_log典型的input.fofn是这样的：
/mnt/data/Analysis_Results/m140913_050931_42139_c100713652400000001823152404301535_s1_p0.1.subreads.fasta
/mnt/data/Analysis_Results/m140913_050931_42139_c100713652400000001823152404301535_s1_p0.2.subreads.fasta
/mnt/data/Analysis_Results/m140913_050931_42139_c100713652400000001823152404301535_s1_p0.3.subreads.fasta目录0-rawreads包括：overlapping原始序列的所有的脚本和数据，它包含不同的job_ *和m_ *目录。
.
例如，如果我们将E. coli的数据分成20个块，目录就会如下：
cns_done job_00011 job_00024 job_00037 job_00050 m_00003 m_00016
da_done job_00012 job_00025 job_00038 job_00051 m_00004 m_00017
job_00000 job_00013 job_00026 job_00039 job_00052 m_00005 m_00018
job_00001 job_00014 job_00027 job_00040 job_00053 m_00006 m_00019
job_00002 job_00015 job_00028 job_00041 job_00054 m_00007 m_00020
job_00003 job_00016 job_00029 job_00042 job_00055 m_00008 preads
job_00004 job_00017 job_00030 job_00043 job_00056 m_00009 prepare_db.sh
job_00005 job_00018 job_00031 job_00044 job_00057 m_00010 raw_reads.db
job_00006 job_00019 job_00032 job_00045 job_00058 m_00011 rdb_build_done
job_00007 job_00020 job_00033 job_00046 job_00059 m_00012 run_jobs.sh
job_00008 job_00021 job_00034 job_00047 las_files m_00013
job_00009 job_00022 job_00035 job_00048 m_00001 m_00014
job_00010 job_00023 job_00036 job_00049 m_00002 m_00015job_ *目录存储了每个daligner的输出工作。
m_ *目录是用于合并工作。
有一些空文件，它们的名字是done结尾。这些文件的时间戳是用来跟踪工作流阶段。您可以修改时间戳和重启fc_run.py来触发为工作流程的一部分做某些做计算。然而，不推荐这样，除非你有一些时间来读fc_run.py的源代码，并且知道如何跟踪工作流内部的依赖性。(例如，如果你在成功组装后touch rdb_build_done，并且重新运行fc_run.py。因为所有中间过程取决于文件，rdb_build_done比任何中间文件的都要新，就会触发fc_run.py，重复整个装配过程。
.
las_files存储了最后比对信息。
如果你不打算重新运行工作，但是想知道如何比对的样子，你可以删除所有job_ *和m_ *目录但保留las_files和preads 目录。
.
这一步的主要输出存储在0-rawreads /preads。0-rawreads/preads内部的out.%04d.fa的是输出的reads的fasta文件。如果这一步成功完成, 哨兵文件cns_done将被创建。
该目录存储着p-read与p-read overlapping的数据。它与0-rawreads大体一致，但是却没有consensus 这一步。主要输出是1-preads_ovl/目录里的 *.las。
这是最终的输出目录。
它包含了组装图和draft contigs的所有信息，细节描述参照Graph output format一节。
下面是fc_ovlp_to_graph.py的使用信息：
usage: fc_ovlp_to_graph.py [-h] [--min_len MIN_LEN] [--min_idt MIN_IDT]
[--lfc]
overlap_file
a example string graph assembler that is desinged for handling diploid genomes
positional arguments:
overlap_file a file that contains the overlap information.
optional arguments:
-h, --help show this help message and exit
--min_len MIN_LEN minimum length of the reads to be considered for
assembling
--min_idt MIN_IDT minimum alignment identity of the reads to be considered
for assembling
--lfc use local flow constraint method rather than best overlap
method to resolve knots in string graphDe-dup example and strategy
TODO
Low coverage assembly
TODO
Assembly layout rule
Applications of the assembly graph
TODO
Reproducibility and replicability
Other TODOs
Incremental overlapping
Pre-processing repeat for overlapping
subread：子read，
full-pass subread：
pre-assembly：
error correction：
proper overlap：
string graph：
contig：
primary contig：
associated contig：
p-read：
compound path：
simple path：
Quiver：
科普时间：
什么是virtualenv？
每个应用可能需要各自拥有一套“独立”的Python运行环境。virtualenv就是用来为一个应用创建一套“隔离”的Python运行环境。
资料：
virtualenv - 廖雪峰
virtualenv -- python虚拟沙盒
Virtualenv入门基础教程
现在默认分支是“master”，其包含Falcon最新的源码。
当前最新的综合版本是v0.3.0，安装指南参照：v0.3+ Integration Installation Guide
上一个版本是v0.2.2，具体参照：v0.2.2 release github repository
标签：
原文地址：http://www.cnblogs.com/leezx/p/5724590.html
