Genome에 존재하는 Variant란?

Variant 즉 다양성은 모든 생물이 가지는 고유한 특성이다. 종 수준에서의 생존력을 높이기 위해서는 개체들에게 다양성을 부여하고 환경에 더 잘 적응해가는 방향으로 진화해야 생존할 수 있기 때문이다. 한 종 밖에 존재하지 않아 멸종 할 지도 모른다는 바나나가 적절한 반대 예시일 것이다.

다른 종은 모르겠지만 인간의 경우 모든 사람이 약 99.9%의 DNA가 일치하며 0.1%가 개체의 unique함을 보장해 준다고 알려져 있다. (여기에 epgenetic한 요소까지 더해지면 일란성 쌍둥이라 할 지라도 다른 표현형을 가질 수 있게 되지만 이 포스팅에서는 다루지 않겠다.)

직접적으로 영향을 주는 변이는 크게 SNPs와 CNVs 그리고 Translocation으로 나눌 수 있다.

Single nucleotide polymorphisms (SNPs)

- 단일 염기 다형성이라고 번역되어지며 하나의 염기가 다를 때를 말한다. 약 300bp마다 하나 꼴로 나타나는것으로 알려져 있으며 표현형에 직접적으로 연관되어 있다. 

- 단순히 염기가 달라졌을 때를 변이(mutation)이라고 한다면 그 변이가 집단내에서 패턴으로 존재할 때 이를 SNPs이라고 부르게 된다. 

- 아주아주 단순한 예를 들어 어떠한 염기가 A일 때 흑발 G일 때 금발 C일 때 백발이라면 해당 염기의 위치 (염색체 3번의 123,245,321) 를 SNPs위치라고 말 할 수 있다는 것이다. 

- 물론 예시처럼 단순한하지는 않다. SNPs을 정의하기 위해서는 수 많은 샘플에서 적어도 약 1%가 넘는 돌연변이를 가진 그룹이 존재해야 한다. (A가 9900명 G가 100명 존재한다면 SNPs으로 정의할 수 있다.)

Copy number variants (CNVs)

- Duplication 또는 deletion으로 DNA의 일정 부분이 삭제되거나 중복으로 존재하여 유전자의 copy수가 바뀌게 될 때를 말한다.

- copy가 바뀐다는 말은 해당 유전자의 번역 과정에서 생겨나는 단백질의 양도 변하게 된다는 말이고 이것이 결과적으로 질병 의 발생 등에 영향을 줄 수 있게 된다.

- 대표적으로는 헌팅턴병이 유전자 말단 부분의 duplication으로 일어난다.

Translocations

- translocations은 DNA복제 과정 중에 DNA영역 일부가 기존에 위치와는 다른 위치로 옮겨지는 것을 말한다. 

- 단순해 보이지만 유전자가 발현되기 위해서는 DNA의 2차구조나 promoter 등에 의한 복잡한 조합이 일어나야 하는데 전체 유전자의 일부만 옮겨지게 되면 발현이 지나치게 되거나 적게 되는 문제가 발생할 수 있다.

Reference -

https://www.genomicseducation.hee.nhs.uk/news/item/289-various-types-of-variant-what-is-genomic-variation/

DESeq2에서 연속적인 값을 condtion으로 받기

DESeq2에서는 2가지 이상의 condition에 대한 차이를 보고 Differentially Expressed Genes를 구할 수도 있지만 시간이나 약물 투여 농도에 따른 DEGs도 계산할 수 있다.

아래의 예시처럼 약물 농도는 0nM, 0.1nM, 1nM을 투여하고 이에 따른 유전자의 변화를 가지고 p-value를 구하게 된다.

Reference -

RSeQC를 사용하여 stranded 데이터 확인하기

지난 글에 이어 RNA-seq에서 데이터가 stranded인지 unstranded인지 구별하는 방법에 대해 더 알아보기로 한다.

2018/08/03 - [bioinformatics] - RNA-seq 라이브러리 종류와 구별법

지난 글에서는 IGV로 직접 reads를 확인하는 방법을 설명했지만 이를 확인시켜주는 프로그램도 당연히 존재한다.

하지만 개인적으로 RNA-seq데이터를 다룰 때는 한 번쯤은 IGV로 데이터가 분석에 사용하기에 적절한 품질을 지녔는지 확인해보는것이 좋다고 생각하기 때문에 직접 확인하는 것을 권장하지만 데이터가 수십, 수백 개가 된다면 일일히 다 확인하는 것은 불가능하다.

그럴때 RSeQC에 포함되어있는 infer_experiment.py 스크립트를 사용하여 stranded인지 여부를 확인해보도록 하자.

#http://rseqc.sourceforge.net/#infer-experiment-py

위의 링크를 타고 들어가면 자세한 설명이 있기때문에 여기선 간단하게 요약만 한다.

필요한 데이터는 Alignement file (BAM,SAM) 파일과 annotation.bed 파일이다.

전자는 당연히 준비가 되어 있을 것이고 후자는 RSeQC 홈페이지에서 아래 링크로 가면 받을 수 있다. (human과 mouse만을 지원한다. 다른 종은 직접 만들거나 다른 데이터베이스에서 찾기 어렵지 않을 것이다.)

http://rseqc.sourceforge.net/#download-gene-models-update-on-08-07-2014

결과해석은 위의 infer-experiment-py 링크에 설명을 참고하면 된다.

infer_experiment.py -r hg19.refseq.bed12 -i Pairend_nonStrandSpecific_36mer_Human_hg19.bam

#Output::

This is PairEnd Data
Fraction of reads failed to determine: 0.0172
Fraction of reads explained by "1++,1--,2+-,2-+": 0.4903
Fraction of reads explained by "1+-,1-+,2++,2--": 0.4925

infer_experiment.py -r hg19.refseq.bed12 -i Pairend_StrandSpecific_51mer_Human_hg19.bam

#Output::

This is PairEnd Data
Fraction of reads failed to determine: 0.0072
Fraction of reads explained by "1++,1--,2+-,2-+": 0.9441
Fraction of reads explained by "1+-,1-+,2++,2--": 0.0487

infer_experiment.py -r hg19.refseq.bed12 -i SingleEnd_StrandSpecific_36mer_Human_hg19.bam

#Output::

This is SingleEnd Data
Fraction of reads failed to determine: 0.0170
Fraction of reads explained by "++,--": 0.9669
Fraction of reads explained by "+-,-+": 0.0161

단순하게 위의 Fraction과 아래의 Fraction이 약 50:50을 이룬다면 annotation되어 있는 유전자의 방향성과 맵핑된 reads의 방향성간의 통일성이 없기 때문에 unstranded 데이터라고 말할 수 있으며 반대로 첫 번째 Fraction과 두 번째 Fraction중에 한 곳으로 몰려있다면 방향성이 통일되어 있다는 뜻이기 때문에 stranded라고 할 수 있다.

위의 Fraction과 아래 Fraction중에 어디에 몰리는지는 stranded 데이터를 만들 때 방법에 차이에서부터 기인한다. library 설명 포스팅을 참고하면 된다.

Reference -

http://rseqc.sourceforge.net/

RNA-seq 라이브러리 종류와 구별법

Susan M. Corely et al., Differentially expressed genes from RNA-Seq and functional enrichment results are affected by the choice of single-end versus paired-end reads and stranded versus non-stranded protocols, BMC Genomics, 2016

Single-end 

- single-end는 우리가 생각할 수 있는 가장 기본적인 형태이다. 단순히 길이가 약 35~150bp 정도가 되는 단일 reads들의 집합으로 존재하며 수 년 전에 생산한 데이터들은 이 형태를 많이 띄고있다. 다른 방법들에 비해 만들기 쉽지만 상대적으로 정확한 위치에 맵핑되지 않는다는 단점이 있기 때문에 Illumina사의 HiSeq등의 기술로 생산하는 short reads에서는 많이 보이지 않는 추세이다.

- 하지만 PacBio, Minion, MiSeq등의 기술은 reads의 길이를 350~수 kb에 해당하게 늘릴 수 있고 길이가 긴 reads는 정확환 위치에 맵핑할 확률이 늘어나기 때문에 특정한 상황에서는 여전히 쓰이고 있다.

Paired-end

- paired-end는 초기 single-end의 단점을 극복하고자 read를 양 쪽에서 읽어 한 쌍의 fragement를 만들어 더 정확한 위치에 맵핑되고자 하는 방식이다. 

Zhernakova, et al., PLoS Genet. 2013

- 위의 사진을 비교해보면 길이가 100bp reads가 있다고 가정할 때, single-end에서는 전체 genome상에 reads와 유사한 서열이 100bp만 있다면 그 곳에 맵핑 될 수 있지만 paired-end같은 경우에는 한 쌍이라서 100bp+약간의 간격+100bp가 모두 일치해야만 되기 때문에 200bp+@가 매치되어야한다. 이는 확률적으로 잘못된 위치에 맵핑될 확률이 single-end보다 낮아진다는 것을 의미한다. (@는 insert size를 말하며 라이브러리 제작시 원하는 만큼 줄 수 있다.)

Single-end와 Paired-end 구별법

- single-end와 paired-end의 구분법은 대부분의 paired-end는 왼쪽과 오른쪽 파일이 각각 존재한다. sample_1.fq와 sample_2.fq가 있다면 각 파일 안에 있는 readID가 같은 것 끼리 짝을 이루는 식이다. 한 파일 양 reads가 동시에 존재하는 경우도 있는데 이때도 파일 안에 readID를 자세히 보면 readID가 동일하고 끝에 1/2 이런식으로 구분되어 있다면 paired-end인 것이다.

Unstranded

- 라이브러리를 제작할 때 시퀀싱 기계가 reads를 읽기 위해서 양 말단에 adapter sequence를 붙이게 되는데 같은 서열의 adapter sequence가 붙게 되면 reads의 앞 뒤를 구분할 수 없게 된다. 

Strand-specific

- Strand-specific은 양 말단의 adapter sequence를 다르게 붙임으로써 reads의 앞 뒤를 구분해 주는 방식을 말하는데, 제작방식은 크게 세 가지가 존재하며, 각 방식에따라 맵핑할 때 옵션을 다르게 줘야 할 수도 있으니 구분할 수 있는것이 좋다.

- 맵핑 프로그램을 보면 fr-strand, rf-strand 등이 존재하는데 앞의 fr은 foward reverse, rf는 reverse foward를 의미하며 각각 앞 뒤의 reads가 정방향인지 역방향인지를 의미한다. 어떤 것을 선택할 지는 paired-end 제작 방식에 따라 달라질 수 있기 때문에 지정해주어야한다. (STAR는 프로그램이 알아서 계산하지만 tophat은 직접 넣어줘야 한다. 아래는 tophat의 메뉴얼중의 일부로 옵션별 설명이 있어서 캡쳐해서 가져와 보았다.)

- 아래는 strand 예측이 왜 중요한지를 보여주는 figure이다. 가끔 genome 상에서 유전자들이 겹치는 부분이 존재한다. (DNA상에서의 방향성이 같고 겹친다면 같은 유전자이기 때문에 당연히도 예시의 두 유전자의 방향성은 반대이다.) Unstranded RNA-seq signal을 보면 어떤 유전자에서 나온 reads인지 구분이 안되지만 stranded RNA-seq이라면 RPD signal처럼 구분할 수 있다. 또한 발현량 분석 프로그램이 유전자의 방향성과 반대의 방향성을 가진 reads들은 유전자의 발현량으로 계산하지 않는 다는 점에서 strand-specific은 더 정확한 유전자 발현량을 계산할 수 있게 해주는 방식이라고 말할 수 있다.

(해당 논문은 unstrand데이터를 stranded데이터로 바꾸는 프로그램을 소개하며 해당 방식을 적용했을 때 민감도와 특이도가 높아지는 것을 보여준다.)

Bo-Hyun You et al., High-confidence coding and noncoding trnascriptome maps, Genome Research, 2017

Unstranded와 Strand-specific 구별법

- 아래 그림을 보면 방향성이 <<< 인 유전자에 있는 reads는 대부분 F1R2를 가지고 있고 >>>에는 F2R1를 가지고 있는 것을 확인할 수 있다. (가끔 by chance로 아닌 reads들이 있지만 5% 미만이다.)

위의 방식을 지원하는 프로그램도 존재한다.

2018/08/06 - [bioinformatics] - RSeQC를 사용하여 stranded 데이터 확인하기

마지막으로 Single-end, Paired-end와 Stranded, Unstranded는 별개의 개념이라는 것을 알아야 한다. Single-end이더라도 Stranded일 수도 있고 Unstranded일 수도 있으며 Paired-end도 마찬가지이다.

Reference -

https://genome.cshlp.org/content/27/6/1050.full?sid=63335fd5-50fe-46bd-b849-c88bbc2a0a84

https://www.biostars.org/p/70833/#70834

http://journals.plos.org/plosgenetics/article?id=10.1371/journal.pgen.1003594

MultiQC 설치 및 실행하기

로그파일의 이름과 무관하게 형식만 맞으면 다 읽는 것으로 보이며 그럼에도 불구하고 실행시간은 오래걸리지 않는다. 

로그파일들의 크기가 대체로 작긴 하지만 폴더내에 있는 파일 중에 매치가되지 않으면 바로바로 점프하는 듯 하다.

설치 방법은 가장 간단한 방법은 PyPI을 사용하는 방법이다.

프로그램이 Python 기반이기 때문에 Python을 설치하고 PyPI를 사용하도록 하자.

2017/08/16 - [programming language/python] - Python 설치 및 실행하기

pip가 설치되어 있다면 아래처럼 입력하고

권한문제가 생긴다면 --user 옵션을 넣도록 하자.

--user 옵션을 사용했다면 실행 파일은 ~/.local/bin/multiqc 에 생겼을 것이다.

설치가 정상적으로 끝났다면 프로젝트 폴더에가서

명령을 수행하면 multiqc_report.html 파일이 생성된다.

html파일을 열면 프로젝트 폴더 안에 있는 대부분의 로그파일들이 올라와있는 것을 확인할 수 있다.

당연하게도 모든 중간작업마다 로그파일을 생성했어야만 정상적으로 파일을 읽고 작동한다.

table이나 plot 복사를 하면 report를 작성하는데 매우 유용하게 사용할 수 있을 것이다.

Reference-

http://multiqc.info/docs/#using-multiqc