macでインフォマティクス

macでインフォマティクス

HTS (NGS) 関連のインフォマティクス情報についてまとめています。

rRNAを除く SortMeRNA

2020 2/5 condaインストール追記

2020 6/16 コマンドが大きく変更したため更新

 

 次世代シーケンシング(NGS)技術を生物群集から直接抽出したRNAに適用すると、コーディング型とノンコーディング型の両方のRNAを特徴づける断片が混在している。これらを区別し、メッセンジャーRNAリボソームRNA(rRNA)のファミリーをさらに分類することは、相互作用環境の遺伝子発現パターンや構成種の系統分類を調べる上で重要なステップである。

 本研究では、メタトランスクリプトームデータからrRNA断片を迅速にフィルタリングするために設計された新しいソフトウェアSortMeRNAを紹介する。このソフトウェアは、大規模なリードセットを扱うことができ、rRNAデータベースに一致する全ての断片を高感度かつ低ランニングタイムでソートすることが可能である。

 

 

出力はfasta、fastq、アライメントのsam、またblastライクな出力も可能である。Illumina, 454, Ion Torrent and PacBioのシーケンスデータに対応している。QIIMEと一緒に使用することで、OTUを検出し系統解析にも利用することができる。

 

 マニュアル

 http://bioinfo.lifl.fr/RNA/sortmerna/code/SortMeRNA-user-manual-v2.1.pdf

FAQ

http://bioinfo.lifl.fr/sortmerna/faqs.php

 

インストール

Github

conda install -c bioconda -y sortmerna

> sortmerna -h

$ sortmerna -h

 

  Program:     SortMeRNA version 2.1b, 03/03/2016

  Copyright:   2012-16 Bonsai Bioinformatics Research Group:

               LIFL, University Lille 1, CNRS UMR 8022, INRIA Nord-Europe

               2014-16 Knight Lab:

               Department of Pediatrics, UCSD, La Jolla,

  Disclaimer:  SortMeRNA comes with ABSOLUTELY NO WARRANTY; without even the

               implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.

               See the GNU Lesser General Public License for more details.

  Contact:     Evguenia Kopylova, jenya.kopylov@gmail.com 

               Laurent Noé, laurent.noe@lifl.fr

               Hélène Touzet, helene.touzet@lifl.fr

 

 

  usage:   ./sortmerna --ref db.fasta,db.idx --reads file.fa --aligned base_name_output [OPTIONS]:

 

  -------------------------------------------------------------------------------------------------------------

  | parameter          value           description                                                    default |

  -------------------------------------------------------------------------------------------------------------

     --ref             STRING,STRING   FASTA reference file, index file                               mandatory

                                         (ex. --ref /path/to/file1.fasta,/path/to/index1)

                                         If passing multiple reference files, separate 

                                         them using the delimiter ':',

                                         (ex. --ref /path/to/file1.fasta,/path/to/index1:/path/to/file2.fasta,path/to/index2)

     --reads           STRING          FASTA/FASTQ reads file                                         mandatory

     --aligned         STRING          aligned reads filepath + base file name                        mandatory

                                         (appropriate extension will be added)

 

   [COMMON OPTIONS]: 

     --other           STRING          rejected reads filepath + base file name

                                         (appropriate extension will be added)

     --fastx           BOOL            output FASTA/FASTQ file                                        off

                                         (for aligned and/or rejected reads)

     --sam             BOOL            output SAM alignment                                           off

                                         (for aligned reads only)

     --SQ              BOOL            add SQ tags to the SAM file                                    off

     --blast           STRING          output alignments in various Blast-like formats                

                                        '0' - pairwise

                                        '1' - tabular (Blast -m 8 format)

                                        '1 cigar' - tabular + column for CIGAR 

                                        '1 cigar qcov' - tabular + columns for CIGAR

                                                         and query coverage

                                        '1 cigar qcov qstrand' - tabular + columns for CIGAR,

                                                                query coverage and strand

     --log             BOOL            output overall statistics                                      off

     --num_alignments  INT             report first INT alignments per read reaching E-value          -1

                                        (--num_alignments 0 signifies all alignments will be output)

       or (default)

     --best            INT             report INT best alignments per read reaching E-value           1

                                         by searching --min_lis INT candidate alignments

                                        (--best 0 signifies all candidate alignments will be searched)

     --min_lis         INT             search all alignments having the first INT longest LIS         2

                                         LIS stands for Longest Increasing Subsequence, it is 

                                         computed using seeds' positions to expand hits into

                                         longer matches prior to Smith-Waterman alignment. 

     --print_all_reads BOOL            output null alignment strings for non-aligned reads            off

                                         to SAM and/or BLAST tabular files

     --paired_in       BOOL            both paired-end reads go in --aligned fasta/q file             off

                                         (interleaved reads only, see Section 4.2.4 of User Manual)

     --paired_out      BOOL            both paired-end reads go in --other fasta/q file               off

                                         (interleaved reads only, see Section 4.2.4 of User Manual)

     --match           INT             SW score (positive integer) for a match                        2

     --mismatch        INT             SW penalty (negative integer) for a mismatch                   -3

     --gap_open        INT             SW penalty (positive integer) for introducing a gap            5

     --gap_ext         INT             SW penalty (positive integer) for extending a gap              2

     -N                INT             SW penalty for ambiguous letters (N's)                         scored as --mismatch

     -F                BOOL            search only the forward strand                                 off

     -R                BOOL            search only the reverse-complementary strand                   off

     -a                INT             number of threads to use                                       1

     -e                DOUBLE          E-value threshold                                              1

     -m                INT             INT Mbytes for loading the reads into memory                   1024

                                        (maximum -m INT is 16384)

     -v                BOOL            verbose                                                        off

 

 

   [OTU PICKING OPTIONS]: 

     --id              DOUBLE          %id similarity threshold (the alignment must                   0.97

                                         still pass the E-value threshold)

     --coverage        DOUBLE          %query coverage threshold (the alignment must                  0.97

                                         still pass the E-value threshold)

     --de_novo_otu     BOOL            FASTA/FASTQ file for reads matching database < %id             off

                                         (set using --id) and < %cov (set using --coverage) 

                                         (alignment must still pass the E-value threshold)

     --otu_map         BOOL            output OTU map (input to QIIME's make_otu_table.py)            off

 

 

   [ADVANCED OPTIONS] (see SortMeRNA user manual for more details): 

    --passes           INT,INT,INT     three intervals at which to place the seed on the read         L,L/2,3

                                         (L is the seed length set in ./indexdb_rna)

    --edges            INT             number (or percent if INT followed by % sign) of               4

                                         nucleotides to add to each edge of the read

                                         prior to SW local alignment 

    --num_seeds        INT             number of seeds matched before searching                       2

                                         for candidate LIS 

    --full_search      BOOL            search for all 0-error and 1-error seed                        off

                                         matches in the index rather than stopping

                                         after finding a 0-error match (<1% gain in

                                         sensitivity with up four-fold decrease in speed)

    --pid              BOOL            add pid to output file names                                   off

 

 

   [HELP]:

     -h                BOOL            help

     --version         BOOL            SortMeRNA version number

 

 

indexdb_rna -h

$ indexdb_rna -h

 

  Program:     SortMeRNA version 2.1b, 03/03/2016

  Copyright:   2012-16 Bonsai Bioinformatics Research Group:

               LIFL, University Lille 1, CNRS UMR 8022, INRIA Nord-Europe

               2014-16 Knight Lab:

               Department of Pediatrics, UCSD, La Jolla,

  Disclaimer:  SortMeRNA comes with ABSOLUTELY NO WARRANTY; without even the

               implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.

               See the GNU Lesser General Public License for more details.

  Contact:     Evguenia Kopylova, jenya.kopylov@gmail.com 

               Laurent Noé, laurent.noe@lifl.fr

               Hélène Touzet, helene.touzet@lifl.fr

 

 

  usage:   ./indexdb_rna --ref db.fasta,db.idx [OPTIONS]:

 

  --------------------------------------------------------------------------------------------------------

  | parameter        value           description                                                 default |

  --------------------------------------------------------------------------------------------------------

     --ref           STRING,STRING   FASTA reference file, index file                            mandatory

                                      (ex. --ref /path/to/file1.fasta,/path/to/index1)

                                       If passing multiple reference sequence files, separate

                                       them by ':',

                                      (ex. --ref /path/to/file1.fasta,/path/to/index1:/path/to/file2.fasta,path/to/index2)

   [OPTIONS]:

     --tmpdir        STRING          directory where to write temporary files

     -m              INT             the amount of memory (in Mbytes) for building the index     3072 

     -L              INT             seed length                                                 18

     --max_pos       INT             maximum number of positions to store for each unique L-mer  10000

                                      (setting --max_pos 0 will store all positions)

     -v              BOOL            verbose

     -h              BOOL            help

 

 

 

 データベース

git clone https://github.com/biocore/sortmerna.git

> ls -lh

$ ls -lh sortmerna/data/rRNA_databases/

total 149672

-rw-r--r--  1 kazu  staff   2.0K  6 16 17:25 README.txt

-rw-r--r--  1 kazu  staff   2.2M  6 16 17:25 rfam-5.8s-database-id98.fasta

-rw-r--r--  1 kazu  staff   8.1M  6 16 17:25 rfam-5s-database-id98.fasta

drwxr-xr-x  5 kazu  staff   160B  6 16 17:25 scripts/

-rw-r--r--  1 kazu  staff   3.7M  6 16 17:25 silva-arc-16s-id95.fasta

-rw-r--r--  1 kazu  staff   734K  6 16 17:25 silva-arc-23s-id98.fasta

-rw-r--r--  1 kazu  staff    19M  6 16 17:25 silva-bac-16s-id90.fasta

-rw-r--r--  1 kazu  staff    12M  6 16 17:25 silva-bac-23s-id98.fasta

-rw-r--r--  1 kazu  staff    13M  6 16 17:25 silva-euk-18s-id95.fasta

-rw-r--r--  1 kazu  staff    14M  6 16 17:25 silva-euk-28s-id98.fasta

-rw-r--r--  1 kazu  staff   591K  6 16 17:25 silva_ids_acc_tax.tar.gz

(データベースにはアーキアと真核生物のrRNA配列もある。SILVAからダウンロードするとarb ファイルで扱いにくいので、最新バージョンではないがこちらの配列セットを様々な処理に使うのもあり)

 

 

ラン

1、indexファイルの準備

解析にはデータベースのrRNA (FASTA) にindexをつける必要がある。ここではSILVAの16S rRNAのindexを作成している。

indexdb_rna --ref silva-arc-16s-id95.fasta,silva-arc-16s-id95.fasta-db
  • -v   verbose
  • --ref   STRING,STRING   FASTA reference file, index file

fastaとinexの間は","で区切る。複数のファイルを指定することもできる(マニュアル参照)。

 

 2、SortMeRNAの実行。

sortmerna --ref silva-arc-16s-id95.fasta,silva-arc-16s-id95.fasta-db \
--reads pair_1.fq --reads pair_2.fq \
--aligned mapped --fastx
  • --fastx     output FASTA/FASTQ file  
  • --aligned      STRING          aligned reads filepath + base file name 

mapped.fastqとnohit.fastqが出力される。

 

--fastx、--otherのほかに、--samや--blastがある。詳細はhelpから確認してください。

引用

SortMeRNA: fast and accurate filtering of ribosomal RNAs in metatranscriptomic data.

Kopylova E, Noé L, Touzet H.

Bioinformatics. 2012 Dec 15;28(24):3211-7.