Identification of novel transcripts in annotated genomes using RNA-Seq

CuffLinks 논문중 하나. CuffLinks의 한 옵션인  -g 옵션의 설명을 위한 논문. RABT(reference annotation based transcript assembly)를 설명하는 논문이다.




INTRODUCTION
RNA-Seq의 탄생으로 좀 더 정교한 transcript abundance estimation이 가능해졌다. 그런데 그것 뿐만 아니라 RNA-Seq을 EST로 봄으로써 genome annotation도 가능함을 이야기 한다. 이런 RNA-Seq으로 annotation하는데 있어서 발현양이 적은 transcript의 경우에는 transcript의 부분만이 커버가 된다는 문제점이 있다. 아래 그림에서와 같이 95%이하로 cover된 transcript가 64.4%를 차지함을 볼 수있다. 이렇듯 naive assembly 방식으로는 제대로된 transcript를 구축하는데 문제가 있을 수 있다. 

기존의 transcript assembly를 아래 3가지로 나눌수 있는데 그 어떠한 방식도 기존의 annotation을 이용하지 않는다.

1. De novo assembly
2. Genome reference based transcript assembly : genome에다가 mapping 하고 assembly 한다. Scripture 프로그램에서 했던 것. cover된 영역 사이의 gap을 채우는 것. 하지만 기존에 존재하던 annotation 파일을 이용하지는 않았다는 한계점이 있단다.
3. RNA-Seq assisted protein coding gene annotation : ab initio gene finding program의 증거로다가 RNA-Seq데이터를 이용하는 것.

그래서 이 논문에서 기존의 annotation을 assembly에 이용하는 방식인 RABT assembly 방식을 소개한다.

METHODS

위 method를 이해하기 위해서는 먼저 cufflinks의 transcript assembly를 봐야한다.
Cufflinks Transcript Assembly
cufflinks의 transcript assembly는 다음과 같은 목적을 갖는다.

1. 모든 fragment(paired-read)는 최소 한개의 assembled transcript와 일치한다
2. 모든 transcript는 read들의 tiling에 의한 것이다
3. transcript의 갯수는 위의 목적 1번을 만족하게 하는 최소한의 transcript 수이다.
4. the resulting RNA-Seq models are identifiable

결국 fragment를 설명하기 위한 최소한의 갯수의 transcript를 찾아내는 것이 cufflinks의 목적.

reference annotation을 assembly algorithm에 이용하기 위해 3가지 접근을 채택했다(위의 그림이 overview).

1. reference transcript로 부터 가상으로 faux-read를 tiled 되게 만들어서 low coverage transcript의 missing 된 부분을 찾을 수 있도록 한다. reference transcript의 15bp마다 405bp 길이의 faux-read를 생성.
2. read와 1번 step의 faux-read 모두를 이용하여 Cufflinks로 parsimonious assembly를 한다. 
3. 2번 step에서 생긴 transfrag들을 reference transcript와 비교해서 다음 다섯가지 모든 조건과 일치하면 제거한다.
• 5' endpoint가 reference transcript에 포함된 경우
• 3' endpoint가 reference transcript의 600bp 이상으로 뻗지 못한 경우(그리고 이 영역에 intron이 없는 경우)
• reference transcript에 없는 intron을 가지지 않은 경우
• a
• a

RESULTS


DISCUSSION

Fast and SNP-tolerant detection of complex variants and splicing in short reads

GSNAP 논문. 아.. BWA-SW도 봐야 하는데.. 뭐 읽어도 읽어도 끝도 없고 모르는게 너무 많다는 생각밖에 안들고.. 에라이..


요즘 short read aligner는 몇가지 고려해야 할 사항이 있는데 speed와 sequence variant 그리고 splicing event. speed는 suffix tree와 Burrows-Wheeler Transform을 활용한 방법들이 많이 나왔다. sequence variant의 경우 SNP이 1000bp 당 하나가 있고 또한 human polymorphisms의 7~8%가 indel이며 이 coding indel 중 25%는 3nt보다 길단다. 이 같은 sequence variant는 read가 길어짐에 따라 더 심각해진다. splicing event를 찾는 방법으로 exon-exon을 이어서 인공적인 sequence를 만들어서 mapping 하는것이 한 방법이 될 수 있다. 아니면 tophat 처럼 exon 주변의 splice site junction을 찾는것. 그러나 이것들은 exon 정보를 미리 알고 있거나 아니면 expression 이 많이 일어나는 exon에만 적용이 가능하다는 한계점이 있다.
뭐 이와 같은 문제점을 고려해서 만든것이 GSNAP(Genomic Short-read Nucleotide Alignment Program). 아래 그림이 GSNAP이 찾을 수 있는 complex variant 의 예. 또한 GSNAP은 single reference sequence 뿐 아니라 dbSNP 같은걸 포함하는 reference, 여기서 표현하는 것을 빌리자면 'space' reference를 이용할 수 있다고 한다.

Overview
alignment는 search problem과 같다고 보고 searching은 generating, filtering, verifying을 포함한다. efficiency는 generating과 filtering에 의존적이다. MAQ과 같은 기존의 프로그램은 read를 먼저 pre-processing 하고 나서 이 read index를 genome에 대해 generating과 filtering 해서 candidate genomic region을 찾는다. genome이 큰 경우에는 genome을 먼저 preprocessing 하는 것이 보다 효율적이다.

TopHat: discovering splice junction with RNA-Seq

항상 느끼는거지만 정석대로 가는게 젤 빠른 길인거 같다. 어떻게든 빨리 가볼려고 인터넷 뒤지고 뭐하고 해도 결국은 첫 시작점으로 돌아간다. 내가 자바스크립트 완벽 가이드를 산 것처럼. 여튼 trinity을 이해할려고 cufflinks를 보게됐고.. cufflinks를 보기 위해 tophat까지 왔다. 그렇다면 bowtie까지 가야 하는건가.. 다행히도 예전에 bowtie 논문을 조금 본 기억이 있는데.. 더 내려가진 않았으면 한다.


논문의 intro에서는 QPALMA 라는 프로그램이랑 비교한다. QPLAMA는 svm을 이용한 machine-learning 방법으로 기존에 알려져 있는 junction에 read가 mapping 되는걸로 training을 해서 새로운 exon 이랑 junction을 찾는단다. 근데 자기껀 아니라는 거다. 결국 자기네 건 training set이 필요 없단다. 그리고 훨 빠르단다. 

위 그림이 전체적인 tophat 알고리즘이다(그림만 보면 컨셉은 심플하다). 


세세하게 보자면 일단 Bowtie를 이용해서 mapping 한다. 그래서 mapping 안된 read를 IUM (initially unmapped) reads 라고 해서 따로 모아 놓는다. mapping 할 때 mismatch 나 뭐 이런 것들은 default 값으로 유지하는 것 같고 mapping redundancy는 10까지 허용. 그 이상은 버림. 10까지 허용하는 건 gene의 multi copy를 위한것. 글고 low complexity¹에 의한 unmapping은 IUM read에 들어가지 않는다.
그담엔 Bowtie의 assembly 모듈을 이용해서 mapping된 read를 assemble 한단다. splice junction에 걸쳐진 read들은 IUM read로 빠졌을 거고 그렇다면 exon의 끝부분에는 read가 몇에 없을거고 게다가 quality도 나쁠거니까(read가 끝에 얼추 정확하게 맞는다는건 exon의 끝부분이 그 read의 끝부분일테니까).. 이런 단순 mapping에 의해서 exon의 끝부분 sequence를 잃을 수 있다는 생각에 default로다가 양쪽으로 45 bp 확장한다. 발현양이 적은건 부분부분 read가 mapping 될거라서 이런 single island가 가까우면 merge 한다(default 6, 허나 mammalian일 경우 70 추천). 
그 담은 각 island 마다 canonical intron² donor 랑 acceptor를 찾는다(GT-AG).  70~20000bp(default) 안에 위치한 island 끼리의 canonical intron donor-acceptor 의 paring을 지어서 이 pair들에다가 IUM read를 맵핑한다. 중요한건 within single island 에서도 GT-AG pair를 찾는다는것. 단 속도를 높이기 위해 모든 single island에서 찾는건 아니고 depth가 높은 것에서만 찾는다.


마지막으로 Velvet+gmap이랑 자기네꺼랑 비교한다. 음 당연히 reference를 참고로 assembly 하는 TopHap이 performance가 좋단다.


음 일단은 parameter를 정리해야 할거 같다. parameter가 곧 알고리즘의 키가 되는 것이기 때문에(사실 항상 어떤 프로그램을 사용하기 전에 논문을 보는 이유가 이거다. 논문을 이해 못하면 프로그램을 제대로 사용 못하기 때문에)..


<possible parameter> 논문 보면서 예상되는 프로그램의 parameter


-D : single island로 assembly 된 것 중에도 alternative splicing에 의해 isoform이 있는 경우가 있기 때문에 single island에서도 canonical intron을 찾는데 speed를 위해 depth가 높은 island에서만 찾는다. 그 값을 정하는게 D
-s : read의 seed look up table을 만들때 사용하는 seed를 생성할 총 영역의 길이(5'쪽부터 얼마의 bp)
-k : junction의 seed라고 할 수 있다. junction에서 앞뒤로 k-mer, 즉 2k-mer를 junction의 seed로 해서 exact match되는 read의 seed index를 찾는다.
-minor isoform expression percentage : algorithm 수행후 찾아진 junction들 중에 그 expression이 주변 exon의 expression보다 15%가 안되면 report 자체를 안하는데 이 비율을 조정하는 값 
----------------------------------reference---------------------------------------


1.low complexity DNA : 항상 이 단어를 보면 그냥 당연스레 sequence complexity 가 낮은거, 그러니까 polyN 같은 걸 말하겠지 하고 그냥 넘어 갔는데. . 최소한 정확한 definition이라도 알아야 할거 같아서 링크를 건다. low complexity 를 찾는데 많이 쓰는 프로그램 중 하나가 repeatmasker 라고 링크 따라 가면 low complexity DNA에 대한 내용이 나온다. 음 대충 보니 poly-purine/ poly-pyrimidine 이랑 high AT/GC content를 이야기 하는거 같네. 


2. canonical intron, donor, acceptor : 위키 RNA splicing에 잘 나와있다. 넘 길어서 담에 봐야지