출처

참조한 글은 문장의 뒤에 링크를 포함한 출처를 표기해두었습니다.

0. 개요

Human genome filtering이란, 주어진 DNA sequencing 데이터에서 사람의 DNA로 추정되는 데이터는 필터링을 해주는 과정을 말한다.

만약 사람의 장내 세균 구성을 살펴보기 위해, 분변샘플로부터 WMS (WMS 개념에 대해서는 이전 글 참고)를 돌려 DNA 데이터들을 얻었다고 했을 때, 이 데이터에는 세균, 바이러스 등의 미생물 DNA 뿐만 아니라 장 상피 세포 안의 DNA 같은, 즉 사람의 DNA도 존재한다.

분석의 목적이 장내 세균의 DNA 데이터 분석이라면, 이러한 human DNA는 제거를 해줘야 추후 왜곡되지 않은 분석을 수행할 수 있을 것이다. 이를 위해 human genome filtering을 수행한다.

1. 수행 과정

이 과정은 보통 다른 생물정보학적 소프트웨어 (Bowtie2, STAR)를 많이 사용한다. 이 글에서는 RNA-Seq의 경우가 아닌 DNA-Seq의 경우를 가정하고, Bowtie2를 사용했을 때를 기반으로 설명한다.

- Human genome dataset 사용

먼저 Human genome dataset을 wget 명령어로 다운로드해준다. 본인은 GRCh38/hg38을 다운로드했다.

이후에는 bowtie2로 이 dataset을 subject로 삼아, raw data를 query로 하여 alignment를 해줄 건데, 이때 다운로드 받은 GRCh38을 그대로 사용하는 것이 아니라, bowtie2 index를 구축하여 genome 내에서 raw data read들의 정렬이 가능하도록 해줘야한다.

bowtie2-build를 사용하면 되고, 명령어를 잘 읽어보고 연구에 맞게 사용하길 권장드린다. 특히 thread 옵션은 속도에 영향을 미치니, 서버의 가용 자원에 맞게 사용해야 한다. 모든 옵션들을 하나씩 읽어보고 천천히 사용을 해보는게 좋은 것 같다고 생각한다. 옵션 하나를 놓쳐서 다시 소프트웨어를 사용하는 일이 비일비재하기 때문이다.

bowtie2-build -f {bowtie2 index 구축하고자 하는 타겟 파일의 경로} --threads {사용할 쓰레드 수} ./{해당 bowtie2 index database의 이름을 결정후, 경로 적기}

bowtie2-build 후, 본인은 다음과 같은 코드를 사용하여 bowtie2를 돌렸다.

"""
작성자 : 김한빈
작성일자 : 22.12.15
기능 : 모든 파일 read들을 GRCH38 reference genome에 bowtie2로 매핑

명령어 : nohup bowtie2 -x {database_path} -1 {input_path}/{file_1[i]} -2 {input_path}/{file_2[i]} -S {file_1[i].split('.')[0][:-2]}.sam --threads 10"
"""

import os 

database_name = "GRCh38"
database_path = f"툴 경로/b_bowtie2/GRCh38_bowtie2/{database_name}"

input_path = "raw data 경로/3_trimmed_fastp"
output_path = "결과를 담아둘 폴더 경로/5_human_reads_filtering"

file_1 = []
file_2 = []
for name in sorted(os.listdir(input_path)):
    if name.split(".")[0][-1] == '1': # fastq paired-end의 첫째 리드 파일을 뜻함
        file_1.append(name)
    if name.split(".")[0][-1] == '2': # fastq paired-end의 둘째 리드 파일을 뜻함
        file_2.append(name)

command = ""
for i in range(len(file_1)):
    check_count += 1
    command = f"nohup bowtie2 -x {database_path} -1 {input_path}/{file_1[i]} -2 {input_path}/{file_2[i]} -S {file_1[i].split('.')[0][:-2]}.sam --threads 10 --no-unal"  # f string을 사용하여 command를 입력

    os.system(command)

이렇게 하고나면 sam file이 나오는데, 각 sam file에 대해서, 다음과 같은 명령어를 쳐보자. 다음 명령어는, bowtie2 align 결과 (sam file)에서, 정렬되지 않은 리드들만 뽑아내라는 뜻의 명령어이다.

우리가 subject로 사람의 genome을 이용했고, query로 넣은 raw data에서 뽑고자 하는 것이 “미생물” DNA data이기 때문에, 정렬된 리드 (사람)가 아닌 정렬되지 않은 리드 (미생물)를 뽑아야 한다.

khb@khb: $samtools view -f 4 test_result5.sam

Imgur

위와 같은 예시 출력 텍스트가 나온다 (본인 노트북 가용데이터 용량이 높지 못해 genome 자체를 데이터베이스화 시키지 못했다..). 이후 이를 파이썬 이용해서 다시 fastq 파일로 바꿔주면 미생물 DNA만 뽑아 내는 과정이 완료된다.

하지만 본 목적이 장내 세균의 DNA 식별이라면, 이렇게 human genome filtering 후의 결과물만을 가지고 연구 데이터로 활용하기는 어렵다. 해당 결과물은 미생물 (사람이외의 바이러스, 세균, 곰팡이 모두 포함함) DNA일 가능성이 높은 데이터이기 때문이다. 따라서 다시 이 결과데이터를 가지고 taxon 식별이 필요한데, 이는 이 다음 글에서 다뤄보려고 한다.

이전글: WMS tech 2- QC와 트리밍