GDAS003-Bioconductor與基因組級(jí)數(shù)據(jù)分析簡(jiǎn)介

title: GDAS003-Bioconductor與基因組級(jí)數(shù)據(jù)分析簡(jiǎn)介
date: 2019-09-03 12:0:00
type: "tags"
tags:

  • Bioconductor
    categories:
  • Genomics Data Analysis Series

前言

本篇筆記主要是介紹了Bioconductor與基因組級(jí)數(shù)據(jù)分析的關(guān)系。

R語(yǔ)言，R包與倉(cāng)庫(kù)(repositories)

學(xué)習(xí)這一系列課程的前提是你已經(jīng)有了一定的R語(yǔ)言基因。并且我們推薦使用RStudio作為R語(yǔ)言的使用環(huán)境。如果你跳過(guò)前面的1x-4x部分，直接進(jìn)入5x，可以先參考一下其它的R語(yǔ)言教程，例如 Try R 。

Why R?

Bioconductor的基礎(chǔ)就是R。使用R有3個(gè)原因：

  • 許多統(tǒng)計(jì)學(xué)家與生物統(tǒng)計(jì)學(xué)家使用R語(yǔ)言，并且他們使用R語(yǔ)言設(shè)計(jì)了許多算法可以幫助我們理解復(fù)雜的實(shí)驗(yàn)數(shù)據(jù)。
  • R語(yǔ)言有著高度的互操作性(interoperable)，促進(jìn)了使用其它語(yǔ)言編寫(xiě)軟件組件的重用。
  • R可移植到在各種計(jì)算機(jī)上的操作系統(tǒng)里(Linux、MacOSX、Windows)，初學(xué)者可以使用這些平臺(tái)中的任一個(gè)來(lái)學(xué)習(xí)。

總之，R的易用性和在統(tǒng)計(jì)學(xué)和“數(shù)據(jù)科學(xué)”中扮演著核心角色，因此當(dāng)生物學(xué)家和統(tǒng)計(jì)學(xué)家面對(duì)基因組級(jí)的實(shí)驗(yàn)數(shù)據(jù)時(shí)，選擇R作為處理工具就是自然而然的事情了。Bioconductor項(xiàng)目開(kāi)始于2001年，自那以后，它就與基因組級(jí)生物學(xué)中不斷增長(zhǎng)的數(shù)據(jù)量和復(fù)雜性保持著同步。

函數(shù)式面向?qū)ο蟪绦蛟O(shè)計(jì)

R語(yǔ)言是一種混合了函數(shù)式編程與面向?qū)ο蟮幕旌闲驼Z(yǔ)言。

在R中進(jìn)行函數(shù)式編程時(shí)，如下所示：

square = function(x) x^2

上面就是一個(gè)典型的函數(shù)式編程。我們定義了一個(gè)函數(shù)square()，它會(huì)計(jì)算輸入的數(shù)據(jù)字的平方，即^2，例如當(dāng)我們輸入3時(shí)，結(jié)果就是在R中，所有的計(jì)算都是通過(guò)函數(shù)進(jìn)行的。

在進(jìn)行面向?qū)ο蟮木幊虝r(shí)，我們所關(guān)注的生就是構(gòu)建數(shù)據(jù)結(jié)構(gòu)，定義能夠操作這些結(jié)構(gòu)化數(shù)據(jù)的方法。這是一種接近于人正常思維的編程思想。最早出現(xiàn)在20世紀(jì)90年代?，F(xiàn)在我們以Bioconductor中的一個(gè)案例說(shuō)明一下：

library(Homo.sapiens)
class(Homo.sapiens)
methods(class=class(Homo.sapiens))

結(jié)果如下所示：

> library(Homo.sapiens)
> class(Homo.sapiens)
[1] "OrganismDb"
attr(,"package")
[1] "OrganismDbi"
> methods(class=class(Homo.sapiens))
 [1] asBED                 asGFF                 cds                   cdsBy                
 [5] coerce<-              columns               dbconn                dbfile               
 [9] disjointExons         distance              exons                 exonsBy              
[13] extractUpstreamSeqs   fiveUTRsByTranscript  genes                 getTxDbIfAvailable   
[17] intronsByTranscript   isActiveSeq           isActiveSeq<-         keys                 
[21] keytypes              mapIds                mapToTranscripts      metadata             
[25] microRNAs             promoters             resources             select               
[29] selectByRanges        selectRangesById      seqinfo               show                 
[33] taxonomyId            threeUTRsByTranscript transcripts           transcriptsBy        
[37] tRNAs                 TxDb                  TxDb<-               
see '?methods' for accessing help and source code

在這個(gè)案例中，我們可以看到：

  1. OrganismDb是一個(gè)類(lèi)(class)；
  2. Homo.sapiens是OrganismDb這個(gè)類(lèi)的一個(gè)實(shí)例(instance)；
  3. 某個(gè)類(lèi)的所有方法都可以應(yīng)用到這個(gè)類(lèi)的所有實(shí)例中去（methods(class=class(Homo.sapiens))就是查看這個(gè)類(lèi)的方法）；
  4. 每個(gè)方法的實(shí)現(xiàn)是通過(guò)R中的函數(shù)進(jìn)行的（我的理解就是，方法其實(shí)本質(zhì)上就是函數(shù)）。并且函數(shù)的運(yùn)行依賴(lài)于它自身的參數(shù)。

其中需要特別注意的是juncture(交界區(qū))的方法，即genes，exons，transcripts，這三個(gè)方法輸出基因組基本組成成分的一些信息。上面我們列出的所有方法都可以應(yīng)用于OrganismDb類(lèi)定義的一些實(shí)例，也就是一些模式生物，例如小鼠(Mus musculus)，酒釀酵母(S. cerevisiae)和秀麗隱桿線(xiàn)蟲(chóng)(C. elegans)。

R包，模塊化與持續(xù)集成(Continuous* *Integration)

這一部分在閱讀的時(shí)候可以跳過(guò)。

R包的結(jié)構(gòu)

我們可以通過(guò)編寫(xiě)R代碼來(lái)用R執(zhí)行面向?qū)ο蟮暮瘮?shù)式編程。一種基本的方法就是創(chuàng)建“腳本”，用于定義數(shù)據(jù)的導(dǎo)入和分析的基本流程。當(dāng)腳本以?xún)H定義函數(shù)和數(shù)據(jù)結(jié)構(gòu)的方式編寫(xiě)時(shí)，我們就可以將這些腳本打包，從而方便地發(fā)布，供其他有類(lèi)似數(shù)據(jù)處理問(wèn)題的用戶(hù)使用。

R軟件的打包協(xié)議規(guī)定了如何將R和其他語(yǔ)言的源代碼與元數(shù)據(jù)以及文檔一起組織起來(lái)，以促進(jìn)簡(jiǎn)便的測(cè)試和發(fā)布的流程。例如，它義些文檔的包的早期版本具有如下的目錄和文件布局：

├── DESCRIPTION  (文本文件，用于說(shuō)明提供的元數(shù)據(jù)，授權(quán)信息)
├── NAMESPACE    (定義輸入與輸出)
├── R            (文件夾，包含了一些R代碼)
├── README.md    (可選，Github上的說(shuō)明內(nèi)容)
├── data         (文件夾，含有示例數(shù)據(jù))
├── man          (文件夾，包的詳細(xì)文檔)
├── tests        (文件夾，正式軟件測(cè)試代碼)
└── vignettes    (文件夾，含有高質(zhì)量的文檔)
    ├── biocOv1.Rmd
    ├── biocOv1.html

打包協(xié)議文檔“Writing R Extensions”提供了完整的詳細(xì)信息。使用R命令 R CMD build [foldername] 將對(duì)程度包文件夾的內(nèi)容進(jìn)行打包，從而創(chuàng)建可使用R CMD install[archivename] 添加到R安裝的文件。不過(guò)現(xiàn)在都使用R Studio來(lái)完成這些任務(wù)。

模塊化與包的相互依賴(lài)

打包協(xié)議可以幫助我們隔離執(zhí)行有限操作集的軟件，并識(shí)別隨時(shí)間固有變化的程序集合的版本。沒(méi)有客觀(guān)的方法來(lái)確定一組操作是否適合打包。一些非常有用的包僅執(zhí)行少量任務(wù)，而其他包則具有非常廣泛的應(yīng)用范圍。重要的是打包的概念允許軟件的模塊化。這在兩個(gè)方面很重要：范圍和時(shí)間。范圍的模塊化對(duì)于允許并行獨(dú)立開(kāi)發(fā)解決不同問(wèn)題的軟件工具非常重要。時(shí)間上的模塊化對(duì)于識(shí)別行為穩(wěn)定的軟件版本是很重要的。

持續(xù)集成：測(cè)試包的正確性和互操作性

下列的圖片就是開(kāi)發(fā)Bioconductor分支的一個(gè)構(gòu)建報(bào)告的快照：

image

在上面的表格中，一共有6列，其中最后1列標(biāo)明了“Installed Pkgs”，其中2000表示用于Linux平臺(tái)，這個(gè)數(shù)據(jù)在不同平臺(tái)之間有所不同，并且通常會(huì)伴隨著開(kāi)發(fā)分支的時(shí)間而增加。

匯總

Bioconductor的核心開(kāi)發(fā)者小組致力于開(kāi)發(fā)數(shù)據(jù)結(jié)構(gòu)，使用戶(hù)能夠方便地處理基因組和基因組規(guī)模的數(shù)據(jù)。用于設(shè)計(jì)的結(jié)構(gòu)是為了支持基因組規(guī)模生物學(xué)實(shí)驗(yàn)的主要目標(biāo)：

  • 解析由微陣列或測(cè)序儀產(chǎn)生的大規(guī)模數(shù)據(jù)集。
  • 對(duì)（相對(duì)）原始數(shù)據(jù)進(jìn)行預(yù)處理，以支持可靠的統(tǒng)計(jì)解釋。
  • 將分析量化與樣本信息數(shù)據(jù)進(jìn)行結(jié)合，用于檢驗(yàn)分子過(guò)程和生物體水平特征（例如生長(zhǎng)、疾病狀態(tài)）之間關(guān)系的假設(shè)。
  • 在本課程里，我們將回顧你可以在自己的研究中用于執(zhí)行這些相關(guān)任務(wù)的無(wú)明和功能。

5x的基因前提和概述

你知道如何使用R來(lái)操作和分析數(shù)據(jù)后，并且還對(duì)統(tǒng)計(jì)建模有不錯(cuò)的理解后。Bioconductor項(xiàng)目就會(huì)證明，在進(jìn)行基因組級(jí)規(guī)模的計(jì)算生物操作時(shí)（但非所有），R是一個(gè)有效的工具。

將Bioconductor與處理基因組數(shù)據(jù)的其它軟件系統(tǒng)的不同之處在于：

  • 使用面向?qū)ο蟮脑O(shè)計(jì)理解統(tǒng)一了基因組實(shí)驗(yàn)中出現(xiàn)的不同數(shù)據(jù)類(lèi)型；
  • 致力于基因組注釋的可互操作結(jié)構(gòu)，從核苷酸到人口規(guī)模；
  • 發(fā)布和開(kāi)發(fā)周期的持續(xù)集成原則，在多個(gè)廣泛使用的計(jì)算平臺(tái)上進(jìn)行日常測(cè)試。我們這個(gè)為期四周的模塊化學(xué)習(xí)目標(biāo)旨在針對(duì)基因組規(guī)模數(shù)據(jù)分析的方方面面建立起對(duì)該系統(tǒng)使用的專(zhuān)業(yè)能力。

525.5x課程主要為了四個(gè)主要部分，每部分花費(fèi)一周時(shí)間來(lái)學(xué)習(xí)，內(nèi)容如下：

  • 動(dòng)機(jī)和技術(shù)：我們?yōu)槭裁匆獧z測(cè)，以及如何檢測(cè)數(shù)據(jù)，以及我們?nèi)绾问褂肦語(yǔ)言來(lái)管理數(shù)據(jù)。
  • 基因組注釋?zhuān)河绕涫且⒁饣蚪M坐標(biāo)中的區(qū)間（ranges，有的譯為“范圍”）在識(shí)別基因組結(jié)構(gòu)方面的作用。
  • 基因組級(jí)數(shù)據(jù)的預(yù)處理概念，重要如何通過(guò)Bioconductor來(lái)實(shí)現(xiàn)。
  • 使用Bioconductor來(lái)做基因組級(jí)數(shù)據(jù)的假設(shè)檢驗(yàn)。

Some of the fundamental concepts that distinguish Bioconductor from other software systems addressing genome-scale data are

本章的不同小節(jié)將簡(jiǎn)要描述一下這些概念，以及說(shuō)明如何進(jìn)行計(jì)算。

動(dòng)機(jī)與技術(shù)

“What we measure and why”這個(gè)視頻（這個(gè)視頻位于Youtube上）給出了基本生物過(guò)程的示意圖，我們可以通過(guò)計(jì)算來(lái)研究這些過(guò)程。我們注意到，對(duì)有機(jī)體的生命過(guò)程至關(guān)重要的所有蛋白質(zhì)的編碼序列都位于有機(jī)體的基因組DNA中。研究生物體之間的差異， 以及生物體內(nèi)的某些變化（例如腫瘤的發(fā)展），通常涉及基因組DNA序列的計(jì)算。

Bioconductor提供了處理許多生物體基因組DNA序列的工具。一種處理序列的基本方法計(jì)算就是給定一個(gè)“參考序列”，然后我們計(jì)算另外一條序列與參考序列之間的差異。

參考序列的獲取

使用Bioconductor很容易處理人類(lèi)(Homo sapiens)的參考序列?，F(xiàn)在我們來(lái)看一下17號(hào)染色體，如下所示：

library(BSgenome.Hsapiens.UCSC.hg19)
Hsapiens$chr17
##   81195210-letter "DNAString" instance
## seq: AAGCTTCTCACCCTGTTCCTGCATAGATAATTGC...GGTGTGGGTGTGGTGTGTGGGTGTGGGTGTGGT

我們注意到：

  • 我們通過(guò)了一個(gè)R包來(lái)提供了序列。
  • 包的名稱(chēng)表明了這個(gè)是來(lái)源于UCSC的參考基因組hg19。
  • 我們使用美元符號(hào) $ 提取了染色體的序列。

表示DNA突變體

與參考序列有單個(gè)偏離的標(biāo)準(zhǔn)表示方法就是使用VCF格式的文件（ Variant Call Format）。VariantAnnotation包就含有這樣一個(gè)案例。我們有兩個(gè)某些DNA突變體的兩個(gè)高質(zhì)量表示方法(high-level representations)，即VCF內(nèi)容中的摘要，以及突變體本身在基因組上的地址，如下所示：

fl <- system.file("extdata", "ex2.vcf", package="VariantAnnotation") 
vcf <- readVcf(fl, "hg19")
vcf
## class: CollapsedVCF 
## dim: 5 3 
## rowRanges(vcf):
##   GRanges with 5 metadata columns: paramRangeID, REF, ALT, QUAL, FILTER
## info(vcf):
##   DataFrame with 6 columns: NS, DP, AF, AA, DB, H2
## info(header(vcf)):
##       Number Type    Description                
##    NS 1      Integer Number of Samples With Data
##    DP 1      Integer Total Depth                
##    AF A      Float   Allele Frequency           
##    AA 1      String  Ancestral Allele           
##    DB 0      Flag    dbSNP membership, build 129
##    H2 0      Flag    HapMap2 membership         
## geno(vcf):
##   SimpleList of length 4: GT, GQ, DP, HQ
## geno(header(vcf)):
##       Number Type    Description      
##    GT 1      String  Genotype         
##    GQ 1      Integer Genotype Quality 
##    DP 1      Integer Read Depth       
##    HQ 2      Integer Haplotype Quality
rowRanges(vcf)
## GRanges object with 5 ranges and 5 metadata columns:
##                  seqnames             ranges strand | paramRangeID
##                     <Rle>          <IRanges>  <Rle> |     <factor>
##        rs6054257       20 [  14370,   14370]      * |         <NA>
##     20:17330_T/A       20 [  17330,   17330]      * |         <NA>
##        rs6040355       20 [1110696, 1110696]      * |         <NA>
##   20:1230237_T/.       20 [1230237, 1230237]      * |         <NA>
##        microsat1       20 [1234567, 1234569]      * |         <NA>
##                             REF                ALT      QUAL      FILTER
##                  <DNAStringSet> <DNAStringSetList> <numeric> <character>
##        rs6054257              G                  A        29        PASS
##     20:17330_T/A              T                  A         3         q10
##        rs6040355              A                G,T        67        PASS
##   20:1230237_T/.              T                           47        PASS
##        microsat1            GTC             G,GTCT        50        PASS
##   -------
##   seqinfo: 1 sequence from hg19 genome

我們需要注意：

  • 演示中的數(shù)據(jù)已經(jīng)包含在了包中，它僅用于說(shuō)明和測(cè)試
  • 變量 vcf 向用戶(hù)簡(jiǎn)潔地展示了信息
  • 使用rowRanges函數(shù)提取了突變體在hg19基因組中的位置，并與標(biāo)簽一同顯示

檢測(cè)基因表達(dá)

我們將以觀(guān)察模式生物釀酒酵母(Sacchomyces Cerevisiae)中基因表達(dá)的檢測(cè)數(shù)值來(lái)結(jié)束這一部分的討論?，F(xiàn)在我們來(lái)看一個(gè)實(shí)驗(yàn)，這個(gè)實(shí)驗(yàn)檢測(cè)了釀酒酵母生殖周期中一系列時(shí)間點(diǎn)上的全基因組mRNA豐度。我們會(huì)使用R包來(lái)管理這些數(shù)據(jù)，并且我們使用特殊的Bioconductor定義 的數(shù)據(jù)結(jié)構(gòu)來(lái)訪(fǎng)問(wèn)有關(guān)的實(shí)驗(yàn)和結(jié)果信息，如下所示：

library(yeastCC)
data(spYCCES)
spYCCES
## ExpressionSet (storageMode: lockedEnvironment)
## assayData: 6178 features, 77 samples 
##   element names: exprs 
## protocolData: none
## phenoData
##   sampleNames: cln3_40 cln3_30 ... elu_390 (77 total)
##   varLabels: syncmeth time
##   varMetadata: labelDescription
## featureData: none
## experimentData: use 'experimentData(object)'
##   pubMedIds: 9843569 
## Annotation:
experimentData(spYCCES)
## Experiment data
##   Experimenter name: Spellman PT 
##   Laboratory: Department of Genetics, Stanford University Medical Center, Stanford, California 94306-5120, USA. 
##   Contact information:  
##   Title: Comprehensive identification of cell cycle-regulated genes of the yeast Saccharomyces cerevisiae by microarray hybridization. 
##   URL:  
##   PMIDs: 9843569 
## 
##   Abstract: A 150 word abstract is available. Use 'abstract' method.

經(jīng)過(guò)一段時(shí)間的學(xué)習(xí)，你就會(huì)在接下來(lái)的幾周內(nèi)成為這方面的高手，你可以將整個(gè)細(xì)胞周期中參與調(diào)控的基因隨時(shí)間變化的過(guò)程繪制出來(lái)。

plot of chunk lkycc2

需要注意：

  • 關(guān)于實(shí)驗(yàn)的元數(shù)據(jù)已經(jīng)綁定到R包的數(shù)據(jù)中（通過(guò) experimentData可以可看PubmedID和摘要）；
  • 我們使用簡(jiǎn)單的語(yǔ)法就能選擇這些復(fù)雜實(shí)驗(yàn)設(shè)計(jì)中的信息；在這個(gè)案例中，我們使用spYCCES[, spYCCES$syncmeth=="alpha"] 就能提取出控制alpha pheromone的基因；
  • R的繪圖工具支持通用的繪圖注釋與增強(qiáng)功能 (plot annotation and enhancemetn)；
  • 統(tǒng)計(jì)模型工具能夠直接幫我們區(qū)分周期與非周期的基因。

結(jié)束語(yǔ)

你即將學(xué)習(xí)一些關(guān)于基因組結(jié)構(gòu)和檢測(cè)它們的分子生物學(xué)技術(shù)的高水平講座。當(dāng)你遇到這些概念時(shí)，請(qǐng)你記住那些與理解結(jié)構(gòu)以及數(shù)據(jù)處理過(guò)程相關(guān)的計(jì)算思路與方法。在Bioconductor中找到相應(yīng)的計(jì)算工具，并熟練掌握這些工具。如果你找不到它們，請(qǐng)通知我們，或者是如果滿(mǎn)足某些需求的工具不存在，你也可以選擇構(gòu)建它們。

參考資料

  1. HarvardX Biomedical Data Science Open Online Training
  2. Bioconductor for genome-scale data -- quick intro
?著作權(quán)歸作者所有,轉(zhuǎn)載或內(nèi)容合作請(qǐng)聯(lián)系作者
【社區(qū)內(nèi)容提示】社區(qū)部分內(nèi)容疑似由AI輔助生成，瀏覽時(shí)請(qǐng)結(jié)合常識(shí)與多方信息審慎甄別。
平臺(tái)聲明：文章內(nèi)容（如有圖片或視頻亦包括在內(nèi)）由作者上傳并發(fā)布，文章內(nèi)容僅代表作者本人觀(guān)點(diǎn)，簡(jiǎn)書(shū)系信息發(fā)布平臺(tái)，僅提供信息存儲(chǔ)服務(wù)。

相關(guān)閱讀更多精彩內(nèi)容

友情鏈接更多精彩內(nèi)容