macでインフォマティクス

macでインフォマティクス

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

(ヒトゲノム向け)ニューラルネットワークを使用したロングリードのSVコーラー NanoVar

 

 構造変異は、多くのヒト疾患の発症に関与しており、ヒト集団の遺伝的変異の大部分を占めている(ref.3,4)。 50 b5を超えるゲノム変化として定義される構造変異(SV)は、遺伝子調節異常または新規遺伝子融合を引き起こす可能性のある遺伝子病変を形成することにより、細胞生理学に機能的に影響を与え、ガン(ref.6,7)、メンデル障害(ref.8、 9)、および複雑な疾患(ref.10)などの疾患の発症を促進する。 SVは、欠失、重複、挿入、逆位、転座などのさまざまなクラスとして存在する。長年にわたり、疾患に関連するSVは、診断、予後、および患者の治療ガイダンスのバイオマーカーとして示されており、クリニックでのシーケンスベースおよび非シーケンスベースの方法でスクリーニングされた。 SVの臨床的影響が明らかになり続けているため、バイオマーカーを促進するために患者のルーチンSVプロファイリングを行うための正確で迅速かつ安価なワークフローに対する明確なニーズがある。現在、包括的なSV検出のためのシーケンスベースの方法には、ロングリードまたは第3世代シーケンシング(3GS)とショートリードまたは第2世代の2つの主要な標準がある。3GSテクノロジーは多くのユーザーが利用できるようになったが、シーケンスエラー率が高くスループットが低いため、2GSテクノロジーに取って代わっていない。 3GSは現在、主に小さなゲノムの研究またはターゲットシーケンスに限定され、最近の研究では哺乳類の全ゲノムシーケンス(WGS)が報告されているが、古い技術と比較してメガベースあたりのシーケンスコストが高くなっている。 SV発見の領域では、多くのグループが、シーケンスエラー率が高いにもかかわらず、3GSアプローチが2GSよりも高いSV検出感度と分解能を提供したことを報告している。これは主に、短い配列では挿入の検出(50-200 bp)が不十分であり、新規配列の挿入やrepetitiveエレメントを含む大きなゲノム変異を解明できないためである。一方、長いリード長(> 1 kb)はマッピングのあいまいさを減らし、リピート配列と複雑なSVを解決し、ショートリードよりもはるかに広範囲のSVを発見する。 SVの検出機能は向上しているが、3GSのメガベースあたりの低スループットと高いシーケンスコストにより、患者の日常的なSV探査で使用する可能性が妨げられている。
これらの問題を克服するために、患者の正確なSV特性評価のために、低デプスのオックスフォードナノポアテクノロジー(ONT)WGSデータを利用する新しいSVコールツールであるNanoVarを開発した。 NanoVarは、すべてのSVクラスに対して、信頼性の高いSV検出およびSV接合性推定のために、ニューラルネットワークベースのアルゴリズムを採用している。フローセルのケミストリの性質、ライブラリー調製に応じて、1から5回のONT MinIONシーケンス実行で達成できる、合計塩基数が4Xまたは12ギガベース(Gb)の最小シーケンスデプスで浅いロングリードのWGSデータを処理するように最適化されている。本論文では、NanoVarのSV検出精度を評価し、シミュレーションデータセットを使用して他のツールと比較した。

 

Githubにシミュレーションデータのダウンロードリンクあり。

 

インストール

ubuntu18.04LTSのpython3.7環境でテストした。

依存

Linux (x86_64 architecture, tested in Ubuntu 14.04, 16.04, 18.04)

  • bedtools >=2.26.0
  • makeblastdb and windowmasker
  • hs-blastn

Github 

#bioconda (link)
conda create -n nanovar -c bioconda -y nanovar python=3.7
source activate nanovar

nanovar -h

$ nanovar -h

usage: nanovar [-h] [-f FILTER_BED] [-l MINLEN] [-p SPLITPCT] [-a MINALIGN]

               [-b BUFFER] [-s SCORE] [-v] [-q] [-t [1-54]] [--force]

               [--mdb MDB] [--wmk WMK] [--hsb HSB]

               [long_reads] [reference_genome] [working_directory]

 

NanoVar is a neural-network-based structural variant (SV) caller that utilizes 

low-depth long-read sequencing data.

 

positional arguments:

  [long_reads]          Path to long reads. Formats: fasta/fa/fa.gzip/fa.gz

                        fastq/fq/fq.gzip/fq.gz

  [reference_genome]    Path to reference genome in FASTA. Genome indexes

                        created will overwrite indexes created by other

                        aligners (e.g. bwa)

  [working_directory]   Path to working directory. Directory will be created

                        if it does not exist

 

optional arguments:

  -h, --help            show this help message and exit

  -f FILTER_BED, --filter_bed FILTER_BED

                        BED file with genomic regions to be excluded. (e.g.

                        telomeres and centromeres). Either specify name of in-

                        built reference genome filter (i.e. hg38, hg19, mm10)

                        or provide FULL path to own BED file. [None]

  -l MINLEN, --minlen MINLEN

                        minimum length of SV to be detected. [25]

  -p SPLITPCT, --splitpct SPLITPCT

                        minimum percentage of unmapped bases within a long

                        read to be considered as a split-read. 0.05<=p<=0.50

                        [0.05]

  -a MINALIGN, --minalign MINALIGN

                        minimum alignment length for single alignment reads.

                        [200]

  -b BUFFER, --buffer BUFFER

                        nucleotide length buffer for SV breakend clustering.

                        [50]

  -s SCORE, --score SCORE

                        score threshold for defining PASS/FAIL SVs in VCF.

                        Default score 2.6 was derived from simulated analysis.

                        [2.6]

  -v, --version         show version and exit

  -q, --quiet           hide verbose

  -t [1-54], --threads [1-54]

                        number of available threads for use, max=54 [2]

  --force               run full pipeline, do not skip any redundant steps

                        (e.g. index generation)

  --mdb MDB             Specify path to 'makeblastdb' executable

  --wmk WMK             Specify path to 'windowmasker' executable

  --hsb HSB             Specify path to 'hs-blastn' executable

 

 

実行方法

ONTのfastqとレファレンス、オプションで除外する領域情報を持つbedファイルを指定する。

nanovar -t 24 -f hg38 long_reads.fq reference.fa working_dir

 

 

  • <long_reads>   Path to long reads. Formats: fasta/fa/fa.gzip/fa.gz fastq/fq/fq.gzip/fq.gz
  • <reference_genome>    Path to reference genome in FASTA. Genome indexes created will overwrite indexes created by other aligners (e.g. bwa)
  • <working_directory>   Path to working directory. Directory will be created if it does not exist
  • -f   BED file with genomic regions to be excluded. (e.g. telomeres and centromeres). Either specify name of in- built reference genome filter (i.e. hg38, hg19, mm10) or provide FULL path to own BED file. [None]

 

引用

NanoVar: Accurate Characterization of Patients’ Genomic Structural Variants Using Low-Depth Nanopore Sequencing
Cheng Yong Tham, Roberto Tirado-Magallanes, Yufen Goh, Melissa J. Fullwood, Bryan T.H. Koh, Wilson Wang, Chin Hin Ng, Wee Joo Chng, Alexandre Thiery, Daniel G. Tenen, Touati Benoukraf

bioRxiv preprint first posted online Jun. 17, 2019

 

関連


 

メタゲノムアセンブラの注意点

2019 11/25 誤字修正

 

 メタゲノムのde novoアセンブラについて少し誤解している人がいたので、注意喚起を兼ねて簡単にまとめておく。

 

 メタゲノムのデータセットは特定の環境の生物の混ぜ物のシーケンシングリードに由来しているため、よく似ているがわずかに異なるclosely relatedなゲノム配列が様々なabundanceで存在している。メタゲノムのアセンブラは、この混合の(バルクの)シーケンシングリードセットから、複数の近縁な株のコンセンサスのバックボーン配列を構築するアセンブラである。この性質により、配列多少の多様性があっても、同種のゲノムはより長い領域に渡って再構成される。配列が高度に断片化していると、ゲノムを分類したり特徴を調べたりすることができなくなるため、このようなメタゲノムアセンブラの一般的性質は、下流解析を効率的に進めるためにも欠かせない特徴と言える。

 しかしながら、この特徴により、異なる生物間で保存性の高い領域も1つの配列にまとめられる可能性がある。このことを理解していないと、ターゲットの配列が1生物からのみ由来すると判断して、解析に影響を及ぼす恐れがある。ここでは、ツールに振り回されないためにも、シミュレーションリードを使ってメタゲノムアセンブラの性質を実感できる簡単なテストを行ってみる。

 

 

異なるクラスタの生物間で長い領域に渡って保存性されている配列として、遺伝子クラスターが挙げられる。ここでは、よく研究されている遺伝子クラスターであるニトロゲナーゼのnifクラスターを例に説明する。下の図は、このnifクラスターの配列をNCBI のGenBnakで公開されている4生物からfetchして、blastnによるアラインメント比較した結果である。ペアワイズの配列比較と作図には、GenBankから自動でblast比較して簡単にゲノム比較の図を作れるEasyFig(紹介)を使っている。

f:id:kazumaxneo:20191126005350p:plain

 

1、上では4つしかないが、以下の表の5生物のnif clusterをfetchした。これらのnif clusterの配列をリファレンスに、ARTのHiseq2500プロファイルでシミュレートして、ペアエンドのリード(150x2)を表のcoverageの深さで発生させた。

f:id:kazumaxneo:20191125212743p:plain

 

ディレクトリ構成

f:id:kazumaxneo:20191125213617p:plain

 

2、リードを結合して1つのファイルにまとめ、 de novoアセンブリを実行する。

cat */paired_1.fq > mix_1.fq
cat */paired_2.fq > mix_2.fq

#spades
spades -1 mix_1.fq -2 mix_2.fq -k auto -o spades

#metaspades
metaspades -1 mix_1.fq -2 mix_2.fq -k auto -o metaspades

 

3、spadesアセンブリのGFAファイルをBandageに読み込んで視覚化した。

f:id:kazumaxneo:20191125212908p:plain

カバレッジが低い生物のnif clusterは断片化しているが、長さとカバレッジを見る限り、カバレッジが十分な生物のnif clusterはほぼ全長が再構成できている。

 

次がmetaspadesのアセンブリ結果。

f:id:kazumaxneo:20191125213023p:plain

metaspadesでは再構成に成功している配列もあるものの、上の交差した配列群のように、共通した配列部分でパスを共有している領域がある。

 

結論

アンプリコンシーケンスで97% OTUなどの基準で比較を行うように、メタゲノムのアセンブラも、株レベルの細かな違いは無視して、近い距離にあるゲノムを1つのコンセンサス配列として表し、下流の様々な解析に適した配列を作成しようとする。そのため、メタゲノムのアセンブラを使うと、保存された遺伝子クラスターがパスを共有したコンセンサス配列として表現される場合があることを示した。では、仮に株レベルの違いまでアセンブラの感度を上げるとどうなるか? 実際のメタゲノムのシーケンシングリードをspadesアセンブラなどでアセンブリすれば分かるが、アセンブリグラフの様々な部位で分岐が生じるため、進化速度が遅い領域を除き、大半の領域では、アセンブリ後の配列は高度に断片化してしまう。これではそもそもアセンブルをする意味が薄くなってしまう。

 メタゲノムのアセンブラは、まじめに考えれば解釈が困難な混合物のシーケンシングデータをなんとかして処理するため、また計算負荷を下げるため、アセンブリグラフをシンプルにし長い配列を出力しようと試みる。しかし、これに対するトレードオフの1つが、よく似た配列のある種誤ったコンセンサス表現である。Nif clusterのようなよく似た長い配列を分離するには、クラスターを構成するリードを回収して精度の高いspadesなどのアセンブラで再アセンブリするか(*1)、倍数体ゲノムからのハプロタイプ再構成のような追加のテクニックが必要になる。ということで、メタゲノムシーケンスの目的が保存された遺伝子クラスターの分析である場合、そのままwholeアセンブリするだけでは良い結果を出せない可能性があることを示した。 

 

 

補足

実際に、これらのnifクラスターはどのくらい似ているのか、追加でall versus allのANI比較とmauveのゲノム比較を行った。
f:id:kazumaxneo:20191125212712p:plain

(ANI比較は4つのnifクラスターしか使用していない)。

ANI比較ではそれぞれのnifクラスターは70%以下の塩基配列相同性になった。

 

mauveによるゲノム比較。シンテシーブロック表現のeasyfigより細かい部分まで視覚化される。

f:id:kazumaxneo:20191125214229p:plain

 

追記

たとえcomplete MAGsが得られたとしても、それをpublic databaseに通常のゲノムとして登録することは避けるべきです。以下を読んでみてください。


 

 

 

*1

コミュニティの複雑性が低ければうまくいくかもしれない。 

GFF3のツールキット GFF3toolkit

 

i5k Workspace @ NAL (HP) でサポートされているGFF3toolkit(https://github.com/NAL-i5K/GFF3toolkit)は、節足動物ゲノムプロジェクトとその研究コミュニティからのGFF3形式の遺伝子アノテーションを処理するためのツールスイートを提供する。 遺伝子アノテーションのGFF3フォーマットを改善するために、品質管理とマージ手順がGFF3toolkitとともに提案されている。 特に、このツールキットは、GFF3ファイルのソート、GFF3形式エラーの検出、2つのGFF3ファイルのマージ、GFF3ファイルからの生物学的シーケンスの生成を行う機能を提供する。

 

 

インストール

ubuntu18.04のpython3.7環境でテストした(docker使用、ホストOS macos10.14)。

依存

Github

pip install gff3tool

#latest version
pip install git+https://github.com/NAL-i5K/GFF3toolkit.git

gff3_QC -h

$ gff3_QC -h

usage: gff3_QC [-h] [-g GFF] [-f FASTA] [-noncg] [-i] [-n ALLOWED_NUM_OF_N]

               [-t [CHECK_N_FEATURE_TYPES [CHECK_N_FEATURE_TYPES ...]]]

               [-o OUTPUT] [-s STATISTIC] [-v]

 

    Testing environment:

    1. Python 2.7

 

    Inputs:

    1. GFF3: Specify the file name with the -g or --gff argument; Please note that this program requires gene/pseudogene and mRNA/pseudogenic_transcript to have an ID attribute in column 9.

    2. fasta file: Specify the file name with the -f or --fasta argument

 

    Outputs:

    1. Error report for the input GFF3 file

* Line_num: Line numbers of the found problematic models in the input GFF3 file.

* Error_code: Error codes for the found problematic models. Please refer to lib/ERROR/ERROR.py to see the full list of Error_code and the corresponding Error_tag.

        * Error_tag: Detail of the found errors for the problematic models. Please refer to lib/ERROR/ERROR.py to see the full list of Error_code and the corresponding Error_tag.

 

    Quick start:

    gff3_QC -g example_file/example.gff3 -f example_file/reference.fa -o test

    or

    gff3_QC --gff example_file/example.gff3 --fasta example_file/reference.fa --output test

 

optional arguments:

  -h, --help            show this help message and exit

  -g GFF, --gff GFF     Genome annotation file, gff3 format

  -f FASTA, --fasta FASTA

                        Genome sequences, fasta format

  -noncg, --noncanonical_gene

                        gff3 file is not formatted in the canonical gene model

                        format.

  -i, --initial_phase   Check whether initial CDS phase is 0 (default: no

                        check)

  -n ALLOWED_NUM_OF_N, --allowed_num_of_n ALLOWED_NUM_OF_N

                        Max number of Ns allowed in a feature, anything more

                        will be reported as an error (default: 0)

  -t [CHECK_N_FEATURE_TYPES [CHECK_N_FEATURE_TYPES ...]], --check_n_feature_types [CHECK_N_FEATURE_TYPES [CHECK_N_FEATURE_TYPES ...]]

                        Count the number of Ns in each feature with the type

                        specified, multiple types may be specified, ex: -t CDS

                        exon (default: "CDS")

  -o OUTPUT, --output OUTPUT

                        output file name (default: report.txt)

  -s STATISTIC, --statistic STATISTIC

                        statistic file name (default: statistic.txt)

  -v, --version         show program's version number and exit

gff3_fix -h

$ gff3_fix -h

usage: gff3_fix [-h] [-qc_r QC_REPORT] [-g GFF] [-og OUTPUT_GFF] [-v]

 

Testing environment:

1. Python 3.*

 

Input:

1. Error report: Error report from gff3_QC.py. Specify the file name with the -qc_r or --qc_report argument;

2. GFF3: Specify the file name with the -g or --gff argument;

 

Output:

1. Corrected GFF3

 

Quick start:

gff3_fix -qc_r error.txt -g example_file/example.gff3 -og corrected.gff3

 

optional arguments:

  -h, --help            show this help message and exit

  -qc_r QC_REPORT, --qc_report QC_REPORT

                        Error report from gff3_QC.py

  -g GFF, --gff GFF     Genome annotation file, gff3 format

  -og OUTPUT_GFF, --output_gff OUTPUT_GFF

                        output gff3 file name

  -v, --version         show program's version number and exit

gff3_sort

$ gff3_sort

usage: gff3_sort [-h] [-g GFF_FILE] [-og OUTPUT_GFF] [-t SORT_TEMPLATE] [-i]

                 [-v] [-r]

 

Sort a GFF3 file according to the order of Scaffold (seqID), coordinates on a Scaffold, and feature relationship based on sequence ontology.

 

Inputs:

1. GFF3 file: Specify the file name with the -g argument

 

Outputs:

1. Sorted GFF3 file: Specify the file name with the -og argument

 

Examples:

1. Specify the input, output file names and options using short arguments:

   gff3_sort -g example_file/example.gff3 -og example_file/example_sorted.gff

2. Specify the input, output file names and options using long arguments:

   gff3_sort --gff_file example_file/example.gff3 --output_gff example_file/example_sorted.gff

 

optional arguments:

  -h, --help            show this help message and exit

  -g GFF_FILE, --gff_file GFF_FILE

                        GFF3 file that you would like to sort.

  -og OUTPUT_GFF, --output_gff OUTPUT_GFF

                        Sorted GFF3 file

  -t SORT_TEMPLATE, --sort_template SORT_TEMPLATE

                        A file that indicates the sorting order of features

                        within a gene model

  -i, --isoform_sort    Sort multi-isoform gene models by feature type

                        (default: False)

  -v, --version         show program's version number and exit

  -r, --reference       Sort scaffold (seqID) by order of appearance in gff3

                        file (default is by number)

gff3_to_fasta -h

$ gff3_to_fasta -h

usage: gff3_to_fasta [-h] [-g GFF] [-f FASTA] [-embf] [-st SEQUENCE_TYPE]

                     [-u [USER_DEFINED [USER_DEFINED ...]]] [-d DEFLINE]

                     [-o OUTPUT_PREFIX] [-noQC] [-v]

 

Extract sequences from specific regions of genome based on gff file.

Testing enviroment:

1. Python 2.7

 

Required inputs:

1. GFF3: specify the file name with the -g argument

2. Fasta file: specify the file name with the -f argument

3. Output prefix: specify with the -o argument

 

Outputs:

1. Fasta formatted sequence file based on the gff3 file.

 

Example command:

gff3_to_fasta -g example_file/example.gff3 -f example_file/reference.fa -st all -d simple -o test_sequences

 

optional arguments:

  -h, --help            show this help message and exit

  -g GFF, --gff GFF     Genome annotation file in GFF3 format

  -f FASTA, --fasta FASTA

                        Genome sequences in FASTA format

  -embf, --embedded_fasta

                        Specify this option if you want to extract sequence from embedded fasta.

  -st SEQUENCE_TYPE, --sequence_type SEQUENCE_TYPE

                        Type of sequences you would like to extract: 

                        "all" - FASTA files for all types of sequences listed below, except user_defined;

                        "gene" - gene sequence for each record;

                        "exon" - exon sequence for each record;

                        "pre_trans" - genomic region of a transcript model (premature transcript);

                        "trans" - spliced transcripts (only exons included);

                        "cds" - coding sequences;

                        "pep" - peptide sequences;

                        "user_defined" - specify parent and child features via the -u argument.

  -u [USER_DEFINED [USER_DEFINED ...]], --user_defined [USER_DEFINED [USER_DEFINED ...]]

                        Specify parent and child features for fasta extraction, format: [parent feature type] [child feature type] (ex: -u mRNA CDS). Required if -st user_defined is given.

  -d DEFLINE, --defline DEFLINE

                        Defline format in the output FASTA file:

                        "simple" - only ID would be shown in the defline;

                        "complete" - complete information of the feature would be shown in the defline.

  -o OUTPUT_PREFIX, --output_prefix OUTPUT_PREFIX

                        Prefix of output file name

  -noQC, --quality_control

                        Specify this option if you do not want to execute quality control for gff file. (default: QC is executed)

  -v, --version         show program's version number and exit

gff3_merge -h

$ gff3_merge -h

usage: gff3_merge [-h] [-g1 GFF_FILE1] [-g2 GFF_FILE2] [-f FASTA]

                  [-u1 USER_DEFINED_FILE1] [-u2 USER_DEFINED_FILE2]

                  [-og OUTPUT_GFF] [-r REPORT_FILE] [-a] [-noAuto] [-v]

 

Merge two gff files of the same genome into one.

 

Testing enviroment:

1. Python 3.*

 

Inputs:

1. GFF3 file 1: Gff with annotations modified relative to the original gff (e.g. output from the Apollo program), specify the file name with the -g1 argument

2. GFF3 file 2: Original/Reference gff, specify the file name with the -g2 argument

3. FASTA: Genomic sequences in the FASTA format with the -f argument

 

Outputs:

1. Merged GFF3: Models from GFF3 file 1 replace Models from GFF3 file 2 based on their replace tag. Specify the output file name with the -og argument

2. Log report for the integration: specify the file name with the -r argument

 

Examples:

1. Specify the input, output file names and options using short arguments:

   gff3_merge -g1 example_file/new_models.gff3 -g2 example_file/reference.gff3 -f example_file/reference.fa -og merged.gff -r merged_report.txt

2. Specify the input, output file names and options using long arguments:

   gff3_merge --gff_file1 example_file/new_models.gff3 --gff_file2 example_file/reference.gff3 --fasta example_file/reference.fa --output_gff merged.gff --report_file merged_report.txt

 

optional arguments:

  -h, --help            show this help message and exit

  -g1 GFF_FILE1, --gff_file1 GFF_FILE1

                        Updated GFF3 file, such as Apollo gff

  -g2 GFF_FILE2, --gff_file2 GFF_FILE2

                        Reference GFF3 file, such as Maker gff or OGS gff

  -f FASTA, --fasta FASTA

                        Genomic sequences in the fasta format

  -u1 USER_DEFINED_FILE1, --user_defined_file1 USER_DEFINED_FILE1

                        File for specifing parent and child features for fasta

                        extraction from updated GFF3 file.

  -u2 USER_DEFINED_FILE2, --user_defined_file2 USER_DEFINED_FILE2

                        File for specifing parent and child features for fasta

                        extraction from reference GFF3 file.

  -og OUTPUT_GFF, --output_gff OUTPUT_GFF

                        The merged GFF3 file

  -r REPORT_FILE, --report_file REPORT_FILE

                        Log file for the integration

  -a, --all             auto-assignment replace tags for all transcript

                        features. (default: Only automatically assign replace

                        tags for the transcript without replace tags)

  -noAuto, --auto_assignment

                        Turn off the auto-assignment of replace tags, if you

                        already have replace tags in your updated gff

                        (default: Automatically assign replace tags and then

                        merge the gff files)

  -v, --version         show program's version number and exit

 

 

実行方法

ここではテストファイルを使う。

git clone https://github.com/NAL-i5K/GFF3toolkit.git
cd GFF3toolkit/

 

gff3_sort - scaffolds順、座標、およびアノテーションのparent-child の関係に従ってGFF3ファイルを並べ替える

gff3_sort -g example_file/example.gff3 -og example-sorted.gff3
  •  -og    Sorted GFF3 file

 

gff3_to_fasta - ゲノムの特定の領域から配列(スプライスされた転写産物、cds、ペプチド)を抽出

gff3_to_fasta -g example_file/example.gff3 -f example_file/reference.fa -st all -d simple -o test_sequences
  • -st    Type of sequences you would like to extract:  
  1. "all" - FASTA files for all types of sequences listed below, except user_defined;
  2. "gene" - gene sequence for each record;
  3. "exon" - exon sequence for each record;
  4. "pre_trans" - genomic region of a transcript model (premature transcript);
  5. "trans" - spliced transcripts (only exons included);
  6. "cds" - coding sequences;
  7. "pep" - peptide sequences;
  8. "user_defined" - specify parent and child features via the -u argument. 

 

gff3_QC - GFFの様々なタイプのエラーを検出 (link)

gff3_QC -g example_file/example.gff3 -f example_file/reference.fa -o error.txt -s statistic.txt

statistic.txt

f:id:kazumaxneo:20191125015900p:plain

error.txt

f:id:kazumaxneo:20191125015845p:plain

 

gff3_fix - QCで検出されたエラーを修正 (link)

gff3_fix -qc_r error.txt -g example_file/example.gff3 -og corrected.gff3
  •  -qc_r   Error report from gff3_QC.py

 

gff3_merge - 2つのGFFをマージ (link)

gff3_merge -g1 example_file/new_models.gff3 -g2 example_file/reference.gff3 -f example_file/reference.fa -og merged.gff -r merged_report.txt

 

 

引用

The GFF3toolkit: QC and Merge Pipeline for Genome Annotation

Methods Mol Biol. 2019;1858:75-87
Chen MM, Lin H, Chiang LM, Childers CP, Poelchau MF

 

関連


 

GenBankから配列やアノテーションを取り出すWebサービス FeatureExtract

 

 イントロン/エクソン構造、プロモーター領域の内容、上流域および下流域における他の遺伝子の位置など、DNA配列の特徴のアノテーションに簡単にアクセスできることは、多くの生物学的問題に取り組むことが非常に有益である。たとえば、遺伝子内のイントロンの配置を考慮すると、相同な遺伝子の系統解析に役立つ。 PCRまたはDNAマイクロアレイを使用してUTR領域を調査するための実験を設計するには、UTR領域の既知のエレメントと染色体上の他の遺伝子の位置とストランドの情報が必要である。このような情報は、GenBankNCBI Human Genome buildsなどのデータベースに豊富に記録され、文書化されている。ただし、この情報にアクセスするには、通常、かなりのバイオインフォマティクススキルとデータ形式の詳細な知識が必要になる。ここでは、GenBankエントリから機能アノテーションを抽出するための非常に柔軟で使いやすいツールを紹介する。このツールは、特定の機能(プロモーターなど)に対応するデータセットの抽出にも役立つ。最も重要なことは、出力データ形式が一貫しており、ユーザーにとって扱いやすく、計算的に解析しやすいことである。 FeatureExtract Webサーバーは、http://www.cbs.dtu.dk/services/FeatureExtract/でアカデミックおよび商用の両方で自由に利用できる。

 

 

使い方

http://www.cbs.dtu.dk/services/FeatureExtract/にアクセスする。

f:id:kazumaxneo:20191124142556p:plain

 

ここではexampleのGenBankファイルを使う。default条件でsubmit。

f:id:kazumaxneo:20191124143718p:plain

 

結果

エントリーの数、intronを含むエントリーの数、トータルサイズ(bp)などがまとめられる。

f:id:kazumaxneo:20191124143828p:plain

 

fasta形式のファイルとしてダウンロードできる。

f:id:kazumaxneo:20191124144059p:plain

 

submit前にカスタマイズすることもできる。

defaultではCDS フィーチャが対象になる。geneやrRNAなどに変更可能。

f:id:kazumaxneo:20191124144951p:plain

 

GenBankのfullバージョンならフィーチャの隣接領域も取り出せる。

f:id:kazumaxneo:20191124150428p:plain

 

隣接領域は小文字になる。

f:id:kazumaxneo:20191124150633p:plain

 

 

補足

SMSのwebツールでもGenBank => fastA変換が利用できます。SMSにはEMBL formatからのアノテーション抽出ルールなどもあります。

https://www.bioinformatics.org/sms2/genbank_fasta.html

 

Sequence format converterというツールもあります。

Sequence format converter

f:id:kazumaxneo:20191124145909p:plain

こちらも様々な出力に対応。

引用

FeatureExtract—extraction of sequence annotation made easy
Rasmus Wernersson
Nucleic Acids Research, Volume 33, Issue suppl_2, 1 July 2005, Pages W567–W569, 

UCSC、NCBI、Ensemblからゲノムをダウンロードする genomepy

2021 10/9 コマンドの修正(バージョンアップ)

 

タイトルの通りのツール。簡単に紹介します。

 

インストール

依存

  • tabix
  • genePredToBed
  • genePredToGtf
  • bedToGenePred
  • gtfToGenePred
  • gff3ToGenePred
conda install -c bioconda -y ucsc-genepredtobed
conda install -c bioconda -y ucsc-genepredtogtf
conda install -c bioconda -y ucsc-gtftogenepred
conda install -c bioconda -y ucsc-bedtogenepred
conda install -c bioconda -y ucsc-gff3togenepred

本体 Github

#bioconda (link)
conda install -c bioconda -y genomepy

#pip
pip install genomepy

genomepy -h

$ genomepy -h

Usage: genomepy [OPTIONS] COMMAND [ARGS]...

 

Options:

  --version   Show the version and exit.

  -h, --help  Show this message and exit.

 

Commands:

  config     manage configuration

  genomes    list available genomes

  install    install genome

  plugin     manage plugins

  providers  list available providers

  search     search for genomes

 

 

実行方法

1、ONにするpluginを確認。

genomepy plugin list

$ genomepy plugin list

plugin              enabled

blacklist           

bowtie2             

bwa                 

gaps                *

gmap                

hisat2              

minimap2            

sizes               *

アライナーはbowtie2、bwa、gmap、hisat2、minimap2、starを指定できる。ONにすると、そのアライナーのindexがgenomeのダウンロード後に作成される。当然、そのアライナーのコマンドが使用できない環境では作成されない。blascklistは以下の論文が詳しい(link)。

 

(optional)、 pluginのbwaをONにする。

genomepy plugin enable bwa

ONになっているか確認する。

> genomepy plugin list 

$ genomepy plugin list

plugin              enabled

blacklist           

bowtie2             

bwa                 *

gaps                *

gmap                

hisat2              

minimap2            

sizes               *

ONになっている。 

 

3、利用できるデータベースを確認する。

genomepy providers

$ genomepy providers

Ensembl

UCSC

NCBI

UCSCNCBIEnsemblのゲノムが利用できる。

 

例えばUCSCからダウンロードできるゲノムを確認する。

genomepy genomes -p UCSC

$ genomepy genomes -p UCSC

UCSC hg38 Human Dec. 2013 (GRCh38/hg38) Genome at UCSC

UCSC hg19 Human Feb. 2009 (GRCh37/hg19) Genome at UCSC

UCSC hg18 Human Mar. 2006 (NCBI36/hg18) Genome at UCSC

UCSC hg17 Human May 2004 (NCBI35/hg17) Genome at UCSC

UCSC hg16 Human July 2003 (NCBI34/hg16) Genome at UCSC

UCSC mm10 Mouse Dec. 2011 (GRCm38/mm10) Genome at UCSC

UCSC mm9 Mouse July 2007 (NCBI37/mm9) Genome at UCSC

UCSC mm8 Mouse Feb. 2006 (NCBI36/mm8) Genome at UCSC

UCSC mm7 Mouse Aug. 2005 (NCBI35/mm7) Genome at UCSC

UCSC anoGam3 A. gambiae Oct. 2006 (AgamP3/anoGam3) Genome at UCSC

UCSC anoGam1 A. gambiae Feb. 2003 (IAGEC MOZ2/anoGam1) Genome at UCSC

UCSC apiMel2 A. mellifera Jan. 2005 (Baylor 2.0/apiMel2) Genome at UCSC

UCSC apiMel1 A. mellifera July 2004 (Baylor 1.2/apiMel1) Genome at UCSC

UCSC xenLae2 African clawed frog Aug. 2016 (Xenopus_laevis_v2/xenLae2) Genome at UCSC

UCSC vicPac2 Alpaca Mar. 2013 (Vicugna_pacos-2.0.1/vicPac2) Genome at UCSC

UCSC vicPac1 Alpaca Jul. 2008 (Broad/vicPac1) Genome at UCSC

UCSC allMis1 American alligator Aug. 2012 (allMis0.2/allMis1) Genome at UCSC

UCSC dasNov3 Armadillo Dec. 2011 (Baylor/dasNov3) Genome at UCSC

UCSC gadMor1 Atlantic cod May 2010 (Genofisk GadMor_May2010/gadMor1) Genome at UCSC

UCSC papAnu4 Baboon Apr. 2017 (Panu_3.0/papAnu4) Genome at UCSC

UCSC papAnu2 Baboon Mar. 2012 (Baylor Panu_2.0/papAnu2) Genome at UCSC

UCSC papHam1 Baboon Nov. 2008 (Baylor Pham_1.0/papHam1) Genome at UCSC

UCSC bisBis1 Bison Oct. 2014 (Bison_UMD1.0/bisBis1) Genome at UCSC

UCSC panPan2 Bonobo Aug. 2015 (MPI-EVA panpan1.1/panPan2) Genome at UCSC

UCSC panPan1 Bonobo May 2012 (Max-Planck/panPan1) Genome at UCSC

UCSC aptMan1 Brown kiwi Jun. 2015 (MPI-EVA AptMant0/aptMan1) Genome at UCSC

UCSC melUnd1 Budgerigar Sep. 2011 (WUSTL v6.3/melUnd1) Genome at UCSC

UCSC otoGar3 Bushbaby Mar. 2011 (Broad/otoGar3) Genome at UCSC

UCSC caePb2 C. brenneri Feb. 2008 (WUGSC 6.0.1/caePb2) Genome at UCSC

UCSC caePb1 C. brenneri Jan. 2007 (WUGSC 4.0/caePb1) Genome at UCSC

UCSC cb3 C. briggsae Jan. 2007 (WUGSC 1.0/cb3) Genome at UCSC

UCSC cb1 C. briggsae July 2002 (WormBase cb25.agp8/cb1) Genome at UCSC

UCSC ce11 C. elegans Feb. 2013 (WBcel235/ce11) Genome at UCSC

UCSC ce10 C. elegans Oct. 2010 (WS220/ce10) Genome at UCSC

UCSC ce6 C. elegans May 2008 (WS190/ce6) Genome at UCSC

UCSC ce4 C. elegans Jan. 2007 (WS170/ce4) Genome at UCSC

UCSC ce2 C. elegans Mar. 2004 (WS120/ce2) Genome at UCSC

UCSC ci3 C. intestinalis Apr. 2011 (Kyoto KH/ci3) Genome at UCSC

UCSC ci2 C. intestinalis Mar. 2005 (JGI 2.1/ci2) Genome at UCSC

UCSC ci1 C. intestinalis Dec. 2002 (JGI 1.0/ci1) Genome at UCSC

UCSC caeJap1 C. japonica Mar. 2008 (WUGSC 3.0.2/caeJap1) Genome at UCSC

UCSC caeRem3 C. remanei May 2007 (WUGSC 15.0.1/caeRem3) Genome at UCSC

UCSC caeRem2 C. remanei Mar. 2006 (WUGSC 1.0/caeRem2) Genome at UCSC

UCSC felCat9 Cat Nov. 2017 (Felis_catus_9.0/felCat9) Genome at UCSC

UCSC felCat8 Cat Nov. 2014 (ICGSC Felis_catus_8.0/felCat8) Genome at UCSC

UCSC felCat5 Cat Sep. 2011 (ICGSC Felis_catus 6.2/felCat5) Genome at UCSC

UCSC felCat4 Cat Dec. 2008 (NHGRI/GTB V17e/felCat4) Genome at UCSC

UCSC felCat3 Cat Mar. 2006 (Broad/felCat3) Genome at UCSC

UCSC galGal6 Chicken Mar. 2018 (GRCg6a/galGal6) Genome at UCSC

UCSC galGal5 Chicken Dec 2015 (Gallus_gallus-5.0/galGal5) Genome at UCSC

UCSC galGal4 Chicken Nov. 2011 (ICGSC Gallus_gallus-4.0/galGal4) Genome at UCSC

UCSC galGal3 Chicken May 2006 (WUGSC 2.1/galGal3) Genome at UCSC

UCSC galGal2 Chicken Feb. 2004 (WUGSC 1.0/galGal2) Genome at UCSC

UCSC panTro6 Chimp Jan. 2018 (Clint_PTRv2/panTro6) Genome at UCSC

UCSC panTro5 Chimp May 2016 (Pan_tro 3.0/panTro5) Genome at UCSC

UCSC panTro4 Chimp Feb. 2011 (CSAC 2.1.4/panTro4) Genome at UCSC

UCSC panTro3 Chimp Oct. 2010 (CGSC 2.1.3/panTro3) Genome at UCSC

UCSC panTro2 Chimp Mar. 2006 (CGSC 2.1/panTro2) Genome at UCSC

UCSC panTro1 Chimp Nov. 2003 (CGSC 1.1/panTro1) Genome at UCSC

UCSC criGriChoV2 Chinese hamster Jun. 2017 (CHOK1S_HZDv1/criGriChoV2) Genome at UCSC

UCSC criGriChoV1 Chinese hamster Aug. 2011 (CHO K1 cell line/criGriChoV1) Genome at UCSC

UCSC criGri1 Chinese hamster Jul. 2013 (C_griseus_v1.0/criGri1) Genome at UCSC

UCSC manPen1 Chinese pangolin Aug 2014 (M_pentadactyla-1.1.1/manPen1) Genome at UCSC

UCSC latCha1 Coelacanth Aug. 2011 (Broad/latCha1) Genome at UCSC

UCSC bosTau9 Cow Apr. 2018 (ARS-UCD1.2/bosTau9) Genome at UCSC

UCSC bosTau8 Cow Jun. 2014 (Bos_taurus_UMD_3.1.1/bosTau8) Genome at UCSC

UCSC bosTau7 Cow Oct. 2011 (Baylor Btau_4.6.1/bosTau7) Genome at UCSC

UCSC bosTau6 Cow Nov. 2009 (Bos_taurus_UMD_3.1/bosTau6) Genome at UCSC

UCSC bosTau4 Cow Oct. 2007 (Baylor 4.0/bosTau4) Genome at UCSC

UCSC bosTau3 Cow Aug. 2006 (Baylor 3.1/bosTau3) Genome at UCSC

UCSC bosTau2 Cow Mar. 2005 (Baylor 2.0/bosTau2) Genome at UCSC

UCSC macFas5 Crab-eating macaque Jun. 2013 (Macaca_fascicularis_5.0/macFas5) Genome at UCSC

UCSC droAna2 D. ananassae Aug. 2005 (Agencourt prelim/droAna2) Genome at UCSC

UCSC droAna1 D. ananassae July 2004 (TIGR/droAna1) Genome at UCSC

UCSC droEre1 D. erecta Aug. 2005 (Agencourt prelim/droEre1) Genome at UCSC

UCSC droGri1 D. grimshawi Aug. 2005 (Agencourt prelim/droGri1) Genome at UCSC

UCSC dm6 D. melanogaster Aug. 2014 (BDGP Release 6 + ISO1 MT/dm6) Genome at UCSC

UCSC dm3 D. melanogaster Apr. 2006 (BDGP R5/dm3) Genome at UCSC

UCSC dm2 D. melanogaster Apr. 2004 (BDGP R4/dm2) Genome at UCSC

UCSC dm1 D. melanogaster Jan. 2003 (BDGP R3/dm1) Genome at UCSC

UCSC droMoj2 D. mojavensis Aug. 2005 (Agencourt prelim/droMoj2) Genome at UCSC

UCSC droMoj1 D. mojavensis Aug. 2004 (Agencourt prelim/droMoj1) Genome at UCSC

UCSC droPer1 D. persimilis Oct. 2005 (Broad/droPer1) Genome at UCSC

UCSC dp3 D. pseudoobscura Nov. 2004 (FlyBase 1.03/dp3) Genome at UCSC

UCSC dp2 D. pseudoobscura Aug. 2003 (Baylor freeze1/dp2) Genome at UCSC

UCSC droSec1 D. sechellia Oct. 2005 (Broad/droSec1) Genome at UCSC

UCSC droSim1 D. simulans Apr. 2005 (WUGSC mosaic 1.0/droSim1) Genome at UCSC

UCSC droVir2 D. virilis Aug. 2005 (Agencourt prelim/droVir2) Genome at UCSC

UCSC droVir1 D. virilis July 2004 (Agencourt prelim/droVir1) Genome at UCSC

UCSC droYak2 D. yakuba Nov. 2005 (WUGSC 7.1/droYak2) Genome at UCSC

UCSC droYak1 D. yakuba Apr. 2004 (WUGSC 1.0/droYak1) Genome at UCSC

UCSC canFam3 Dog Sep. 2011 (Broad CanFam3.1/canFam3) Genome at UCSC

UCSC canFam2 Dog May 2005 (Broad/canFam2) Genome at UCSC

UCSC canFam1 Dog July 2004 (Broad/canFam1) Genome at UCSC

UCSC turTru2 Dolphin Oct. 2011 (Baylor Ttru_1.4/turTru2) Genome at UCSC

UCSC eboVir3 Ebola virus Sierra Leone 2014 (G3683/KM034562.1/eboVir3) Genome at UCSC

UCSC loxAfr3 Elephant Jul. 2009 (Broad/loxAfr3) Genome at UCSC

UCSC calMil1 Elephant shark Dec. 2013 (Callorhinchus_milii-6.1.3/calMil1) Genome at UCSC

UCSC musFur1 Ferret  Apr. 2011 (MusPutFur1.0/musFur1) Genome at UCSC

UCSC fr3 Fugu Oct. 2011 (FUGU5/fr3) Genome at UCSC

UCSC fr2 Fugu Oct. 2004 (JGI 4.0/fr2) Genome at UCSC

UCSC fr1 Fugu Aug. 2002 (JGI 3.0/fr1) Genome at UCSC

UCSC thaSir1 Garter snake Jun. 2015 (Thamnophis_sirtalis-6.0/thaSir1) Genome at UCSC

UCSC nomLeu3 Gibbon Oct. 2012 (GGSC Nleu3.0/nomLeu3) Genome at UCSC

UCSC nomLeu2 Gibbon Jun. 2011 (GGSC Nleu1.1/nomLeu2) Genome at UCSC

UCSC nomLeu1 Gibbon Jan. 2010 (GGSC Nleu1.0/nomLeu1) Genome at UCSC

UCSC aquChr2 Golden eagle Oct. 2014 (aquChr-1.0.2/aquChr2) Genome at UCSC

UCSC rhiRox1 Golden snub-nosed monkey Oct. 2014 (Rrox_v1/rhiRox1) Genome at UCSC

UCSC gorGor5 Gorilla Mar. 2016 (GSMRT3/gorGor5) Genome at UCSC

UCSC gorGor4 Gorilla Dec 2014 (gorGor4.1/gorGor4) Genome at UCSC

UCSC gorGor3 Gorilla May 2011 (gorGor3.1/gorGor3) Genome at UCSC

UCSC chlSab2 Green monkey Mar. 2014 (Chlorocebus_sabeus 1.1/chlSab2) Genome at UCSC

UCSC cavPor3 Guinea pig Feb. 2008 (Broad/cavPor3) Genome at UCSC

UCSC eriEur2 Hedgehog May 2012 (EriEur2.0/eriEur2) Genome at UCSC

UCSC eriEur1 Hedgehog June 2006 (Broad/eriEur1) Genome at UCSC

UCSC equCab3 Horse Jan. 2018 (EquCab3.0/equCab3) Genome at UCSC

UCSC equCab2 Horse Sep. 2007 (Broad/equCab2) Genome at UCSC

UCSC equCab1 Horse Jan. 2007 (Broad/equCab1) Genome at UCSC

UCSC dipOrd1 Kangaroo rat Jul. 2008 (Broad/dipOrd1) Genome at UCSC

UCSC petMar3 Lamprey Dec. 2017 (Pmar_germline 1.0/petMar3) Genome at UCSC

UCSC petMar2 Lamprey Sep. 2010 (WUGSC 7.0/petMar2) Genome at UCSC

UCSC petMar1 Lamprey Mar. 2007 (WUGSC 3.0/petMar1) Genome at UCSC

UCSC braFlo1 Lancelet Mar. 2006 (JGI 1.0/braFlo1) Genome at UCSC

UCSC anoCar2 Lizard May 2010 (Broad AnoCar2.0/anoCar2) Genome at UCSC

UCSC anoCar1 Lizard Feb. 2007 (Broad/anoCar1) Genome at UCSC

UCSC galVar1 Malayan flying lemur Jun. 2014 (G_variegatus-3.0.2/galVar1) Genome at UCSC

UCSC triMan1 Manatee Oct. 2011 (Broad v1.0/triMan1) Genome at UCSC

UCSC calJac3 Marmoset March 2009 (WUGSC 3.2/calJac3) Genome at UCSC

UCSC calJac1 Marmoset June 2007 (WUGSC 2.0.2/calJac1) Genome at UCSC

UCSC oryLat2 Medaka Oct. 2005 (NIG/UT MEDAKA1/oryLat2) Genome at UCSC

UCSC geoFor1 Medium ground finch Apr. 2012 (GeoFor_1.0/geoFor1) Genome at UCSC

UCSC pteVam1 Megabat Jul. 2008 (Broad/pteVam1) Genome at UCSC

UCSC myoLuc2 Microbat Jul. 2010 (Broad Institute Myoluc2.0/myoLuc2) Genome at UCSC

UCSC balAcu1 Minke whale Oct. 2013 (BalAcu1.0/balAcu1) Genome at UCSC

UCSC micMur2 Mouse lemur May 2015 (Mouse lemur/micMur2) Genome at UCSC

UCSC micMur1 Mouse lemur Jul. 2007 (Broad/micMur1) Genome at UCSC

UCSC hetGla2 Naked mole-rat Jan. 2012 (Broad HetGla_female_1.0/hetGla2) Genome at UCSC

UCSC hetGla1 Naked mole-rat Jul. 2011 (BGI HetGla_1.0/hetGla1) Genome at UCSC

UCSC oreNil2 Nile tilapia Jan. 2011 (Broad oreNil1.1/oreNil2) Genome at UCSC

UCSC monDom5 Opossum Oct. 2006 (Broad/monDom5) Genome at UCSC

UCSC monDom4 Opossum Jan. 2006 (Broad/monDom4) Genome at UCSC

UCSC monDom1 Opossum Oct. 2004 (Broad prelim/monDom1) Genome at UCSC

UCSC ponAbe3 Orangutan Jan. 2018 (Susie_PABv2/ponAbe3) Genome at UCSC

UCSC ponAbe2 Orangutan July 2007 (WUGSC 2.0.2/ponAbe2) Genome at UCSC

UCSC priPac1 P. pacificus Feb. 2007 (WUGSC 5.0/priPac1) Genome at UCSC

UCSC chrPic1 Painted turtle Dec. 2011 (v3.0.1/chrPic1) Genome at UCSC

UCSC ailMel1 Panda Dec. 2009 (BGI-Shenzhen 1.0/ailMel1) Genome at UCSC

UCSC susScr11 Pig Feb. 2017 (Sscrofa11.1/susScr11) Genome at UCSC

UCSC susScr3 Pig Aug. 2011 (SGSC Sscrofa10.2/susScr3) Genome at UCSC

UCSC susScr2 Pig Nov. 2009 (SGSC Sscrofa9.2/susScr2) Genome at UCSC

UCSC ochPri3 Pika May 2012 (OchPri3.0/ochPri3) Genome at UCSC

UCSC ochPri2 Pika Jul. 2008 (Broad/ochPri2) Genome at UCSC

UCSC ornAna2 Platypus Feb. 2007 (ASM227v2/ornAna2) Genome at UCSC

UCSC ornAna1 Platypus Mar. 2007 (WUGSC 5.0.1/ornAna1) Genome at UCSC

UCSC nasLar1 Proboscis monkey Nov. 2014 (Charlie1.0/nasLar1) Genome at UCSC

UCSC oryCun2 Rabbit Apr. 2009 (Broad/oryCun2) Genome at UCSC

UCSC rn6 Rat Jul. 2014 (RGSC 6.0/rn6) Genome at UCSC

UCSC rn5 Rat Mar. 2012 (RGSC 5.0/rn5) Genome at UCSC

UCSC rn4 Rat Nov. 2004 (Baylor 3.4/rn4) Genome at UCSC

UCSC rn3 Rat June 2003 (Baylor 3.1/rn3) Genome at UCSC

UCSC rheMac10 Rhesus Feb. 2019 (Mmul_10/rheMac10) Genome at UCSC

UCSC rheMac8 Rhesus Nov. 2015 (BCM Mmul_8.0.1/rheMac8) Genome at UCSC

UCSC rheMac3 Rhesus Oct. 2010 (BGI CR_1.0/rheMac3) Genome at UCSC

UCSC rheMac2 Rhesus Jan. 2006 (MGSC Merged 1.0/rheMac2) Genome at UCSC

UCSC proCap1 Rock hyrax Jul. 2008 (Broad/proCap1) Genome at UCSC

UCSC sacCer3 S. cerevisiae Apr. 2011 (SacCer_Apr2011/sacCer3) Genome at UCSC

UCSC sacCer2 S. cerevisiae June 2008 (SGD/sacCer2) Genome at UCSC

UCSC sacCer1 S. cerevisiae Oct. 2003 (SGD/sacCer1) Genome at UCSC

UCSC strPur2 S. purpuratus Sep. 2006 (Baylor 2.1/strPur2) Genome at UCSC

UCSC strPur1 S. purpuratus Apr. 2005 (Baylor 1.1/strPur1) Genome at UCSC

UCSC aplCal1 Sea hare Sept. 2008 (Broad 2.0/aplCal1) Genome at UCSC

UCSC oviAri4 Sheep Nov. 2015 (Oar_v4.0/oviAri4) Genome at UCSC

UCSC oviAri3 Sheep Aug. 2012 (ISGC Oar_v3.1/oviAri3) Genome at UCSC

UCSC oviAri1 Sheep Feb. 2010 (ISGC Ovis_aries_1.0/oviAri1) Genome at UCSC

UCSC sorAra2 Shrew Aug. 2008 (Broad/sorAra2) Genome at UCSC

UCSC sorAra1 Shrew June 2006 (Broad/sorAra1) Genome at UCSC

UCSC choHof1 Sloth Jul. 2008 (Broad/choHof1) Genome at UCSC

UCSC speTri2 Squirrel Nov. 2011 (Broad/speTri2) Genome at UCSC

UCSC saiBol1 Squirrel monkey Oct. 2011 (Broad/saiBol1) Genome at UCSC

UCSC gasAcu1 Stickleback Feb. 2006 (Broad/gasAcu1) Genome at UCSC

UCSC tarSyr2 Tarsier Sep. 2013 (Tarsius_syrichta-2.0.1/tarSyr2) Genome at UCSC

UCSC tarSyr1 Tarsier Aug. 2008 (Broad/tarSyr1) Genome at UCSC

UCSC sarHar1 Tasmanian devil Feb. 2011 (WTSI Devil_ref v7.0/sarHar1) Genome at UCSC

UCSC echTel2 Tenrec Nov. 2012 (Broad/echTel2) Genome at UCSC

UCSC echTel1 Tenrec July 2005 (Broad/echTel1) Genome at UCSC

UCSC tetNig2 Tetraodon Mar. 2007 (Genoscope 8.0/tetNig2) Genome at UCSC

UCSC tetNig1 Tetraodon Feb. 2004 (Genoscope 7/tetNig1) Genome at UCSC

UCSC nanPar1 Tibetan frog Mar. 2015 (BGI_ZX_2015/nanPar1) Genome at UCSC

UCSC tupBel1 Tree shrew Dec. 2006 (Broad/tupBel1) Genome at UCSC

UCSC melGal5 Turkey Nov. 2014 (Turkey_5.0/melGal5) Genome at UCSC

UCSC melGal1 Turkey Dec. 2009 (TGC Turkey_2.01/melGal1) Genome at UCSC

UCSC macEug2 Wallaby Sep. 2009 (TWGS Meug_1.1/macEug2) Genome at UCSC

UCSC cerSim1 White rhinoceros May 2012 (CerSimSim1.0/cerSim1) Genome at UCSC

UCSC xenTro9 X. tropicalis Jul. 2016 (Xenopus_tropicalis_v9.1/xenTro9) Genome at UCSC

UCSC xenTro7 X. tropicalis Sep. 2012 (JGI 7.0/xenTro7) Genome at UCSC

UCSC xenTro3 X. tropicalis Nov. 2009 (JGI 4.2/xenTro3) Genome at UCSC

UCSC xenTro2 X. tropicalis Aug. 2005 (JGI 4.1/xenTro2) Genome at UCSC

UCSC xenTro1 X. tropicalis Oct. 2004 (JGI 3.0/xenTro1) Genome at UCSC

UCSC taeGut2 Zebra finch Feb. 2013 (WashU taeGut324/taeGut2) Genome at UCSC

UCSC taeGut1 Zebra finch Jul. 2008 (WUGSC 3.2.4/taeGut1) Genome at UCSC

UCSC danRer11 Zebrafish May 2017 (GRCz11/danRer11) Genome at UCSC

UCSC danRer10 Zebrafish Sep. 2014 (GRCz10/danRer10) Genome at UCSC

UCSC danRer7 Zebrafish Jul. 2010 (Zv9/danRer7) Genome at UCSC

UCSC danRer6 Zebrafish Dec. 2008 (Zv8/danRer6) Genome at UCSC

UCSC danRer5 Zebrafish July 2007 (Zv7/danRer5) Genome at UCSC

UCSC danRer4 Zebrafish Mar. 2006 (Zv6/danRer4) Genome at UCSC

UCSC danRer3 Zebrafish May 2005 (Zv5/danRer3) Genome at UCSC

NCBIEnsembl多いので保存した方がよい。多いので検索にも数分かかる。

#Ensembl
genomepy genomes -p Ensembl > list

#NCBI
genomepy genomes -p NCBI > list

 

4、ダウンロード

例えば以下のようなリストの、

$ genomepy genomes -p Ensembl

 

Ensembl Red5_PS1_1.69.0 actinidia_chinensis

Ensembl AMTR1.0 amborella_trichopoda

Ensembl ASM9120v1 cyanidioschyzon_merolae

Ensembl AUK_PRJEB4211_v1 coffea_canephora

Ensembl CcrdV1 cynara_cardunculus

(以下略)

cyanidioschyzon merolaeのゲノムアセンブリASM9120v1をダウンロードしてみる。

以下のように打つ。

genomepy install ASM9120v1 -p Ensembl

-gをつけるとパス指定。

 

追記

マウスのmm10。デフォルトではEnsemblだったが、ーpでUCSCに変更。

genomepy install mm10 -p UCSC

 

デフォルトでは~/.local/share/genomes/に保存される。開いてみる。

open ~/.local/share/genomes/

f:id:kazumaxneo:20191123183432p:plain

ゲノムのfasta以外に、bwaのpluginをonにしたため、bwaアライナーのindexファイルも作成される。他のアライナーはパスを通してないためindexが作成されなかった。

 

optional

生物名で検索してダウンロードする例。

UCSCから利用できるヒトゲノムバージョンをサーチ。

genomepy search -p UCSC human

#データベース指定なしなら(少し時間がかかる)
genomepy search human

$ genomepy search -p UCSC human

UCSC hg38 Human Dec. 2013 (GRCh38/hg38) Genome at UCSC

UCSC hg19 Human Feb. 2009 (GRCh37/hg19) Genome at UCSC

UCSC hg18 Human Mar. 2006 (NCBI36/hg18) Genome at UCSC

UCSC hg17 Human May 2004 (NCBI35/hg17) Genome at UCSC

UCSC hg16 Human July 2003 (NCBI34/hg16) Genome at UCSC

最新のhg38以外に4つ利用できる。

 

今度はアノテーションも含めてダウンロードしてみる。ソフトフィルタリングもONにする。

genomepy install -g UCSC_human_genome -m soft --annotation hg38 UCSC
  • -g <TEXT>      genome directory
  • -m <TEXT>     mask (hard or soft)
  • --annotation   download annotation

ダウンロードされた。

f:id:kazumaxneo:20191123191757p:plain

READMEを開く。

f:id:kazumaxneo:20191123191830p:plain

関連


 

 

 

 

ノンスペシャリストのための堅牢なphylogenetic analysis webサービス Phylogeny.fr

 

 系統解析は、生物学の多くの研究分野の中心であり、通常、相同な配列の同定、それらのマルチプルアライメント、系統樹再構築、および推定されたツリーのグラフィカルな表現を伴う。 Phylogeny.frプラットフォームは、これらのタスクを自動的に実行するプログラムを透過的にチェーンする。主にphylogenyの経験のない生物学者向けに設計されているが、専門家のニーズにも対応できる。最初のものは、phylogenyパイプラインにチェーンされた最新のツールを見つけて、シンプルかつ堅牢な方法でデータを分析し、専門家は洗練された分析を簡単に構築して実行できる。 Phylogeny.frには3つの主要なモードがある。 「ワンクリック」モードは、非専門家を対象とし、精度と速度を備えたすぐに使用できるパイプラインチェーンプログラムを提供する: MUSCLE for multiple alignment、PhyML for tree building、TreeDyn for tree rendering。すべてのパラメータはほとんどの研究に適合するように設定されており、ユーザーは入力シーケンスを入力するだけですぐにプリントできるツリーを取得できる。 「詳細」モードでは同じパイプラインが使用されるが、ユーザーが各プログラムのパラメーターをカスタマイズできる。 「アラカルト」モードは、ユーザーが特定のニーズに合わせてツールの幅広い選択肢から必要なステップを選択および設定することにより、独自のパイプラインを構築できるため、柔軟性と洗練性が向上する。phylogenetic analysisの前に、ユーザーは、一般的なデータベースまたは特殊なデータベースでBLASTを実行することにより、クエリシーケンスの近傍を収集することもできる。次に、ガイドツリーは、phylogeneticパイプラインの入力として使用される隣接シーケンスを選択するのに役立つ。 Phylogeny.frは、http://www.phylogeny.fr/で入手できる。

 

使用されているソフトウェアのダウンロード

http://www.phylogeny.fr/downloads.cgi

Documentation(FASTA、EMBL、またはNEXUSなどのフォーマットの解説、各ツールのフォーマット対応状況のテーブルもあり)

http://www.phylogeny.fr/documentation.cgi

 

使い方

http://www.phylogeny.fr/phylogeny.cgiにアクセスする。

f:id:kazumaxneo:20191123022635p:plain

 

ここではone clickモードの流れを確認する。

f:id:kazumaxneo:20191123022834p:plain

 

FASTA、EMBL、またはNEXUSフォーマットの配列ファイルをアップロードする。

f:id:kazumaxneo:20191123022818p:plain 

テストデータをロードした。

f:id:kazumaxneo:20191123023047p:plain

alignのタブではマルチプルシーケンスアラインメントの結果を確認でき、結果をfasta format、Phylip format、Clustal formatのいずれかでダウンロードできる。

 

Gblocksプログラムを使用して、アラインメントがpoorな領域や分岐した領域を排除する。

f:id:kazumaxneo:20191123024552p:plain

 

PhyMLを使ったphylogenetic reconstruction。Newick formatでダウンロードできる。

f:id:kazumaxneo:20191123025405p:plain

 

系統樹レンダリング

f:id:kazumaxneo:20191123023019p:plain

 

ツリーはその場で編集できる。

f:id:kazumaxneo:20191123025900p:plain

 

引用
Phylogeny.fr: robust phylogenetic analysis for the non-specialist

Dereeper A, Guignon V, Blanc G, Audic S, Buffet S, Chevenet F, Dufayard JF, Guindon S, Lefort V, Lescot M, Claverie JM, Gascuel O

Nucleic Acids Res. 2008 Jul 1;36(Web Server issue):W465-9

 

関連


 

ラージデータにも対応したマルチプルシーケンスアラインメントツール Kalign 3

 

 マルチプルシーケンスアラインメント(MSA)は、生物学的シーケンス解析の重要なタスクのままである。 MSAプログラムは、Consistency-based のメソッドとprogressive なメソッドに分けることができる。後者は、ペアワイズシーケンス距離を推定し、ガイドツリーを構築し、ガイドツリーの順序に従ってシーケンスをアラインする。Consistency-basedのメソッドは、progressive なメソッドと比較してより正確である傾向があるが、桁違いに遅いため、数千のシーケンスを整列する場合には実用的ではない。 Kalign(Lassmann et al、2008)は、一連の一般的なベンチマークデータセット上の他のアラインメントプログラムと比較して、精度と速度のバランスが優れたプログレッシブアライメント手法である(例:Sievers et al、2011を参照)。十分に熟成したにもかかわらず、Kalignは、今日頻繁に遭遇する何万ものシーケンスを処理するようには設計されていない。特に、元のKalignプログラムは、算術平均(UPGMA)アルゴリズムを使用した非加重ペアグループ法を使用して、ガイドツリーを構築し、2次時間の複雑さをもたらす。より最近のアライメントプログラムは、ガイドツリーを構築するためのヒューリスティックを実装することにより、このハードルを克服している(Blackshields et al、2010; Katoh and Toh、2006)。

 ここでは、SIMD(single instruction, multiple data)加速バージョンのGene Myersのビット並列アルゴリズム(Myers、1999)を導入して、ペアワイズシーケンス距離を推定し、Blackshields et al(2010)によって導入されたシーケンス埋め込み戦略を採用して ガイドツリーの構築を高速化するKalignの新しいバージョンを紹介する。

 Kalign2(Muth and Manber、1996)で使用されている高速文字列照合アルゴリズムを、Gene Myersの近似文字列照合アルゴリズムの新しい実装に置き換えまた。このアルゴリズムは、ビット並列命令を使用して、2つの文字列間の正確な編集距離を計算する。標準実装では、クエリの最大長はコンピューターワード(64ビットアーキテクチャでは64文字)のサイズに相当する。ただし、アルゴリズムは、最新のすべてのコンピューターで使用可能なAVXおよびAVX2命令を含むSIMD命令を使用して、さらに並列化できる。これらの命令を使用すると、長さ256のシーケンスを比較することが可能になる。AVX命令を使用したジーンマイヤーズアルゴリズムの実装はかなり単純だが、AVX命令セットにない操作があり、個別に実装する必要があった。アルゴリズムスタンドアロン実装はKalignと一緒に配布され、ダウンストリームの採用と開発を促進する。(以下略)

 

 

インストール

macos10.14でテストした。

GIthub

git clone https://github.com/TimoLassmann/kalign.git 
cd kalign
./autogen.sh
./configure
make
make check
make install

#homebrew
brew install brewsci/bio/kalign

#2019 11月現在、kalign2はcondaで導入可能
conda install -c bioconda kalign2

> kalign -h

$ kalign -h

 

Kalign (3.1.1)

 

Copyright (C) 2006,2019 Timo Lassmann

 

This program comes with ABSOLUTELY NO WARRANTY; for details type:

`kalign -showw'.

This is free software, and you are welcome to redistribute it

under certain conditions; consult the COPYING file for details.

 

Please cite:

  Lassmann, Timo.

  "Kalign 3: multiple sequence alignment of large data sets."

  Bioinformatics (2019) 

  https://doi.org/10.1093/bioinformatics/btz795

 

 

Usage: kalign  -i <seq file> -o <out aln> 

 

Options:

 

   --format           : Output format. [Fasta]

   --reformat         : Reformat existing alignment. [NA]

   --version (-V/-v)  : Prints version. [NA]

 

 

 

実行方法

配列を指定する。

#fasta formatで出力
kalign -i input.fa -f fasta -o out.fa

#clustal formatで出力
kalign -i input.fa -f clu -o out.clu

 

引用

Kalign 3: multiple sequence alignment of large datasets
Timo Lassmann
Bioinformatics, btz795, Published: 26 October 2019