首页 / 资讯中心 / 文章详情

用TBtools轻松完成拟南芥-水稻共线性分析,告别命令行

用TBtools轻松完成拟南芥-水稻共线性分析,告别命令行 ★ FEATURED ARTICLE
做基因组共线性分析很多人的第一反应就是又要打开命令行解压缩文件、跑BLAST、跑MCScanX、再写脚本画图。这一连串操作不要说初学者就连一些有几年经验的人也会觉得烦。我见过太多人数据已经下载好了论文里急需一张共线性图结果卡在命令行环境配置上整整一周。后来我换了个思路直接用TBtools的MCScanX插件把原来要敲十几条命令的流程压缩成几次鼠标点击拟南芥和水稻这两个模式物种之间的共线性分析5分钟就能跑出核心结果。这篇文章就是把我自己从0到1跑通这个流程的经验完整写出来包括数据下载、参数设置、结果解读、报错排查。目标是让没接触过命令行的人也能顺利出图。适合谁看一是准备做基因家族分析的同学二是想在论文里加共线性证据的科研党三是对生物信息学工具感兴趣但被命令行劝退的人。文章不会讲太深的算法原理但会把关键参数、文件格式、容易踩的坑都说清楚照着做基本能复现。1. 共线性分析的基本概念与应用场景1.1 什么是共线性分析为什么重要共线性英文叫Synteny或Collinearity通俗讲就是两个基因组之间基因在染色体上的排列顺序保持一致。你可以把它想成两栋建筑虽然建造年代不同、外观也不同但承重墙的位置、管线的走向却出奇地一致。这种一致性不是巧合而是因为两栋楼都沿用了同一张原始设计图。在基因组里这张设计图就是共同祖先的基因组。拟南芥和水稻在进化上分道扬镳了大约1.5亿年这期间基因组经历了大量的基因丢失、串联重复、倒位、易位等事件但仍然能在染色体上找到许多区域基因的排列顺序和方向保持高度一致。共线性分析就是系统性地把这些保守区域找出来形成一张基因组之间的对应关系图。从分析角度看共线性分析和简单的同源基因搜索不一样。同源基因搜索BLAST只看单个基因之间的相似性而共线性分析看的是一连串基因在染色体顺序上的呼应。后者证据更强更能说明两个基因组区域确实源自一个共同的祖先片段。这也是为什么共线性区块经常被称为同源区块或上古保守片段。1.2 拟南芥-水稻共线性分析的典型价值拟南芥是双子叶植物的模式物种水稻是单子叶植物的模式物种。这两个物种的基因组都已经组装得非常完善注释质量也很高做它们之间的共线性分析至少有三个层面的价值。第一只从进化角度看拟南芥-水稻的比较能揭示单双子叶植物分化后基因组结构的演化规律。比如哪些染色体片段更保守、哪些区域发生了频繁重排这些都能通过共线性区块的大小和密度直观看出来。第二从基因家族分析角度看这是最普遍的需求。比如你在水稻里鉴定了一个NBS-LRR抗病基因家族想看看拟南芥中对应的直系同源基因分布在哪些染色体区域共线性分析就是最直接的手段。通过共线性区块中的基因对可以区分哪些同源基因来自全基因组复制WGD哪些来自串联复制哪些是物种特异性的。第三从功能迁移角度看拟南芥作为研究最深入的植物大量基因有明确的功能注释。如果水稻某个基因与拟南芥中一个已知功能基因位于共线性区块内就可以基于共线性邻位关系推测这个水稻基因可能参与类似的生物学过程。这种方法在候选基因预测中非常实用。1.3 适合人群与前置知识清单这个流程的门槛其实很低但也不是完全零基础。我盘点了一下需要用到的前置知识大概有四块。第一必须知道什么是蛋白序列和基因注释文件。蛋白序列是fasta格式以大于号开头后面跟基因ID再下面是氨基酸序列。基因注释文件常见的是GFF3格式用于记录基因在染色体上的位置。学的不是特别精通没关系能识别文件后缀和大致结构就够用。第二要对染色体有概念。拟南芥有5条染色体水稻有12条染色体这个数字在后续参数检查和结果解读中会反复出现。第三最好能理解BLAST的基本原理。你不需要自己跑命令但至少要知道BLASTP是蛋白序列比对E-value越小结果越可信。第四电脑上要装好Java环境。TBtools是Java程序没有Java环境软件打不开。此外Windows用户建议把压缩软件、文本编辑器比如Notepad准备好后面会用到解压和查看结果文件。2. 工具选型与原理MCScanX在TBtools中的角色2.1 MCScanX算法核心逻辑MCScanX全称是Multiple Collineation Scan前身是多重共线性扫描工具。它检测共线性区块的核心思路可以拆成三步来理解。第一步接收BLASTP比对结果。两个物种分别提供蛋白序列BLASTP跑完之后会得到一个列表告诉你每条基因在其他物种中找到哪些同源基因以及这些同源基因的相似性得分。这一步等于是先把可能的同源基因对找出来。第二步把同源基因对映射到染色体上。这一步由注释文件GFF完成。GFF文件提供了每个基因所在的染色体ID、起始位置、终止位置和方向。MCScanX把BLAST结果中的基因对放到真实的染色体坐标上形成一个候选锚点集合。第三步用动态规划算法搜索共线性区块。从这些候选锚点中寻找那些在两条染色体上连续出现的同源基因对。动态规划的思想本质上就是如果某个基因对在两条染色体上的上下游基因也都能形成同源对就认为这段区域属于同一个共线性区块。算法会沿染色体顺序逐步延伸直到连续的同源基因对中断。参数中match_score控制同源基因匹配的得分gap_penalty控制允许的空位惩罚max_gaps限制在一个共线性区块内允许存在的最大非匹配基因数。这些参数都会影响最终检测到的区块边界。默认参数是在拟南芥、水稻这类基因组上反复调试出来的大部分情况下直接沿用即可。2.2 原生MCScanX与TBtools插件的对比原生MCScanX必须以命令行方式运行。它要求输入文件符合严格的命名规则假设你的文件前缀是 ath_osa那么你必须有 ath_osa.blast 和 ath_osa.gff 两个文件放在同一目录下然后执行makeblastdb -in ath.pep.fasta -dbtype prot blastp -query osa.pep.fasta -db ath.pep.fasta -out ath_osa.blast -outfmt 6 -evalue 1e-5 -num_threads 8 MCScanX ath_osa看起来没几行但操作中到处都是细节BLAST结果是不是标准outfmt 6格式GFF文件里有没有多余的空格染色体ID是否统一MCScanX是否在你的环境变量里。任何一环出问题程序都会报一些非常不友好的错误比如failed to parse line或者直接exit。我当时排查一个路径问题花了好几个小时最后发现只是Win系统下的反斜杠和正斜杠混用导致的。TBtools的MCScanX插件把这些坑都处理掉了。它通过图形界面接收你的fasta和gff路径自动在后面拼装makeblastdb、blastp、MCScanX这些命令行工具并自动检查格式、统一ID映射关系。对大多数人来说效果就是点几下鼠标结果文件就出来了。而且TBtools在可视化层面也做了深度集成结果可以直接导入画图组件。2.3 TBtools环境准备与版本选择TBtools目前有1.x和2.x两个大版本系列。1.x系列非常经典各种教程和插件资源丰富菜单结构相对稳定。2.x系列界面变化较大功能和性能有升级但有些老教程里的菜单位置对不上。我的建议是如果你是第一次使用直接下载官网或GitHub上的最新版本如果你需要在固定环境下复现操作选择一个版本后就不要频繁切换。安装方面Windows下解压后双击 exe 或 jar 文件即可macOS 和 Linux 下需要配置Java运行环境建议安装OpenJDK 17或更高版本。打开TBtools后如果看到Java Version相关提示按照提示操作就行。内存设置方面共线性分析涉及的序列数据不算特别大默认内存通常够用如果分析特别大的植物基因组可以在启动脚本中把-Xmx参数调大比如 -Xmx4g。比较重要的一点TBtools会保留上次运行的日志如果某一步报错可以在日志窗口里看到具体调用命令。这个功能在排查问题时非常有用我后面讲报错时会专门提到。3. 数据下载与预处理拟南芥和水稻数据准备3.1 基因组数据下载做拟南芥和水稻共线性分析数据源很多但为了省事我推荐统一从Ensembl Plants下载。访问 plants.ensembl.org搜索 Arabidopsis thaliana然后在基因组版本下拉框里选择 TAIR10搜索 Oryza sativa版本选择 MSU 7.0 或 RGAP 7.0。每个物种都需要下载两类文件。第一类是蛋白序列文件通常命名为 xxx.pep.all.fa.gz 或 xxx_pep.fasta.gz第二类是基因注释文件通常命名为 xxx.gff3.gz。注意不要下载成CDS序列文件因为MCScanX的输入要求是蛋白序列用CDS核酸序列会导致BLASTP阶段无法识别。下载完成后先解压然后用文本编辑器打开看几行。蛋白序列文件里的基因ID格式大致是这样的AT1G01010.1 MASS... Os01g0100100 MASS...GFF3文件里则记录了基因的详细信息。重点看第1列的染色体ID和第9列的注释属性。比如拟南芥GFF是Chr1开头水稻GFF可能是Chr1也有的版本是1开头。如果染色体ID对不上后面BLAST没问题但MCScanX建立位置关系时会出错。3.2 文件格式与ID一致性检查这一步是整个流程里最容易被忽视、却最容易导致失败的环节确保fasta的基因ID和gff的基因ID能对应上。以拟南芥为例TAIR10版本的蛋白序号可能是 AT1G01010.1而GFF3中的第9列写的是 IDAT1G01010.1;ParentAT1G01010.1;biotypeprotein_coding。这里ID可以对应上。但Ensembl下载的新版本有时会在ID后面加版本号比如 AT1G01010.10而gff里可能是 AT1G01010不一致的情况比比皆是。水稻的情况更特殊。MSU 7.0 注释里的基因ID可能是 LOC_Os01g01000而蛋白fasta里的ID可能是 Os01g0100100。两者通过gff里的Dbxref或Name属性对应。这需要你在下载后做一个简单检查。如果发现ID不一致有几种处理办法。最简单的用文本编辑器打开fasta和gff把带版本号的后缀统一去掉。但文件大时手动改不现实我通常用TBtools自带的File Tool或Sequence Toolkit里的ID处理功能批量删除或替换ID中的点号后缀。如果涉及两个不同命名体系的转换可以用TBtools的GFF Tools菜单下的ID映射功能选择基于gff本身建立的ID关系进行转换。3.3 目录整理与命名规范在正式开始分析前我强烈建议花两分钟整理好目录。这一步看着不起眼但能帮你少踩一半的坑。创建的工作目录路径不要出现中文、空格和特殊字符比如 D:\synteny_ath_osa 就比 D:\基因组共线性分析\data 稳妥得多。原因前面说过TBtools虽然界面是图形化但底层调用的BLAST和MCScanX还是命令行时代的程序对路径中的空格和中文支持不佳。目录内放四个文件就行文件建议命名作用拟南芥蛋白序列ath.pep.fastaBLASTP的query和subject拟南芥注释GFFath.gff3基因位置信息水稻蛋白序列osa.pep.fastaBLASTP的query和subject水稻注释GFFosa.gff3基因位置信息命名简洁统一后续在TBtools里选择输入时会非常方便。文件名不要太长前缀尽量简短这样后面生成的 .collinearity 文件名也不会因为太长而影响操作。4. 核心实操5分钟跑通共线性分析4.1 One Step MCScanX操作详解TBtools中的 One Step MCScanX 是把BLASTP和MCScanX这两个步骤封装在一起的集成功能界面上一共就几个输入框非常适合第一次使用。我在常用版本上验证过的操作路径是这样的打开TBtools后在顶部搜索框搜MCScanX选择Commonly Used下面的 One Step MCScanX。弹出的界面大致有四个重要区域CDS/Protein File A填入拟南芥蛋白序列 ath.pep.fastaGFF/GTF File A填入拟南芥注释 ath.gff3CDS/Protein File B填入水稻蛋白序列 osa.pep.fastaGFF/GTF File B填入水稻注释 osa.gff3数据填好后下面还有几个参数项。E-value默认一般是1e-5对于拟南芥和水稻这个亲缘关系1e-5是完全够用的。Threads这个参数代表线程数建议填你电脑CPU核心数减一比如8核机器填7跑起来会快一些。点击Start后TBtools会弹出一个运行窗口显示命令执行过程。整个过程一般会经历两个阶段第一阶段是BLASTP比对可能会持续几分钟第二阶段是MCScanX运行速度很快几秒到几十秒就完成。看到Finished字样后说明核心分析已经结束可以关闭窗口去找结果文件了。4.2 分步模式BLASTP MCScanX如果你想把每一步看得更清楚或者需要自己对BLAST结果做中间过滤可以选择分步模式。分步模式的操作是先运行TBtools中的BLASTP工具把两个蛋白序列文件传进去得到BLAST比对结果然后运行MCScanX工具将BLAST结果和GFF文件一起输入。使用分步模式的场景通常是你已有现成的BLAST结果文件不想重新跑一遍或者你想尝试不同的E-value阈值看几组参数下的BLAST输出对最终共线性结果的影响。这种情况下BLAST输出格式需要是MCScanX能识别的标准tabular格式也就是BLAST默认outfmt 6列顺序是query id、subject id、identity、alignment length、mismatches、gap opens、q.start、q.end、s.start、s.end、evalue、bit score。需要注意分步模式中BLAST结果文件和GFF文件的文件名前缀要一致。假设BLAST结果命名为 ath_osa.blast那么GFF文件也要命名为 ath_osa.gffMCScanX才能正确读取。TBtools插件界面上一般会提示文件命名要求跟着界面提示操作即可。4.3 参数选择与调整经验关于MCScanX插件中的参数我多说几句经验。默认参数在拟南芥和水稻这个组合上是经过大量验证的绝大多数情况下不需要动。但如果遇到以下特殊场景可以有针对性地调整。第一个场景亲缘关系较远的物种。拟南芥和水稻虽然分属单双子叶但共线性信号仍然比较强。如果你以后拿真菌和植物去做比较E-value可以放宽到1e-3同时可以降低min align score让算法更容易识别弱同源关系。第二个场景目标区域基因密度极其不均匀。比如某条染色体上有一段非常保守的基因簇而邻接区域是大量重复序列。如果max_gaps设得太小算法会在重复区域把共线性区块截断导致检测结果碎片化。这时可以适当增大max_gaps但不要超过50否则会引入很多假阳性。第三个场景想得到更“严苛”的共线性结论。如果你的论文需要强调高置信共线性区块可以把E-value收紧到1e-10并调整match_score到较高的值比如70或80。这样得到的区块数量会下降但每个区块的可信度更高。调整参数的原则是一切以生物学问题为准。不要为了追求区块数量而无限放宽参数也不要因为默认参数结果少就急着改先看看结果文件里有哪些内容再决定怎么调。5. 结果解读与可视化进阶5.1 .collinearity文件的深入解读运行完成后输出目录里会生成以你设置的输出前缀命名的核心结果文件常见的有两个以.collinearity为后缀的文件以及以.gff为后缀的处理后注释文件。.collinearity文件是文本格式可以用任何文本编辑器打开。文件开头通常有一段以## Alignment开头的区块声明每个区块代表一段共线性区域。区块内部每行记录一对共线基因的对应关系格式类似## Alignment 0: Chr1: 1234 - 5678 chr1: 4321 - 8765 0-1: AT1G01010.1 Os01g0100100 1e-100 2-3: AT1G01020.1 Os01g0100200 1e-87每一行的前两个数字代表这对基因分别在该物种基因列表中的索引冒号后是基因ID最后是这个基因对的E-value。一个区块中基因对越多说明这一段共线性关系越强。解读结果时除了看单个基因对的对应关系还需要关注区块的整体特征。比如区块内基因对的顺序是否完全一致有没有局部翻转两个物种间共线性区块的染色体分布是否偏好某几条染色体。拟南芥和水稻之间的共线性区块通常会覆盖大多数染色体但覆盖程度不同这些差异本身就很有生物学意义。5.2 点图与圆图的绘制拿到.collinearity文件后最关心的就是怎么画一张能放进论文的图。TBtools提供了两个主流的可视化组件Dual Synteny Plot和Circle Gene View。Dual Synteny Plot适合展示两个物种之间的共线性对应关系本质上是点图的增强版。横坐标是一个物种的染色体纵坐标是另一个物种的染色体每个点代表一个同源基因对。正相关的点连成斜线就代表共线性区域反向区域则呈负斜率。TBtools会在共线性区块背景上画斜线或色块让对应关系非常直观。操作上只需要加载.collinearity文件和对应的.gff文件然后可以选染色体范围、点颜色、连线颜色等。Circle Gene View适合展示一个物种内部或两个物种之间的整体共线性关系它把所有染色体放在一个圆环上共线性区块通过内部连线表示。这个图尤其适合展示全基因组复制事件WGD或大规模染色体重排。驱动Circle Gene View需要提供染色体长度信息以及共线性区块的连接信息。TBtools可以从.gff文件自动提取染色体长度连接信息从.collinearity自动解析所以作图过程非常简单。5.3 多物种共线性的扩展玩法当你掌握了两个物种的分析后多物种共线性分析只是把同样的流程重复几遍。比如你测定了水稻、玉米、高粱、二穗短柄草等禾本科物种的基因组可以把它们作为一组输入两两运行MCScanX后把多次结果整合到同一个圆图中。TBtools在一些版本中提供了Multiple Synteny Plot或Advanced Circos功能可以一次加载多个物种的共线性结果用不同颜色的线条区分不同物种之间的连接。我做禾本科多物种比较时就习惯把水稻作为参考物种分别跑水稻-玉米、水稻-高粱、水稻-短柄草的共线性分析然后把结果合在一张圆图里。不过要提前提醒一句多物种圆图对内存和显存的要求会明显上升而且连线太多会造成视觉混乱。建议先把每个物种的染色体内区段划分好只显示你关注的区域而不是把所有染色体全画出来否则图面会被密密麻麻的连线占满。6. 常见问题与排错指南6.1 高频报错速查表实际操作中大家遇到的问题其实高度集中我整理了一份速查表方便你直接对照。现象可能原因解决办法选择fasta文件后提示无法解析文件不是标准fasta格式或存在空行用文本编辑器重新保存为UTF-8无BOM格式BLASTP运行后没有任何输出E-value设置过严或两物种序列中没有相似片段放宽E-value到1e-3测试MCScanX运行中直接崩溃GFF文件中染色体ID与fasta不一致检查并统一ID结果中有很多区块但基本是乱的输入序列并非蛋白序列而是CDS确认fasta序列由氨基酸字母组成插件显示Java OutOfMemory数据量太大或内存设置不足调整TBtools启动脚本的-Xmx参数圆图生成后大量区域空白染色体长度设置错误或gff文件缺少部分染色体检查.gff中染色体数量和长度6.2 实战踩坑记我第一次做拟南芥和水稻共线性分析的时候踩过一个让我印象极深的坑。当时我下载的拟南芥注释文件是Araport11版本蛋白序列是TAIR10版本两者ID体系基本兼容但有一小部分基因在Araport11中被合并或拆分导致BLAST比对结果本来很漂亮MCScanX却一直在某个位置崩溃。后来我同时换成TAIR10版本的注释和蛋白序列问题立刻消失。所以我要强调下载数据时最好在同一来源、同一页面中同时下载同类文件不要混搭不同版本的注释和蛋白序列。版本不一致是很多诡异报错的根源。另一个坑是水稻染色体的命名。部分水稻GFF文件把不同亚基因组定义为chr01到chr12但蛋白fasta里却是chr1、chr2这种没补零的写法。如果没统一MCScanX会把chr01和chr1当成两个不同的染色体导致共线性区块无法正确关联结果里可能出现大量染色体级错配。遇到这种情况可以用TBtools的序列工具批量替换把补零的版本统一成不带零的形式。6.3 可选如何从TBtools过渡到命令行版本虽然这篇文章的主题是告别命令行但如果你想进一步做大项目或者批量处理了解命令行版本的MCScanX仍然有价值。TBtools插件能帮你完成大部分常规分析但当数据量上升到几十个物种时图形界面的操作效率会下降。从TBtools过渡到命令行推荐路径是先在TBtools里查看生成的日志复制出当时实际调用的BLAST和MCScanX命令。TBtools的日志窗口中会显示底层命令的详细内容包括参数和路径。你只需要在终端里手动运行这些命令就能复现同样的分析过程并且可以在此基础上去写循环脚本。在这个阶段我建议系统地学习一下文件格式处理和Shell脚本基础会大幅提升你的分析效率。比如批量给多个物种运行BLAST、统一文件命名、整合结果这些需求如果靠点击图形界面逐次操作会非常费时间。我个人学习命令行版本的经验是不要一开始就背命令而是先从TBtools生成的实际命令开始理解每条命令的作用然后慢慢替换自己的数据。这样既不会在初期被命令行吓退又能逐步掌握自动化分析的能力。最后的实操心得跑过几次共线性分析之后我最大的体会是工具的发展真的把很多高不可攀的分析拉到了普通实验室都能操作的层面。TBtools的MCScanX插件不能说多么完美但它至少让一堆不会写脚本的研究者也能在短时间内拿到可靠结果这对科研效率的提升是实打实的。最后再分享一个小技巧分析完成后一定把运行日志保存一份连同输入文件的版本信息、修改日期一起记录在项目文件夹里。论文写方法部分时需要准确写出数据版本和分析参数没有这些记录你会非常被动。我自己就因为当时偷懒没留日志后来补材料时不得不重新跑了一遍分析。省那两分钟往往后面要花两个小时来还。
阅读完成 · 觉得有帮助?
咨询建站