데이터분석/Bioinformatics

[한주현님강의] DNA 분석 파이프라인(4)

묘걍 2023. 10. 25. 15:55

Genome Analysis Tutorial - Duplication 리드 marking

🎨Duplication Read

- NGS에서 서열을 읽는 과정에서 필연적으로 PCR과 같이 DNA를 증폭시키는 과정이 있다

더보기

PCR은 DNA 증폭 방법이다. 증폭이라는 말은 '많이 만들어내는 것'이다. PCR을 사용하면 작은 양의 DNA를 가지고 그 DNA의 많은 복사본을 만들 수 있다.

생각해보면, 복사기를 사용해서 한 장의 문서를 여러 장으로 복사하는 것과 비슷하다. PCR은 DNA 복사기와 같은 역할을 하는데, 이를 통해 DNA의 매우 작은 부분을 많이 복사하여 연구나 검사에 사용할 수 있다

더보기

간단히 설명하면, PCR 검사는 바이러스의 조그만 부분을 크게 확대해서 볼 수 있게 만드는 도구라고 생각하면 된다. 코로나 바이러스가 우리 몸 안에 조금 있을 때, 이 작은 양을 찾아내기는 어렵다. 그래서 PCR을 사용하여 그 작은 양을 많이 만들어서 (증폭시켜서) 잘 볼 수 있게 하는 것이다.

즉, PCR 검사를 받으면, 만약 우리 몸 안에 코로나 바이러스가 있으면, 그 바이러스의 RNA를 많이 만들어서 바이러스 유무를 확인하는 것이다

- 이 과정에서 동일한 리드가 나오게 된다. (말 그대로 Duplicate된 것)

- 동일한 리드들은 추후에 variant calling(변이 찾는 과정)에서 문제가 될 수 있다.

생물정보학에서의 "duplication read"는 DNA 시퀀싱 과정에서 얻은 데이터에서 중복으로 나타나는 시퀀스이다. 동일한 DNA 조각이 여러 번 읽힌 것이다. 이는 다양한 원인으로 발생할 수 있는데, 그 중 하나는 시퀀싱 기술의 한계나 실험 과정에서의 오류 때문이다. 이런 중복 읽기는 데이터 분석 시 정확성을 해칠 수 있기 때문에 주의가 필요하다.

 

👀문제 상황 예시

출처: 강사님 화면

- 위 예시의 경우 기준 서열은 G인데 특정 리드에서 T로 되어 있음 (표시 안 된 곳은 G)

- T가 4개, G가 4개인 셈이니까 T로 calling 될 것 같다(??)

- 만약 T로 나온 리드들이 PCR 과정에 의한 아티펙트(artifact)라면

"아티펙트"는 실험하면서 생기는 원하지 않는 실수나 오류 같은 것이다. 그러니까, 우리가 원하는 결과와는 다른, 예기치 않은 결과가 나오는 것이다.

PCR은 DNA를 많이 복사해내는 방법이다. 그런데 이 과정에서 아티펙트가 생긴다는 것은, DNA를 복사할 때 생기는 작은 실수나 오류 때문에 DNA의 글자(시퀀스)가 조금 바뀌는 경우가 있다는 걸 의미한다. 이런 실수로 바뀐 DNA 글자는 진짜 DNA의 변화가 아니라, 실험 도중에 생긴 오류일 가능성이 크다.

- 실제로 생물 안에서는 이렇게 T가 많지 않은데 리드가 많이 복제돼서 생긴 인위적인 요소인 것

출처: 강사님 화면

- 아티펙트를 제거하는 과정을 Mark Duplicate라고 한다

- Duplication을 제거하면 4개에서 하나로 줄어든다

- Duplicate를 Mark해주는 툴로는 피카드와 샘툴즈가 있다

피카드(Picard):
피카드는 DNA 시퀀싱 데이터를 처리하고 분석하는 데 사용되는 도구 모음이다. 여기에는 여러 기능이 있는데, 그 중 하나는 "듀플리케이션 마킹"이다. 이 기능은 같은 DNA 조각이 여러 번 읽혀서 중복되는 데이터를 찾아 표시해주는 역할을 한다.

샘툴즈(Samtools):
샘툴즈는 DNA 시퀀싱 데이터를 다루는데 필요한 여러 도구들을 포함한 프로그램이다. 이 도구를 사용하면 DNA 데이터를 보기, 정렬하기, 변이 찾기 등 다양한 작업을 할 수 있다. 샘툴즈 자체에는 직접적인 듀플리케이션 마킹 기능은 없지만, 시퀀싱 데이터를 처리하는 다양한 과정에서 중요한 역할을 한다.

🧩 Samtools Command

fixmate

- Aligner가 reference 서열에 붙인 리드의 짝을 교정하는 과정

- Aligner는 DNA나 RNA 서열을 기준이 되는 서열(reference)에 맞춰서 정렬하는 도구이다. 예를 들면, 우리가 새로 얻은 DNA 서열이 어디에 위치하는지, 어떤 위치와 가장 잘 맞는지 찾아주는 도구라고 생각하면 된다.

- 우리가 DNA 정보를 얻을 때, 대부분의 경우 "리드"라는 작은 조각들로 얻게 된다. 이 리드들은 대부분 쌍으로 이루어져 있다. 이제, Aligner 도구를 사용해서 이 리드들을 기준 서열에 붙여보면(정렬하면), 때로는 이 쌍 중 하나는 제대로 붙지만 다른 하나는 잘못된 위치에 붙을 수 있다. fixmate는 바로 이렇게 잘못 붙은 리드의 짝을 바로잡아주는 작업을 하는 도구이다. 즉, 두 리드가 제대로 된 위치에 붙어 있는지 확인하고, 잘못된 경우 교정해 주는 것이다.

sort

- BAM 파일에서 리드들의 순서를 정렬하는 과정

- Mark Duplication을 하려면 fixmate 뿐만 아니라 순서도 정렬해야 한다

더보기

BAM 파일에서 리드들의 순서를 정렬하는 것은 여러 가지 이유 때문이다:

  1. 효율적인 접근: 특정 위치에 대한 리드를 빠르게 접근하려면 해당 리드들이 정렬된 상태여야 한다. 예를 들어, 특정 유전자 또는 지역에 대한 변이를 빠르게 찾아볼 수 있다.
  2. 인덱싱: BAM 파일을 정렬하면 해당 파일의 인덱스를 생성할 수 있다. 인덱싱된 BAM 파일은 특정 크로모솜 또는 위치 범위에 빠르게 접근할 수 있게 해준다.
  3. 변이 탐지: 대부분의 변이 탐지 도구는 리드가 정렬된 BAM 파일을 입력으로 사용한다. 리드들이 참조 서열에 대해 올바르게 정렬되어 있어야 정확한 변이 탐지가 가능하다.
  4. 시각화: 리드의 정렬 상태는 시각화 도구(예: IGV, UCSC Genome Browser)에서 유전체 데이터를 탐색할 때 중요하다. 정렬된 리드를 사용하면 해당 영역의 시퀀싱 깊이나 패턴을 더욱 명확하게 볼 수 있다.
  5. 중복 제거: 중복된 리드(duplicates)는 PCR 중복 또는 기타 원인으로 인해 발생할 수 있다. 정렬된 BAM 파일을 사용하면 이러한 중복 리드를 쉽게 식별하고 제거할 수 있다.
  6. 다른 분석과 연계: 많은 생물정보학 도구와 분석 파이프라인은 정렬된 BAM 파일을 기본 입력으로 사용한다. 따라서 리드를 정렬하는 것은 다양한 분석 작업에 필수적인 전처리 단계다.

이러한 이유들로 인해 BAM 파일에서 리드들을 정렬하는 것은 생물정보학 분석에서 흔히 수행되는 중요한 단계다.

markdup

- BAM 파일에서의 duplication이 있는 리드를 표시하는 과정

- (마크듑이라고 발음)

index

- BAM 파일의 index(색인)인 BAI 파일을 생성하는 과정

 

🎮실습

실습 시작 위치

https://github.com/KennethJHan/GenomeAnalysisTutorial

 

GitHub - KennethJHan/GenomeAnalysisTutorial: Genome Analysis Pipeline for Lecture

Genome Analysis Pipeline for Lecture. Contribute to KennethJHan/GenomeAnalysisTutorial development by creating an account on GitHub.

github.com

여기서 4번 참고

 

🧩 fixmate

$SAMTOOLS fixmate -m sample.mapped.bam sample.fixmate.bam
  • -m 옵션: fixmate 기능의 특별한 옵션. mate 즉 짝 정보를 고려하여 추가적인 처리를 수행한다. 두 리드 쌍의 거리(distance)나 플래그(flag) 정보를 조정할 때 사용된다.
  • sample.mapped.bam 파일에 fixmate 처리를 -m 옵션과 함께 수행한 후 그 결과를 smaple.fixmate.bam 파일로 저장하라는 뜻이다

위 명령어를 입력하고 ll로 확인해보면

위 파일이 생성된 것을 알 수 있다.

 

🧩 sort

$SAMTOOLS sort -o sample.fixmate.sorted.bam sample.fixmate.bam
  • -o: 정렬된 결과를 어느 파일에 저장할 것인지 지정한다. -o 다음에 오는 파일 이름(sample.fixmate.sorted.bam)은 정렬된 결과가 저장될 파일 이름이다
  • smaple.fixmate.bam은 입력 파일이다
  • smaple.fixmate.bam 파일을 정렬하고 그 정렬된 결과를 sample.fixmate.sorted.bam파일로 저장하라는 의미이다

완료 후 ll로 확인해보면

 

🧩 markdup

$SAMTOOLS markdup sample.fixmate.sorted.bam sample.markdup.bam
  • sample.fixmate.sorted.bam 파일의 중복된 리드를 표시하고, 그 결과를 smaple.markdup.bam 파일로 저장하라는 의미이다

완료 후 ll -tr로 확인해보면

 

👉🏻 실제로 BAM 파일에 Duplication이 마킹 되었는지 확인해보기

🎮 BAM이 어떻게 생겼는지 확인해보기

이 파일을 확인해보기

더보기
  1. $SAMTOOLS view: SAMTOOLS의 view 기능은 BAM/SAM 파일을 보거나 변환하기 위해 사용된다. 이 명령어를 사용하면 BAM 파일의 내용을 텍스트 형식으로 볼 수 있다.
  2. h: 이 옵션은 헤더 정보를 포함하여 결과를 출력하라는 의미다. 헤더는 BAM/SAM 파일의 상단에 위치하며, 리드 데이터와는 별개로 파일에 대한 여러 정보와 설정을 담고 있다.
  3. sample.markdup.bam: view 명령어로 열고자 하는 BAM 파일의 이름이다.
  4. |: 파이프라고 부르는 이 기호는 앞선 명령어의 출력을 다음 명령어의 입력으로 전달한다.
  5. less: 텍스트 파일의 내용을 페이지 단위로 보여주는 프로그램이다. 큰 파일을 한 번에 볼 때 유용하다.
  6. S: less 프로그램의 옵션으로, 이 옵션을 사용하면 오버플로우 문자가 줄바꿈 없이 화면에 표시된다. 즉, 긴 줄이 여러 줄로 나누어져 표시되는 것을 방지하고, 사용자가 수평으로 스크롤하여 내용을 볼 수 있게 해준다.

따라서, $SAMTOOLS view -h sample.markdup.bam | less -S 명령어 전체의 의미는 "sample.markdup.bam 파일의 내용을 헤더와 함께 보여주되, 그 결과를 페이지 단위로 볼 수 있게 less 프로그램을 사용하고, 긴 줄은 줄바꿈 없이 그대로 표시하라"는 것이다.

🎨 samflag

- BAM 리드에 대한 정보들이 들어 있다

- SAM 파일과 그 압축된 형태인 BAM 파일은 시퀀싱 데이터를 표현하는 형식이다. 각 리드에 대한 정보는 줄 단위로 표현되며, 여러 칼럼(column)으로 구성된다. 그 중에서 samflag 또는 종종 FLAG라고 불리는 칼럼은 해당 리드의 여러 속성을 2진수(비트) 플래그 형식으로 나타낸다.

- samflag에서 각 비트는 특정 의미를 가진다

더보기
  1. 0x1: 템플릿에 존재하는 리드가 쌍을 이루는 경우이다.
  2. 0x2: 각 리드 쌍이 올바르게 정렬된 경우이다.
  3. 0x4: 리드가 참조 서열에 정렬되지 않은 경우이다.
  4. 0x8: 짝(paired) 리드가 참조 서열에 정렬되지 않은 경우이다.
  5. 0x10: 리드가 역방향에 정렬된 경우이다.
  6. 0x20: 짝(paired) 리드가 역방향에 정렬된 경우이다.
  7. 0x40: 첫 번째 리드인 경우이다 (리드 쌍 중).
  8. 0x80: 두 번째 리드인 경우이다 (리드 쌍 중).
  9. 0x100: '보조 정렬'으로 판명된 리드이다 (주 정렬이 아닌 다른 위치에도 정렬 가능한 경우).
  10. 0x200: 리드가 플랫폼 또는 파이프라인에 의해 제외된 경우이다.
  11. 0x400: PCR 또는 기타 중복으로 판명된 리드이다.
  12. 0x800: 서열 결정이나 정렬이 실패한 리드이다.

- 이 플래그 값들은 여러 비트가 동시에 설정될 수 있으므로, 실제 samflag 값은 위의 값들의 합으로 나타난다. 주어진 samflag 값을 해석하기 위해서는 각 비트의 설정 여부를 확인하여 어떤 속성이 해당 리드에 적용되었는지를 판별할 수 있다.

https://broadinstitute.github.io/picard/explain-flags.html

 

Explain SAM Flags

 

broadinstitute.github.io

- 각각에 해당하는 숫자가 이진법으로 올라간다

BAM 파일 내의 Duplication 마킹 여부를 파악하려면 SAM flag의 특정 비트를 확인하면 된다. BAM은 SAM 파일의 이진 버전이다. 각 레코드에는 다양한 정보를 나타내는 SAM flag가 있다.

Duplication 마킹 여부는 SAM flag의 1024 (0x400) 비트로 확인할 수 있다. 예를 들어 SAM flag 값이 1089라면, 그 값은 64 + 1024 + 1의 합으로 이루어진다. 이때 1024 비트가 설정되어 있으면 해당 리드는 중복으로 마킹된 것이다.

따라서 BAM 파일의 SAM flag를 분석함으로써 중복 마킹 여부를 알 수 있다. samtools view와 같은 도구를 사용하면 이를 쉽게 확인할 수 있다.

🎮 위 리드 플래그 중 1024가 포함된 리드만 골라내보자

더보기
  1. $SAMTOOLS view: SAMTOOLS의 view 기능은 BAM/SAM 파일을 보거나 변환하기 위해 사용된다.
  2. h: 이 옵션은 헤더 정보를 포함하여 결과를 출력하라는 의미다.
  3. f 1024: f 옵션은 특정 플래그를 가진 리드만 출력하도록 필터링하는 옵션이다. 여기서 1024는 SAM flag 값으로, 이 값은 "PCR 또는 기타 중복으로 판명된 리드"를 의미한다. 따라서 f 1024는 중복된 리드만 출력하라는 의미다.
  4. sample.markdup.bam: view 명령어로 확인하고자 하는 BAM 파일의 이름이다.
  5. |: 파이프라고 부르는 이 기호는 앞선 명령어의 출력을 다음 명령어의 입력으로 전달한다.
  6. less: 텍스트 파일의 내용을 페이지 단위로 보여주는 프로그램이다.
  7. S: less의 옵션으로, 긴 줄을 줄바꿈 없이 그대로 화면에 표시하도록 한다.

결국, $SAMTOOLS view -h -f 1024 sample.markdup.bam | less -S 명령어 전체의 의미는 "sample.markdup.bam 파일에서 PCR 또는 기타 중복으로 판명된 리드만 헤더와 함께 출력하고, 그 결과를 페이지 단위로 볼 수 있게 less를 사용하되, 긴 줄은 줄바꿈 없이 그대로 표시하라"는 것이다.

만약 대문자 F를 쓰게 되면

포함되지 않은 것을 보여준다!

 

소문자 f로 했을 때 다음과 같은 결과가 나온다

이 중 첫 번째 리드인 1153을 위 싸이트에 넣어보면

해당하는 것에 체크가 되고 summary가 나온다

 

🎮 -h 말고 -H를 넣는다면

헤더만 보여준다

🎮 옵션을 안 준다면

리드만 나온다

 

🎮 wc -l

라인을 카운트해준다

전체 리드 갯수

🎮 Duplicate가 있는 라인 수

 

🧩 index

$SAMTOOLS index sample.markdup.bam
  • SAMTOOLS의 index 기능은 BAM 파일에 대한 인덱스를 생성한다. 인덱스는 BAM파일 내의 특정 위치에 빠르게 접근하기 위한 목적으로 사용된다
  • smaple.markdup.bam파일에 대한 인덱스를 생성하라는 것이다.
  • 인덱스 파일은 보통 bai확장자를 가진다

ll -tr로 확인해보면

 

🎨 tview

$SAMTOOLS tview sample.markdup.bam

 

tview
SAMTOOLS의 tview 기능은 BAM 파일에 정렬된 리드들을 텍스트 기반의 그래픽 인터페이스로 시각화해주는 도구다. 즉, tview를 사용하면 터미널 창에서 참조 서열과 그에 정렬된 리드들을 직접 보면서 데이터를 탐색할 수 있다. 이를 통해 특정 영역의 시퀀싱 결과를 확인하거나 문제점이 있는 지역을 직접 관찰할 수 있다.

입력하면 다음과 같이 나온다

🎮 포지션 검색하기

더보기

포지션 (Position)

포지션은 참조 서열에 대한 특정 위치를 가리킨다. 예를 들어, 사람의 유전체는 수십억 개의 서열로 이루어져 있으며, 이 서열의 각 부분에는 숫자로 된 위치 번호가 있다. 이 숫자로 된 위치를 '포지션'이라고 부른다.

/키를 누르면 검색창이 나온다

(강사님은 Goto만 나오는데 나는 윗칸과 디폴트값이 왜 들어가있는지 모르겠다..)

chr21:5012650를 입력한다

리드가 쌓인 것을 눈으로 볼 수 있게 해준다

더보기

리드가 쌓여 있다

시퀀싱 데이터를 분석할 때, 우리는 많은 양의 DNA 조각(리드)을 얻게 된다. 이 리드들을 참조 서열에 '정렬'하면, 같은 위치의 서열에 여러 리드가 겹쳐서 쌓이게 된다. 이렇게 여러 리드가 같은 위치에 겹쳐 있는 상태를 "리드가 쌓여 있다"고 표현한다. 이는 해당 영역의 DNA가 여러 번 읽혔다는 것을 의미하며, 이 정보는 해당 영역의 시퀀싱 품질을 알려주는 중요한 지표가 될 수 있다.

🎮 도움말

shift + ?키

q키 누르면 사라짐

🎮 베이스마다 컬러 입히기

c키

🎮 이동하기

방향키

🎮 읽는 방법

  1. 첫 번째 줄: 이 줄은 포지션(position)을 나타낸다. 이것은 참조 서열의 특정 위치를 나타내는 숫자다. 이 예시에서는 5012661부터 5012831까지의 위치를 나타낸다.
  2. 두 번째 줄: 이 줄은 참조 서열(reference sequence)을 표시한다. 여기서는 대부분 'N'으로 나와 있어 해당 영역의 정확한 염기서열이 무엇인지 알 수 없다. 그러나 마지막 부분에서는 염기서열 (CCTGGACGTC...)을 확인할 수 있다.
  3. 세 번째 줄 이후: 여기에는 각 리드(read)의 서열이 표시된다. 리드는 개인의 DNA 샘플에서 읽은 짧은 DNA 조각을 의미한다. 각 리드는 참조 서열에 '정렬'되어 있다.
    • 대문자는 해당 위치에서 리드와 참조 서열이 일치함을 나타낸다.
    • 소문자는 리드와 참조 서열 사이에 불일치가 있다는 것을 나타낸다.

즉, SAMTOOLS tview를 사용하면 참조 서열에 정렬된 리드의 시각적인 표현을 볼 수 있다. 이를 통해 특정 영역에서 리드와 참조 서열 간의 일치 및 불일치를 쉽게 확인할 수 있다.

🎮 나가기

q키

터미널로 빠져나감

 

 

Genome Analysis Tutorial - Variant calling

출처: 강사님 화면

👉🏻복습

- Raw 파일인 FASTQ가 있다

- 각각의 리드들을 BWA2 라는 툴을 가지고 기준 서열(Reference Sequence)에 매핑을 한다

- 리드에서 나오는 아티펙트를 제거하는 과정이 Mark Duplicate

- 쌓인 리드 내에서 기준 서열과 다른 것을 찾아내는 Variant Calling

- 그러면 최종적으로 VCF 파일이 나온다

 

🎨 GATK

- Genomic Analysis Tool Kit

GATK(유전체 분석 도구킷, Genome Analysis Toolkit)은 유전체 분석에 사용되는 소프트웨어 도구 모음이다. GATK는 DNA 또는 RNA 시퀀싱 데이터를 분석하여 유전체 정보를 추출하고 유전체 변이를 탐지하는 데 도움을 주는 강력한 도구 이다.

1. 유전체 분석: GATK는 생물학적 시료에서 얻은 DNA 또는 RNA의 염기서열(시퀀스)를 분석한다. 이것은 생물체의 유전 정보를 이해하고, 유전체 내에서 중요한 부분을 찾는 것을 의미한다.
2. 유전체 변이: 유전체 변이는 생물체의 유전 정보에서 발생하는 다양한 변화를 가리킨다. 이 변화에는 단일 염기의 변화(일명 SNP)나 염기의 삽입, 삭제(인델) 등이 포함된다. GATK는 이러한 변이를 탐지하고 정확하게 기록한다.
3. 도구 모음: GATK는 다양한 유전체 분석 작업을 수행하는 도구들의 모음이다. 예를 들어, 데이터 정제, SNP/인델 탐지, 유전체 변이 호출 등을 처리하는 다양한 도구가 포함되어 있다.
4. 정확성: GATK는 고도로 정확한 결과를 얻을 수 있는 도구로 알려져 있습니다. 이는 유전체 분석에서 오류를 최소화하고 신뢰성 있는 결과를 얻기 위해 중요하다.
5. 연구 및 의료 응용: GATK는 연구 및 의료 분야에서 널리 사용되며, 질병 연구, 진단, 치료 계획 등 다양한 응용 분야에서 활용된다.

요약하면, GATK는 유전체 분석에 사용되는 도구 모음으로, DNA 또는 RNA 데이터를 분석하여 유전체 정보와 변이를 탐지하는 데 도움을 주는 정확하고 강력한 도구 이다.

🧩 GATK Command

BaseRecalibrator

- Base Quality Score를 재계산하여 변이를 calling할 때 정확도를 높이는 과정

Base quality score(베이스 품질 점수)는 DNA 시퀀싱 데이터에서 각 염기(아데닌, 티민, 구아닌, 사이토신)의 품질을 나타내는 값이다. 이 품질 점수는 해당 염기가 시퀀싱 시에 얼마나 정확하게 기록되었는지를 나타내며, 높은 품질 점수는 더 정확한 염기 시퀀스를 나타낸다.

Base quality score는 일반적으로 0에서 40 사이의 숫자로 표현되며, 더 높은 점수는 더 정확한 염기 시퀀스를 의미한다. 예를 들어, 품질 점수가 40인 염기는 매우 정확하게 시퀀싱된 것을 나타내며, 데이터의 신뢰성이 높다.

BaseRecalibrator는 이러한 베이스 품질 점수를 재 계산하고 보정하여, 초기 시퀀싱 데이터에서 발생한 오차를 보상하고 더 정확한 시퀀스 정보를 얻을 수 있도록 도와준다. 이는 유전체 분석의 정확성을 향상시키는 중요한 단계 중 하나이다.

- DNA 시퀀싱 데이터의 정확성을 향상시키는 역할

- 시퀀싱 데이터에서 발생하는 오류를 식별하고, 그 오류를 보정하는 역할

- 이렇게 하면 더 정확한 유전체 분석 결과를 얻을 수 있게 도와준다.

- 데이터 품질을 향상시키는 중요한 단계 중 하나

 

ApplyBQSR

- 계산한 base quality score를 BAM파일에 반영

- DNA 시퀀싱 데이터의 정확성을 개선하는 역할

- BaseRecalibrator가 식별한 오류를 고친 정보를 시퀀싱 데이터에 적용하여 더 정확한 결과를 얻도록 도움

 

HaplotypeCaller

- BAM파일로부터 gvcf(Genomic VCF)를 생성하는 과정

     - gvcf는 genomic 포지션의 전체가 나온다

gVCF 파일은 유전체 분석에서 사용되는 특별한 종류의 파일이다. 이 파일은 DNA 시퀀싱 데이터에서 발생하는 다양한 유전체 정보를 저장하는데 사용된다.

여기에는 몇 가지 중요한 내용이 포함됩니다:
- 염기 정보: gVCF 파일은 DNA의 각 염기(아데닌, 티민, 구아닌, 사이토신)에 대한 정보를 담고 있다. 이 정보는 어떤 염기가 어느 위치에 있는지를 설명한다.
- 품질 정보: 파일에는 각 염기의 품질에 대한 정보도 포함되어 있다. 이것은 각 염기의 정확성을 나타내며, 높은 품질은 더 신뢰할 수 있는 정보를 의미한다.
- 유전체 변이 정보: gVCF 파일은 유전체 내에서 발생한 다양한 변이(예: 돌연변이)에 대한 정보를 담고 있다. 이것은 DNA의 어떤 부분이 바뀌었는지를 나타내며, 이러한 정보는 질병 연구나 다양한 유전체 연구에 사용된다.

- 유전체 변이를 탐지하는데 도움을 준다

- 유전체 변이 탐지: DNA 시퀀싱 데이터를 분석하여 다양한 종류의 유전체 변이를 식별. 이러한 변이에는 단일 염기 다형성(SNP), 인델(염기 삽입 또는 삭제), 대립 염기 서열 등이 포함된다
- 해플로타입 결정: 이 도구는 시퀀싱 데이터에서 유용한 정보를 활용해 해플로타입이라불리는 유전체 단위의 염기 서열을 결정한다. 이를 통해 정확한 변이 정보를 추출할 수 있다.
- 고품질 결과: HaplotypeCaller는 높은 품질의 결과를 생성하여, 정확한 유전체 변이를 찾는 데 도움을 준다. 이는 연구, 진단, 치료, 계획 등 다양한 응용 분야에서 유용하게 활용된다.

 

GenotypeGVCFs

- gvcf파일로부터 변이만 있는 position을 추출하여 VCF를 생성

- 다양한 개체의 유전체 정보를 결합하여 각 염기의 상태를 판단하는 도구

- 유전체 변이 데이터 수집: 먼저, 다양한 개체(예: 사람, 동물)에서 얻은 유전체 변이 데이터(GVCF 파일)를 수집한다. 이 파일들은 각 개체의 유전체 정보에 대한 단서를 제공한다
- 유전체 상태 결정: 그 다음, GenotypeGVCFs는 이러한 데이터를 결합하여 각 염기(아데닌, 티민, 구아닌, 사이토신)의 상태를 판단한다. 이는 각 염기가 어떤 형태인지(예: A, T, G, C)를 결정하는 것을 의미한다.
- 유전체 정보 제공: 마지막으로, GenotypeGVCFs는 결정된 유전체 정보를 제공한다. 이 정보는 어떤 개체에서 어떤 염기가 어떤 상태에 있는지에 대한 정보를 제공하며, 이것은 유전체 분석 결과를 이해하고 다양한 연구 및 응용 분야에서 활용할 수 있다.

 

https://github.com/KennethJHan/GenomeAnalysisTutorial

 

GitHub - KennethJHan/GenomeAnalysisTutorial: Genome Analysis Pipeline for Lecture

Genome Analysis Pipeline for Lecture. Contribute to KennethJHan/GenomeAnalysisTutorial development by creating an account on GitHub.

github.com

5번 Variant Calling 부분 따라하면 됨

시작 위치는 당연히 data디렉토리

 

🎮 BaseRecalibrator

java -jar $GATK4 BaseRecalibrator -I sample.markdup.bam -R ../resource/reference/hg38.chr21.fa --known-sites ../resource/knownsites/hg38_v0_Homo_sapiens_assembly38.known_indels.chr21.vcf.gz -L chr21 -O sample.recal_data.table
  • BaseRecalibrator 도구를 사용하여 DNA 시퀀싱 데이터를 정정하고, recalibration 데이터 테이블을 생성하는 작업을 수행
  • java -jar $GATK4: Java로 GATK 버전 4를 실행
  • BaseRecalibrator: GATK에서 사용하는 도구 중 하나로, 시퀀싱 데이터의 염기별 품질 점수를 재 계산하고 보정
  • I sample.markdup.bam: 입력 파일로 sample.markdup.bam을 사용합니다. 이 파일은 시퀀싱된 DNA 데이터가 정렬된 형태로 들어있는 BAM 형식의 파일
  • -known-sites (생략)vcf.gz: 유전체 내에서 이미 알려진 인델(염기 삽입) 변이에 대한 정보를 제공하는 VCF 파일을 지정. 이 정보를 사용하여 시퀀싱 데이터의 품질 점수를 보정할 때 참고
  • L chr21: 작업을 chr21(염기서열 정보가 있는 특정 염기)에 대해서만 수행하도록 제한. 이렇게 하면 작업이 더 빨리 완료되고, 메모리 사용이 줄어든다
  • O sample.recal_data.table: 출력 파일의 이름을 sample.recal_data.table로 지정. 이 파일은 품질 보정 정보가 저장된다

을 통해 확인해보면

table 파일이 생겼다

 

🎮 ApplyBQSR

java -jar $GATK4 ApplyBQSR -R ../resource/reference/hg38.chr21.fa -I sample.markdup.bam --bqsr-recal-file sample.recal_data.table -L chr21 -O sample.recal.bam
  • DNA 시퀀싱 데이터에 대한 Base Quality Score Recalibration (BQSR)를 적용하고 수정된 BAM 파일을 생성하는 작업을 수행
  • java -jar $GATK4: Java로 GATK 버전 4를 실행
  • ApplyBQSR: GATK에서 사용하는 도구 중 하나로, Base Quality Score Recalibration을 적용하는 역할을 한다
  • R ../resource/reference/hg38.chr21.fa: 참조 유전체 시퀀스 파일을 지정한다. 이 파일은 시퀀싱 데이터의 비교 기준이 되는 유전체 정보를 담고 있다.
  • I sample.markdup.bam: 입력 파일로 sample.markdup.bam을 사용한다. 이 파일은 시퀀싱된 DNA 데이터가 정렬된 형태로 들어있는 BAM 파일이다.
  • -bqsr-recal-file sample.recal_data.table: 이 옵션은 앞서 생성한 Recalibration 데이터 테이블인 sample.recal_data.table을 지정한다. 이 테이블은 품질 보정에 사용되는 정보를 담고 있다.
  • L chr21: 작업을 chr21(염기서열 정보가 있는 특정 염기)에 대해서만 수행하도록 제한한다. 이렇게 하면 작업이 더 빨리 완료되고, 메모리 사용이 줄어든다.
  • O sample.recal.bam: 출력 파일의 이름을 sample.recal.bam로 지정한다. 이 파일은 BQSR를 적용한 수정된 BAM 파일이다.

로 확인해보면

recal.bam 파일이 생겼다

 

🎮 HaplotypeCaller

java -jar $GATK4 HaplotypeCaller -R ../resource/reference/hg38.chr21.fa -I sample.recal.bam -L chr21 -O sample.g.vcf -ERC GVCF
  • DNA 시퀀싱 데이터에서 유전체 변이를 탐지하고, 그 정보를 gVCF(Genomic Variant Call Format) 파일에 저장하는 작업을 수행
  • java -jar $GATK4: Java로 GATK 버전 4를 실행
  • HaplotypeCaller: GATK에서 사용하는 도구 중 하나로, DNA 시퀀싱 데이터에서 유전체 변이를 탐지하고 호출하는 역할을 한다.
  • R ../resource/reference/hg38.chr21.fa: 참조 유전체 시퀀스 파일을 지정한다. 이 파일은 시퀀싱 데이터의 비교 기준이 되는 유전체 정보를 담고 있다.
  • I sample.recal.bam: 입력 파일로, Base Quality Score Recalibration(BQSR)가 적용된 sample.recal.bam을 사용한다. 이 파일은 수정된 BAM 파일로, 품질 보정이 적용되었다.
  • L chr21: 작업을 chr21(염기서열 정보가 있는 특정 염기)에 대해서만 수행하도록 제한한다. 이렇게 하면 작업이 더 빨리 완료되고, 메모리 사용이 줄어든다.
  • O sample.g.vcf: 출력 파일의 이름을 sample.g.vcf로 지정한다. 이 파일은 gVCF 형식으로 유전체 변이 정보를 저장한다.
  • ERC GVCF: 이 옵션은 HaplotypeCaller를 GVCF 모드로 설정한다. GVCF 모드에서는 모든 가능한 변이와 해당 위치에서의 품질 정보를 출력한다.

3분 걸려서 완료

로 확인해보면

gvcf 파일이 생성됨

로 뭐가 들었는지 확인해보면

모든 genomic position들이 나온 걸 볼 수 있다

쭉 보다보면 이렇게 <NON_REF> 말고 변이가 있는 포지션이 있다.

  • <NON_REF>: 이 표현은 참조 유전체(Reference Genome)와 비교하여 유전체 변이가 없는 부분을 나타낸다. 즉, 해당 염기 위치에서는 변이가 없고, 유전체가 참조 상태로 유지되는 것을 의미한다. 이는 대부분의 염기 위치에 해당하며, 대부분의 DNA가 참조 유전체와 일치하는 부분이다.
  • *, <NON_REF>: 이 표현은 참조 유전체와 비교하여 유전체 변이가 감지되지만, 정확한 변이 유형이 아직 정해지지 않았음을 나타낸다. 이 경우, 해당 위치에서는 변이가 발생했지만 정확한 유전체 변이 유형(Genotype)이 아직 정해지지 않았다. 이 표현은 보통 다양한 변이 유형을 탐지하기 위한 초기 단계에서 사용된다.

🎮 GenotypeGVCFs

java -jar $GATK4 GenotypeGVCFs -R ../resource/reference/hg38.chr21.fa -V sample.g.vcf -L chr21 -O sample.vcf

 

로 확인해보면

VCF 파일이 생성되었다

로 파일을 확인해보면

- '#'이 붙은 것이 헤더

- 나머지가 내용

  • CHROM:
    • 의미: 이 열은 유전체 상의 염기 서열 위치를 나타낸다.
    • 예시: "chr21" 또는 "1"과 같이 염기 서열이 위치한 염색체를 나타낸다.
  • POS (Position):
    • 의미: 이 열은 변이가 발생한 염기 서열 내의 위치를 나타낸다.
    • 예시: "1000"과 같이 위치를 숫자로 표현한다.
  • ID:
    • 의미: 이 열은 변이에 대한 고유한 식별자 또는 참조되는 ID를 포함할 수 있습니다. 없는 경우 "."으로 표시된다.
    • 예시: "rs12345" 또는 "." (없음).
  • REF (Reference Allele):
    • 의미: 이 열은 참조 서열에 해당하는 염기를 나타낸다.
    • 예시: "A" (참조 서열에 있는 염기).
  • ALT (Alternate Allele):
    • 의미: 이 열은 변이로 인해 참조 서열에서 다른 염기로 바뀐 염기를 나타낸다.
    • 예시: "T" (참조 서열의 "A"가 "T"로 대체된 경우).
  • QUAL (Quality):
    • 의미: 이 열은 변이의 품질 점수를 나타냅니다. 높은 값은 높은 신뢰도를 나타낸다.
    • 예시: "50.0" (품질 점수).
  • FILTER:
    • 의미: 이 열은 변이의 필터링 상태를 나타냅니다. 품질 점수 및 다른 통계 정보를 기반으로 필터링될 수 있다.
    • 예시: "PASS" (합격) 또는 "LowQual" (품질이 낮음).
  • INFO:
    • 의미: 이 열은 변이에 대한 추가 정보를 포함하는 필드다. 다양한 유전체 정보 및 변이 유형에 대한 정보를 제공한다.
    • 예시: "AC=2;AN=4" (다양한 정보 필드).
  • FORMAT:
    • 의미: 이 열은 각 샘플(개체)에 대한 정보 형식을 정의한다. 각 샘플의 염기 서열 및 유전체 정보를 나타낸다.
    • 예시: "GT:DP:GQ:AD:PL" (Genotype 정보, 커버리지, 품질, 대립 염기, 가능성 점수).
  • sample:
    • 의미: 이 열은 각 샘플(개체)에 대한 실제 유전체 정보를 포함한다. Genotype 및 다른 유전체 정보가 포함된다.
    • 예시: "0/1:30:99:10,20:500,600,700" (유전형 정보 및 다른 정보).

 

🧩 INDEL

  • Insertion (인서션):
    • 정의: Insertion은 DNA 서열에 새로운 염기가 추가되는 변이 유형을 나타낸다. 이것은 원래 서열에 새로운 염기가 삽입되어, 서열이 길어지는 것을 의미한다.
    • 예시: 원래 서열이 "ATCG"일 때, Insertion으로 인해 "ATCCG"와 같이 새로운 "C" 염기가 추가되는 것
  • Deletion (딜리션):
    • 정의: Deletion은 DNA 서열에서 하나 이상의 염기가 삭제되는 변이 유형을 나타낸다. 이것은 원래 서열에서 염기가 빠져, 서열이 짧아지는 것을 의미한다.
    • 예시: 원래 서열이 "ATCG"일 때, Deletion으로 인해 "ATG"와 같이 "C" 염기가 삭제되는 것

 

🔗 VCF 파일 Specification

https://samtools.github.io/hts-specs/VCFv4.2.pdf

 

 

🔗 gVCF, VCF 설명

https://korbillgates.tistory.com/136

 

[GATK] VCF란, gVCF란, VCF와 gVCF의 차이점 - Variant Calling Format 설명, genomic VCF

안녕하세요 한주현입니다. 오늘은 GATK의 결과 파일 중 gVCF와 VCF의 차이점과 gVCF가 무엇인지 대해서 알아보겠습니다. 요새 들어 GATK에 대해서 포스팅을 많이 쓰고 있네요 ㅎㅎ;; 워낙 GATK 안에 많

korbillgates.tistory.com

 

 

 

 

Genome Analysis Tutorial - 마무리

출처: 강사님 화면

  • DNA 추출
  • 랜덤하게 자름
  • 시퀀서에 넣음
  • 거기서 나오는 raw data인 FASTQ
  • FASTQ를 BWA를 사용하여 Reference 서열에 Mapping
  • 그 결과 BAM 파일이 나옴
  • BAM파일에서 Variant calling
  • 그 결과 최종적으로 VCF 파일이 나옴

👀 앞으로

- 위 과정을 통해 나온 VCF 내의 변이 정보들

- 그 변이들이 생물학적으로 어떤 의미를 가지고 있는지 알아야해

- 우리가 직접 알아낼 수도 있겠지만 이미 DataBase에 정보가 다 들어있다

🧩 Annotation

- 변이들에 그 정보를 붙이는 과정

- 툴로 진행함

- 어떤 이펙트가 있는지

- 병적으로 아무런 의미가 없다는 디나잉 (?? 이건 잘 알아들은건지 모르겠다)

- 병을 일으키는 페소제닉 등으로 구분

더보기
  • 이펙트(Effect): 이펙트는 변이가 유전체에서 어떠한 영향을 미칠 수 있는 것을 설명하는 용어이다. 이것은 해당 변이가 유전자를 어떻게 영향을 주는지, 유전자의 기능을 변경하는지 등을 나타낸다. 예를 들어, 이펙트가 "논실레섬"인 경우, 해당 변이는 단백질 생산에 영향을 미치지 않으며 기능적으로 중요하지 않을 수 있다.
  • 디나잉(Denying): "디나잉"이라는 용어는 일반적으로 사용되지 않는 용어이며, 여기서 제공된 문맥에서는 어떠한 변이도 아무런 의미가 없거나 중요하지 않다는 것을 나타내는 것으로 이해될 수 있다. 이것은 유전체 데이터에서 중요하지 않은 변이를 나타낼 때 사용될 수 있다.
  • 페소제닉(Pathogenic): "페소제닉"은 변이가 유전체에서 병을 일으킬 수 있는 것을 의미한다. 이것은 해당 변이가 유전자의 기능을 손상시키거나 유전체의 다른 부분과 상호 작용하여 질병을 유발할 수 있는 것을 나타낸다. 페소제닉 변이는 유전병 또는 다른 질병의 원인으로 작용할 수 있다.

VCF 파일을 어노테이션할 때, 변이가 병적인 의미가 없거나 기능적으로 중요하지 않다면 "병적으로 무의미 (Pathogenicly Non-Significant)" 또는 "기능적으로 중요하지 않음 (Functionally Insignificant)" 등의 용어를 사용하여 해당 변이가 병적인 영향을 미치지 않음을 나타낼 수 있다. 이러한 용어를 사용하여 해당 변이가 질병과 관련이 없거나 기능적으로 중요하지 않다는 것을 명확히할 수 있다. 이것은 어노테이션에서 해당 변이의 중요성을 강조하지 않는 것을 의미한다.

 

 

 

출처: 강사님 화면

- 우리가 사용한 데이터는 아주 작은 데이터

- 현실 세계의 데이터는 이정도 됨

- Exome 영역만 Sequencing한 것

Exome 영역은 유전체 내의 중요한 부분 중 하나로, 유전체 내의 단백질을 코딩하는 염기서열을 포함하는 영역이다. 유전체는 모든 유전정보를 담고 있지만, 단백질을 만들어내는 유전자(단백질 코딩 유전자)는 상대적으로 적은 비율을 차지한다. Exome은 이 단백질 코딩 유전자들을 포함하는 영역으로 제한된 영역을 가리킨다.
Exome 시퀀싱은 유전체 전체를 시퀀싱하는 비용과 시간을 줄이기 위해 주로 사용되며, 단백질 코딩 영역을 중점적으로 조사한다. 이러한 시퀀싱은 유전자 변이, 돌연변이, 질병과 관련된 유전자 등을 식별하고 연구하는 데 매우 유용하다. Exome 시퀀싱은 질병 연구, 진단, 유전자 패널 테스트 등 다양한 응용 분야에서 사용된다.
Exome 영역은 유전체의 일부이며, 크게는 인간의 DNA 중에서 약 1~2% 정도를 차지한다. 그러나 이 영역은 유전학적 연구 및 질병 연구에서 중요한 정보를 제공하며, 단백질 코딩 영역의 변이를 조사하여 유용한 유전체 정보를 얻는 데 도움이 된다.

- Genomic 전체 영역을 시퀀싱

 

 

 

 

 

 

 

출처: https://youtu.be/BBtFdJCvEoY?si=g_F_hQOkiueEGxzZ

https://youtu.be/Y5z_guJA9Z8?si=B1aCECa7DpsF5mcN

https://youtu.be/h0h2_Oj-xro?si=joeJleHkzW3ipXP_