모티프 검색 — 트라이로 프로모터 서열에서 전사인자 자리 찾기

이 토픽을 마치면

교과서에서 배운 트라이(Trie) 와 이진 탐색 을 엮어서, 프로모터 서열에서 여러 전사인자 결합 자리를 한 번에 찾는 도구를 직접 만들 수 있습니다. Aho-Corasick 같은 실전 다중 패턴 매칭 알고리즘의 기초 자료구조를 코드로 이해합니다.

이 글은 교육용 일반 예제 입니다. 실전 모티프 검색은 위치 가중치 행렬(PWM), 배경 모델 등 훨씬 정교합니다.

"500개 모티프 × 100kb 서열" — 순진한 검색의 함정

여러분이 마우스 프로모터 100kb 서열에서 알려진 전사인자 결합 부위 500개(각각 6~12bp)를 모두 찾으려 한다고 합시다.

순진한 접근:

python

def naive_search(sequence: str, motifs: list[str]) -> list[tuple[int, str]]:
    hits = []
    for motif in motifs:
        for i in range(len(sequence) - len(motif) + 1):
            if sequence[i:i+len(motif)] == motif:
                hits.append((i, motif))
    return hits

이 접근의 문제는 순회를 모티프마다 반복한다는 것입니다. 시간 복잡도 O(모티프 수 × 서열 길이) = O(500 × 100,000) = 5 × 10⁷. 파이썬으로 몇 초.

만약 서열이 30억 bp짜리 사람 유전체라면? 5 × 10⁵ × 3 × 10⁹ = 1.5 × 10¹⁵. 우주가 끝나기 전에 안 끝납니다.

진짜 접근은 트라이입니다. 500개 모티프를 하나의 트리로 병합하고, 서열을 한 번만 순회하면서 트리를 따라가며 매칭을 검출합니다. 시간 복잡도가 O(서열 길이 + 매칭 수) 로 압축됩니다.

블랙박스에서 부품으로 — 트라이 열어보기

부품 1: 트라이 자료구조

트라이는 각 노드가 문자 하나를 대표하고, 루트에서 잎으로 가는 경로가 하나의 문자열을 나타내는 트리입니다.

예시 모티프 집합 {"TATA", "TATT", "GCGC"} 를 트라이로 만들면:

text

root
        /  \
       T    G
       |    |
       A    C
       |    |
       T    G
      / \   |
     A*  T* C*

* 는 완결된 모티프의 끝을 표시합니다. 세 모티프가 앞부분(TA)을 공유하면 그 부분이 한 번만 트리에 저장됩니다.

파이썬 구현:

python

class TrieNode:
    def __init__(self) -> None:
        self.children: dict[str, TrieNode] = {}
        self.matches: list[str] = []


class Trie:
    def __init__(self) -> None:
        self.root = TrieNode()
    
    def insert(self, pattern: str) -> None:
        node = self.root
        for char in pattern:
            if char not in node.children:
                node.children[char] = TrieNode()
            node = node.children[char]
        node.matches.append(pattern)

부품 2: 트라이 순회 검색

이제 서열의 각 위치에서 트라이를 따라 내려가며 매칭을 확인합니다.

python

def search(trie: Trie, sequence: str) -> list[tuple[int, str]]:
    hits: list[tuple[int, str]] = []
    
    for i in range(len(sequence)):
        node = trie.root
        j = i
        while j < len(sequence) and sequence[j] in node.children:
            node = node.children[sequence[j]]
            for match in node.matches:
                hits.append((i, match))
            j += 1
    
    return hits

각 위치 i에서 시작해 트리를 최대 깊이(가장 긴 모티프 길이)만큼 따라갑니다. 500개 모티프의 최대 길이가 12라면, 각 위치당 최대 12번 조회. 총 O(서열 길이 × 최대 모티프 길이).

실전 최적화(Aho-Corasick)는 여기서 failure link 를 추가해 재시작 위치를 미리 계산합니다. 이 튜토리얼에서는 기본 트라이로 개념만 확립합니다.

이진 탐색과의 대비

트라이의 접근을 이해하려면 대안적 접근인 정렬된 배열 + 이진 탐색 과 비교하면 명확해집니다.

정렬된 모티프 배열 접근:

python

def sorted_search(sequence: str, sorted_motifs: list[str], max_motif_len: int) -> list[tuple[int, str]]:
    from bisect import bisect_left
    
    hits = []
    for i in range(len(sequence)):
        for length in range(1, max_motif_len + 1):
            substring = sequence[i:i+length]
            idx = bisect_left(sorted_motifs, substring)
            if idx < len(sorted_motifs) and sorted_motifs[idx] == substring:
                hits.append((i, substring))
    return hits

이 접근의 시간복잡도: O(서열 길이 × 최대 모티프 길이 × log 모티프 수).

트라이 대비 log 인자만 있지만, 실제로는 트라이가 더 빠릅니다. 이유는 캐시 지역성. 트라이 순회는 관련 노드들이 연속된 메모리 접근인 반면, 이진 탐색은 배열 여기저기를 랜덤 접근합니다.

이진 탐색은 모티프 집합이 자주 바뀌지 않고 로드 시간이 중요할 때 트라이보다 유리합니다. 트라이 구축이 느리지만 검색은 빠릅니다. 이진 탐색은 정렬만 하면 되므로 초기 로드가 빠릅니다. 트레이드오프입니다.

실전 사용 예시

간단한 시나리오로 확인해 봅시다.

python

# 알려진 전사인자 결합 모티프 (예시)
motifs = [
    "TATAAA",       # TATA box
    "CAAT",         # CAAT box
    "GGGCGG",       # GC box
    "CACGTG",       # E-box
    "TGACTCA"       # AP-1
]

# 가상 프로모터 서열
promoter = "GCTATAAACCAATGGGCGGATGCACGTGCCCTGACTCAAG"

# 트라이 구축
trie = Trie()
for motif in motifs:
    trie.insert(motif)

# 검색
hits = search(trie, promoter)
for pos, motif in sorted(hits):
    print(f"위치 {pos}: {motif}")

# 출력:
# 위치 2: TATAAA
# 위치 9: CAAT
# 위치 12: GGGCGG
# 위치 21: CACGTG
# 위치 28: TGACTCA

한 번의 서열 순회로 5개 모티프가 모두 검출됐습니다.

페이딩 — 여러분이 채워야 할 두 빈칸

빈칸 1: IUPAC 코드 지원

실제 전사인자 결합 자리는 가변성을 가집니다. IUPAC 코드(R = A/G, Y = C/T, W = A/T 등)로 표현됩니다.

python

IUPAC = {
    "A": {"A"}, "C": {"C"}, "G": {"G"}, "T": {"T"},
    "R": {"A", "G"}, "Y": {"C", "T"}, "S": {"C", "G"},
    "W": {"A", "T"}, "K": {"G", "T"}, "M": {"A", "C"},
    "B": {"C", "G", "T"}, "D": {"A", "G", "T"},
    "H": {"A", "C", "T"}, "V": {"A", "C", "G"},
    "N": {"A", "C", "G", "T"}
}


def search_with_iupac(trie: Trie, sequence: str) -> list[tuple[int, str]]:
    """
    모티프에 IUPAC 코드가 포함될 수 있음.
    예: "TATA**W**A" 는 TATAAA 또는 TATATA 모두 매칭.
    """
    # TODO: 트라이 순회 시 각 노드의 children 키가 IUPAC 코드면 확장 세트로 매칭
    pass

힌트: 트리 구축은 그대로 두고, 검색에서 sequence[j] 와 각 자식 노드 키를 IUPAC 확장 세트로 매칭.

빈칸 2: 역상보 자동 검색

DNA는 이중가닥이므로 모티프가 역상보(reverse complement) 로도 나타날 수 있습니다.

python

def reverse_complement(seq: str) -> str:
    complement = {"A": "T", "T": "A", "G": "C", "C": "G"}
    return "".join(complement.get(b, b) for b in reversed(seq))


def search_both_strands(motifs: list[str], sequence: str) -> list[tuple[int, str, str]]:
    """
    각 모티프와 그 역상보를 모두 트라이에 넣고 검색.
    반환: (위치, 모티프, "forward" or "reverse")
    """
    # TODO: 트라이 구축 시 각 모티프 정방향과 역상보 모두 삽입
    # 각 엔트리의 원래 모티프와 방향을 기록
    pass

힌트: TrieNode.matches 를 list[tuple[str, str]] (원래 모티프, 방향)로 확장.

성찰 — 실전 모티프 검색과의 차이

Position Weight Matrix (PWM): 실전에서 대부분의 전사인자 결합 자리는 엄격한 문자열 매칭이 아닌 확률적 매칭입니다. 각 위치의 각 염기에 확률이 할당되어, 후보 서열의 점수를 계산합니다. JASPAR 같은 DB의 모티프가 이 형식입니다.

배경 모델: 실전은 무작위 매칭 확률을 배경 모델로 계산해 통계적으로 유의미한 매칭만 보고합니다. FIMO(MEME suite)가 표준 도구.

Aho-Corasick: 여러분의 트라이에 failure link 를 추가하면 O(서열 길이 + 매칭 수) 로 진짜 선형 시간이 됩니다. 실전 다중 패턴 매칭의 표준.

전장 유전체 규모: 사람 유전체(3B bp) 검색에는 접미사 배열, FM-index 같은 색인 자료구조가 대안입니다. BWA, Bowtie 같은 정렬기의 심장부.

확장 프로젝트

1. JASPAR 데이터 로드: 실제 JASPAR DB에서 인간 전사인자 100개 PWM을 다운로드해 여러분의 도구로 프로모터 스캔.

2. 시각화: 검색 결과를 프로모터 지도로 그리기. matplotlib 나 IGV 스타일 트랙.

3. Aho-Corasick 확장: failure link를 추가해 실제 선형 시간 매칭기 완성.

이 편의 부품 지도

[F] 트라이(Trie): 여러 문자열을 트리로 병합. 노드 = 문자, 경로 = 문자열.
[F] 이진 탐색: 정렬 배열 + bisect 로 대안 접근. 트라이와의 트레이드오프 이해.
[W] 파일 I/O: FASTA 파싱 등 (완성 스크립트로 제공).

[F] = 여러분이 직접 구현 / [W] = 완성 코드로 제공.