Skip to content

Annotating_Genomic_Ranges #19

New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

Sign up for GitHub

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Open
lixin856 wants to merge 2 commits into
base: master
Choose a base branch
from
Open

Annotating_Genomic_Ranges #19

Changes from all commits
Commits
Show all changes
2 commits
Select commit Hold shift + click to select a range
321ef9a
Merge pull request #1 from openbiox/master
lixin856 Jun 19, 2019
77295a2
注释基因组区域
lixin856 Jul 6, 2019
File filter

Filter by extension

Filter by extension
Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
198 changes: 99 additions & 99 deletions BiocWorkflow_todo/05annotation_en/Annotating_Genomic_Ranges.Rmd
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -12,7 +12,7 @@ vignette: >
%\VignetteEngine{knitr::rmarkdown}
---

# Version Info
# 版本信息
```{r, echo=FALSE, results="hide", warning=FALSE}
suppressPackageStartupMessages({
library('annotation')
Expand All @@ -26,66 +26,66 @@ library('annotation')
**Package version**: `r packageVersion("annotation")`
</p>

# Background
# 背景

Bioconductor can import diverse sequence-related file types, including fasta,
fastq, BAM, VCF, gff, bed, and wig files, among others. Packages support common
and advanced sequence manipulation operations such as trimming, transformation,
and alignment. Domain-specific analyses include quality assessment, ChIP-seq,
differential expression, RNA-seq, and other approaches. Bioconductor includes
an interface to the Sequence Read Archive (via the
[SRAdb](/packages/release/bioc/html/SRAdb.html) package).
Bioconductor可以导入各种序列相关的文件类型，
包括 fasta、fastq、BAM、VCF、gff、bed 和 wig 文件等。
Bioconductor 包支持常规和高级的序列操作，例如 trimming，转换和对齐；
也支持特定领域的分析，包括质量评估、ChIP-seq、差异表达、RNA-seq等。
Bioconductor包含序列读取存档的接口（通过
[SRAdb](/ packages / release / bioc / html / SRAdb.html)包）。

This workflow walks through the annotation of a generic set of ranges with
Bioconductor packages. The ranges can be any user-defined region of interest or
can be from a public file.

# Data Preparation
此工作流程将介绍用 Bioconductor 包进行通用区域的注释，
区域可以是任何用户定义的感兴趣区域，也可以来自公共文件。


# 数据准备

## Human hg19

As a first step, data are put into a GRanges object so we can take
advantage of overlap operations and store identifiers as
metadata columns.
第一步，将数据放入 GRanges 对象中，以便利用
overlap operations，以及存储标识符作为元数据列。


第一组区域是来自 dbSNP 的 VCF 文件。该文件可以从 NCBI 的 ftp 站点下载
[ftp://ftp.ncbi.nlm.nih.gov/snp/](ftp://ftp.ncbi.nlm.nih.gov/snp/），
然后使用 VariantAnnotation 包中的 readVcf() 导入；或者作为 AnnotationHub
中预解析的 VCF 对象使用。

The first set of ranges are variants from a dbSNP Variant Call Format (VCF)
file. This file can be downloaded from the ftp site at NCBI
[ftp://ftp.ncbi.nlm.nih.gov/snp/](ftp://ftp.ncbi.nlm.nih.gov/snp/) and
imported with readVcf() from the VariantAnnotation package. Alternatively,
the file is available as a pre-parsed VCF object in the AnnotationHub.

```{r, echo=FALSE}
suppressPackageStartupMessages(library(annotation))
suppressPackageStartupMessages(library(annotation)) #导入"annotation"包，屏蔽导入时的信息输出。（译者注：如果之前尚未安装 annotation 包，可参考 https://bioconductor.org/install/ 安装）
```


The Hub returns a VcfFile object with a reference to the file on disk.
返回一个引用磁盘上文件的 VcfFile 对象--hub。

```{r}
hub <- AnnotationHub()
```

Query the Hub for clinvar VCF files build GRCh37:
查询 hub，返回 clinvar VCF 文件 GRCh37 版本的 sourceurl（vcf文件的原始URL）：

```{r}
mcols(query(hub, "clinvar.vcf", "GRCh37"))[,"sourceurl", drop=FALSE]
```

Retrieve one of the files:
检索其中一个文件：

```{r}
fl <- query(hub, "clinvar.vcf", "GRCh37")[[1]]
```

Read the data into a VCF object:
将数据读入 VCF 对象：

```{r}
vcf <- readVcf(fl, "hg19")
dim(vcf)
```

Overlap operations require that seqlevels and the genome of the objects match.
Here the VCF seqlevels are modified to match the TxDb.
Overlap operations 要求 seqlevels 和对象的基因组匹配。
这里修改 VCF seqlevels 以匹配 TxDb。


```{r}
Expand All @@ -95,64 +95,64 @@ seqlevels(vcf)
seqlevels(vcf) <- paste0("chr", seqlevels(vcf))
```

For this example we'll annotate chromosomes 3 and 18:
这里我们以注释染色体 3 和 18 为例：

```{r}
seqlevels(vcf, pruning.mode="coarse") <- c("chr3", "chr18")
seqlevels(txdb_hg19) <- c("chr3", "chr18")
```

Sanity check to confirm we have matching seqlevels.
检查确认 seqlevels 是匹配的。
```{r}
intersect(seqlevels(txdb_hg19), seqlevels(vcf))
```

The genomes already match so no change is needed.
基因组已经匹配，因而无需进行任何更改。

```{r}
unique(genome(txdb_hg19))
unique(genome(vcf))
```

使用'rowRanges()'提取 VCF 对象中的 GRanges。

The GRanges in a VCF object is extracted with 'rowRanges()'.
```{r}
gr_hg19 <- rowRanges(vcf)
```

## Mouse mm10

The second set of ranges is a user-defined region of chromosome 4 in mouse.
The idea here is that any region, known or unknown, can be annotated with
the following steps.
第二组注释区域是用户自定义的小鼠染色体 4 的区域。
这里是想说明，任何已知或未知的区域都可以通过以下
步骤进行注释。

Load the TxDb package and keep only the standard chromosomes.
加载 TxDb 包且只保留标准染色体。
```{r}
txdb_mm10 <- keepStandardChromosomes(TxDb.Mmusculus.UCSC.mm10.ensGene)
```
We are creating the GRanges from scratch and can specify the seqlevels
(chromosome names) to match the TxDb.

我们从头开始创建 GRanges 并指定 seqlevels（染色体名称）以匹配 TxDb。
```{r}
head(seqlevels(txdb_mm10))
gr_mm10 <- GRanges("chr4", IRanges(c(4000000, 107889000), width=1000))
```

Now assign the genome.
分配基因组。
```{r}
unique(genome(txdb_mm10))
genome(gr_mm10) <- "mm10"
```

# Location in and Around Genes
# 区域内及周围的基因

locateVariants() in the VariantAnnotation package annotates ranges
with transcript, exon, cds and gene ID's from a TxDb. Various
extractions are performed on the TxDb (exonsBy(), transcripts(),
cdsBy(), etc.) and the result is overlapped with the ranges. An appropriate
GRangesList can also be supplied as the annotation. Different variants
such as 'coding', 'fiveUTR', 'threeUTR', 'spliceSite', 'intron',
'promoter', and 'intergenic' can be searched for by passing the appropriate
constructor as the 'region' argument. See ?locateVariants for details.
VariantAnnotation 包中的 locateVariants() 使用来自 TxDb
的 transcript，exon，cds 和 gene ID 来注释基因区域。可以
对 TxDb 执行各种提取（exonsBy()，transcript()，cdsBy()等
），注释结果与提取范围对应。适当的 GRangesList 也可以提供
作为注释。通过将适当的构造函数作为 “region” 的参数传递，
可以搜索不同的变体，例如 'coding'，'fiveUTR'，'threeUTR'
，'spliceSite'，'intron'，'promoter' 和 'intergenic'。详
细信息，请参阅 ?locateVariants。

```{r}
loc_hg19 <- locateVariants(gr_hg19, txdb_hg19, AllVariants())
Expand All @@ -161,10 +161,10 @@ loc_mm10 <- locateVariants(gr_mm10, txdb_mm10, AllVariants())
table(loc_mm10$LOCATION)
```

# Annotate by ID
# 通过 ID 注释

The ID's returned from locateVariants() can be used in select() to map
to ID's in other annotation packages.
从 locateVariants() 返回的 ID 可以在 select() 中用于
映射到其他注释包中的 ID。


```{r}
Expand All @@ -173,25 +173,25 @@ keys <- na.omit(unique(loc_hg19$GENEID))
head(select(org.Hs.eg.db, keys, cols, keytype="ENTREZID"))
```

The 'keytype' argument specifies that the mouse TxDb contains
Ensembl instead of Entrez gene id's.
参数 'keytype' 指定小鼠 TxDb 包含 Ensembl 而不是
Entrez gene id。


```{r}
keys <- unique(loc_mm10$GENEID)
head(select(org.Mm.eg.db, keys, cols, keytype="ENSEMBL"))
```

# Annotate by Position
# 通过位置注释

Files stored in the AnnotationHub have been pre-processed into
ranged-based R objects such as a GRanges, GAlignments and VCF.
The positions in our GRanges can be overlapped with the ranges in
the AnnotationHub files. This allows for easy subsetting
of multiple files, resulting in only the ranges of interest.
存储在 AnnotationHub 中的文件已经预处理为基于
范围的 R 对象，例如 GRanges，GAlignments 和
VCF。GRanges 中的位置可以与 AnnotationHub 文
件中的范围重叠，这对多个文件取子集很容易，从
而仅产生感兴趣的范围。

Create a 'hub' from AnnotationHub and filter the files based
on organism and genome build.
从 AnnotationHub 创建一个 'hub' ，并根据生物
和基因组版本过滤文件。

```{r}
hub <- AnnotationHub()
Expand All @@ -200,29 +200,29 @@ hub_hg19 <- subset(hub,
length(hub_hg19)
```

Iterate over the first 3 files and extract ranges that overlap 'gr_hg19'.
迭代前 3 个文件并提取与 'gr_hg19' 重叠的范围。
```{r, echo=FALSE}
ov_hg19 <- lapply(1:3, function(i) subsetByOverlaps(hub_hg19[[i]], gr_hg19))
```
```{r}
ov_hg19 <- lapply(1:3, function(i) subsetByOverlaps(hub_hg19[[i]], gr_hg19))
```

Inspect the results.
检查结果。
```{r}
names(ov_hg19) <- names(hub_hg19)[1:3]
lapply(ov_hg19, head, n=3)
```
Annotating the mouse ranges in the same fashion is left
as an exercise.
作为练习，以相同的方式注释小鼠范围。


# Annotating Variants
# 注释变体

<h4 id=amino-acid-coding-changes">氨基酸编码变化</h4>
对于落在编码区中的一组 dbSNP 变体，可以计算氨基酸变化。
一个变体-转录本的匹配输出一行，这导致每个变体可能输出
多行。

<h4 id=amino-acid-coding-changes">Amino acid coding changes</h4>
For the set of dbSNP variants that fall in coding regions, amino
acid changes can be computed. The output contains one line for each
variant-transcript match which can result in multiple lines for
each variant.


```{r}
Expand All @@ -233,34 +233,34 @@ head(predictCoding(vcf, txdb_hg19, Hsapiens), 3)
sessionInfo()
```

# Exercises

Exercise 1: VCF header and reading data subsets.

VCF files can be large and it's often the case that only a subset of
variables or genomic positions are of interest. The scanVcfHeader()
function in the VariantAnnotation package retrieves header information
from a VCF file. Based on the information returned from scanVcfHeader()
a ScanVcfParam() object can be created to read in a subset of data from
a VCF file.
* Use scanVcfHeader() to inspect the header information in the
'chr22.vcf.gz' file in VariantAnnotation package.
* Select a few 'info' or 'geno' variables and create a ScanVcfParam object.
* Use the ScanVcfParam object as the 'param' argument to readVcf()
to read in a subset of data.
Note that the header() accessor operates on VCF objects in the R
workspace. Try header(vcf) on the dbSNP file from AnnotationHub.

Exercise 2: Annotate the mouse ranges in 'gr_mm10' with AnnotationHub files.
* Create a new 'hub' and filter on organism.
* Isolate the files for the appropriate genome build and perform overlaps.
# 练习

练习 1 ：VCF 头部和读取数据子集。

VCF 文件可能很大，而且通常情况下只有一部分变量或
基因组位置是有意义的。 VariantAnnotation 包中的
scanVcfHeader() 函数从 VCF 文件中检索标头信息。根
据 scanVcfHeader() 返回的信息，可以创建 ScanVcfParam()
对象以从 VCF 文件中读取数据子集。

* 使用 scanVcfHeader() 检查 VariantAnnotation 包中 'chr22.vcf.gz'
文件中的标题信息。
* 选择一些 'info' 或 'geno' 变量并创建一个 ScanVcfParam 对象。
* 使用 ScanVcfParam 对象作为 readvcf() 的 'param' 参数来读入
数据子集。
请注意，header() 访问器对 R 工作空间中的 VCF 对象进行操作。
在 AnnotationHub 的 dbSNP 文件中尝试 header(vcf)。

练习 2 ：使用 AnnotationHub 文件在 'gr_mm10' 中注释小鼠区域。
* 创建一个新的 “hub” 并过滤 organism。
* 从适当基因组版本提取出文件并执行 overlaps。

Exercise 3: Annotate a gene range from Saccharomyces Scerevisiae.
* Load TxDb.Scerevisiae.UCSC.sacCer3.sgdGene and extract the
gene ranges. (Hint: use transcriptsBy() and range()).
* Isolate the range for gene "YBL086C".
* Create a new 'hub' from AnnotationHub and filter by organism.
(You should see >= 39 files.)
* Select the files for 'sacCer3' and perform overlaps.
练习 3 ：注释来自 Saccharomyces Scerevisiae 的基因区域。
* 加载 TxDb.Scerevisiae.UCSC.sacCer3.sgdGene 并提取基因
范围。 （提示：使用 transcriptsBy() 和 range()）。
* 分离基因 “YBL086C” 的区域。
* 从 AnnotationHub 创建一个新的 “hub” 并按 organism 过滤。
（你应该看到 >= 39 个文件。）
* 选择 'sacCer3' 的文件并执行 overlaps。

<p class="back_to_top">[ <a href="#top">Back to top</a> ]</p>