蛋白工厂（rRNA）和物流系统（tRNA）的识别鉴定-Barrnap，Infernal和tRNA-scan-SE
文章目录
- rRNA和tRNA：蛋白工厂和物流系统
- tRNAscan-SE-鉴定tRNA
- trnascan-se安装
- usage
- Barrnap-鉴定rRNA
- installtion
- Usage
- 选项
- 与 RNAmmer 的比较
- Infernal-包括rRNA在内的多种RNA鉴定
- Infernal installtion and database download
- 使用
- 输出
- 提取序列
- Reference
rRNA和tRNA：蛋白工厂和物流系统
rRNA和tRNA是生命分子中两种最基本的构件，前者是蛋白质加工的工厂，后者是负责运输原料（氨基酸）的物流系统，二者共同作用完成蛋白质的合成。在原核生物中，上述两种基因与其他单拷贝基因不同，往往具备多个拷贝基因。
rRNA在细菌中一般1-15个，在古菌中一般1-4 个拷贝rrnDB；tRNA至少18/20个，且与rRNA的数量呈现正相关。这个很好理解，蛋白加工工厂越多，负责提供货源的卡车自然也越多。
tRNA拷贝数与rRNA拷贝数呈现正相关关系rrnDB
16S rRNA在基因组中的拷贝数情况，多数在1-8个拷贝
如何从原核生物基因组序列中识别和鉴定上述元件呢，这里介绍2个最为常用的软件，tRNAscan-SE和Infernal。
tRNAscan-SE-鉴定tRNA
tRNAscan-SE是识别鉴定tRNA最经典的软件，如果待鉴定基因组比较少，可以直接使用tRNAscan-SE在线版，上传序列后，指定相应的序列源（真核？细菌？古菌？线粒体。。。）和搜索模式运行即可，具体不表。如果基因组比较多，则需要安装本地版，批量运行搜索。
trnascan-se安装
# 1.conda install
conda create -n trna -c bioconda trnascan-se=2
# or 2.install from source code
wget http://trna.ucsc.edu/software/trnascan-se-2.0.7.tar.gz
tar -zxvf trnascan-se-2.0.7.tar.gz
cd tRNAscan-SE-2.0/
./configure && sudo make install
# 加入环境变量
usage
- 基本用法
[Neptuneyut]$ tRNAscan-SEtRNAscan-SE 2.0.9 (July 2021)FATAL: No sequence file(s) specified.Usage: tRNAscan-SE [-options] <FASTA file(s)>Scan a sequence file for tRNAs-- default: use Infernal & tRNA covariance modelswith eukaryotic sequences(use -B, -A, -M, -O or -G to scan other types of sequences)Basic Options-E : search for eukaryotic tRNAs (default)，真核-B : search for bacterial tRNAs，细菌-A : search for archaeal tRNAs，古菌-M <model> : search for mitochondrial tRNAs，线粒体options: mammal, vert，哺乳动物或脊椎动物-O : search for other organellar tRNAs，其他细胞器tRNA-G : use general tRNA model (cytoslic tRNAs from all 3 domains included)，三域生物的细胞质tRNA-L : search using the legacy method (tRNAscan, EufindtRNA, and COVE)，传统搜索模式use with -E, -B, -A, -O, or -G-I : search using Infernal (default)，Infernal 搜索模式use with -E, -B, -A, -O, or -G-o <file> : save final results in <file>-f <file> : save tRNA secondary structures to <file>，输出tRNA二级结构-m <file> : save statistics summary for run in <file>(speed, # tRNAs found in each part of search, etc)-H : show both primary and secondary structure components tocovariance model bit scores-q : quiet mode (credits & run option selections suppressed)-h : print full list (long) of available options- 搜索大肠杆菌（Ecoli str.K-12 substr. MG1655）示例：
首先，我们在NCBI中查到Ecoli str.K-12 substr. MG1655的基因组信息如下：基因组4.64M，4285个蛋白，22个核糖体以及86个tRNA
下载Ecoli str.K-12 substr. MG1655基因组序列，使用tRNA对大肠杆菌搜索，速度非常块（不到30s）
(trna) [u@h@Ecoli_k12]$ time tRNAscan-SE -B -I -o Ecoli_k12_trnascan_output.txt -f Ecoli_k12_trnascan_secondary_structure.txt -m Ecoli_k12_trnascan_summary.txt Ecoli_k12.fastatRNAscan-SE v.2.0.9 (July 2021) - scan sequences for transfer RNAs
Copyright (C) 2020 Patricia Chan and Todd LoweUniversity of California Santa Cruz
Freely distributed under the GNU General Public License (GPLv3)------------------------------------------------------------
Sequence file(s) to search: Ecoli_k12.fasta
Search Mode: Bacterial
Results written to: Ecoli_k12_trnascan_output.txt
Output format: Tabular
Searching with: Infernal First Pass->Infernal
Isotype-specific model scan: Yes
Covariance model: /home/chenyy/anaconda3/envs/trna/lib/tRNAscan-SE/models/TRNAinf-bact.cm/home/chenyy/anaconda3/envs/trna/lib/tRNAscan-SE/models/TRNAinf-bact-SeC.cm
Infernal first pass cutoff score: 10Temporary directory: /tmp
tRNA secondary structurepredictions saved to: Ecoli_k12_trnascan_secondary_structure.txt
Search statistics saved in: Ecoli_k12_trnascan_summary.txt
------------------------------------------------------------Status: Phase I: Searching for tRNAs with HMM-enabled Infernal
Status: Phase II: Infernal verification of candidate tRNAs detected with first-pass scanreal 0m27.328s
user 0m54.958s
sys 1m20.211s
- 结果
- summary
(trna) [u@h@Ecoli_k12]$ cat Ecoli_k12_trnascan_summary.txttRNAscan-SE v.2.0.9 (July 2021) scan results (on host super-AS-2023US-TR4)
Started: 2021年 10月 30日 星期六 17:35:15 CST
命令行参数及模型使用情况：
------------------------------------------------------------
Sequence file(s) to search: Ecoli_k12.fasta
Search Mode: Bacterial
Results written to: Ecoli_k12_trnascan_output.txt
Output format: Tabular
Searching with: Infernal First Pass->Infernal
Isotype-specific model scan: Yes
Covariance model: /home/anaconda3/envs/trna/lib/tRNAscan-SE/models/TRNAinf-bact.cm/home/anaconda3/envs/trna/lib/tRNAscan-SE/models/TRNAinf-bact-SeC.cm
Infernal first pass cutoff score: 10Temporary directory: /tmp
tRNA secondary structurepredictions saved to: Ecoli_k12_trnascan_secondary_structure.txt
Search statistics saved in: Ecoli_k12_trnascan_summary.txt
------------------------------------------------------------First-pass Stats:第一轮搜索结果统计
---------------
Sequences read: 1，参考基因组序列数量
Seqs w/at least 1 hit: 1
Bases read: 4641652 (x2 for both strands)
Bases in tRNAs: 7219
tRNAs predicted: 91，预测tRNA数量
Av. tRNA length: 79，tRNA平均长度
Script CPU time: 0.22 s
Scan CPU time: 18.00 s
Scan speed: 515.7 Kbp/secFirst pass search(es) ended: 2021年 10月 30日 星期六 17:35:17 CSTInfernal Stats:第二轮Infernal搜索统计
-----------
Candidate tRNAs read: 91
Infernal-confirmed tRNAs: 89，Infernal确认tRNA数量
Bases scanned by Infernal: 9039
% seq scanned by Infernal: 0.1 %
Script CPU time: 0.20 s
Infernal CPU time: 36.47 s
Scan speed: 247.8 bp/secInfernal analysis of tRNAs ended: 2021年 10月 30日 星期六 17:35:42 CSTOverall scan speed: 169125.6 bp/sectRNAs decoding Standard 20 AA: 86，编码20种基本氨基酸的tRNA数量，我们看到和NCBI结果一致
Selenocysteine tRNAs (TCA): 1，硒代半胱氨酸1个
Possible suppressor tRNAs (CTA,TTA,TCA): 0
tRNAs with undetermined/unknown isotypes: 1，不确定/未知的转运rna
Predicted pseudogenes: 1，假基因-------
Total tRNAs: 89tRNAs with introns: 0|Isotype / Anticodon Counts:tRNA密码子使用情况，后续可以用于密码使用偏好性研究Ala : 5 AGC: GGC: 2 CGC: TGC: 3 丙氨酸有5个tRNA,2个是GGC编码，3个是TGC
Gly : 6 ACC: GCC: 4 CCC: 1 TCC: 1
Pro : 3 AGG: GGG: 1 CGG: 1 TGG: 1
Thr : 5 AGT: GGT: 2 CGT: 2 TGT: 1
Val : 7 AAC: GAC: 2 CAC: TAC: 5
Ser : 5 AGA: GGA: 2 CGA: 1 TGA: 1 ACT: GCT: 1
Arg : 7 ACG: 4 GCG: CCG: 1 TCG: CCT: 1 TCT: 1
Leu : 8 AAG: GAG: 1 CAG: 4 TAG: 1 CAA: 1 TAA: 1
Phe : 2 AAA: GAA: 2
Asn : 4 ATT: GTT: 4
Lys : 6 CTT: TTT: 6
Asp : 3 ATC: GTC: 3
Glu : 4 CTC: TTC: 4
His : 1 ATG: GTG: 1
Gln : 4 CTG: 2 TTG: 2
Ile : 5 AAT: GAT: 3 CAT: 2 TAT:
Met/fMet: 6 CAT: 6
Tyr : 3 ATA: GTA: 3
Supres : 0 CTA: TTA: TCA:
Cys : 1 ACA: GCA: 1
Trp : 1 CCA: 1
SelCys : 1 TCA: 1
- output:预测tRNA基本信息
(trna) [u@h@Ecoli_k12]$ cat Ecoli_k12_trnascan_output.txt
Sequence tRNA Bounds tRNA Anti Intron Bounds Inf
Name tRNA # Begin End Type Codon Begin End Score Note
-------- ------ ----- ------ ---- ----- ----- ---- ------ ------
NC_000913.3 1 225381 225457 Ile GAT 0 0 75.8
NC_000913.3 2 225500 225575 Ala TGC 0 0 78.1
NC_000913.3 3 228928 229004 Asp GTC 0 0 85.5
NC_000913.3 4 236931 237007 Asp GTC 0 0 85.5
NC_000913.3 5 262871 262946 Thr CGT 0 0 85.1
NC_000913.3 6 297184 297254 Thr CGT 0 0 44.2
NC_000913.3 7 564723 564799 Arg TCT 0 0 86.6
NC_000913.3 8 586007 586101 Undet NNN 0 0 31.7 pseudo，得分较低，假基因
NC_000913.3 9 780554 780629 Lys TTT 0 0 94.9
NC_000913.3 10 780765 780840 Val TAC 0 0 84.6
NC_000913.3 11 780843 780918 Lys TTT 0 0 94.9
NC_000913.3 12 781068 781143 Val TAC 0 0 84.6
NC_000913.3 13 781147 781222 Lys TTT 0 0 94.9
NC_000913.3 14 781369 781444 Lys TTT 0 0 94.9
NC_000913.3 15 781577 781652 Lys TTT 0 0 94.9
NC_000913.3 16 1746435 1746511 Val GAC 0 0 82.5
NC_000913.3 17 1746516 1746592 Val GAC 0 0 85.9
NC_000913.3 18 2044549 2044624 Asn GTT 0 0 79.7
NC_000913.3 19 2059851 2059926 Asn GTT 0 0 79.7
NC_000913.3 20 2062260 2062335 Asn GTT 0 0 79.7
NC_000913.3 21 2286211 2286287 Pro GGG 0 0 67.3
NC_000913.3 22 2466309 2466383 Arg CCT 0 0 62.5
NC_000913.3 23 2520931 2521006 Val TAC 0 0 84.6
NC_000913.3 24 2521051 2521126 Val TAC 0 0 84.6
NC_000913.3 25 2521173 2521248 Val TAC 0 0 84.6
NC_000913.3 26 2521253 2521328 Lys TTT 0 0 94.9
NC_000913.3 27 2755859 2755955 Undet NNN 0 0 29.6 得分较低，未知类型
NC_000913.3 28 2947387 2947463 fMet CAT 0 0 77.5
NC_000913.3 29 2947497 2947573 fMet CAT 0 0 77.5
NC_000913.3 30 2947607 2947683 fMet CAT 0 0 77.5
NC_000913.3 31 3110366 3110441 Phe GAA 0 0 75.3
NC_000913.3 32 3215598 3215673 Ile2 CAT 0 0 85.1
NC_000913.3 33 3836222 3836316 SeC TCA 0 0 142.5
NC_000913.3 34 3943435 3943510 Glu TTC 0 0 66.6
NC_000913.3 35 3946872 3946948 Asp GTC 0 0 85.5
NC_000913.3 36 3946957 3947032 Trp CCA 0 0 79.9
NC_000913.3 37 3982375 3982451 Arg CCG 0 0 82.2
NC_000913.3 38 3982510 3982585 His GTG 0 0 78.7
NC_000913.3 39 3982606 3982692 Leu CAG 0 0 67.1
NC_000913.3 40 3982735 3982811 Pro TGG 0 0 82.4
NC_000913.3 41 4037141 4037217 Ile GAT 0 0 75.8
NC_000913.3 42 4037260 4037335 Ala TGC 0 0 78.1
NC_000913.3 43 4168372 4168447 Glu TTC 0 0 66.6
NC_000913.3 44 4175388 4175463 Thr TGT 0 0 81.3
NC_000913.3 45 4175472 4175556 Tyr GTA 0 0 73.0
NC_000913.3 46 4175673 4175747 Gly TCC 0 0 70.6
NC_000913.3 47 4175754 4175829 Thr GGT 0 0 84.7
NC_000913.3 48 4209774 4209849 Glu TTC 0 0 66.6
NC_000913.3 49 4392360 4392435 Gly GCC 0 0 88.4
NC_000913.3 50 4392472 4392547 Gly GCC 0 0 88.4
NC_000913.3 51 4392583 4392658 Gly GCC 0 0 88.4
NC_000913.3 52 4496405 4496489 Leu CAA 0 0 79.3
NC_000913.3 53 4606401 4606315 Leu CAG 0 0 67.1
NC_000913.3 54 4606286 4606200 Leu CAG 0 0 65.0
NC_000913.3 55 4606165 4606079 Leu CAG 0 0 67.1
NC_000913.3 56 4362626 4362551 Phe GAA 0 0 75.3
NC_000913.3 57 3708692 3708616 Pro CGG 0 0 76.0
NC_000913.3 58 3427152 3427076 Ile GAT 0 0 75.8
NC_000913.3 59 3427033 3426958 Ala TGC 0 0 78.1
NC_000913.3 60 3423655 3423580 Thr GGT 0 0 79.8
NC_000913.3 61 3322158 3322072 Leu GAG 0 0 64.5
NC_000913.3 62 3318289 3318213 fMet CAT 0 0 75.5
NC_000913.3 63 2999057 2998984 Gly CCC 0 0 81.1
NC_000913.3 64 2818645 2818553 Ser GCT 0 0 83.4
NC_000913.3 65 2818549 2818473 Arg ACG 0 0 85.1
NC_000913.3 66 2818274 2818198 Arg ACG 0 0 85.1
NC_000913.3 67 2818135 2818059 Arg ACG 0 0 85.1
NC_000913.3 68 2817860 2817784 Arg ACG 0 0 85.1
NC_000913.3 69 2785837 2785762 Ile2 CAT 0 0 84.7
NC_000913.3 70 2729444 2729369 Glu TTC 0 0 66.6
NC_000913.3 71 2518231 2518156 Ala GGC 0 0 75.0
NC_000913.3 72 2518116 2518041 Ala GGC 0 0 75.0
NC_000913.3 73 2058102 2058027 Asn GTT 0 0 79.7
NC_000913.3 74 2043557 2043468 Ser CGA 0 0 90.4
NC_000913.3 75 1992117 1992042 Gly GCC 0 0 88.4
NC_000913.3 76 1991987 1991914 Cys GCA 0 0 63.0
NC_000913.3 77 1991901 1991815 Leu TAA 0 0 71.7
NC_000913.3 78 1287622 1287538 Tyr GTA 0 0 72.2
NC_000913.3 79 1287328 1287244 Tyr GTA 0 0 72.2
NC_000913.3 80 1097652 1097565 Ser GGA 0 0 76.9
NC_000913.3 81 1031712 1031625 Ser TGA 0 0 83.5
NC_000913.3 82 925971 925884 Ser GGA 0 0 76.9
NC_000913.3 83 697133 697057 Met CAT 0 0 84.4
NC_000913.3 84 697047 696963 Leu TAG 0 0 74.1
NC_000913.3 85 696939 696865 Gln TTG 0 0 70.4
NC_000913.3 86 696830 696756 Gln TTG 0 0 70.4
NC_000913.3 87 696740 696664 Met CAT 0 0 84.4
NC_000913.3 88 696616 696542 Gln CTG 0 0 70.8
NC_000913.3 89 696504 696430 Gln CTG 0 0 70.8
- tRNA二级结构
(trna) [u@h@Ecoli_k12]$ head Ecoli_k12_trnascan_secondary_structure.txt
NC_000913.3.trna1 (225381-225457) Length: 77 bp
Type: Ile Anticodon: GAT at 35-37 (225415-225417) Score: 75.8* | * | * | * | * | * | * | *
Seq: AGGCTTGTAGCTCAGGTGGTtAGAGCGCACCCCTGATAAGGGTGAGGtCGGTGGTTCAAGTCCACTCAGGCCTACCA
Str: >>>>>>>..>>>>.........<<<<.>>>>>.......<<<<<.....>>>>>.......<<<<<<<<<<<<....NC_000913.3.trna2 (225500-225575) Length: 76 bp
Type: Ala Anticodon: TGC at 34-36 (225533-225535) Score: 78.1* | * | * | * | * | * | * | *
Seq: GGGGCTATAGCTCAGCTGGGAGAGCGCCTGCTTTGCACGCAGGAGGtCTGCGGTTCGATCCCGCATAGCTCCACCA
Barrnap-鉴定rRNA
Barrnap(BAsic Rapid Ribosomal RNA Predictor)，Barrnap预测了基因组中核糖体RNA基因的位置。它支持细菌（5S,23S,16S）、古细菌（5S,5.8S,23S,16S）、多细胞动物线粒体（12S,16S）和真核生物（5S,5.8S,28S,18S）。
它将FASTA DNA序列作为输入，并将GFF3作为输出。它使用HMMER 3.1中的新工具nhmmer，用于RNA：DNA类型的HMM搜索。支持多线程，可以期望在更多的CPU上实现大致的线性加速。
installtion
conda install -c bioconda -c conda-forge barrnap
Usage
$ barrnap --quiet examples/small.fna
##gff-version 3
P.marinus barrnap:0.8 rRNA 353314 354793 0 + . Name=16S_rRNA;product=16S ribosomal RNA
P.marinus barrnap:0.8 rRNA 355464 358334 0 + . Name=23S_rRNA;product=23S ribosomal RNA
P.marinus barrnap:0.8 rRNA 358433 358536 7.5e-07 + . Name=5S_rRNA;product=5S ribosomal RNA$ barrnap -q -k mito examples/mitochondria.fna
##gff-version 3
AF346967.1 barrnap:0.8 rRNA 643 1610 . + . Name=12S_rRNA;product=12S ribosomal RNA
AF346967.1 barrnap:0.8 rRNA 1672 3228 . + . Name=16S_rRNA;product=16S ribosomal RNA$ barrnap -o rrna.fa < contigs.fa > rrna.gff
$ head -n 3 rrna.fa
>16S_rRNA::gi|329138943|tpg|BK006945.2|:455935-456864(-)
ACGGTCGGGGGCATCAGTATTCAATTGTCAGAGGTGAAATTCTTGGATT
TATTGAAGACTAACTACTGCGAAAGCATTTGCCAAGGACGTTTTCATTA
选项
General
--help show help and exit
--version print version in form barrnap X.Y and exit
--citation print a citation and exit
Search
--kingdom is the database to use: Bacteria:bac, Archaea:arc, Eukaryota:euk, Metazoan Mitochondria:mito
--threads is how many CPUs to assign to nhmmer search
--evalue is the cut-off for nhmmer reporting, before further scrutiny
--lencutoff is the proportion of the full length that qualifies as partial match
--reject will not include hits below this proportion of the expected length
Output
--quiet will not print any messages to stderr
--incseq will include the full input sequences in the output GFF
--outseq creates a FASTA file with the hit sequences
与 RNAmmer 的比较
Barrnap被设计为RNAmmer的替代品。它的动机是我想消除Prokka对RNAmmer的依赖，因为RNAmmer受制于免费的学术注册许可，而且RNAmmer对传统的HMMER 2.x的依赖与大多数人现在使用的HMMER 3.x相冲突。
RNAmmer比Barrnap更复杂，也更准确，因为它在glocal对齐模式下使用HMMER 2.x，而NHMMER 3.x目前只支持局部对齐（Sean Eddy预计glocal在2014年就会被支持，但在2018年仍然没有支持）。
在实践中，Barrnap会在几秒钟内找到所有典型的rRNA基因（在细菌中），但可能会把端点弄错几个碱基，并可能会错过奇怪的rRNAs。它使用的HMM模型来自Rfam、Silva和RefSeq。
Infernal-包括rRNA在内的多种RNA鉴定
Infernal(INFERence of RNA ALignment)是用于搜索DNA序列数据库中的RNA结构和序列相似性。它是一种基于被称为协方差概率模型(CMs)，这是随机上下文自由语法（SCFG）的一种专门类型。CM就像一个序列剖面，并有一个特定位置的评分系统，用于碱基替换、修补、插入和删除。但它的得分是序列一致性和RNA二级结构一致性的组合，所以在许多情况下，它更有能力识别保守二级结构的RNA同源物，而不是保守一级序列。
与其他仅基于序列比对和数据库搜索工具相比，Infernal明显更准确，更能够检测到远程同源物，因为它对序列和结构进行建模。但是建立结构模型需要很高的计算成本，所以CM同源搜索的速度很慢（以前的版本）。在1.1版本中，典型的同源搜索速度是以前的100倍。这要归功于来自HMMER3软件包的加速HMM方法的加入。
Infernal主要包括以下常用程序：
Infernal installtion and database download
1.安装Infernal软件套装
sudo apt-get install infernal infernal-doc
# or conda
conda install -c bioconda infernal=1.1.2
2.下载Rfam数据库
Rfam数据库是一个RNA家族的集合，每个家族由一个CM和用于建立该CM的多序列排列代表。截至2020年12月，目前的版本是14.4，其中包括3446个家族。 Rfam网站允许使用cmscan对Rfam CM数据库进行基于网络的搜索，而用户可以上传查询序列。 另外，你也可以通过在本地运行cmscan来进行同样的搜索，如本例所示。 通过用你的序列数据集搜索所有的Rfam，你将用一个命令为你的数据集注释大多数已知类型的结构性RNAs。
# download Rfam 14.4
wget ftp://ftp.ebi.ac.uk/pub/databases/Rfam/14.4/Rfam.cm.gz
gunzip Rfam.cm.gz
# Then, as in the previous example, you’ll need to runcmpresson this CM database:
cmpress Rfam.cm
# CM database: /mnt/nfs/database/Rfam/Rfam.cm
使用
# 下一步是运行cmscan。 为了重现Rfam搜索的执行方式（Nawrockiet al., 2015），需要几个命令行选项。下面将对每个选项进行解释。完整的命令是（分成两行，以便在页面上显示）。
# quick start
time cmscan --rfam --cpu 100 --tblout ecoli.tblout --fmt 2 /mnt/nfs/database/Rfam/Rfam.cm Ecoli_k12/Ecoli_k12.fasta > ecoli.cmscan
# full usage
cmscan --rfam --cut_ga --nohmmonly --tblout mrum-genome.tblout --fmt 2 \
--clanin tutorial/Rfam.14.3.clanin Rfam.cm tutorial/mrum-genome.fa --cpu 10 > mrum-genome.cmscan
大肠杆菌100核用时2min,ecoli.tblout第二列记录着tRNA/rRNA/gene type,可以根据该列过滤，ecoli.cmscan是所有信息
参数意义：
--cpu <n> ：用于多线程的并行CPU数量
--rfam： 指定过滤管道在快速模式下运行，使用与Rfamsearches和其他大于20Gb的序列数据库相同的严格过滤器（见第4节）。
--cut_ga：指定使用特殊的RfamGA(收集)阈值来确定报告哪些命中。这些阈值存储在thermam .cmfile中。每个模型都有自己的GA比特评分阈值，这个阈值是由Rfam管理器确定的，当比特评分达到或超过该值时，所有命中都被认为与模型是真实的同构。这些决定是基于对Rfam使用的大型Rfamseq数据库的观察结果做出的(Nawrocki等人，2015)
--nohmmonly：所有的模型，甚至那些零基对的模型，都是在CM模式下运行的（而不是HMM模式）。这就保证了所有的GAcutoffs（在CM模式下为每个模型确定的）都是有效的。
--tblout：指定应该创建一个表格式的输出文件，见第6节。
--fmt 2：表格输出文件将采用格式2，其中包括对重叠命中的注释。 见第61页对该格式的完整描述。
--clanin： Clan information这个文件列出了哪些模型属于同一个族。 族是同源的模型组，因此预计这些模型的一些命中率会重叠。 例如，LSUrRNAarchaea和LSUrRNAbacteria模型都在同一个族中。
当cmscan命令运行完毕后，文件mrum-genome.cmscan将包含程序的标准输出。这个文件将与我们在前面的cmscan例子中看到的类似。**文件mrum-genome.tblout也已经被创建，它是所有命中的表格，每一个命中有一行表示。看一下这个文件。前两行是注释行（前缀为#字符），文件中27列数据的标签。后面的每一行都有27个以空格分隔的标记。这些标记的具体含义将在第6节中详细描述。下面我将包括文件的前24行，但删除了第3-5、7-9和13-16列（用…代替），这样文本就可以放在这一页。
这个表格格式包括目标模型名称、序列名称（在第3栏，上面省略了，以节省空间）、宗族名称（clan name）、序列坐标、位点得分、E值等等。 因为使用了–fmt 2选项，这个文件包括了哪些命中与其他命中重叠的信息，从标有 "OLP "的那一列开始，以 "wfrct2 "结束。 在 "OLP "列中带有 “"字符的命题不与任何其他命题重叠。 那些带有 "ˆ “的命题确实与其他至少一个命题重叠，但这些重叠的命题没有一个有更好的分数（出现在列表的更高位置）。那些带有”="的命题也与至少一个其他的命题重叠，这些命题的索引在 "anyidx "列中给出。关于这些栏目的更详细的解释，请参见第61页。 这是在Methanobrevibacter ruminantium基因组中的LSU rRNA的两个拷贝。 3号和4号是LSUrRNAbacteri模型，与1号和2号几乎完全重叠（1号是从762872到765862的序列位置，3号是从762874到765862的序列位置）。这种重叠并不奇怪，因为细菌和古细菌的LSU rRNA模型非常相似，所以对相同的子序列给予了高分。此外，命中5是LSUrRNAeukarya的，也与命中1和3相重叠。因为这三个LSU模型由于其同源性，预计都会产生重叠的命中，所以Rfam将它们归入同一族，注意所有三个命中的 "族名 "列中的 "CL00112 "值。 这个氏族信息是在我们提供给ocmscan的fam.14.3.claninfoinput文件中，通过使用–claninoption提供的。"lp "列表明，命题1是三个重叠命题中得分最高的，因为它含有 "ˆ "字符。 命中3和5在 “lp “栏中都有”=”，表明有另一个与这些命中重叠的模型的命中，并且有更好的得分。如果你使用这些结果为Methanobrevibacter ruminantium基因组产生注释，你可能想忽略任何有更高得分重叠的命中。 要做到这一点，你只需删除 "OLP “栏中带有”="的任何命中。 另外，你也可以通过在ocmscan中增加 "oskip "选项，使这些命题不被打印到表格输出文件中。您也可以通过–oclanoption来修改重叠标注行为，该选项将重叠标注限制在同一族内模型的命中。 不在同一族内的模型的重叠点击不会被标记为重叠，而是在 "OLP "字段中标记为 "”。
输出
(gtdbtk) [u@h@rRNA_and_tRNA_out]$ head -n40 ecoli.cmscan
# cmscan :: search sequence(s) against a CM database
# INFERNAL 1.1.3 (Nov 2019)
# Copyright (C) 2019 Howard Hughes Medical Institute.
# Freely distributed under the BSD open source license.
# - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
# query sequence file: ../Genomes/ISO_Y784-1.fa
# target CM database: /mnt/nfs/database/Rfam/Rfam.cm
# tabular output of hits: ecoli.tblout
# tabular output format: 2
# Rfam pipeline mode: on [strict filtering]
# number of worker threads: 100 [--cpu]
# - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -Query: 4 [L=379216]
Hit scores:rank E-value score bias modelname start end mdl trunc gc description---- --------- ------ ----- ------------ ------ ------ --- ----- ---- -----------(1) ! 7.7e-31 141.0 0.0 T-box 103732 103948 + cm no 0.41 T-box leader(2) ! 1.4e-30 140.0 0.0 T-box 54103 53866 - cm no 0.44 T-box leader(3) ! 6.4e-28 129.7 0.0 T-box 53823 53579 - cm no 0.47 T-box leader(4) ! 3.3e-21 96.3 0.0 SAM 287143 287037 - cm no 0.46 SAM riboswitch (S box leader)(5) ! 1.6e-17 80.9 13.5 cspA 147215 146861 - cm no 0.34 cspA thermoregulator(6) ! 2.5e-14 97.0 0.0 RsaE 85440 85320 - cm no 0.36 RNA Staph. aureus E (RoxS)(7) ! 1.5e-08 55.0 0.0 tRNA 20433 20505 + cm no 0.64 tRNA------ inclusion threshold ------(8) ? 0.14 31.5 0.0 MIR394 222554 222650 + cm no 0.49 microRNA MIR394(9) ? 0.21 33.8 0.0 BsrC 264593 264673 + cm no 0.41 BsrC(10) ? 0.24 28.6 0.0 mir-14 167868 167813 - cm no 0.34 microRNA mir-14(11) ? 0.25 34.8 0.8 Phe_leader 146681 146822 + cm no 0.31 Phenylalanine leader peptide提取序列
Because the programs for CM fetching (cmfetch) and sequence fetching (esl-sfetch) can fetch any number of profiles or (sub)sequences by names/accessions given in a list,and these programs can also read these lists from astdin pipe, you can craft incantations that generate subsets of queries or targets on the fly. For example, you can extractand align all hits found by cmsearch with an E-value below the inclusion threshold as follows (using the\character twice below to split up the final command onto three lines):
/home/chenyy/anaconda3/envs/kofam/bin/esl-sfetch --index ISO_Y438.fa # esl-sfetch构建index
grep -v ˆ# 2782Genomes_CMSCAN_for_tRNA_and_rRNA/ISO_Y438.tblout| grep ! | awk '{ printf("%s/%s-%s %s %s %s\n", $1, $8, $9, $8, $9, $1); }'|/home/chenyy/anaconda3/envs/kofam/bin/esl-sfetch -Cf links/ISO_Y438.fa # 官方示例
cmsearch --tblout tRNA5.mrum-genome.tbl tRNA5.cm mrum-genome.faesl-sfetch --index mrum-genome.facat tRNA5.mrum-genome.tbl | grep -v ˆ# | grep ! \| awk ’{ printf(‘‘%s/%s-%s %s %s %s\n’’, $1, $8, $9, $8, $9, $1); }’ \| esl-sfetch -Cf mrum-genome.fa - | ../src/cmalign tRNA5.cm -
#In retrieving subsequences listed in a file (-C -f, or just -Cf), each line of the file
# is in GDF format: <newname> <from> <to> <source seqname>, space/tab delimited.# 批量提取核糖体RNA
for type in LSU_rRNA_bacteria LSU_rRNA_archaea SSU_rRNA_bacteria SSU_rRNA_archaea 5_8S_rRNA
dofor file in *tbloutdo grep -v "ˆ#" $file| grep ! | awk '{ printf("%s %s %s %s\n", $2, $8, $9, $1); }'|\/home/chenyy/anaconda3/envs/kofam/bin/esl-sfetch -Cf \/mnt/nfs/yutao/ISOs_and_MAGs/2782rRNA_and_tRNA_out/links/${file%.tblout}.fa -|\bioawk -c fastx -v type=$type -v file=${file%.tblout} '$name==type{print ">"file"_"$name"\n"$seq}'done > /mnt/nfs/yutao/ISOs_and_MAGs/2782rRNA_and_tRNA_out/Extract_RNA_sequence/${type}.fa
done第一个命令使用CM文件tRNA5.c.cm和序列文件mrum-genome.fa（这些在教程中使用）进行搜索，并将表格式输出保存为tRNA5.mrum-genome.tbl。 第二条命令对基因组文件进行索引，为快速（子）序列检索做准备。在第三条命令中，我们只从tRNA5.mrum-genome.tbl中提取了那些不以a#开头的线（这些是注释线），同时也包括a！（这些是E值低于包含阈值的命中），使用了前两条ogrep命令。 这个输出通过hawk重新格式化成 "GDF "格式，即atesl-sfetchexpects: 。然后这些行被输送到esl-sfetch（使用’-‘参数）以检索每个命中（仅包括每个命中的子序列，而不是完整的目标序列）。esl-sfetch将输出一个FASTA文件，最后被输送到ocmalign，同样使用’-'参数。实际打印到屏幕上的输出结果将是所有包含的tRNA命中的多重对齐。
Reference
tRNAscan-SE
Infernal
Infernal mannual
蛋白工厂（rRNA）和物流系统（tRNA）的识别鉴定-Barrnap，Infernal和tRNA-scan-SE相关推荐
- 【智能制造】索菲亚家居智能工厂与物流系统建设
家具制造业的智能化转型是必然趋势,索菲亚家居股份有限公司作为行业领军企业,经过了信息化管理.自动化与信息化齐头并进两大发展阶段,于近两年继续大力推进信息化智能化建设,加速向智能制造转型.索菲亚家居的智 ...
- 案例！智慧物流系统的规划思路
导语 大家好,我是智能仓储物流技术研习社的社长,你的老朋友,老K. 知识星球 * 原创电子书 * 深海社区 * 微信群 来源:意尔康 01 [智慧物流,驱动因素] 随着科技高速发展,智慧物流自动化仓储 ...
- 智能工厂物流系统总体规划导向与逻辑
导语 在价值链运营环境下,物流已经成为智能工厂的核心要素,工厂规划和运营管理必须要具备"流动思维"和"供应链交付思维"."大交付.大物流.小生产&qu ...
- 智能工厂物流系统规划步骤与关键要素
导读:智能工厂物流规划是一个系统规划的过程,需要遵循一定的规划步骤.智能工厂物流规划的一般步骤,包括需求梳理.概念设计.初步规划.详细规划.方案验证.实施落地等六个阶段. 消费者和客户需求拉动交付,从 ...
- [渝粤教育] 西南科技大学 物流系统规划与设计 在线考试复习资料
物流系统规划与设计--在线考试复习资料 一.单选题 1.TDC型物流中心的管理重点是( ). A.储位的管理 B.客户订单管理 C.产品货物管理 D.销售点管理 2.一般来讲,批量大.价值低.运距短的 ...
- agv系统介绍_AGV物流系统工作流程及模块介绍
原标题:AGV物流系统工作流程及模块介绍 众所周知,企业内部物流体系是-个复杂系统,很多企业都逐渐开始重视物流提高自己竞争力. 因此能降低零部件库存,降低周转箱数量,平衡物料接收,提高装货卸货效率的A ...
- 五步实现“中国制造”精益物流系统
作者:Smarttomato log96 2007-10-19 内容导航: 建立PFEP和零部件供应超市 第1页: 建立PFEP和零部件供应超市 第2页: 设计物料上线路线 完善需求信号 <s ...
- 考虑人机协同的智能工厂多AGV物流调度仿真研究
摘要 智能工厂的多AGV物流调度作业,对多AGV的调度优化有较高要求,存在上料种类多.频次高.干涉强和动线复杂等问题.鉴于此,从人机协同视角,对智能工厂多AGV物流调度进行仿真研究.分析仿真需求,在空 ...
- 蒙牛六期的高度自动化物流系统
蒙牛乳业集团成立于1999年,总部设在内蒙古呼和浩特市和林格尔县盛乐经济园区,企业总资产达76亿元,乳制品年生产能力500万吨.随着生产规模的不断扩大,自2002年起,蒙牛乳业集团开始采用自动化立体仓 ...
- 【智能物流】PPT干货，智能物流系统
卷烟厂物流自动化系统智能解决方案分析 精智工厂 物流自动化系统迅速发展起来的新型自动化立体仓库系统,集光.机.电子一体的系统工程,把物流.信息流用计算机和现代信息技术集成在一起,涉及激光导航.红外通讯 ...
最新文章
- Jedis无法远程连接阿里云服务器的redis问题
- linux内核线程创建销毁机制
- Kali安装之后必做20件事 第二版(2017-07-07不断更新)
- C++(一)——存储持续性、作用域、链接性
- printf格式化字符串用法
- Nginx+keepalived从入门到集群搭建（手把手教学，建议收藏）
- 卷积神经网络迁移学习
- 安装MySQL-python报错 error: command 'gcc' failed with exit status 1解决方法
- Selenium2Library+ride学习笔记
- 唱好铁血丹心谐音正规_长沙正规的音乐高考培训学校
- Android xml manifest属性详解
- phpstom怎样导出数据库_用phpStorm的数据库工具来管理你的数据库
- oracle通信通道的文件结尾_oracle里执行full join 报通信通道的文件结尾问题
- 百度蜘蛛ip地址大全，百度搜索引擎蜘蛛的IP地址段
- python垃圾分类程序_如何利用Python进行垃圾分类
- 【捞】明朝灭亡的经济原因
- 如何用安装启动盘启动计算机,用u盘启动电脑进入系统安装 如何进入启动u盘安装系统...
- STM32之DS1682
- 美国网络再次“瘫痪”，华为意外“出头”，网络服务器世界第一
- 免费下论文及查重投稿的10来个方法