全基因組重測序數據分析

1、全基因組重測序數據分析1. 簡介(Introduction)通過高通量測序識別發(fā)現de novo的somatic和germ line 突變，結構變異-SNV，包括重排突變（deletioin, duplication 以及copy number variation）以及SNP的座位；針對重排突變和SNP的功能性進行綜合分析；我們將分析基因功能（包括miRNA），重組率（Recombination）情況，雜合性缺失（LOH）以及進化選擇與mutation之間的關系；以及這些關系將怎樣使得在disease（cancer）genome中的mutation產生對應的易感機制和功能。我們將在基因組學以及

2、比較基因組學，群體遺傳學綜合層面上深入探索疾病基因組和癌癥基因組。實驗設計與樣本（1）Case-Control 對照組設計 ；（2）家庭成員組設計：父母-子女組（4人、3人組或多人）；初級數據分析1數據量產出： 總堿基數量、Total Mapping Reads、Uniquely Mapping Reads統(tǒng)計，測序深度分析。2一致性序列組裝：與參考基因組序列（Reference genome sequence）的比對分析，利用貝葉斯統(tǒng)計模型檢測出每個堿基位點的最大可能性基因型，并組裝出該個體基因組的一致序列。3SNP檢測及在基因組中的分布：提取全基因組中所有多態(tài)性位點，結合質量值、測序深度、

3、重復性等因素作進一步的過濾篩選，最終得到可信度高的SNP數據集。并根據參考基因組信息對檢測到的變異進行注釋。4InDel檢測及在基因組的分布: 在進行mapping的過程中，進行容gap的比對并檢測可信的short InDel。在檢測過程中，gap的長度為15個堿基。對于每個InDel的檢測，至少需要3個Paired-End序列的支持。5Structure Variation檢測及在基因組中的分布: 能夠檢測到的結構變異類型主要有：插入、缺失、復制、倒位、易位等。根據測序個體序列與參考基因組序列比對分析結果，檢測全基因組水平的結構變異并對檢測到的變異進行注釋。高級數據分析1.測序短序列匹配（R

4、ead Mapping）（1）屏蔽掉Y染色體上假體染色體區(qū)域（pseudo-autosomal region）, 將Read與參考序列NCBI36進行匹配（包括所有染色體，未定位的contig，以及線粒體序列mtDNA（將用校正的劍橋參考序列做替代）)。采用標準序列匹配處理對原始序列文件進行基因組匹配， 將Read與參考基因組進行初始匹配；給出匹配的平均質量得分分布；（2）堿基質量得分的校準。我們采用堿基質量校準算法對每個Read中每個堿基的質量進行評分，并校準一些顯著性誤差，包括來自測序循環(huán)和雙核苷酸結構導致的誤差。（3）測序誤差率估計。 pseudoautosomal contigs，sh

5、ort repeat regions（包括segmental duplication，simple repeat sequence-通過tandem repeat識別算法識別）將被過濾；2. SNP Calling 計算 （SNP Calling）我們可以采用整合多種SNP探測算法的結果，綜合地，更準確地識別出SNP。通過對多種算法各自識別的SNP進行一致性分析，保留具有高度一致性的SNP作為最終SNP結果。這些具有高度一致性的SNP同時具有非常高的可信度。在分析中使用到的SNP識別算法包括基于貝葉斯和基因型似然值計算的方法，以及使用連鎖不平衡LD或推斷技術用于優(yōu)化SNP識別檢出的準確性。統(tǒng)計

6、SNV的等位基因頻率在全基因組上的分布稀有等位基因數目在不同類別的SNV中的比率分布（a）；SNV的類別主要考慮：（1）無義（nonsense）,（2）化學結構中非同義，（3）所有非同義，（4）保守的非同義，（5）非編碼，（6）同義，等類型SNV； 另外，針對保守性的討論，我們將分析非編碼區(qū)域SNV的保守型情況及其分布（圖a, b）3. 短插入/缺失探測（Short Insertion /Deletion （Indel）Call）(1). 計算全基因組的indel變異和基因型檢出值的過程計算過程主要包含3步：（1）潛在的indel的探測；（2）通過局部重匹配計算基因型的似然值；（3）基于LD連

7、鎖不平衡的基因型推斷和檢出識別。Indel在X，Y染色體上沒有檢出值得出。(2). Indel 過濾處理4. 融合基因的發(fā)現（Fusion gene Discovery）選擇注釋的基因信息來自于當前最新版本的Ensemble Gene數據庫，RefSeq數據庫和Vega Gene數據庫。下面圖例給出的是融合基因的形成，即來自不同染色體的各自外顯子經過重組形成融合基因的模式圖。5.結構變異（Structure Variation）結構變異（Structure VariationSV）是基因組變異的一類主要來源，主要由大片段序列（一般1kb）的拷貝數變異（copy number variation

8、, CNV）以及非平衡倒位（unbalance inversion）事件構成。目前主要一些基因組研究探測識別的SV大約有20,000個（DGV數據庫）。在某些區(qū)域上，甚至SV形成的速率要大于SNP的速率，并與疾病臨床表型具有很大關聯(lián)。我們不僅可以通過測序方式識別公共的SV，也可以識別全新的SV。全新的SV的生成一般在germ line和突變機制方面都具有所報道。然而，當前對SV的精確解析需要更好的算法實現。同時，我們也需要對SV的形成機制要有更重要的認知，尤其是SV否起始于祖先基因組座位的插入或缺失，而不簡單的根據等位基因頻率或則與參考基因組序列比對判斷。SV的功能性也結合群體遺傳學和進化生物

9、學結合起來，我們綜合的考察SV的形成機制類別。SV形成機制分析，包括以下幾種可能存在的主要機制的識別發(fā)現：（A）同源性介導的直系同源序列區(qū)段重組（NAHR）；（B）與DNA雙鏈斷裂修復或復制叉停頓修復相關的非同源重組（NHR）；（C）通過擴展和壓縮機制形成可變數量的串聯(lián)重復序列（VNTR）；（D）轉座元件插入（一般主要是長短間隔序列元件LINE/SINE或者伴隨TEI相關事件的兩者的組合）。結構變異探測和擴增子（Amplicon）的探測與識別分析:如下圖所示6. 測序深度分析測序深度分析就是指根據基因組框內覆蓋度深度與期望覆蓋度深度進行關聯(lián)，并識別出SV。我們也將采用不同算法識別原始測序數據中

10、的缺失片段（deletion）和重復片段（duplication）。7. SV探測識別結果的整合與FDR推斷(可選步驟)(1). PCR或者芯片方式驗證SV(2). 計算FDR-錯誤發(fā)現率（配合驗證試驗由客戶指定）(3)篩選SV檢出結果用于SV的合并和后續(xù)分析：我們通過不同方式探測識別SV的目的極大程度的檢出SV，并且降低其FDR（=10%）。通過下屬篩選方法決定后續(xù)分析所使用到的SV集合。每種SV探測識別算法得到的SV的FDR要求小于10%，并將各自符合條件的SV合并；對于FDR大于10% 的算法計算識別的SV結果，如果有PCR和芯片平臺驗證數據，同樣可以納入后續(xù)SV分析中。最后，針對不同算

11、法得到的SV，整合處理根據breakpoint斷點左右重合覆蓋度的置信區(qū)間來評定；8. 變異屬性分析(1) neutral coalescent分析測序數據可以探測到低頻率的變異體（MAF=5%）。根據來自群體遺傳學理論（neutral coalescent理論）的期望值可以計算低頻度變異的分布。我們用不同等位基因頻率下每Mb變異數目與neutral coalescent 選擇下的期望值比值，即每Mb 基因組windows內的theta觀測值，來刻畫和反映自然純化選擇與種群（cancer cell-line可以特定的認為是可以區(qū)分的種群）增長速率。該分布分別考察SNP（藍色線），Indel（紅

12、色線），具有基因型的大片段缺失（黑色線），以及外顯子區(qū)域上的 SNP（綠色線）在不同等位基因頻率區(qū)間上的theta情況（參見下圖）。(2). 全新變異體(novel variant)的等位基因頻率和數量分布分析對象包括全新預測的SNP，indel，large deletion, 以及外顯子SNP在每個等位基因頻率類別下的數目比率（fraction）（參見下圖）；全新預測是指預測分析結果與dbSNP（當前版本129）以及deletion數據庫dbVar（2010年6月份版本）和已經發(fā)表的有關indels研究的基因組數據經過比較后識別確定的全新的SNP，indel以及deletion。dbSNP包

13、含SNP和indels; dbVAR包含有deletion,duplication,以及mobile element insertion。dbRIP以及其他基因組學研究（JC Ventrer 以及Watson 基因組，炎黃計劃亞洲人基因組）結果提供的short indels和large deletion。(3). 變異體的大小分布以及新穎性分布計算SNP，Deletion，以及Insertion 大小分布；計算SNP，Deletion，以及Insertion中屬于全新預測結果的數目占已有各自參考數據庫數目的比例（相對于dbSNP數據庫；dbSNP包含SNP和indels;dbVAR包含有del

14、etion,duplication,以及mobile element insertion。dbRIP以及其他基因組學研究（JC Ventrer 以及Watson 基因組，炎黃計劃亞洲人基因組）結果提供的short indels和large deletion）其中，可以給出LINE，Alu的特征位置。(4). 結構變異SV的斷點聯(lián)結點(BreakPoint Junction)分析根據SV不同檢出結果經過一些列篩選步驟構建所有結構變異SV的斷點聯(lián)結點數據庫，保留長度大于等于50bp的SV；分析斷點聯(lián)結點處具有homology或者microhomology的SV；并將同一染色體，起始和終止位置坐標下

15、的不同SV進行去冗余處理。分析識別SV 的斷點聯(lián)結點（Breakpoint）: 將Breakpoint按照可能形成的方式可以分類為以下幾類：（a）非等位基因同源重組型（non-allelic homologous recombination-NAHR）;（b）非同源重組（nonhomologous recombination-NHR），包括nonhomologous end-joining (NHEJ)和fork stalling /template switching（FoSTeS/MMBIR）；（c）可變串聯(lián)重復（VNTR）（d）轉座插入元件（TEI）。圖 CSV形成偏好性分析分析SV形成

16、機制與斷裂點臨近區(qū)域序列的關系，包括染色質界標（端粒，中心粒），重組高發(fā)熱點區(qū)域，重復序列以及含量，短DNA motif和微同源區(qū)域（microhomology region）。9.突變率估計針對以家庭成員為單位的測序方案，我們主要探測de novo的突變（DNM）；通過采用不同的方法/算法，我們給出每個家庭一份推斷的DNM報表；(1) 根據基因型推斷結果，分別對每人每堿基位置上的de novo突變進行綜合度量；(2) 采用貝葉斯方法計算家庭組設計中DNM的后驗概率10. SNP，SNV功能分析與注釋(1). 祖先等位基因的注釋通過將人類（NCBI36），黑猩猩（chimpanzee2.1），

17、猩猩（PPYG2）以及恒河猴（MMUL1）4種基因組進行基因組比對，發(fā)現保守的序列區(qū)域，計算祖先等位基因；以及duplication/deletion事件的進化分析。(2). 分析基因結構序列上不同區(qū)域的多樣性（Diversity）與分歧進化（divergence）根據基因型分析結果計算基因結構序列上的多樣性程度，即雜合度(heterozygosity); 雜合度指標可以說明選擇效應的存在以及局部變異的結構分布特征模式。我們將考慮基因5UTR上游200bp ，5UTR ，第一個外顯子，第一個內含子，中間外顯子，中間內含子，最末外顯子和內含子，以及3UTR及其下游200bp區(qū)域左右考察的范圍(參

18、見下圖a)。分析編碼轉錄本的起始/終止位置臨近區(qū)域的多樣性和進化分歧度（參見下圖b）。(3). 疾病變異體探測將樣本測序中分析得到SV與HGMD疾病變異體數據進行比對，得到交叉記錄的錯義和無義的SNP；通過將HGMD疾病關聯(lián)突變與CUI（疾病概念分類標識數據庫）比對獲得HGMD中所有SV的疾病表型，并獲得HGMD與測序數據分析得到的SV的疾病表型；并通過Fisher檢驗和Bonferroni多重假設檢驗校正計算樣本SV所富集的疾病表型。(4). 拷貝數變異CNV所含基因的功能注釋將CNV是否覆蓋區(qū)段重復SD區(qū)域分類為2大類，每類CNV的所含基因的功能富集情況計算，顯著性在橫軸表示；各種顯著性功

19、能在縱軸表示。(5). 變異的功能性分析與注釋（a）. SNP, Indels以及大的結構變異SV的功能注釋;（b）. 對包含翻譯起始注釋信息的轉錄本編碼區(qū)上的SNP分類為：同義SNP，非同義SNP和無義SNP（引入終止子），干擾終止子的SNP，以及干擾剪接位點的SNP；為了降低假陽性，我們采用嚴格的篩選方式過濾來自indels的錯誤；（c）.對錯義編碼區(qū)突變的功能性分析: 通過信息學分析算法評估相對于生殖系變異的體細胞突變對蛋白質的結構和功能的影響效應。(6). SNV，SNP與miRNA研究之間的關聯(lián)分析miRNA是起重要的調控作用的小分子，我們將對miRNA的pri-mRNA，pre-m

20、iRNA以及miRNA靶基因序列進行分析，識別潛在的SNP功能位點。據文獻研究提供證據表明Human pre-miRNA的二級結構中存在不同位置上的SNP，我們將通過熱力學穩(wěn)定性分析方法評估SNP對pre-miRNA結構的影響；另外，我們也將對miRNA-Target靶基因相互作用位點做分析，評估對SNP對靶基因靶向性的影響。(7). SNV，SNP與GWAS研究之間的關聯(lián)分析分析GWAS研究中得到的易感基因在基因組上不同坐標上的OR值分布情況； 將當前已知的GWAS研究成果與SNP進行比較；根據LD連鎖不平衡將SNP與易感基因的關系進行深入討論;直接與間接關聯(lián)方法可以分別識別與表型相關的SN

21、P，對于不易獲得（missing）和定位的SNP，通過LD連鎖不平衡推斷疾病易感基因突變座位。(8) 生物學通路（代謝通路，信號通路）分析生物學通路（Biological pathway），包括代謝通路和信號轉導通路是生物功能的重要組成部分，我們將各種形式的突變、變異，包括SNV和SNP，的對應基因放到生物學通路中進行綜合分析，考察功能性突變對pathway的影響程度和影響的規(guī)律。通過GSEA（配合芯片表達譜數據），KS檢驗，超幾何分布檢驗等方法對變異基因在某些pathway的富集程度進行排序，識別發(fā)生功能改變的潛在通路。(9). 蛋白質-蛋白質相互作用（PPI）網絡分析蛋白質相互作用也是生物

22、分子功能增益和缺失的重要途徑，因此我們針對蛋白質相互作用網絡中的突變的蛋白及其收到影響的網絡節(jié)點蛋白進行系統(tǒng)分析，并對收到影響的網絡子結構進行功能注釋分析和聚類富分析。我們采用網絡分析算法對由于各種突變所受到影響的子網絡（subnetwork）進行功能富集度的分析；(10). 順式基因調控網絡模塊（CRM）分析(a) 啟動子序列分析包括動子區(qū)域上的Motif預測，并與已知轉錄因子數據庫TRANSFAC和JASPAR中的TFBS結合位點進行比對；啟動子區(qū)域上保守性分析，分析突變位置和保守性區(qū)域的關聯(lián)；(b) 計算全基因組保守性。確定TFBS的保守性以及mutation位置的保守性；（11）重排（

23、arrangements）與突變（mutation）的全基因組統(tǒng)計（a）. 體細胞(somatic)和生殖系（germline）重排（arrangements）體細胞突變是相對于germ line 突變的一類需要重要分析的內容，我們針對Case-control設計的測序方案可以分別分析突變的情況，包括SNV，indel，以及CNV；如果僅在tumor/disease(Case組)出現而不在normal（對照組）出現的突變我們可以認為是somatic體細胞突變。將somatic mutation 與dbSNP數據庫比對可以發(fā)現潛在的全新的突變和有記錄的突變位置。然后，將突變分別比對到基因區(qū)域和非基因區(qū)域?；騾^(qū)域具體包括：內含子區(qū)，UTR，剪接位點區(qū)和外顯子區(qū)。其中外顯子區(qū)分別統(tǒng)計：同義（synonymous），缺失（deletion），閱讀框移位（frameshift），插入（insertion）,錯義（missense）,無義（nonsense）以及非編碼蛋白外顯子（non-protein coding exon）等不同類型。綜合不同方面分析的結果，并按照突變分類給出各重排(arrangements)類型：SNV，CNV的數目統(tǒng)計數據表（參見下圖） 。對每