Upgrade to Pro — share decks privately, control downloads, hide ads and more …

Finding the Hamming Distance: Counting Point Mutations

onouyek
January 07, 2022

Finding the Hamming Distance: Counting Point Mutations

onouyek

January 07, 2022
Tweet

More Decks by onouyek

Other Decks in Programming

Transcript

  1. 【第6回】ゼロから始めるゲノム解析
    (Python編)
    Finding the Hamming Distance: Counting Point Mutations
    @onouyek

    View full-size slide

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

    View full-size slide

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

    View full-size slide

  4. 環境構築 - 必要パッケージ群のインストール
    # 公開されているレポジトリからファイル群を取得
    $ 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

    View full-size slide

  5. 本日のお題 http://rosalind.info/problems/hamm/
    DNA配列sとtのハミング距離を求めよ

    View full-size slide

  6. 前提知識編

    View full-size slide

  7. ハミング距離
    文字列sからtに変換するのに必要な編集(置換)の数
    AAACCCGGGTTT
    || |
    CGACGATATGTC
    AAACCCGGGTTT
    |||||||||||
    AACCCGGGTTTA
    例1:sとtで共通のものは12塩基のうち3塩基なのでハミング距離は 9
    例2:sとtの配列はほぼ共通( 92%)だが開始位置から比べるのでハミング距離は 4
    AAACCCGGGTTT
    || || || ||
    AACCCGGGTTTA
    s
    t
    s
    t
    s
    t

    View full-size slide

  8. 位置引数の複数指定
    位置引数は順番があるのでargparseで指定した順番にして渡す
    $ ./hamm.py -h
    usage: hamm.py [-h] str str
    Hamming distance
    positional arguments:
    str Sequence 1
    str Sequence 2
    optional arguments:
    -h, --help show this help message and exit

    View full-size slide

  9. abs():絶対値を求める
    >>> seq1, seq2 = 'AC', 'ACGT'
    >>> len(seq1) - len(seq2)
    -2
    2つの配列が異なる場合、差をとると負になる場合がある
    >>> distance = abs(len(seq1) - len(seq2))
    2
    abs()で絶対値を求めることができる

    View full-size slide

  10. min():最小値を求める
    range()の引数に求めた最小値を渡せばインデックスを使って文字を比較できる
    >>> seq1, seq2 = 'AC', 'ACGT'
    >>> min(len(seq1), len(seq2))
    2
    2つの配列を比べる場合、短い方の長さだけ調べれば良い
    min()で最小値が求められる
    >>> for i in range(min(len(seq1), len(seq2))):
    ... print(seq1[i], seq2[i])
    ...
    A A
    C C

    View full-size slide

  11. Solution 1: min()
    ①seq1とseq2の長さの差の絶対値をdistanceに代入
    ②seq1とseq2の短い方の長さでイテレーション
    ③seq1とseq2の同じ位置の文字を比較して異なればdistanceをインクリメント



    View full-size slide

  12. Solution 2: zip()
    min()で短い方の長さを求めなくても、zip()で短い方に合わせて2つの
    リストからそれぞれの要素を取り出せる
    >>> list(zip('AC', 'ACGT'))
    [('A', 'A'), ('C', 'C')]

    View full-size slide

  13. Solution 3: zip_longest()
    zip_longest()で長い方に合わせて2つのリストからそれぞれの要素を取り出せる
    短い方は要素がなくなるとNoneが返されるのでabs()で初期値を求める必要がなくなる
    >>> from itertools import zip_longest
    >>> list(zip_longest('AC', 'ACGT'))
    [('A', 'A'), ('C', 'C'), (None, 'G'), (None, 'T')]

    View full-size slide

  14. Solution 4: リスト内包表記
    リスト内包表記で1行で書くことができる
    >>> seq1, seq2 = 'GAGCCTACTAACGGGAT', 'CATCGTAATGACGGCCT'
    >>> [1 if c1 != c2 else 0 for c1, c2 in zip_longest(seq1, seq2)]
    [1, 0, 1, 0, 1, 0, 0, 1, 0, 1, 0, 0, 0, 0, 1, 1, 0]
    demo

    View full-size slide

  15. Solution 5: filter()
    filter()とlambda関数を使って2つの配列の異なる要素のイテレーターが得られる
    >>> seq1, seq2 = 'AC', 'ACGT'
    >>> list(filter(lambda t: t[0] != t[1], zip_longest(seq1, seq2)))
    [(None, 'G'), (None, 'T')]

    View full-size slide

  16. Solution 6: map()
    map()とlambda関数を使って2つの配列の要素を比較したbool値のイテレーターが得られる
    map(lambda t: t[0] != t[1], zip_longest('AC', 'ACGT')))
    map(lambda t: t[0] != t[1], [('A', 'A'), ('C', 'C'), (None, 'G'), (None, 'T')])
    lambda('A', 'A'): 'A' != 'A' -> False
    lambda('C', 'C'): 'C' != 'C' -> False
    lambda(None, 'G'): None != 'G' -> True
    lambda(None, 'T'): None != 'T' -> True

    View full-size slide

  17. Solution 7: starmap()
    starmap()を使うことで、関数にタプルを展開して渡してくれる
    def not_same(t):
    return t[0] != t[1]
    starmap(not_same, zip_longest('AC', 'ACGT'))
    starmap(lambda a, b: a != b, [('A', 'A'), ('C', 'C'), (None, 'G'), (None, 'T')])
    *('A', 'A')
    lambda 'A', 'A': 'A' != 'A' -> False
    demo

    View full-size slide

  18. operatorモジュール:関数形式の標準演算子
    >>> import operator
    >>> operator.ne('A', 'A')
    False
    >>> operator.ne('A', 'T')
    True
    >>> 1 + 2
    3
    >>> operator.add(1, 2)
    3
    >>> 'AC' + 'GT'
    'ACGT'
    >>> operator.concat('AC', 'GT')
    'ACGT'

    View full-size slide

  19. 追加課題
    プログラムを拡張して、2つ以上の入力配列を扱えるようにしてください。各配列のペア
    間のハミング距離を表示するようにしてください。

    View full-size slide

  20. 本日学んだこと
    ● abs()
    ● min()
    ● zip()
    ● zip_longest()
    ● filter()
    ● map()
    ● starmap()
    ● operatorモジュール

    View full-size slide