% \ VignetteIndexEntry{07。注释和范围练习}%\VignetteEngine{knitr::knitr} \documentclass{article} < >= options(max.print=1000) BiocStyle::latex() library(knitr) opts_chunk$set(cache=TRUE, tidy=FALSE) @ < >= suppressPackageStartupMessages({library(ShortRead) library(VariantAnnotation) library(parallel);options(mc.cores=detectCores())库(ggplot2)库(rnaseqdata . hnrnpc .bam. ch14)库(org.Hs.eg.db)库(TxDb.Hsapiens.UCSC.hg19.knownGene)库(BSgenome.Hsapiens.UCSC.hg19)库(AnnotationHub)库(rtracklayer)}) @ \title{实用:注释}\作者{Martin Morgan (mtmorgan@fhcrc.org)} \日期{2014年2月27-28日}\newcommand{\Hsap}{\emph{H. cn}~智人}}\ newcommand {\ Dmel} {\ emph {D。~melanogaster}} \usepackage{Exercise} \begin{document} \maketitle \tableofcontents \section{基因注释}\分段{数据包}有机体级(' org')包包含一个中心标识符(例如，Entrez基因id)和其他标识符(例如，GenBank或Uniprot登录号，RefSeq id等)之间的映射。org包的名称始终是\texttt{org的形式。 ． .db}(例如\Biocannopkg{org.Sc.sgd.db}) where \texttt{ }是生物体的两个字母缩写(例如\texttt{Sc}代表\emph{Saccharomyces cerevisiae})和\texttt{ }是描述中心标识符类型的缩写(小写)(例如\texttt{sgd}表示\emph{Saccharomyces}基因组数据库分配的基因标识符，\texttt{eg}表示Entrez基因标识符)。\Biocpkg{AnnotationDbi}包(org包只是一种“。db”注释包)中的“如何使用“。db”注释包”是一个关键参考。' .db'和大多数其他\Bioconductor{}注释包每6个月更新一次。注释包通常包含一个以包本身命名的对象。这些对象统称为\Rclass{AnnotationDb}对象，更具体的类命名为\Rclass{OrgDb}、\Rclass{ChipDb}或\Rclass{TranscriptDb}对象。可以应用于这些对象的方法包括\Rfunction{cols}、\Rfunction{keys}、\Rfunction{keytypes}和\Rfunction{select}。表~\ref{tab:select-ops}中总结了检索注释的常用操作。\begin{table} \居中\标题{检索和操作注释的常用操作。}{选项卡:select-ops} \ \标签开始{表格}{微光}类别和功能和描述\ \ \线\ noalign {\ smallskip}发现& \ Rfunction{列}&列表的列可以返回\ \ & \ Rfunction {keytypes} &列列表可以用作键\ \ & \ Rfunction{键}&列表值,可以预期给定keytype \ \ & \ Rfunction{选择}&匹配检索注释\ Rcode}{键,\ Rcode {keytype}和\ Rcode{列}\ \操作& \ Rfunction {setdiff}, \ Rfunction{联盟},\Rfunction{intersect} &对集合\\ & \Rfunction{duplicate}， \Rfunction{unique} &标记或删除重复的\\ & \Rfunction{\%in\%}， \Rfunction{match} &查找匹配\\ & \Rfunction{any}， \Rfunction{all} &是否有任何\Rcode{TRUE}? Are all? \\ & \Rfunction{merge} & Combine two different \Robject{data.frames} based on shared keys \\ \Rclass{GRanges*} & \Rfunction{transcripts}, \Rfunction{exons}, \Rfunction{cds} & Features (transcripts, exons, coding sequence) as \Rclass{GRanges}. \\ & \Rfunction{transcriptsBy} , \Rfunction{exonsBy} & Features group by gene, transcript, etc., as \Rclass{GRangesList}.\\ & \Rfunction{cdsBy}\\ \hline \end{tabular} \end{table} \begin{Exercise} This exercise illustrates basic use of the `select' interface to annotation packages. \begin{enumerate} \item What is the name of the org package for \emph{Homo sapiens}? Load it. Display the \Rclass{OrgDb} object for the \Biocpkg{org.Hs.eg.db} package. Use the \Rfunction{columns} method to discover which sorts of annotations can be extracted from it. \item Use the \Rfunction{keys} method to extract ENSEMBL identifiers and then pass those keys in to the \Rfunction{select} method in such a way that you extract the SYMBOL (gene symbol) and GENENAME information for each. Use the following ENSEMBL ids. \end{enumerate} < >= ensid <- c("ENSG00000130720"， "ENSG00000103257"， "ENSG00000156414"， "ENSG00000144644"， "ENSG00000159307"， "ENSG00000144485") @ \end{练习}\begin{解决方案}\Rclass{OrgDb}对象命名为\Rcode{org.Hs.eg.db}。<>= library(org.Hs.eg.db) keytypes(org.Hs.eg.db) columns(org. hs .db) cols <- c("SYMBOL"， "GENENAME") select(org.Hs.eg.db, keys=ensids, columns=cols, keytype="ENSEMBL") @ \end{解决方案}\分段{互联网资源}select \Bioconductor{}支持基于web查询的软件包的简短摘要在表~\ref{tab:webannoservices}中。\begin{table} \居中\标题{查询基于web的注释服务的选定包。} \label{tab:webannoservices} \begin{tabular}{ll} Package & Description\\ hline\noalign{\smallskip} \Biocpkg{AnnotationHub} & Ensembl, Encode, dbSNP, UCSC数据对象\\ \Biocpkg{biomaRt} & \url{http://biomart.org}， Ensembl和其他注释\\ \Biocpkg{PSICQUIC} & \url{https://code.google.com/p/psicquic.org}，蛋白质相互作用\\ \Biocpkg{uniprot。ws} & \url{http://uniprot.org}，蛋白标注\\ \Biocpkg{KEGGREST} & \url{http://www.genome.jp/kegg}， KEGG通路\\ \Biocpkg{SRAdb} & \url{http://www.ncbi.nlm.nih.gov/sra}，测序实验。\\ \Biocpkg{rtracklayer} & \url{http://genome.ucsc.edu}，基因组跟踪。\\ \Biocpkg{GEOquery} & \url{http://www.ncbi.nlm.nih.gov/geo/}，阵列和其他数据\\ \Biocpkg{ArrayExpress} & \url{http://www.ebi.ac.uk/arrayexpress/}，阵列和其他数据\\ \hline \end{tabular} \end{table} \段{Using biomaRt} \Biocpkg{biomaRt}包提供对在线\href{http://www.biomart.org}{biomart}资源的访问。它由几个数据库资源组成，称为“marts”。每个集市允许访问多个数据集;\Biocpkg{biomaRt}包提供了用于mart和数据集发现的方法，以及用于检索数据的标准方法\Rfunction{getBM}。\begin{练习}\警告{此练习需要INTERNET访问}\begin{enumerate} \item加载\Biocpkg{biomaRt}包并列出可用的marts。选择\emph{ensembl} mart并列出该mart的数据集。设置一个mart来使用\emph{ensembl} mart和\emph{hsapiens\_gene\_ensembl}数据集。 \item A \Biocpkg{biomaRt} dataset can be accessed via \Rfunction{getBM}. In addition to the mart to be accessed, this function takes filters and attributes as arguments. Use \Rfunction{filterOptions} and \Rfunction{listAttributes} to discover values for these arguments. Call \Rfunction{getBM} using filters and attributes of your choosing. \end{enumerate} \end{Exercise} \begin{Solution} <>= ##需要上网!!library(biomaRt) head(listMarts()， 3) ##列出marts head(listDatasets(useMart("ensembl"))， 3) ## mart数据集ensembl <- ##完全指定的mart useMart("ensembl"， dataset = "hsapiens_gene_ensembl") head(listFilters(ensembl)， 3) ## filters myFilter <- "chromosome_name" head(filterOptions(myFilter, ensembl)， 3) ##返回值myValues <- c("21"， "22") head(listAttributes(ensembl)，3) ## attributes myAttributes <- c("ensembl_gene_id"，"chromosome_name") ##组装和查询集市res <- getBM(attributes = myAttributes, filters = myFilter, values = myValues, mart = ensembl) @使用\Rcode{head(res)}查看结果。\end{解决方案}\begin{练习}作为一个可选练习，注释在DESeq2实验室中差异表达的基因，例如，找到与五个差异表达最多的基因相关的\texttt{GENENAME}。这些在生物学上有意义吗?你能将注释结果与“top table”结果合并，以提供统计和生物学上有信息的摘要吗?\end{练习}\section{基因组注释}有多种可用来表示大基因组的包和类。其中包括:\begin{description} \item [\Rpackage{TxDb。*}]用于转录本和其他基因组/坐标注释。\item [\Biocpkg{BSgenome}]用于全基因组表示。 See \Rfunction{available.packages} for pre-packaged genomes, and the vignette `How to forge a BSgenome data package' in the \item [\Biocannopkg{Homo.sapiens}] For integrating \Rpackage{TxDb*} and \Rpackage{org.*} packages. \item [\Rpackage{SNPlocs.*}] For model organism SNP locations derived from dbSNP. \item [\Rfunction{FaFile}] (\Biocpkg{Rsamtools}) for accessing indexed FASTA files. \item [\Rpackage{SIFT.*}, \Rpackage{PolyPhen}, \Rpackage{ensemblVEP}] Variant effect scores. \end{description} \subsection{Transcript annotation packages} Genome-centric packages are very useful for annotations involving genomic coordinates. It is straight-forward, for instance, to discover the coordinates of coding sequences in regions of interest, and from these retrieve corresponding DNA or protein coding sequences. Other examples of the types of operations that are easy to perform with genome-centric annotations include defining regions of interest for counting aligned reads in RNA-seq experiments and retrieving DNA sequences underlying regions of interest in ChIP-seq analysis, e.g., for motif characterization. \subsection{\Rpackage{rtracklayer}} \label{subsec:rtracklayer} The \Biocpkg{rtracklayer} package allows us to query the UCSC genome browser, as well as providing \Rfunction{import} and \Rfunction{export} functions for common annotation file formats like GFF, GTF, and BED. \begin{Exercise} \warning{This exericse requires INTERNET ACCESS} Here we use \Biocpkg{rtracklayer} to retrieve estrogen receptor binding sites identified across cell lines in the ENCODE project. We focus on binding sites in the vicinity of a particularly interesting region of interest. \begin{enumerate} \item Define our region of interest by creating a \Rclass{GRanges} instance with appropriate genomic coordinates. Our region corresponds to 10Mb up- and down-stream of a particular gene. \item Create a session for the UCSC genome browser \item Query the UCSC genome browser for ENCODE estrogen receptor ERalpha$_a$ transcription marks; identifying the appropriate track, table, and transcription factor requires biological knowledge and detective work. \item Visualize the location of the binding sites and their scores; annotate the mid-point of the region of interest. \end{enumerate} \end{Exercise} \begin{Solution} Define the region of interest <>= roi <- GRanges("chr10"， IRanges(92106877, 112106876, names="ENSG00000099194")) @ \noindent创建会话<->= library(rtracklayer) session <- browserSession() @ \noindent在我们感兴趣的区域查询UCSC中的特定轨道，表和转录因子<>= trackName <- "wgEncodeRegTfbsClusteredV2" tableName <- "wgEncodeRegTfbsClusteredV2" trFactor <- "ERalpha_a" ucscTable <- getttable (ucscTableQuery(session, track=trackName, range=roi, table=tableName, name=trFactor)) @ \noindent显示结果<>= plot(score ~ chromStart, ucscTable, pch="+") abline(v=start(roi) + (end(roi) - start(roi) + 1) / 2, col="blue") @ \end{解决方案}\end{文档}