使用基因功能分類GO(Gene Ontology)分析(二):GO term的結構

GO term定義相當嚴謹,每個GO term間都存在一定的關係,可以使用DAG(Directed Acyclic Graph),即“有向無環”圖表示:

screenshot.png

某個圖中的一個圈圈就是一個GO term,也是一個節點(node),這些點的關係可以用parent和child來解釋,而整組網絡至少會有一個root,和leaf,在此圖來看,root就是A,而leaf為D, E, F,每一個點在網絡中的特性可以尤其所謂的depth 和 height表示,depth為這個節點距離root需要跨過幾個節點,而leaf則是這個節點到整個圖的最後一個時,要經過幾個節點。

另一個重點就是節點和節點的關係,有至少六種關係存在,分別是:regulates, negatively regulates, positively regulates, contribute to , colocalize with, not等等。

一個GO identifier由GO:nnnnnn組成,即GO字串加上七位數字,每個gene都有獨一無二的GO identifier, 但有多個形容此GO identified特性的GO term。

目前有需多使用GO來分析的網路工具UniGene, LocusLink(已不使用), Swiss-Prot, GeneMAPP, KEGG, PathDB. Reactome、DAVID。(遠不止這些…數繁不及備載)

文獻閱讀:RNAseq在臨床應用的可能和挑戰(下)

原始論文:Byron, S. A., Van Keuren-Jensen, K. R., Engelthaler, D. M., Carpten, J. D., & Craig, D. W. (2016). Translating RNA sequencing into clinical diagnostics: opportunities and challenges. Nature Reviews Genetics, 17(5), 257–271. doi:10.1038/nrg.2016.10

 

論文的最後,談到這些新科技如何整合進臨床之中

第一個重點就是“新的檢測(這邊指得當然就是RNAseq)”,要能變成現行臨床工作的一部份,需要符合許多特性,可以由三個層面一步步來看,分別是Analytical validity (檢測本身的效力)、Clinical validity(臨床上的檢測能力)、Clinical utility(對於最終的臨床決策的貢獻度)。

第一步:Analytical validity

相對於實驗室裡的檢測,處理的檢體都是在control condition下,所以要能讓每一次的實驗都能reproducible的技術門檻就比較低,對比於臨床簡體的複雜性,還要能在此狀況下維持檢驗的穩定性、感測極限等等,就是第一關的重要門檻,要能在有良好的檢測標準來提供檢測sensitivity、specificity的驗證,再來是每次重複的結果應該要相似的(reproducibility),不會因為微小的變數造成數值的大波動。(robustness),如前exosome的研究便是因為其容易波動而造成無法定量。

以RNAseq來看,撇開前製定序所花的不算,光後端的分析方式、參數的差異,就會造成結果有所不同,在2013年,Genetic European Variant in Disease consortium為了解決這問題設計了一組實驗來處理這些technical reproducibility和feasible的問題,收集了465個人的lymphoblastoid cell,在七個定序中心分別執行,最後總結了一些在後端上可以用來驗證定序結果的,像是GC content、fragment size、transcript length、percentage of reads mapped to annotated exon (這部分可以參考所謂的annotation QC)。除此之外,由FDA領頭的SEQC studyABRF study也有大量研究為了解決和提出相關標準。

Su, Z. et al. A comprehensive assessment of RNA-seq accuracy, reproducibility and information content by the Sequencing Quality Control Consortium. Nat. Biotechnol. 32, 903–914 (2014).

Li, S. et al. Multi-platform assessment of transcriptome profiling using RNA-seq in the ABRF next-generation sequencing study. Nat. Biotechnol. 32, 915–925 (2014).

 

第二步: Clinical validity

當檢測技術能在臨床檢體上保持者可重複性時,接下來檢測就要達成可以實際區分biological difference的區別力。

第三步: Clinical utility

當前兩部都達成後,再來要看這個檢驗的結果能如何改變治療,如同companion diagnostic等等,通常此時就需要setup 臨床的試驗,來看實際在臨床場域中的效果!

screenshot.png

 

除此之外,FDA也有針對RNAseq為主的檢測項目做出一些法規的架構,可以由此深入看美國FDA對於NGS技術延生的檢測所抱持的想法

US Department of Health & Human Services. Center for Devices and Radiological Health. FDA notification and medical device reporting for laboratory developed tests (LDTs) — draft guidance. [online]  (2014).

US Department of Health & Human Services. Center for Devices and Radiological Health. Framework for regulatory oversight of laboratory developed tests (LDTs) — draft guidance. [online]  (2014)

US Department of Health & Human Services. Optimizing FDA’ s regulatory oversight of next generation sequencing diagnostic tests — preliminary discussion paper. [online]  (2014).

Evans, B. J., Burke, W. & Jarvik, G. P. The FDA and genomic tests—getting regulation right. N. Engl. J. Med. 372, 2258–2264 (2015).

 

文獻閱讀:RNAseq在臨床應用的可能和挑戰(中)

原始論文:Byron, S. A., Van Keuren-Jensen, K. R., Engelthaler, D. M., Carpten, J. D., & Craig, D. W. (2016). Translating RNA sequencing into clinical diagnostics: opportunities and challenges. Nature Reviews Genetics, 17(5), 257–271. doi:10.1038/nrg.2016.10

文獻閱讀:RNAseq在臨床應用的可能和挑戰(中)

承襲上篇閱讀的內容,這邊談到論文的第二部分,即如何應用RNA sequencing於感染性疾病:

 可將RNAsequencing用來作為RNA-based相關病原體的診斷使用

目前臨床上重要的RNA virus如愛滋病(HIV)、伊波拉病毒(Ebola)、登革熱(Dengue)、肝炎(Hepatitis) 、流行性感冒等在醫院裡面都是使用qRT-PCR assay作為檢驗工具,這些都有可能在未來可直接用RNAsequencing來取代,除此之外,也可以用來追蹤某個疫情的爆發,像2014年的伊波拉病毒大爆發,便是使用amplicon sequencing的技術來追蹤

用來診斷傳統培養方式難以確診的細菌感染

一部份的感染疾病,其感染部位的不可性,再加上感染的細菌量少,用培養的方式通常要很久,且不一定正確,如臨床上腦膜炎的感染,確診的話會需要抽脊髓液且不一定能培養出東西,要是能從病人血液中檢驗到少許病原體的RNA便能解決,畢竟抽脊髓液本身在臨床處置上是很侵入式的,另一方面,如Mycobacterium tuberculosis的診斷用傳統方式會有問題,使用 RNAseq來提早確診外(這類細菌培養的時間都非常長,甚至幾週),也能從其transcriptome來判斷其抗藥性。

但目前來說,相對於傳統的 qt-PCR assay,使用RNAsequence在傳統的實驗中會有很大的分佈統計上的問題!

文獻閱讀:RNAseq在臨床應用的可能和挑戰(上)

原始論文:Byron, S. A., Van Keuren-Jensen, K. R., Engelthaler, D. M., Carpten, J. D., & Craig, D. W. (2016). Translating RNA sequencing into clinical diagnostics: opportunities and challenges. Nature Reviews Genetics, 17(5), 257–271. doi:10.1038/nrg.2016.10

screenshot.png

篇文章Translating RNA sequencing into clinical diagnostics: opportunities and challenges,鉅細彌遺地把RNAsequencing的應用與相關的疾病相連結,並且闡述其在臨床上的潛力,其主要分三大部分來談論RNAseq,第一部分先講利用RNAseq可以探討用來看哪些東西(轉錄變異transcript variant、non-codingRNA等),第二部份則是針對RNAseq在感染疾病的診斷和追蹤的應用,第三部份則是談論將RNAseq技術實際整合進目前的臨床情境中可能面臨的validation和法規問題。

 第一部分介紹從RNAseq可以用來看各式aberrant transcription

  • mRNA expression profile是最基本用RNAsequencing想看的,而其在臨床上實際的應用相對於DNAsequencing 可以有更多的空間,因為其可以看到動態的基因表現狀況,可用來監測、診斷疾病的狀態,如OncotypeDx 21-gene expression assay,其可以用在乳癌病患之基因檢測,用來預測癌細胞切除手術後復發的機率,或是化療及放射線治療的反應、AlloMap則是用來看心臟移植病患的急性細胞排斥,這項基因檢測服務非常完整,其臨床試驗結果還發表在NEJM上,主要使用11個基因作為signature來預測風險,收受檢體後2-3個工作天就有報告,除此之外,其他基因檢測如免疫系統的監測Adaptive Biotechnology,其主要看T-cell receptor的變異,不過也自己研發用來分析RNAseq的分析工具immunoSEQ。(雖然上述公司所發展的基因檢測商品並非直接使用RNAseq)
  • gene fusion所造成的transcripts異常,可用來區分癌細胞的subtype,不過最常見得還是在血液疾病,其中急性骨髓白血病(AML)最具有代表性,其在t(8;21)(q22;q22) translocation造成的AML1-ETO嵌合,在攝護腺癌細胞中也發現有TMPRSS2-ERG的嵌合現象,跟疾病預後相關。
  • alternative transcripts(由splicing或是structural variant造成)在發育疾病、神經退化疾病或是癌症中有角色存在,如AR-V7在頑固型攝護腺癌、EGFRvIII在glioblastoma、乳癌的BRCA1/BRCA2。
  • 特殊的RNA種類隨者技術的改善和推陳出新,開始有人研究nc-RNA或是extracelluar RNA在疾病中可用來診斷的機會。
    • extracellular RNA 是在biofluid中的胞外RNA,其被體內細胞分泌出後由extracellular vesicles 或是RNA-binding protein及lipoprotein。去量測exRNA的好處是其可以由血液中獲取(及最近很火俗稱的blood biopsy),但其lack of tissue specificity是目前的問題所在,因為所有的細胞都可能分泌出exRNA,而目前有科學家使用抽取特定組織間液的方式來解決這問題。而exRNA可以用來檢測或是追蹤癌症進展,ExoDx Lung 是今天2016年一月發表的商業基因檢測商品,其可以由血液檢體中檢測肺癌病人是否有EML4-ALK fusion transcripts或是T790M。

第二部分則介紹將RNAseq使用在感染疾病的追蹤、診斷上,如RNA-virus相關的疾病登革熱、A肝、D肝、SARS等等,在感染疾病的追蹤上則可以利用RNAseq來做metagenomic的應用,看是哪種type的感染或是其可能的傳播路徑。

第三部份則是介紹實際將RNAseq導入臨床應用,撇開本身實驗及分析RNAseq結果的困難度外,最重要的是相關法規,如執行RNAseq的醫院實驗室要怎麼管理和通過驗證,都算有其難度。

如何組織每次的分析專案檔案

每一次的資料分析,都會產生大量各式檔案,有原始資料、前處理過的、分析的程式碼,所以多花點心思在每個步驟都稍微整理回顧一下,長久來說會減少很多“災難”,讓自己能有reproducible research skill. (其實就算是wet實驗也都是如此,只是dry的分析,幾乎都在電腦上,所以好的習慣變得很重要)。

論壇上其實有不少人也有這樣的問題,故整理一些不錯的文章,雖然每一個人的做法些許不同,但重點是要建立自己一套流程。

screenshot.png

第一個重點是以project為單位來整理,不要用時間或是分散式的,所有跟某一次實驗或是分析相關的就用一組架構,架構裡最主要就分三類:資料、分析結果、分析程式碼,也是根據這三類來分,另外,可以用subproject來往下。而可以多多利用command line的技巧來減緩創建資料夾的麻煩(使用wild card),而在README文件要有良好的documentation,且可以用純文字撰寫,但可以用markdown的風格撰寫,可以在使用工具轉換成pdf、doc,使用像是pandoc的工具。

主要多花時間整理自己電腦裡的資料,整體分析的工作流程會進步更快,雖然一開始在資料量少的時候會比較煩!之後再分享用git來做版本控制!

Ensembl API(二):安裝 API

假如不害怕英文的話,這邊有很詳細的官方安裝簡介還有影片

screenshot.png

Ensembl資料主要是用MySQL所搭建的關聯資料庫,而Ensembl API則是寫好的各種Perl package來讀取裡頭儲存之資料,在上一篇有介紹過,Ensembl資料庫有三大主要資料庫:Core、Compara和Variation,每個資料庫都有不同的API需要安裝,當然,Ensembl團隊有準備好可以直接全部API安裝的方式,這邊的API不只是單純可以用來Query Ensembl資料庫,也可以用來搭建自己資料庫使用Ensembl的Schema來管理和運作,這部分又更進階了些。

 

舉例如何安裝Ensembl API

  • 下載他們的壓縮檔

           mkdir src
           cd src
           wget  ftp://ftp.ensembl.org/pub/ensembl-api.tar.gz
           wget  http://bioperl.org/DIST/BioPerl-1.6.1.tar.gz

  • 解壓縮

            tar zvxf ftp://ftp.ensembl.org/pub/ensembl-api.tar.gz 
            tar zvxf  http://bioperl.org/DIST/BioPerl-1.6.1.tar.gz

  • 重設環境變數
    PERL5LIB=${PERL5LIB}:${HOME}/src/bioperl-1.6.1
    PERL5LIB=${PERL5LIB}:${HOME}/src/ensembl/modules
    PERL5LIB=${PERL5LIB}:${HOME}/src/ensembl-compara/modules
    PERL5LIB=${PERL5LIB}:${HOME}/src/ensembl-variation/modules
    PERL5LIB=${PERL5LIB}:${HOME}/src/ensembl-funcgen/modules
    export PERL5LIB
    

另外一種安裝方式,可以使用git來安裝,這種方式比較scalable,是可以來練習一下!

Ensembl API(一):簡介

Ensembl是一個由European Bioinformatics Institute和Wellcome Trust Sanger Institute 的科學家所推動的計畫,從1999年便開始參與人類基因體計畫的執行,其致力於統整基因註釋和定序資料的整合,讓其他科學家可以輕易地使用網路來獲取需要的資料,這計畫下面主要分成8個小組(Genebuild小組產生基因定序資料、 Core Softward小組發展核心資料庫的api、 Compara小組負責進行各物種間比較計算、 Variation小組處理基因變異的資料、Regulation小組處理整體工作準則的制定、Web team小組建構網站和網路應用、Production小組打造一些分之的軟體和版本的維持、 Outreach推展ensembl計畫)約莫50人在進行。

為何要學習如何使用Ensembl API

Ensembl資料庫跟NCBI以及UCSC並列目前三大資料庫,相對於NCBI和UCSC的資料庫,Ensembl在處理個版本和註釋的處理比較一致,另外,理解Ensembl API可以讓我們將原本的資料處理workflow更好的整合在一起

Ensembl 資料庫的基本架構

Ensembl的資料庫主要是以MySQL來儲存資料的,所以基本上也可以用MySQL的語法來query裡頭的資料,基本上,在Ensembl的架構下有很多不同種類的資料庫,最主要的為Core database,儲存基因體的序列、基因註釋等等,而每一種物種都有一個獨立的Core Database,除了Core database,還有Variation database裡頭置放SNP、strains等等資料、Compara database則是放由ensembl團隊將各物種core database裡的基因序列比較後的資訊。

Ensembl資料庫幾乎每三個月都會更新一次,而直接使用Ensembl API就可以設定固定版本的database來獲取其資料,且ensembl資料庫裡放的所有資料都由實驗結果整理而成,有一套嚴謹的歸類方式。

Ensembl API主要以Perl所撰寫

Ensembl API是以物件導向的方式來撰寫(object-oriented perl),裡頭的documentation非常仔細,大觀念就是其分為Data object,和Data adaptor兩大類,Data object就是把資料連起來的集合,而Data adaptor就是取的這個方法,之後再仔細介紹其撰寫的細節。

RNAseq: Annotation Based Quality Control

當我們將定序出來的reads使用參考基因序列對應回位置後,我們便可以進一步計算每個基因上、exon、transcripts上對應的reads數量,來定量基因的表現,除此之外我們還可以用這些資訊來做進一步的實驗或是定序品質檢驗。

 

理論上來說,計算每個基因上所對應的reads數,可以很簡單的定量出基因表現,但在真核生物中,同樣的基因會有不同的isoforms,這樣便會出現一些問題,除此之外,在製備定序使用的library其protocol也會有影響,如rRNA有無去除乾淨等等,這些可以在alignment後對bam檔做annotation-based的quality control

 

雖然我們已經對fastq檔案做過初步的fastqc,但許多重要的資訊單純的fastq檔案無法提供,因為在fastq檔案裡我們只能對每個base-pair其定序出來的品質來檢查,如有沒有被rRNA汙染、或是製備protocol的檔案有沒有造成一些bias等等。

 

以下列出一些重要的品質指標(annotation-based quality metrics)

 

定序深度是否過深Saturation of sequencing depth

定序的深度其實會決定整個RNAseq實驗後續的基因表現分析、RNA剪接或是尋找特殊transcripts,所以到底實驗設計想要回答的問題需要定序多少深度會適合,而這決定於目標轉錄體的複雜性。

定序片段的在特定功能基因區分布是否有問題Read distribution between different genomic features

 可以看定序片段(reads)在不同的genomic features上的分佈,像是外顯子區域(exonic region)、內含子區域(intronic regions)、intergenic regions(基因間隔區),看分佈是否有特殊的地方,可能會發現新穎的基因或是isoforms,當然也可能會發現有無被污染,尤其是ribosomal RNA的比例,這可以知道是否在製備library時有清除乾淨。

transcripts上是否有均一地被定序到Coverage uniformity along transcripts

  不同實驗的步驟可能會造成transcripts上某一部份location被定序到的比例不同,比如說要是在製備library的時候,有poly-A capture的步驟,則其3’ ends就會有比較高被定序到的現象.

 

可以用來做annotation-based quality control的工具

 

  • RNA-SeQC
    • 為Java撰寫的程式,採用GTF格式的annotation來分析,另外,其需要bam所使用的reference fasta且要有.fai index和sequence dictionary file(,dict),其會輸出html格式的檔案,裡面包含mean coverage、針對不同location transcripts coverage的計算(3’或是5’ bias),會將基因的表現和GC content分成三個程度:low、medium、high。這篇文章有進一步介紹其功能
  • Qualimap
    • 為Java撰寫的程式,裡頭會使用到R/bioconductor的package,除了使用GTF格式,也採用BED格式作為annotation,其對於saturation和biotype 分佈的分析還不錯。其對於測序深度的分析會包含不同genomic features的深度分析、且會提供每增加一百萬的定序深度會提高發現多少新的features的機會,在biotype分佈分析,會將定序片段(reads)分佈在蛋白質轉譯區域(protein coding region)、偽基因(pseudogens)、rRNA、miRNA等地方顯示出來。
  • RseQC
    • 為python撰寫的程式,其特色是在計算不同genomic features上定序片段分佈時,會把上游和下游的transcripts都一起展現出來,另外一點特色則是其對splice junction的分析會標示出是否為novel
  • Picard’s CollectRNASeqMetrics
    • Picard工具裡其中一個Java 程式,其特色為會分析出strand specific的程度

使用samtools來看特定區域的alignments狀況

samtools的功能非常強大,也有很多用法這篇有概略介紹而這邊仔細講解一下比較細緻關於alignment檔案的extraction,就是檢視特定區域的alignment狀況。

看特定區域的alignment狀況

在position-sorted和indexed過的bam file,我們可以擷取特定區域的alignment狀況使用samtools view這指令。

首先,先將想要使用的bam file用samtools index過

$samtools index sample.bam

接下來可以看指定的區域

$samtools view sample.bam 12:7788426-7796061

是將這段區域儲存成另一個bam 檔

$samtools view sample.bam 12:7788426-7796061 > Nanog_region.bam

或是有使用BED file儲存多的region,也可直接用samtools來看這些區域

$samtools view -L cancer_exons.bed sample.bam

Cufflinks: Transcripts Assembly Tools

reference:
Trapnell, C., Williams, B. a, Pertea, G., Mortazavi, A., Kwan, G., van Baren, M. J., … Pachter, L. (2010).,Transcript assembly and quantification by RNA-Seq reveals unannotated transcripts and isoform switching during cell differentiation. Nature Biotechnology, 28(5), 511–515. http://doi.org/10.1038/nbt.1621

cufflinks是一款用c++撰寫的程式,可以用來處理RNAseq已經align玩的資料,將其alignmenta組合成transcipts,並且找尋不同isoforms的可能,操作上在command line下使用。可以在這邊的網站下載以及使用說明書。

原理

第一步:在paired-end的設計下,其會把一對的reads一起看,將其排列成overlapped的bundles,程式會先找出incompatible的fragment,因為這些特別的點代表者一種spliced isoforms的可能,如下圖所示,黃色藍色和紅色分別都是被標示的incompatible fragments,整個資訊會用node和edge的方式串接起來

screenshot.png

第二步,會使用這種graph產生出isoforms的可能組合

screenshot

軟體使用的注意事項:

  1. 假如是paired-end reads實驗設計,要把reads identified的數值要改成一樣的。
  2. samtools必須安裝好,且路徑設計在PATH變數下

語法:
screenshot

連結

–overlap–radius:default 50 base-pair
 這是用來找出比較長距離的gene models

-g  reference.gtf
可以給定已知的gene model reference

 

輸出:

總共會輸出四個檔案:transcripts.gtf   isoforms.fpkm_tracking genes.fpkm_tracking

小建議:

假如要有多組sample要做cufflinks的話,建議先把bam檔合成單一個檔案,因為當alignment數越多,使用其演算法能找出的獨特isoforms較多!