[바이오파이썬] 03. Sequence 객체

안녕하세요

한주현입니다.

이번 강좌에서는 바이오파이썬의 기본이 되는 Sequence 객체에 대해서 알아보고자 합니다.

이전 강좌 (02. 바이오파이썬으로 할 수 있는 일 들) 에서 잠깐 보여드렸듯이

Sequence 객체는 다음과 같이 만들 수 있습니다.

>>> from Bio.Seq import Seq
>>> my_seq = Seq("AGTACACTGGT")
>>> my_seq
Seq('AGTACACTGGT', Alphabet())
>>> my_seq.alphabet
Alphabet()
>>> print(my_seq)
AGTACACTGGT
>>> print(my_seq.lower())
agtacactggt
>>> len(my_seq)
11 
>>>

위의 예제에서 보는 것 처럼 Seq 을 통해 Sequence 객체를 만들었습니다.

그런데 Seq 객체를 print 하거나 .lower()  또는 len() 처럼 파이썬 String 자료형과 같은 형태가 보입니다.

이전 강좌에서 말씀 드렸 듯 Seq객체의 method는 다음과 같으며 파이썬의 String 자료형과 같은 부분도 있으면서,

Seq객체에만 있는 메소드 (complement, reverse_complement 등) 들도 있습니다.

>>> from Bio.Seq import Seq
>>> my_seq = Seq("AGTACACTGGT")
>>> dir(my_seq)
['__add__', '__class__', '__contains__', '__delattr__', '__dict__', '__dir__', '__doc__', '__eq__', '__format__',

 '__ge__', '__getattribute__', '__getitem__', '__gt__', '__hash__', '__init__', '__init_subclass__', '__le__', 

'__len__', '__lt__', '__module__', '__ne__', '__new__', '__radd__', '__reduce__', '__reduce_ex__', '__repr__', 

'__setattr__', '__sizeof__', '__str__', '__subclasshook__', '__weakref__', '_data', 

'_get_seq_str_and_check_alphabet', 'alphabet', 'back_transcribe', 'complement', 'count', 'endswith', 'find', 

'lower', 'lstrip', 'reverse_complement', 'rfind', 'rsplit', 'rstrip', 'split', 'startswith', 'strip', 'tomutable',

'tostring', 'transcribe', 'translate', 'ungap', 'upper']

Sequence 와 Alphabet

Sequence는 파이썬의 String과 비슷하면서도 추가적인 메소드 들이 있다고 말씀드렸습니다.

그리고 아래의 예제에서 보듯 Alphabet이란게 자꾸 나오는데 이것은 무엇을 의미하는 것 일까요?

>>> from Bio.Seq import Seq
>>> my_seq = Seq("AGTACACTGGT")
>>> my_seq
Seq('AGTACACTGGT', Alphabet())
>>> my_seq.alphabet 

Alphabet()

Alphabet 객체는 Sequence 객체가 어떠한 객체를 담고 있는지 정의해줍니다.

예를 들어 Sequence 객체 안에 담고 있는 서열이 DNA 이면 Alphabet에 DNA라고 선언을 해주고,

RNA라면 RNA라고 선언을 해주는 것 입니다.

선언을 하는 방식은 IUPAC alphabet을 따릅니다. (http://www.chem.qmul.ac.uk/iupac/)

기본 IUPACProtein class를 예로 들자면 Alanine 은 'A', Arginine 은 'R' 인 것 이지요.

우리가 평소에 보던 것과 같습니다.

또한 추가적으로 Extended IUPACProtein class 를 소개해드리겠습니다.

U	Sec	selenocysteinie
O	Pyl	pyrrolysine
B	Asx	asparagine 또는 aspartic acid
Z	Glx	glutamine 또는 glutamic acid
J	Xle	leucine isoleucine
X	Xxx	unknown 아미노산

IUPAC의 메소드에는 다음과 같습니다.

ambiguous_dna	IUPACAmbiguousDNA
ambiguous_rna	IUPACAmbiguousRNA
extended_dna	ExtendedIUPACDNA
extended_protein	ExtendedIUPACProtein
protein	IUPACProtein
unambiguous_dna	IUPACUnambiguousDNA
unambiguous_rna	IUPACUnambiguousRNA
	IUPACData

Sequence 객체를 만들 때는,

1. 기본 alphabet

>>> from Bio.Seq import Seq
>>> my_seq = Seq("AGTACACTGGT")
>>> my_seq
Seq('AGTACACTGGT', Alphabet())
>>> my_seq.alphabet
Alphabet

2. IUPAC 지정 

>>> from Bio.Seq import Seq

>>> from Bio.Alphabet import IUPAC
>>> my_seq = Seq('AGTACACTGGT', IUPAC.unambiguous_dna)
>>> my_seq
Seq('AGTACACTGGT', IUPACUnambiguousDNA())
>>> my_seq.alphabet
IUPACUnambiguousDNA()

파이썬의 String 과 비슷한 Sequence

앞서 파이썬의 자료구조 중 String과 Sequence가 비슷하다고 말씀드렸습니다.

아래의 예시를 보시죠.

Sequence 객체의 iterable, len()

>>> from Bio.Seq import Seq
>>> from Bio.Alphabet import IUPAC
>>> my_seq = Seq("GATCG", IUPAC.unambiguous_dna)
>>> for index, letter in enumerate(my_seq):
...     print(index, letter)
... 
0 G
1 A
2 T
3 C
4 G
>>> print(len(my_seq)) 

5

String의 iterable 함과 길이를 반환하는 함수인 len() 가 Seq 객체에서도 확인하실 수 있습니다.

Indexing 도 가능합니다.

Sequence 객체의 indexing 과 slicing

>>> print(my_seq[0])
G
>>> print(my_seq[2])
T
>>> print(my_seq[-1])
G
>>> print(my_seq[1:3])
AT

Sequence 객체의 count 메소드

>>> "AAAA".count("AA")
2

 >>> Seq("AAAA").count("AA")

2

-> Count 메소드는 String 메소드의 count와 같은 결과를 보여주며, 상황에 따라 위의 예제의 결과값으로 오버랩 되도록 3이 나와야 할 때는 별도로 스크립트를 작성해주어야 합니다.

GC contents 계산 - 1 (Count 메소드를 사용)

>>> from Bio.Seq import Seq

>>> from Bio.Alphabet import IUPAC
>>> my_seq = Seq('GATCGATGGGCCTATATAGGATCGAAAATCGC', IUPAC.unambiguous_dna)
>>> len(my_seq)
32
>>> my_seq.count("G")
9
>>> 100 * float(my_seq.count("G") + my_seq.count("C")) / len(my_seq)
46.875

-> GC(%)의 계산은 이렇게 count와 전체 길이를 계산하여 계산 할 수 있습니다.

GC contents 계산 - 2 (Bio.SeqUtils.GC() 사용)

>>> from Bio.Seq import Seq

>>> from Bio.Alphabet import IUPAC

>>> from Bio.SeqUtils import GC

>>> my_seq = Seq('GATCGATGGGCCTATATAGGATCGAAAATCGC', IUPAC.unambiguous_dna)

>>> GC(my_seq)

46.875

Sequence 객체를 String으로 변환

Sequence 객체를 String 으로 바꾸는 것은 매우 쉽습니다.

파이썬 내장함수인 str() 을 사용하시면 됩니다.

>>> from Bio.Seq import Seq
>>> from Bio.Alphabet import IUPAC
>>> from Bio.SeqUtils import GC
>>> my_seq = Seq('GATCGATGGGCCTATATAGGATCGAAAATCGC', IUPAC.unambiguous_dna)
>>> str_seq = str(my_seq)
>>> print(str_seq)
GATCGATGGGCCTATATAGGATCGAAAATCGC

Sequence 객체 더하기

Sequence 객체끼리 더하여 문자열의 합을 만들 수 있습니다.

>>> from Bio.Seq import Seq
>>> from Bio.Alphabet import generic_dna
>>> my_seq1 = Seq("ACGT", generic_dna)
>>> my_seq2 = Seq("AACC", generic_dna)
>>> my_seq3 = my_seq1 + my_seq2
>>> my_seq3 

Seq('ACGTAACC', DNAAlphabet())

Sequence에는 alphabet 객체가 있다고 말씀드렸습니다.

만약 alphabet 객체가 다른 종류끼리 더하면 어떻게 될까요? 오류가 나려나요? 아니면 그냥 생각없이 더해지려나요?

>>> from Bio.Alphabet import IUPAC
>>> from Bio.Seq import Seq
>>> protein_seq = Seq("EVRNAK", IUPAC.protein)
>>> dna_seq = Seq("ACGT", IUPAC.unambiguous_dna)
>>> protein_seq + dna_seq
Traceback (most recent call last):
  File "<stdin>", line 1, in <module>
  File "/anaconda/lib/python3.6/site-packages/Bio/Seq.py", line 292, in __add__
    self.alphabet, other.alphabet)) 

TypeError: Incompatible alphabets IUPACProtein() and IUPACUnambiguousDNA()

-> 바이오파이썬은 생각없이 Sequence객체를 더하지는 않는군요.

보시는 것 처럼 IUPAC.protein + IUPAC.unambiguous_dna 를 더하니 오류가 납니다.

그러면 generic_nucleotide 와 IUPAC.unambiguous_dna를 더해보면 어떨까요?

>>> from Bio.Seq import Seq

>>> from Bio.Alphabet import generic_nucleotide
>>> from Bio.Alphabet import IUPAC
>>> nuc_seq = Seq("GATCGATGC", generic_nucleotide)
>>> dna_seq = Seq("ACGT", IUPAC.unambiguous_dna)
>>> nuc_seq
Seq("GATCGATGC", NucleotideAlphabet())
>>> dna_seq
Seq("ACGT", IUPACUnambiguousDNA())
>>> nuc_seq + dna_seq
Seq("GATCGATGCACGT", NucleotideAlphabet())

-> generic_nucleotide와 IUPAC.unambiguous_dna 는 잘 더해짐을 확인하실 수 있습니다.

대소문자 변환

Sequence 객체에서 대소문자 변환은 String 구조와 마찬가지로 upper, lower 메소드를 사용하시면 됩니다.

>>> from Bio.Seq import Seq
>>> from Bio.Alphabet import generic_dna
>>> dna_seq = Seq("acgtACGT", generic_dna)
>>> dna_seq
Seq(’acgtACGT’, DNAAlphabet())
>>> dna_seq.upper()
Seq(’ACGTACGT’, DNAAlphabet())
>>> dna_seq.lower() 

Seq(’acgtacgt’, DNAAlphabet())

Sequence의 상보, 역상보 서열 (complement, reverse complement)

Sequence객체의 nucleotide sequence에서 complement나 reverse complement  를 한 줄이면 구할 수 있습니다.

메소드를 사용하면 되는데요 아래의 예제를 보시죠.

>>> from Bio.Seq import Seq

>>> from Bio.Alphabet import IUPAC
>>> my_seq = Seq("GATCGATGGGCCTATATAGGATCGAAAATCGC", IUPAC.unambiguous_dna)
>>> my_seq
Seq("GATCGATGGGCCTATATAGGATCGAAAATCGC", IUPACUnambiguousDNA())
>>> my_seq.complement()
Seq("CTAGCTACCCGGATATATCCTAGCTTTTAGCG", IUPACUnambiguousDNA())
>>> my_seq.reverse_complement()
Seq("GCGATTTTCGATCCTATATAGGCCCATCGATC", IUPACUnambiguousDNA())

Sequence객체에 complement, reverse_complement 메소드를 사용하니 바로 상응하는 결과가 나왔습니다.

>>> my_seq.complement()[::-1] 

Seq("GCGATTTTCGATCCTATATAGGCCCATCGATC", IUPACUnambiguousDNA())

Sequence객체는 슬라이싱이 된다고 말씀드렸듯이 [::-1] 을 하여 complement 서열을 reverse_compelment 로 바꿀 수도 있습니다.

단백질 서열의 상보 또는 역상보 서열은 어떻게 될까요?

>>> from Bio.Seq import Seq
>>> from Bio.Alphabet import IUPAC
>>> protein_seq = Seq("EVRNAK", IUPAC.protein)
>>> protein_seq.complement()
Traceback (most recent call last):
  File "<stdin>", line 1, in <module>
  File "/anaconda/lib/python3.6/site-packages/Bio/Seq.py", line 760, in complement
    raise ValueError("Proteins do not have complements!") 

ValueError: Proteins do not have complements!

네, 상식적으로 생각하였을 때 단백질 서열에는 상보, 역상보 서열이 없지요.

그래서 오류가 나게 됩니다.

전사 (Transcription)

전사(Transcription) 에 대하여 말씀드리기 앞서 우선 strand에 대해 이야기 해보겠습니다.

5' - ATGGCCATTGTAATGGGCCGCTGAAAGGGTGCCCGATAG - 3'  (DNA coding strand, +1 strand)

3' - TACCGGTAACATTACCCGGCGACTTTCCCACGGGCTATC - 5'  (DNA template strand, -1 strand)

-- 전사 --> 

5' - AUGGCCAUUGUAAUGGGCCGCUGAAAGGGUGCCCGAUAG - 3'  (Single stranded messenger RNA)

우리는 coding strand를 입력값으로 받아서 template strand 와  mRNA를 만들어 보겠습니다.

>>> from Bio.Seq import Seq
>>> from Bio.Alphabet import IUPAC
>>> coding_dna = Seq("ATGGCCATTGTAATGGGCCGCTGAAAGGGTGCCCGATAG", IUPAC.unambiguous_dna)
>>> coding_dna
Seq('ATGGCCATTGTAATGGGCCGCTGAAAGGGTGCCCGATAG', IUPACUnambiguousDNA())
>>> template_dna = coding_dna.reverse_complement()
>>> template_dna
Seq('CTATCGGGCACCCTTTCAGCGGCCCATTACAATGGCCAT', IUPACUnambiguousDNA())
>>> messenger_rna = coding_dna.transcribe()
>>> messenger_rna
Seq('AUGGCCAUUGUAAUGGGCCGCUGAAAGGGUGCCCGAUAG', IUPACUnambiguousRNA()) 

번역 (Translation)

전사에 이어 번역을 해보겠습니다.

Sequence 객체에서 translate 메소드를 사용하시면 됩니다.

>>> messenger_rna
Seq('AUGGCCAUUGUAAUGGGCCGCUGAAAGGGUGCCCGAUAG', IUPACUnambiguousRNA())
>>> messenger_rna.translate() 

Seq('MAIVMGR*KGAR*', HasStopCodon(IUPACProtein(), '*'))

번역 테이블 (Translation Table)

앞서 Sequence 메소드 중 translate를 사용하여 번역과정을 진행하였습니다.

이처럼 번역과정에는 각 코돈 테이블에 맞는 아미노산이 번역됩니다.

그리고 코돈 테이블에는 사람만 있는 것이 아니라 미토콘드리아 등 여러 코돈 테이블이 있습니다.

여러 종에 대한 코돈 테이블은 다음 링크를 참조하시면 됩니다.

https://www.ncbi.nlm.nih.gov/Taxonomy/Utils/wprintgc.cgi

아래의 예는 사람과 미토콘드리아에 관한 번역 테이블입니다.

>>> from Bio.Data import CodonTable
>>> standard_table = CodonTable.unambiguous_dna_by_name["Standard"]
>>> mito_table = CodonTable.unambiguous_dna_by_name["Vertebrate Mitochondrial"]
>>> print(standard_table)
Table 1 Standard, SGC0

  |  T      |  C      |  A      |  G      |
--+---------+---------+---------+---------+--
T | TTT F   | TCT S   | TAT Y   | TGT C   | T
T | TTC F   | TCC S   | TAC Y   | TGC C   | C
T | TTA L   | TCA S   | TAA Stop| TGA Stop| A
T | TTG L(s)| TCG S   | TAG Stop| TGG W   | G
--+---------+---------+---------+---------+--
C | CTT L   | CCT P   | CAT H   | CGT R   | T
C | CTC L   | CCC P   | CAC H   | CGC R   | C
C | CTA L   | CCA P   | CAA Q   | CGA R   | A
C | CTG L(s)| CCG P   | CAG Q   | CGG R   | G
--+---------+---------+---------+---------+--
A | ATT I   | ACT T   | AAT N   | AGT S   | T
A | ATC I   | ACC T   | AAC N   | AGC S   | C
A | ATA I   | ACA T   | AAA K   | AGA R   | A
A | ATG M(s)| ACG T   | AAG K   | AGG R   | G
--+---------+---------+---------+---------+--
G | GTT V   | GCT A   | GAT D   | GGT G   | T
G | GTC V   | GCC A   | GAC D   | GGC G   | C
G | GTA V   | GCA A   | GAA E   | GGA G   | A
G | GTG V   | GCG A   | GAG E   | GGG G   | G
--+---------+---------+---------+---------+--
>>> print(mito_table)
Table 2 Vertebrate Mitochondrial, SGC1

  |  T      |  C      |  A      |  G      |
--+---------+---------+---------+---------+--
T | TTT F   | TCT S   | TAT Y   | TGT C   | T
T | TTC F   | TCC S   | TAC Y   | TGC C   | C
T | TTA L   | TCA S   | TAA Stop| TGA W   | A
T | TTG L   | TCG S   | TAG Stop| TGG W   | G
--+---------+---------+---------+---------+--
C | CTT L   | CCT P   | CAT H   | CGT R   | T
C | CTC L   | CCC P   | CAC H   | CGC R   | C
C | CTA L   | CCA P   | CAA Q   | CGA R   | A
C | CTG L   | CCG P   | CAG Q   | CGG R   | G
--+---------+---------+---------+---------+--
A | ATT I(s)| ACT T   | AAT N   | AGT S   | T
A | ATC I(s)| ACC T   | AAC N   | AGC S   | C
A | ATA M(s)| ACA T   | AAA K   | AGA Stop| A
A | ATG M(s)| ACG T   | AAG K   | AGG Stop| G
--+---------+---------+---------+---------+--
G | GTT V   | GCT A   | GAT D   | GGT G   | T
G | GTC V   | GCC A   | GAC D   | GGC G   | C
G | GTA V   | GCA A   | GAA E   | GGA G   | A
G | GTG V(s)| GCG A   | GAG E   | GGG G   | G 

--+---------+---------+---------+---------+--

코돈테이블에 대해서는 다음 메소드도 가능합니다.

>>> mito_table.stop_codons

['TAA', 'TAG', 'AGA', 'AGG']
>>> mito_table.start_codons
['ATT', 'ATC', 'ATA', 'ATG', 'GTG']

Sequence 비교하기

Sequence 객체는 서로 비교가 가능합니다.

>>> from Bio.Seq import Seq
>>> from Bio.Alphabet import IUPAC
>>> seq1 = Seq("ACGT", IUPAC.unambiguous_dna)
>>> seq2 = Seq("ACGT", IUPAC.ambiguous_dna)
>>> str(seq1) == str(seq2)
True
>>> seq1 == seq2
True
>>> seq1 == "ACGT" 

True

각 객체끼리 또는 문자열 데이터와도 비교 가능합니다.

만약 Sequence객체에서 비교를 할 수없는 객체끼리의 비교는 어떨까요?

예를 들면, DNA와 Protein 서열의 비교를 해봅시다.

>>> from Bio.Seq import Seq
>>> from Bio.Alphabet import generic_dna, generic_protein
>>> dna_seq = Seq("ACGT", generic_dna)
>>> ptn_seq = Seq("ACGT", generic_protein)
>>> dna_seq == ptn_seq
  BiopythonWarning: Incompatible alphabets DNAAlphabet() and ProteinAlphabet()
  BiopythonWarning) 

True

True가 결과로 나오기는 합니다만 Warning 사인이 나오는군요.

MutableSeq 객체

파이썬에서 String 객체는 immutable 입니다.

>>> s = "ACGT"
>>> s[1]
'C'
>>> s[1] = "G"
Traceback (most recent call last):
  File "<stdin>", line 1, in <module> 

TypeError: 'str' object does not support item assignment

Sequence 객체도 String 과 같이 immutable 입니다.

>>> from Bio.Seq import Seq
>>> from Bio.Alphabet import IUPAC
>>> my_seq = Seq("GCCATTGTTAATGGGCCGCTGAAAGGGTGCCCGA", IUPAC.unambiguous_dna)
>>> my_seq[5] = "G"
Traceback (most recent call last):
  File "<stdin>", line 1, in <module> 

TypeError: 'Seq' object does not support item assignment

그런데 immutable인 Sequence객체를 mutable로 바꿀 수 있습니다.

>>> mutable_seq = my_seq.tomutable()
>>> mutable_seq
MutableSeq('GCCATTGTTAATGGGCCGCTGAAAGGGTGCCCGA', IUPACUnambiguousDNA())
>>> mutable_seq[5] = "C"
>>> mutable_seq
MutableSeq('GCCATCGTTAATGGGCCGCTGAAAGGGTGCCCGA', IUPACUnambiguousDNA())
>>> mutable_seq.remove("T")
>>> new_seq = mutable_seq.toseq()
>>> new_seq
Seq('GCCACGTTAATGGGCCGCTGAAAGGGTGCCCGA', IUPACUnambiguousDNA()) 

>>> 

Sequence 객체에서 tomutable 메소드를 사용하면 MutableSeq 이 됩니다.

그 반대면 toseq 메소드를 사용하면 다시 Sequence 객체가 됩니다.

UnknownSeq 객체

Unknown Sequence 란 길이가 있고 염기가 있는 것은 알지만 정확히 어떠한 염기서열이 있는지 모르는 경우로 보통 'N' 으로 채우게 됩니다.

>>> from Bio.Seq import UnknownSeq
>>> from Bio.Alphabet import IUPAC
>>> unk_dna = UnknownSeq(20, alphabet=IUPAC.ambiguous_dna)
>>> unk_dna
UnknownSeq(20, alphabet = IUPACAmbiguousDNA(), character = 'N')
>>> print(unk_dna) 

NNNNNNNNNNNNNNNNNNNN

String 자료형으로 바로 작업하기

이번 Sequence 객체 강좌를 마치면서 마지막으로 String 자료형을 바로 '역상보서열', '전사', '번역' 하는 방법을 소개하며 마치겠습니다.

>>> from Bio.Seq import reverse_complement, transcribe, translate
>>> my_string = "GCTGTTATGGGTCGTTGGAAGGGTGGTCGTGCTGCTGGTTAG"
>>> reverse_complement(my_string)
'CTAACCAGCAGCACGACCACCCTTCCAACGACCCATAACAGC'
>>> transcribe(my_string)
'GCUGUUAUGGGUCGUUGGAAGGGUGGUCGUGCUGCUGGUUAG'
>>> translate(my_string) 

'AVMGRWKGGRAAG*'

다음 시간은 Sequence annotation 객체에 대하여 말씀드리겠습니다.

그럼 다음 시간에 만나요~

--- 2019.03.13

안녕하세요 여러분~ 한주현입니다!

제가 BJ퍼블릭 출판사와 함께 

"바이오파이썬으로 만나는 생물정보학" 

을 출판하였습니다.

여러 예제들과 함께 준비하였습니다 ㅎㅎ 

자세한 내용은 아래 링크로 ㅎㅎ

생물정보학자 블로그 링크

https://korbillgates.tistory.com/211

출판사 링크

https://bjpublic.tistory.com/324

감사합니다!!

https://bjpublic.tistory.com/324