Help us understand the problem. What is going on with this article?

Bioconductor解説 Biostrings編 その3

過去記事

前アカウントで書いた記事です。

Bioconductor 解説 Biostrings 編 その1

Bioconductor 解説 Biostrings 編 その2

もう前記事書いたの2年前ですね。。。責任をもって終わらせます。

はじめに

今回は、前回紹介したBString、DNAString、RNAString、AAStringオブジェクト(XStringオブジェクト)をリスト化したオブジェクトである BStringSet, DNAStringSet, RNAStringSet, AAStringSetオブジェクト(XStringSetオブジェクト)を紹介したいと思います。
そして次の回で前回予告したこれらオブジェクトを引数として要求する諸関数についてみていきたいと思います。

BStringSet, DNAStringSet, RNAStringSet, AAStringSet

一つの配列情報を記録するXString オブジェクトですが、時として一つでなく複数の配列を一度に処理する場面もありえます。そういった場合に一括処理するオブジェクトがこのXStringSetオブジェクトです。

オブジェクト(インスタンス)を作成するコンストラクタは以下のとおりです。

BStringSetコンストラクタ
BStringSet(x=c( "ABCDEFGHIJK", "LMNOPQRSTUVWXYZ" ) ) 
DNAStringSetコンストラクタ
DNAStringSet( x=c( "ATCGGCAAT", "TTGTACCC" ) ) 
RNAStringSetコンストラクタ
RNAStringSet( x=c( "AUCGGCAAU", "UUGUACCC" ) ) 
AAStringSetコンストラクタ
AAStringSet( x=c( "MARKSLEMSIR", "SRRRIKSMLE" ) )

出力されるインスタンスは以下の通りです。

(出力)BStringSetコンストラクタ
  A BStringSet instance of length 2
    width seq
[1]    11 ABCDEFGHIJK
[2]    15 LMNOPQRSTUVWXYZ
(出力)DNAStringSetコンストラクタ
  A DNAStringSet instance of length 2
    width seq
[1]     9 ATCGGCAAT
[2]     8 TTGTACCC
(出力)RNAStringSetコンストラクタ
  A RNAStringSet instance of length 2
    width seq
[1]     9 AUCGGCAAU
[2]     8 UUGUACCC
(出力)AAStringSetコンストラクタ
  A AAStringSet instance of length 2
    width seq
[1]    11 MARKSLEMSIR
[2]    10 SRRRIKSMLE

引数は以下のようになっています。

x:{文字列型を要素としてもつベクトル}

:配列データのベクトル

start:{整数}

:読み込み開始地点(負の数で末尾からの参照)

end:{整数}

:読み込み終了地点

width:{整数}

:読み込む文字の数

use.names:{論理値}

:TRUEの場合、ベクトルの要素名をそのまま使用する

ほとんどXStringオブジェクトの場合と同じです。
ただxがベクトル型であることと、use.namesで配列名を使用できる点が異なっています。

引数について

まず引数xに関して詳しくみていきましょう。
さっきも言いましたが、XStringオブジェクトでは引数xは文字列型(例:"xxx")で、XStringSetオブジェクトでは引数xは文字列ベクトル(例:c( "xxx", "yyy" ))となっています。

配列の数に関してはメモリがある限りはどれだけでも格納できそうです。

試しに1千万個の配列を格納してみます。

BStringSet( x=rep("A", 1E+7) )
  A BStringSet instance of length 10000000
           width seq
       [1]     1 A
       [2]     1 A
       [3]     1 A
       [4]     1 A
       [5]     1 A
       ...   ... ...
 [9999996]     1 A
 [9999997]     1 A
 [9999998]     1 A
 [9999999]     1 A
[10000000]     1 A

1千億個の配列だと私のPCのメモリ限界でした。

BStringSet( x=rep("A", 1E+10) )
Error: vector memory exhausted (limit reached?)

ちなみにベクトル内の要素が複数である必要はありません。
以下は単一の配列を入力した例です。

DNAStringSetコンストラクタ
DNAStringSet( x=c( "ATCGGCAAT") ) 
  A DNAStringSet instance of length 1
    width seq
[1]     9 ATCGGCAAT

このような性質からある程度はXStringオブジェクトの代替として使うことができます。たとえばXStringオブジェクトでは配列に名前をつけて管理できるようにはなっていないので、名前をつけた状態で管理するためにXStringSetオブジェクトを使う例が考えられます。

名前入りで登録する例
DNAStringSet( x=c( my.seq="ATCGGCAAT") ) 
  A DNAStringSet instance of length 1
    width seq       names               
[1]     9 ATCGGCAAT my.seq

逆に名前を消して登録する場合は、use.names=FALSEにします。

名前入りで登録する例
DNAStringSet( x=c( my.seq="ATCGGCAAT"), use.names=FALSE ) 
  A DNAStringSet instance of length 1
    width seq
[1]     9 ATCGGCAAT

ちなみに後付けで命名したい場合はnames()を使えばできます。

後で配列名を登録する例
my.seq <- DNAStringSet( x=c( "ATCGGCAAT", "CTT" ) ) 
message("配列の名前を表示します。")
names(my.seq) #名前がない状態

message("命名します。")
names(my.seq) <- c( "new.seq1", "new.seq2" ) # 名前を文字列ベクトルで登録

message("配列の名前を表示します。")
names(my.seq) #名前がある状態
配列の名前を表示します。
NULL
命名します
配列の名前を表示します。
[1] "new.seq1" "new.seq2"

この名前の登録の仕方は、ベクトル型やリスト型で命名する時と同じですね。
Rになれない人は、関数(<オブジェクト>) <- "<名前>"という形に違和感を覚えるかもしれません。Rで要素名や列名を記載する際はこのような記法でSetter的な動きをします。

XStringSetのメソッド・その他操作に関して

length()に気を付ける

一部のXStringSetオブジェクトのメソッドは、XStringオブジェクトのメソッドとよく似ています。ただ微妙に仕様が違う点に注意してください。その一例としてlength()があります。

XStringオブジェクトにlength()を適用した場合は、配列長を返します。

length( DNAString("ATCGGCAAT") )
[1] 9

しかし、XStringSetオブジェクトにlength()を適用した場合は、配列の数を返します。

length( DNAStringSet( x=c( A="ATCGGCAAT", B="CTT" ) ) )
[1] 2

XStringSetオブジェクトの各配列長を参照する場合は、nchar()もしくはwidth()を利用します。

nchar( DNAStringSet( x=c( A="ATCGGCAAT", B="CTT" ) ) )
[1] 9 3

ちなみにXStringオブジェクトの場合も、nchar()で配列長を参照できます。ややこしいので「配列長の参照はnchar()だけを使う」というように決めておくほうが混乱しなくてすむかもしれませんね。

nchar( DNAString("ATCGGCAAT") )
[1] 9

subseq()で一括切り出し

XStringSetを使うことの一番わかりやすい利点は一括で配列を処理できる点にあります。
次の例ではsubseq()を使って1番目から3番目までの配列を切り出しています。

ori.seq <- DNAStringSet( x=c( A="ATCGGCAAT", B="CTT" ) )
subseq( x=ori.seq, start=1, end=3 )
  A DNAStringSet instance of length 2
    width seq   names               
[1]     3 ATC   A
[2]     3 CTT   B

この処理は高速ですので、たとえ一千万個程度の配列でも一秒程度で処理できます。

numerous.seq <- BStringSet( x=rep("ATTTCGCGTTTGCGATGA", 1E+7) )
subseq( x=numerous.seq, start=1, end=15 )
  A BStringSet instance of length 10000000
           width seq
       [1]    15 ATTTCGCGTTTGCGA
       [2]    15 ATTTCGCGTTTGCGA
       [3]    15 ATTTCGCGTTTGCGA
       [4]    15 ATTTCGCGTTTGCGA
       [5]    15 ATTTCGCGTTTGCGA
       ...   ... ...
 [9999996]    15 ATTTCGCGTTTGCGA
 [9999997]    15 ATTTCGCGTTTGCGA
 [9999998]    15 ATTTCGCGTTTGCGA
 [9999999]    15 ATTTCGCGTTTGCGA
[10000000]    15 ATTTCGCGTTTGCGA

配列の抽出・変換

XStringSetオブジェクトにはいった配列の一部を取り出したい場合は、[<インデックスのベクトル>]で取り出し可能です。取り出された配列は同じくXStringSetオブジェクトとして生成されます。

抽出(DNAStringSetからDNAStringSetへの変換)
ori.seq <- DNAStringSet( x=c( A="ATCGGCAAT", B="CTT", C="AATCGGGTGT" ) )
sub.seq <- ori.seq[ c(1, 3) ] # 1番目と3番目の配列を取り出す。
sub.seq
  A DNAStringSet instance of length 2
    width seq        names               
[1]     9 ATCGGCAAT  A
[2]    10 AATCGGGTGT C

配列をXStringオブジェクトとして取り出したい場合は、[[<インデックス>]]とすればできます。
ちなみにこの場合で取り出せる配列は一つだけという点には注意してください。

抽出(DNAStringSetからDNAStringへの変換)
ori.seq <- DNAStringSet( x=c( A="ATCGGCAAT", B="CTT", C="AATCGGGTGT" ) )
one.seq <- ori.seq[[ 1 ]] # 1番目の配列を取り出す。
one.seq
  9-letter "DNAString" instance
seq: ATCGGCAAT

文字列ベクトルとして取り出したい場合は、as.character()が使えます。

抽出(DNAStringSetから文字列ベクトルへの変換)
ori.seq <- DNAStringSet( x=c( A="ATCGGCAAT", B="CTT", C="AATCGGGTGT" ) )
ori.vec <- as.character( ori.seq )
ori.vec
           A            B            C 
 "ATCGGCAAT"        "CTT" "AATCGGGTGT" 

配列のマトリックスとして表示したい場合は、as.matrix()が使えます。しかしこの場合はXStringSetオブジェクト内の全ての配列長が同じである必要があります。

抽出(DNAStringSetから配列のマトリックスへの変換)
ori2.seq <- DNAStringSet( x=c( A="ATC", B="CTT", C="AAT" ) )
ori.mat <- as.matrix( ori2.seq )
ori.mat
  [,1] [,2] [,3]
A "A"  "T"  "C" 
B "C"  "T"  "T" 
C "A"  "A"  "T" 

このように手軽に変換できるので、処理方法に策定に行き詰まったら一度変換して処理することも考えてみてください。ただ文字列ベクトル等の基本型で処理した場合はやや処理速度が遅くなると思うので、大規模な配列を処理する場合はまず専用のパッケージがBioconductor内に存在しないか探してみるのが良いかと思います。

さいごに

XStringSetオブジェクトを使うことで複数の配列を能率的に処理することができます。ぜひ利用してみてください。
次回はXStringオブジェクトとXStringSetオブジェクトを引数とする諸関数を紹介していきたいと思います。

KazukiNakamae
プラチナバイオ株式会社 主任研究員 広島大学 研究員兼任 専門分野:ゲノム編集のツール開発・解析
https://scholar.google.co.jp/citations?user=Wje3pk8AAAAJ
Why not register and get more from Qiita?
  1. We will deliver articles that match you
    By following users and tags, you can catch up information on technical fields that you are interested in as a whole
  2. you can read useful information later efficiently
    By "stocking" the articles you like, you can search right away