参考基因组

 参考基因组

rockygao · 2019年12月09日 · 201 次阅读

 目

在学习进行人的 RNA-SEQ 数据分析之前，需要做许多准备工作：

  • 分析环境的搭建以及各个分析软件的安装；
  • 参考基因组的选择；
  • 待分析数据准备；
  • 分析流程

本文章主要介绍参考基因组选择、不同数据中参考基因下组载，以及参考基因组格式说明。

参考基因组及注释文件

  • 参考基因组文件
    通常是 fasta 或者 fasta 格式：以>开头的行标注染色体信息，后续行为该条染色体的碱基信息。
  • 基因注释文件
    通常有GTF(General Transfer Format)和GFF(general feature format)两种，其中 GFF 又可分为 GFF3 和 GFF2。这些格式大同小异，主要作用就是 注释基因组 。上述基因组文件中只有每条染色体的碱基序列，并没有各个基因信息；而 GTF/GFF 文件中包含的是每个转录本的信息。
    1
2
3
4
5
6
7
8
9
10
11
12
13

    # FA文件
>1 dna:chromosome chromosome:GRCh38:1:1:248956422:1 REF
NNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNN
NNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNN
NNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNN
NNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNN
NNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNN# GTF文件
chr1    hg38_ncbiRefSeq stop_codon      67093005        67093007        0.000000        -       .       gene_id "XM_011541469.1"; transcript_id "XM_011541469.1";
chr1    hg38_ncbiRefSeq CDS     67093008        67093604        0.000000        -       0       gene_id "XM_011541469.1"; transcript_id "XM_011541469.1";
chr1    hg38_ncbiRefSeq exon    67092176        67093604        0.000000        -       .       gene_id "XM_011541469.1"; transcript_id "XM_011541469.1";
chr1    hg38_ncbiRefSeq CDS     67095235        67095421        0.000000        -       1       gene_id "XM_011541469.1"; transcript_id "XM_011541469.1";

本流程所用参考基因组为 UCSC，下载方法如下：

1
2
3
4

# 参考基因组文件，UCSC下载
wget -c ftp://hgdownload.soe.UCSC.edu/goldenPath/hg38/bigZips/hg38.fa.gz
# -c 断点续传，避免因为网络中断导致下载失败。
# GEF文件下载见下文

1 参考基因组数据来源

参考基因组主要有三个数据库来源 ENSEMBL、NCBI、UCSC，还有一些特定数据库。各个数据库中格式不尽相同。

1.1 ENSEMBL 参考基因组

ENSMBL列出了各个不同物种的参考基因组和参考转录组，以表格的样式展示，各个文件整理的非常清楚，点击即可到达 FTP 下载页面。

常用物种数据界面 http://asia.ensembl.org/info/data/ftp/index.html
FTP 下载网址:ftp://ftp.ensembl.org/pub

Ensembl 提供的参考基因组有 2 种组装形式和 3 种重复序列处理方式, 分别是primary, toplevel和unmasked (dna)、oft-masked (dna_sm)和masked (dna_rm)。一般选择dna.primary或dna_sm.primary为后缀的文件。各个格式的含义见下文。

1
2
3
4
5

# 参考基因组 ensemb 下载
wget -c ftp://ftp.ensembl.org/pub/release-96/fasta/homo_sapiens/dna/Homo_sapiens.GRCh38.dna.primary_assembly.fa.gz# gff文件 ensemb 下载 
wget ftp://ftp.ensembl.org/pub/release-96/gtf/homo_sapiens/Homo_sapiens.GRCh38.96.gtf.gz

1.2 NCBI

NCBI 数据库好像没有整个基因组的 FA 文件，只有单个染色体的。下载页面如下：

ftp://ftp.ncbi.nlm.nih.gov/genomes/Homo_sapiens

NCBI 只有 gff3 格式下载

1

wget ftp://ftp.ncbi.nlm.nih.gov/genomes/Homo_sapiens/GFF/ref_GRCh38.p12_top_level.gff3.gz

1.3 UCSC

UCSC 下载参考基因组比较简单,如下：

1

wget -c ftp://hgdownload.soe.UCSC.edu/goldenPath/hg38/bigZips/hg38.fa.gz

UCSC 下载基因组注释文件比较麻烦，没有直接的 FTP 下载链接,有两种方式可以下载

  • 第一种方法,界面版

这种方式得到的 GTF 文件只有 gene_id 和 transcript_id，而没有其他信息,如:基因名，外显子。
首先，需要先进入http://genome.UCSC.edu/cgi-bin/hgTables然后进行选择，最后点击get output即可下载。
具体参数如下：
clade: Mammal
genome: human
assembly: Dec. 2013 (GRCh38/hg38)
group: Genes and Gene Predictions
track: NCBI RefSeq
table: UCSC RefSeq (refGene) 或者 RefSeq All (ncbiRefSeq)
当 table 种选择 UCSC RefSeq (refGene) 时与第二种方式基本得到 GTF 文件基本相似。但是这种方式得到的 GTF 文件最后一列只有 gene_id 和 transcript_id，而没有其他信息,如:基因名，外显子。而第二张方式的信息比较全面。

NCBI 中提供各个 gene_id 或 transcript_id 的详细信息文件，下载后然后手动匹配。

1

wget ftp://ftp.ncbi.nlm.nih.gov/gene/DATA/gene2refseq.gz
  • 第二种方式，FTP 下载

UCSC 提供一种非冗余的格式GenePred table format，文件较小，但是需要用genePredToGtf工具转换。
GenePred 格式说明参见下面：

https://genome.UCSC.edu/FAQ/FAQformat.html#format9
http://genomewiki.UCSC.edu/index.php/Genes_in_gtf_or_gff_format

1
2
3
4
5
6
7
8
9

# 下载转换工具
wget http://hgdownload.soe.UCSC.edu/admin/exe/linux.x86_64.v369/genePredToGtf#下载GenePred并转换
wget http://hgdownload.soe.UCSC.edu/goldenPath/hg38/database/refGene.txt.gz
gzip -d refGene.txt.gz
cut -f 2- refGene.txt > refGene.input
./genePredToGtf file refGene.input hg38refGene.gtf
cat hg38refGene.gtf  | sort -k1,1 -k4,4n > hg38refGene.gtf.sorted

1.4 illumina 给的参考基因组

illumina 官网也给了参考基因组的相关文件。

iGenomes

1.5 ENSEMBL、NCBI、UCSC 参考基因组版本对应关系

NCBI    UCSC    ENSEMBL
GRCh36  hg18   release_52.
GRCh37  hg19   release_59/61/64/68/69/75.
GRCh38  hg38   release_76/77/78/80/81/82.

2. 参考基因组格式说明

FA　基因序列
GTF　基因组注释,详细说明见GTF 文件 或者UCSC 对各种格式说明
GFF/GFF3 基因组注释，与 GTF 类似

2.1 参考基因组后缀说明

Primary or Toplevel

  • 参考基因组优先选择的 primary 版本，因为 toplevel 版本会包含 haplotype 信息，多余的信息会增加比对工具的工作，选择 primary 就可以
  • Primary assembly contains all toplevel sequence regions excluding haplotypes and patches. This file is best used for performing sequence similarity searches where patch and haplotype sequences would confuse analysis.

unmasked (dna)、soft-masked (dna_sm) 和 masked (dna_rm，Repeats masked)

  • Masked 基因组是指所有重复区和低复杂区被 N 代替的基因组序列，比对时就不会有 reads 比对到这些区域。一般不推荐用 masked 的基因组，因为它造成了信息的丢失，由此带来的一个问题是 uniquely 比对到 masked 基因组上的 reads 实际上可能不是 unique 的。而且 masked 基因组还会带来比对错误，使得在允许错配的情况下，本来来自重复区的 reads 比对到基因组的其它位置。另外检测重复区和低复杂区的软件不可能是完美的，这就造成遮盖住的重复序列和低复杂区并不一定是 100% 准确和敏感的。在NGS处理过程中，避免使用masked (dna_rm)
  • soft-masked 基因组是指把所有重复区和低复杂区的序列用小写字母标出的基因组，由于主要的比对软件，比如 BWA、bowtie2 等都忽略这些 soft-mask，直接把小写字母当做大写字母比对，所以使用 soft-masked 基因组的比对效果和使用 unmasked 基因组的比对效果是相同的。
  • 当软件不识别大小写时，unmasked等同于soft-masked

2.2 gff/gtf 格式

GFF3

1
2
3
4
5
6
7
8
9
10
11
12
13
14

一个物种的基因组测序完成后，需要对这些数据进行解读，首先要先找到这些序列中转录起始位点、基因、外显子、内含子等组成元件在染色体中的位置信息(即注释)后才能再进行深入的分析。gff/gtf是贮存这些注释信息的两种文件格式。
GFF(general feature format)：这种格式主要是用来注释基因组。 现大部分利用的是第三版，即gff3。
GTF(gene transfer format)：主要是用来对基因进行注释。当前所广泛使用的gtf格式为第二版，即gtf 。
##### GFF3
GFF3允许使用#作为注释符号 ，除去注释外，主体部分共有9列。 GFF3中每一列的含义：seqid source type start end score strand strand attributes
1. seqid ：序列的id。(The name of the sequence where the feature is located.)
2. source：注释的来源，一般指明产生此gff3文件的软件或方法(e.g. Augustus or RepeatMasker)。如果未知，则用点（.）代替。
3. type： 类型，此处不受约束，但为下游分析方便，建议使用gene，repeat_region，exon，CDS，或SO对应编号等。
4. start：起始位置，从1开始计数（区别于bed文件从0开始计数）。
5. end：终止位置。
6. score：得分，注释信息可能性说明，可以是序列相似性比对时的E-values值或者基因预测是的P-values值。”.”表示为空。(indicates the confidence of the source on the annotated feature)
7. strand：“＋”表示正链，“－”表示负链，“.”表示不需要指定正负链，“?” 表示未知. 
8. phase ：步进。仅对编码蛋白质的CDS有效，本列指定下一个密码子开始的位置。可以是0、1或2，表示到达下一个密码子需要跳过碱基个数。
9. attributes：属性。一个包含众多属性的列表，格式为“标签＝值”（tag=value），不同属性之间以分号相隔。

GTF

gtf 文件也是由 9 列组成,其中每一列含义：seqname source，feature，start，end，score，strand，frame，attributes
第 9 列必须以 gene_id 以及 transcript_id 开头

1
2
3
4
5
6
7
8
9
10
11

1. seqname: 序列的名字。通常格式染色体ID或是contig ID。
2. source：注释的来源。通常是预测软件名或是公共数据库。
3. start：起始位置，从1开始计数。
4. end：终止位置。
5. feature ：基因结构.根据所使用软件不同，feature types必须注明。CDS，start_codon，stop_codon是一定要含有的类型。
6. score ：这一列的值表示对该类型存在性和其坐标的可信度，不是必须的，可以用点“.”代替。
7. strand：链的正向与负向，分别用加号+和减号-表示。
8. frame：密码子偏移，可以是0、1或2。
9. attributes：必须要有以下两个值：gene_id value: 表示转录本在基因组上的基因座的唯一的ID。gene_id与value值用空格分开，如果值为空，则表示没有对应的基因。transcript_id value: 预测的转录本的唯一ID。transcript_id与value值用空格分开，空表示没有转录本。

2.3 两种文件的差异比较

GFF与GTF差异

2.4 两种文件格式之间的转换

直接使用 Cufflinks 里面的工具 gffread

1
2
3
4

#gff2gtf
gffread my.gff3 -T -o my.gtf
#gtf2gff
gffread merged.gtf -o- > merged.gff3

参考来源

  1. 基因组各种版本对应关系
  2. NGS 基础 – 参考基因组和基因注释文件
  3. gff/gtf 格式
  4. 从 UCSC 下载基因组的 GTF 文件

Published by

风君子

独自遨游何稽首 揭天掀地慰生平