A Bayesian deconvolution strategy for immunoprecipitation-based DNA methylome analysis

MBD-Seq 을 분석하기 위해서 결정한 MEDIPS R package에서 많이 언급하는 논문이다. 일명 Batman이라는 프로그램을 소개한 논문. 그리고 보건원에서 발표할 다음 논문에서도 이용한는 프로그램이다. immunoprecipitation-based DNA methylation profiling에서 absolute methylation level을 표현하기 위한 algorithm이 소개된다.




서론서는 immunoprecipitated-based DNA methylation profiling이 가격대비 성능이 젤 좋다고..

A Bayesian deconvolution strategy for immunoprecipitation-based DNA methylation analysis

BATMAN(Bayesian tool for methylation analysis) 라는 프로그램은 만든 논문. 그 계보가 MEDIPS 라는 프로그램으로 이어진다. 


이 논문 intro에서도 나온다 MeDIP은 absolute 한 methylation 정도를 판단하기에는 CpG density의 영향을 받는다고. 해서 자기네가 이러한 문제를 해결하고자 Batman을 만들었다고.


첨엔 MeDIP-chip 실험의 feature에 대해 살펴보고 이에 맞는 수식 설명힌다.

C_cp(coupling factor) 라는 것을 정의, 이는 probe p와 CpG dinucleotide c와의 값으로 probe p에 붙은 c를 포함한 DNA 조각의 fraction. 이건 어떠쓰는 용어 정의냐? 이건 probe의 CpG density 를 정의한 C_tot를 구할때 사용한다. C_tot의 정의가 곧 sum of Ccp. probe 주변에 CpG 가 가까이 그리고  많이 있다면 당연히 probe에 붙은 DNA 조각 중에 그 CpG를 갖는 조각들이 많을거고  또 그런식으로 probe에 CpG가 많다면 그 Ccp 값의 총 합이 높을 것이므로(아.. 설명이 참..안되네). 여튼 Ctot를 probe의 CpG density로 사용하였고 그래서 array signal과 그래프로 그려보니(아래 그림 b)..

CpG poor region이 methylation이 많이 되어 있고 CpG density가 높은 곳은 methylation이 안되어 있는걸로 보임. 그리고 CpG poor region에서 CpG density와 array signal이 linear 한 관계를 보인다고(아.. 난 왜 못느끼겟지..). 
그래서 이 signal이 methylated CpG에 의해서만 생긴다고 가정하면, 곧 signal 값은 methylation에 CpG density 만큼 weighted 한 값. 그러므로 array signal의 확률 분포를 표현하자면 아래와 같이 된다.

여기서 A는 array 결과, m은 methylation 상태, G(x|u,a²) 는 가우시안 확률밀도 함수(u는 평균을 a는 표준 편차 의미..) 라고 표현이 되고 결국 우리가 얻어야 하는 값은 f(m|A). 이는 Bayesian inference를 통해서 구한다(이거 고딩때 배운거, 근데 이상하게 어렵게 느껴지긴 한다).

컴퓨터 파워를 낮추기 위해 각 CpG 를 모델링 하기보다는 이를 50- 또는 100 bp로 묶어서 그 안에 들어가 있는 CpG는 같은 level의 methylation level이라고 가정. 그러고 나서 inference를 위해 nested sampling을 사용했다는데.. 원.. 뭔소린지 모르겠다(여기 참조, 사실 이게 bioinformatics를 한다고 하는 나의 가장 큰 약점이다.사실 이 때문에 창피해서 bioinfo한다고 잘 이야기 못하는 듯).

다음은 MeDIP-Seq 에 대한 feature 설명과 BATMAN을 통한 분석 에 대한 설명. 일단은 MeDIP-Seq을 하고 read 만들어서 Maq으로 genome에 mapping하고 보니

위 그림 a와 같이 low mapping quality를 갖는 window가 보이더라(내가 이해하기론 Maqs의 mapping quality가 10 이상인 read mapping 만 이용했는데, genome을 50 kb로 window를 만들어서 각 window별로 read에 의해 얼마나 커버 됐는지에 따라 frequency를  그림으로 표현했을때). 이런 것이 알고 보니 duplication/structural variation이더라. 그래서 이런 부분은 mask 하고 분석했다. 두번째 그림은 read 양별로 얼마나 methylation 부위를 커버하는지. 보면 특정 read 수가 넘어가면 뭐 cover 하는 양이 크게 변화가 없는걸로 보인다.

BATMAN으로 MeDIP-Seq을 분석할려면 두가지에서 변형이 필요한다. 1.linear 모델대신 2차의 polynomial model 사용.2.Gaussian error model 대신 rectified Gaussian model 사용. 100bp resolution으로 결과를 뽑아 냈단다. 그리고 repeat이 있는 부위의 methylation은 under-estimate 될거란 생각에 버린다.

Scripture (Ab initio reconstrunction of cell type-specific trascriptomes in mouse reveals the conserved multi-exonic structure of lincRNAs)

얼마전에 Kobic에 김남신 박사님의 수업자료(?)를 봤는데.. 허걱했다.. 정말 RNA-seq관련 논문은 죄다 꾀고 있는듯이 보였다.. 헐.. 정말 박사라는 이름에 걸맞는 느낌이랄까. 
여튼 자극 받아 scripture 논문도 볼 예정이다. cufflinks랑 비슷한 계열인거 같아서 안볼려 했더만..  

Transcript assembly and quantification by RNA-Seq reveals unannotated transcripts and isoform switching during cell differentiation

작년 이맘때쯤 나온 논문. RNA-seq에 대해서.. 많이들 쓰는거 같다. 프로그램 이름은 Cufflinks.앞전에 trinity는 assembly-first 방식인데 이건 mapping-first 방식이다. 살짝 그림으로 봤을 때 이 프로그램 자체가 mapping 부터 하는거 같진 않고 TopHat 같은 프로그램의 결과를 가지고 뭔가 하는거 같다 (재밌는건 이 논문의 저자 중에 Barbara J Wold가 들어가 있는데 이 저자의 lab에서 ERANGE라는 RNA-seq 분석 프로그램으로 논문을 냈다.. 프로그램을 또 만들다니.. 물론 이사람이 주가 아니라 이건 CBCB가 주가 된거 같지만.. 난 과학자다. 뭐 이런건가.. 예전에 내가 만든 프로그램 따윈 과학의 발전을 위해선 중요치 않다 이건가..). 
어짜피 Trinity setting 하고 나서 cufflinks도 해야 해서.. 그리고 RPKM이 아닌 FPKM 개념이 여기서 나온거 같단 생각이 들어서(trinity output manual 보면 butterfly 결과 fasta 파일의 accession에 rel_FPKM 이 나오는데 이 부분이 아직 fancy 하지 않다고 cufflinks 처럼 발전 시키겠다고 하는걸로 보아)


intro에서 이야기 하는 것이 전에 썻던 개념(ERANGE의 RPKM을 의미하는 것 같은데)은 gene 단위. 그러니까 trascript 단위의 개념이 아니라 gene 단위로 expression을 측정 했기 때문에 promoter isoform 같은건 측정 불가 였다. 또 annotated gene model을 사용했다는 약점이 있다고 말한다. 이것이 예전 프로그램 ERANGE와 Cufflinks의 차이점. 그리고 이것이 Cufflinks의 보안점이다.


아래 그림이 Cufflinks의 algorithm을 보여주는데..

설명하자면.. 일단을 TopHat으로 mapping 한다. 그 결과를 가지고 각 bundle of fragements(=reads) alignment (overlapped 된 mapping read들) 에서 overlap graph를 생성해서 가능한 모든 길을 찾는다(여기서 말하길 Dilworth's Theorem에 따라 incompatible 한 리드가 곧 서로 다른 transcript isoform이 있다는 걸 의미한다고). 이 각각의 길이 transcript isoform. 그리고 나서 각 transcript의 abundance를 측정하는데.. paired-end read 의 길이에 대한 확률분표를 적용한다(그림 d에서 보자면 보라색 paired-end read가 하나 있는데 이는 red trascript와 blue transcript에서 나옹 가능성이 있는 read이다. 그런데 red transcript에서 나왔다면 평균적인 library 길이에서 extremely 하게 벗어난것이기 때문에 이런 것을 수식에 적용하기 위해 길이에 대한 확률 분포 사용). 그런 담에 각각의 transcript isoform의 abundance의 likelihood 값이 최대가 되도록 해서 abundance를 구한단다.


자세한건 supplementary에 있는데 만만찮다. 에라이..

Full-length transcriptome assembly from RNA-Seq data without a reference genome

de novo transcriptome assembler 인 Trinity를 소개하는 논문. 일단은 full-length transcript의 de novo assembler로 Trinity의 방법을 소개하고 fission yeast, mouse, whitefly에 대해 테스트 해본다(참고로 whitefly는 reference genome이 없단다). 또한 다른 프로그램들 (abyss, trans-abyss, oases, SOAPdenovo)와도 비교 분석한다. 

<RNA-seq 분석(특히 full-length transcript reconstruction)의 challenge(이 리스트는 아마도 다음 논문에서 따온듯 하다)>

1. transcript 마다 expression 양이 다르다.

2.한 transcript에서 나온 read들이라도 coverage가 고르지 못하다

3.많이 발현되는 transcript의 시퀀싱 error가 있는 read가 적게 발현되는 transcript의 read보다 많을 수 있다(이는 아마도 유사한 sequence일때 문제가 되지 않을까 싶네).

4.가까이에 있는 transcript 들이 ovelap되서 잘못된 chimeric transcript로 예측할 수 있다.

5.alternative splicing을 고려한 data structure가 필요.

6.서로 다른 transcript 사이에 유사한 sequence가 있을때 생기는 ambiguity.

<RNA-seq 분석에 대한 두가지 접근법>
1.Mapping-first approach : Scripture, Cufflinks 프로그램이 그 예. 일단 reference에 mapping and merging.
2.Assembly-first method : ABySS, SOAPdenovo, Oases. 일단 asseble 먼저 하고 reference 있음 그담에 mapping.

<Trinity 소개>
아래 그림이 Trinity의 대략적 개요인데.. 이걸 보면 처음에.. 뭐지 이거.. 3개의 species에 대해서 test를 해봤다는건가 이렇게 생각이 드는데. 나름 위트있게 각 모듈의 이름과 그 과정을 잘 표현한 그림이다. Trinity는 3개의 모듈(어떻게 보면 프로세스)로 되어 있는데 1.Inchworm, 2.Chrysalis, 3.Butterfly 로 되어있다(그래서 프로그램 이름이 Trinity인거 같기도). 
각 단계별로 보자면
1.Inchworm : greedy k-mer-based approach 방식. 음 뭐냐면 read를 k-mer로 profiling해서 사전을 만들고 가장 frequency가 높은 k-mer를 seed로 삼아서 이 seed를 양방향으로 overlap 되는 것끼리(k-1 mer overlap) 연결시켜 contig를 만드는데 만약 k-1 overlap이 k-mer가 많을 시에는  k-mer 중 frequency가 높은 것을 연결한다. seed 선정부터 contig extend 까지 계속적으로 반복한다.
2.Chrysalis : clustering 후 complete de Bruijn graph를 생성. Inchworm 단계에서 생성된 contig를 clustering하는데 k-1 mer 이상 겹치는 contig들 또는 read에 의해 바로 연결되는 contig들을 grouping. 그 뒤 각 cluster(component)마다 de Bruijn graph를 생성한다.이때 node가 k-1mer이고 edge가 k mer이다. 각 edge는 read frequency로 weight를 구한다. 그 후 각 read들을 가장 많은 k-mer가 공통이 되는 component에 할당한다.
3.Butterfly : 생성된 component들에 read를 대입해서 가능한 모든 full-length transcript를 찾는다. 이는 두 단계로 구성되어 있는데 첫번째로 graph simplification 를 수행한다.  graph simplification은 Chrysalis에서 만든 graph의 node를 merging 하고 error로 판단되는 edge를 pruning 하는 것으로 이루어진다. 두번째로 plausible path scoring 작업을 한다. 이는 실제 read를 simplified 된 graph에 대입함으로써 path를 찾게 되어 full-length transcript를 생성하는 것이다.

더 자세한 내용은 논문 참조.