【第13回】ゼロから始めるゲノム解析 （Python編） Location Restriction Sites: Using, Testing, and Sharing Code @kimoton

本勉強会の概要・目的 書籍名 対象者/目的 Mastering Python for Bioinformatics Python・バイオインフォ知識ほぼゼロの人 を対象に、正しいPythonのコーディング手 法について学ぶ 頻度 毎週〜隔週開催予定 登壇者 募集中！

Rosalindとは ● 問題解決を通じてバイオインフォマティク ス、プログラミング、およびアルゴリズムを 学習するためのプラットフォーム ● 大学やハッカソン、就職の面接にも 600回 以上の採用実績あり 参考：https://qiita.com/_kimoton/items/d534d0fa9b83dd7dc412 概要

環境構築 - 必要パッケージ群のインストール # 公開されているレポジトリからファイル群を取得 $ git clone https://github.com/kyclark/biofx_python $ cd biofx_python # requirements.txt に記載のパッケージをインストール $ pip3 install -r requirements.txt # pylintの設定ファイルをホームディレクトリに移動 $ cp pylintrc ~/.pylintrc # mypyの設定ファイルをホームディレクトリに移動 $ cp mypy.ini ~/.mypy.ini

本日のお題 4～12の長さを持つすべての回文配列の位置と長さを出力せよ https://rosalind.info/problems/revp/

● 回文配列を見つける方法 ● 汎用的な関数をシェアするためのモジュールの作成方法 ● PYTHONPATH環境変数について 本日学ぶこと

前提知識編

回文配列 https://en.wikipedia.org/wiki/Palindromic_sequence 回文配列とは、二本鎖DNA分子またはRNA分子内の核酸配列であり、片方の一本鎖における特定の方向での塩 基配列の読み取りが、もう片方の相補鎖における同じ方向での読み取りと一致する配列 AGATTCGAATCT TCTAAGCTTAGA 1. Complement 2. Reverse Complement 5’ 3’ 3’ 5’ >>> from Bio import Seq >>> seq = 'GCATGC' >>> Seq.reverse_complement(seq) == seq True

回文配列は制限酵素の認識領域 https://en.wikipedia.org/wiki/Palindromic_sequence 制限酵素は制限部位として知られるDNAの特定の回文配列を認識し、その中で切断することが知られている 制限酵素はもともとバクテリアが ファージなどの感染から自身を 守るために外来DNAを切断する 役割

問題のIN/OUT TCAATGCATGCGGGTCTATATGCAT $ ./revp.py tests/inputs/1.fa 5 4 7 4 17 4 18 4 21 4 4 6 6 6 20 6 与えられた配列に含まれている回文配列の位置と長さを出力することが求められている

（振り返り）k-merの抽出 配列から全てのk-mer（k文字の部分配列）を抽出する関数は以下のように書ける >>> def find_kmers(seq, k): ... n = len(seq) - k + 1 ... return [] if n < 1 else [seq[i:i + k] for i in range(n)] ...

部分配列を抽出 find_kmers関数を使うと、指定したkの範囲それぞれについて部分配列を抽出できる # k=4～k=12のkmerがいくつ含まれているか出力 >>> for k in range(4, 13): ... print(k, len(find_kmers(seq, k))) ... 4 22 5 21 6 20 7 19 8 18 9 17 10 16 11 15 12 14

部分配列の逆相補鎖配列を抽出 取得した部分配列に対してreverse_complementメソッドを適用することで逆相補鎖配列を抽出する # find_kmersは文字列のリストで返す >>> kmers = find_kmers(seq, 12) >>> kmers ['TCAATGCATGCG', 'CAATGCATGCGG', 'AATGCATGCGGG', 'ATGCATGCGGGT', 'TGCATGCGGGTC', 'GCATGCGGGTCT', 'CATGCGGGTCTA', 'ATGCGGGTCTAT', 'TGCGGGTCTATA', 'GCGGGTCTATAT', 'CGGGTCTATATG', 'GGGTCTATATGC', 'GGTCTATATGCA', 'GTCTATATGCAT'] # map関数で全ての要素に関数を適用 >>> revc = list(map(Seq.reverse_complement, kmers))

モジュールの活用 再利用可能な処理に関しては、関数化して分離しておく。 使用する際は 「from ファイル名 import 関数名」のように使う >>> from common import find_kmers >>> kmers = find_kmers(seq, 12) >>> kmers ['TCAATGCATGCG', 'CAATGCATGCGG', 'AATGCATGCGGG', 'ATGCATGCGGGT', 'TGCATGCGGGTC', 'GCATGCGGGTCT', 'CATGCGGGTCTA', 'ATGCGGGTCTAT', 'TGCGGGTCTATA', 'GCGGGTCTATAT', 'CGGGTCTATATG', 'GGGTCTATATGC', 'GGTCTATATGCA', 'GTCTATATGCAT'] common.py

モジュールの呼び先 - PYTHONPATH 環境変数PYTHONPATHに検索対象のパスを指定することで、デフォルトの検索パスに加えて指定したパスを検索 するように変更することができる # 例：$HOME/.local/binを追加 (bash) $ export PYTHONPATH="$HOME/.local/bin:$PYTHONPATH" https://docs.python.org/ja/3/using/cmdline.html#envvar-PYTHONPATH

解法編

解法１zip()関数とenumerate()関数を使用 1 指定した範囲の部分配列を抽出 1 2 部分配列を逆相補鎖に変換 2 3 部分配列とその逆相補鎖が一致して いれば出力 3

解法 2 .eq()関数を使用 1 指定した範囲の部分配列を抽出 1 2 部分配列を逆相補鎖に変換 2 3 部分配列とその逆相補鎖が一致して いれば出力 3

解法 3 revp()関数を作成する 1 与えられた文字列からk文字の回文配 列を抽出する関数を作成 2 指定したkの範囲で回文配列の検索を 行い出力 2 1

追加課題 プログラムを拡張して探索する回文配列の長さ（ k）をコマンドラインパラメータとして渡せるようにし て下さい。