ノート/ノート
訪問者数 174      最終更新 2017-04-10 (月) 09:11:37

ノート/breseq/チュートリアルを試す1 > ノート/breseq/チュートリアルを試す2 > ノート/breseq/UT-Austinの進化データを試す
ノート/breseq/Multi-Isolate1 > ノート/breseq/Multi-Isolate2

breseq Multi-Isolate その1 2017-04-01

DDBJから、multi-isolateのプロジェクトを探して、breseq処理をしてみる

出発点は、multi-isolateだと変異の同定が難しいのではないか、と言う疑問。

DDBJのBioProjectのプロパティでmonoisolateでないものを探す

DDBJのデータベースの関係は、BioProject, BioSample, DRA へのデータ登録の7枚目のスライドに図示されている。

この中の、BioProjectデータベースには、プロジェクトを登録されているが、 そのプロパティの中に、Project typeとして、Sample Scope が Monoisolate であるという情報がある。欲しいのはヘテロな環境でのNGSデータなので、Monoisolateではないものを探す。

データベースの記述によると、Multiisolateというのがあるらしい。この種類のデータを実験材料にしたい。

それで、プロジェクトを1つずつチェックしていくのも手ではあるが、面倒なので、 データベースのデータをダウンロードして、ローカルで検索することにする。

データは、

XMLなのだが、面倒ならテキストファイルとしてエディタで検索してもよい。

      <ProjectType>
        <ProjectTypeSubmission>
          <Target capture="eWhole" material="eGenome" sample_scope="eMonoisolate">
            <Organism species="139" taxID="224326">
              <OrganismName>Borrelia burgdorferi B31</OrganismName>

というような顔をしているので、 ねらい目は<Target>の中のsample_scopeがeMultiisolateであるかどうか、である。 また、更に解析対象物を知りたいので、同じ<Target>の中の<OrganismName>を取り出しておく。

pythonにxmlパーザーがあるので、ちょっとしたプログラムを書けば済む。たとえば、

#!/usr/bin/env python
# coding: UTF-8
import xml.etree.ElementTree as ET
#tree = ET.parse("bioproject.xml")
tree = ET.parse("ddbj_core_bioproject.xml")
root = tree.getroot()
for u in root.iter('Package'):
    str = ''
    for v in u.iter('ArchiveID'):
        str += v.get('accession')
    for v in u.iter('Target'):
        str += ',  ' + v.get('sample_scope')
        #text = ''
        for w in v.iter('OrganismName'):
            text = (w.text).strip()
        str += ',  ' + text
    print(str )

この出力をファイルに取ったうえで、sample_scopeのところがMultiisolateである行 だけをgrepなどで拾った。プログラム内でフィルタしてもいいのだが、全体が見たかった ので2段階処理にしてみた。

結果は、fileこのようになった。

ここから、プロジェクトを選ぼうと思う。 登録されているプロジェクトの1つに、PRJDB3087があって、 eMultiisolateであり、 Escherichia coliを対象としているので、おもしろそうだ。


PRJDB3087にあるデータをマップしてみる

DRA Searchにより、キーワードに PRJDB3087 を入れて検索すると3つの結果が得られ、その中の DRA002503.experiment.xml を見ると、3つのAccession  DRX021052 、DRX021053 、DRX021054 がある。 それぞれのリンクをクリックしてページを開くと、

が得られる。

これらの6つのデータそれぞれについてクオリティチェックをした結果と、2つずつのペアリードデータを併せてbreseqでマップした結果を、以下に置く。手順は

fastqc *.bz2
bunzip2 DRR022997_1.fastq.b2z 他のファイルも同じ
breseq -j 16 -o output -r REL606.gbk DRR022997_1.fastq DRR022997_2.fastq

breseqによるマッピング処理について、予め以下の点に疑問がある。

Trimmingをかけてみる

Trimmingの効果を見るため、trimmaticで下記の処理をする

DRR22997の結果
TrimmomaticPE: Started with arguments:
 -threads 16 -phred33 ../DRR022997_1.fastq ../DRR022997_2.fastq
  DRR022997_1_trim20.fastq DRR022997_1_unpaired.fastq
  DRR022997_2_trim20.fastq DRR022997_2_unpaired.fastq
  TRAILING:20 MINLEN:36 SLIDINGWINDOW:4:15
Input Read Pairs: 2496778 Both Surviving: 2456608 (98.39%)
  Forward Only Surviving: 34113 (1.37%) Reverse Only Surviving: 4714 (0.19%)
  Dropped: 1343 (0.05%)
TrimmomaticPE: Completed successfully
DRR22998の結果
TrimmomaticPE: Started with arguments:
 -threads 16 -phred33 ../DRR022998_1.fastq ../DRR022998_2.fastq
  DRR022998_1_trim20.fastq DRR022998_1_unpaired.fastq
  DRR022998_2_trim20.fastq DRR022998_2_unpaired.fastq
  TRAILING:20 MINLEN:36 SLIDINGWINDOW:4:15
Input Read Pairs: 1151142 Both Surviving: 1137953 (98.85%)
  Forward Only Surviving: 12861 (1.12%) Reverse Only Surviving: 137 (0.01%)
  Dropped: 191 (0.02%)
TrimmomaticPE: Completed successfully
DRR22999の結果
TrimmomaticPE: Started with arguments:
 -threads 16 -phred33 ../DRR022999_1.fastq ../DRR022999_2.fastq
  DRR022999_1_trim20.fastq DRR022999_1_unpaired.fastq
  DRR022999_2_trim20.fastq DRR022999_2_unpaired.fastq
  TRAILING:20 MINLEN:36 SLIDINGWINDOW:4:15
Input Read Pairs: 1262796 Both Surviving: 1249828 (98.97%)
  Forward Only Surviving: 12617 (1.00%) Reverse Only Surviving: 171 (0.01%)
  Dropped: 180 (0.01%)
TrimmomaticPE: Completed successfully

改めてFastQCでクオリティをチェック。上記ログのDroppedがそれほど多くない割には、FastQCの形はまあきれいになっていると言えるか。頭(左)の側の乱れはまだかなり大きいと感じる。
このトリムしたデータを、breseqでマップしてみる。


メモ

ncbiのBioProject内でST131をサーチ  36エントリ。かなりたくさんmultiisolateがある。

以下、いくつか目を引いたもの。

Global diversity and evolutionary history of Escherichia coli sequence type 131 (ST131)  これも面白いかな?
表のNumber of Linksのエントリー部分をクリックすると、データへも行けるみたい。
Links from BioProject 138このプロジェクトデータ
SRX1460182: Whole genome sequencing of E. coli ST131  1つ目のデータ例
SRR2970783
SRA-Toolkit(のfastq-dump)を使えとさ。

PRJDB4303  Comparative genomics of ESBL-producing Escherichia coli ST131 isolates

PRJDB3868  Whole genome sequencing of Escherichia coli ST131

PRJNA311313  Escherichia coli JJ1887 Genome sequencing PRJNA307507  Escherichia coli Genome sequencing and assembly Comparative sequencing of ExPEC ST131


添付ファイル: fileacc_multi.list 50件 [詳細]

トップ   編集 凍結 差分 バックアップ 添付 複製 名前変更 リロード   新規 一覧 単語検索 最終更新   ヘルプ   最終更新のRSS
Last-modified: 2017-04-10 (月) 09:11:37 (195d)