国际千人基因组计划数据库（二）

在国际千人基因组计划数据库（一）中，我们介绍了如何下载1000G数据库中特定样本的列表。样本列表文件名：igsr_samples.tsv，共216个样本。

在collection中勾选“1000 Genomes phase 3 release”后样本数没有变化，说明样本来自“Phase 3”：

更多详尽的操作说明可查看IGSR网站的FAQ：

https://www.internationalgenome.org/help

关于数据的下载和查看：

下面仍旧以镰刀型贫血症致病基因HBB的rs334致病位点（详情链接：OMIM）为例。对这个位点，还是要先查询一下dbSNP数据库，结果如下

突变频率：

染色体坐标：

https://www.ncbi.nlm.nih.gov/snp/rs334?horizontal_tab=true#variant_details

因此，rs334的等位基因突变频率（T>A）为0.0274（137/5008）；GRCh37版的染色体坐标为chr11:5248232；GRCh38的为chr11:5227002。

1000G数据库的FTP

在IGSR官网搜索“FTP”：

https://www.internationalgenome.org/home

1000G的FTP地址：

http://ftp.1000genomes.ebi.ac.uk/vol1/ftp/

---

FTP与HTTP

ftp：文件传输协议。HTTP：超文本传输协议（交换文本，图形，声音，动画等的规则集）。

二者主要区别在于ftp需要用户名和密码才能访问（匿名登陆除外），http协议一般是任何人都可以访问的。

ftp在浏览器的地址栏中的格式是：ftp://用户名:密码@网站地址/文件目录/

---

FTP中的内容太多，初次使用时很难理清楚。倒不如参考一下其它数据库如何调用其中的数据，比如Ensembl。

参考Ensembl使用的1000G数据源

因此先看一下Ensembl的工具集中“Data Slicer”（对在线数据库中的BAM或VCF文件进行切割）如何调用1000G的数据。

Ensembl GRCh37

https://grch37.ensembl.org/info/docs/tools/index.html（有很多其它生信工具）

严格按照提示输入染色体区域（GRCh37）：

发现“Data Slicer”调用了1000G的如下FTP站点：

ftp://ftp.1000genomes.ebi.ac.uk/vol1/ftp/release/20130502/ALL.chr11.phase3_shapeit2_mvncall_integrated_v5b.20130502.genotypes.vcf.gz

对下载后的文件（上图红色框），使用bash命令查看有多少个样本：

grep CHROM 11.5248232-5248233.ALL.chr11.phase3_shapeit2_mvncall_integrated_v5b.20130502.genotypes.vcf  | tr "\t" "\n" | tail -n +10 | sort -u | wc -l
# 2504

注意：对于二倍体生物，任何一个变异位点，等位基因的个数（number of alleles）（野生+突变）是样本个数的2倍。

故样本个数为2504时，等位基因个数为5008。如下图：

AN：5008（Total number of alleles in called genotypes）

AC：137（突变的等位基因个数：Total number of alternate alleles in called genotypes）

AF：137/5008 ≈ 0.02735（与上述dbSNP数据库收录的1000G数据一致）。Estimated allele frequency in the range (0,1)。

人群的缩写：Africans (AFR), Europeans (EUR), Admixed Americans (AMR), East Asians (EAS) and South Asians (SAS)。

发现非洲人群的等位基因频率（AFR_AF）特别高（~0.1），是所有人群平均值的~3.7倍（0.1/0.02735）。至于原因，不必再赘述。相关链接：蛋白质二级结构、结构域及蛋白修饰预测；生命科学研究的"上帝视角"—组学(Omics)。

Ensembl GRCh38（过程同上）

使用的FTP文件是GRCh38版：

显然，除了坐标系，GRCh38版的数据与上述GRCh37的完全相同，只是网址不同。可能是Ensembl做了坐标系的转换。所在Ensembl的站点如下（发现了一堆重量级VCF文件）：

tabix切割FTP上的VCF文件

总结Ensembl使用的1000G数据库的数据：

GRCh37：

ftp://ftp.1000genomes.ebi.ac.uk/vol1/ftp/release/20130502/ALL.chr11.phase3_shapeit2_mvncall_integrated_v5b.20130502.genotypes.vcf.gz

http://ftp.1000genomes.ebi.ac.uk/vol1/ftp/release/20130502/ALL.chr11.phase3_shapeit2_mvncall_integrated_v5b.20130502.genotypes.vcf.gz

GRCh38：

ftp://ftp.ensembl.org/pub/data_files/homo_sapiens/GRCh38/variation_genotype/ALL.chr11.phase3_shapeit2_mvncall_integrated_v3plus_nounphased.rsID.genotypes.GRCh38_dbSNP.vcf.gz

http://ftp.ensembl.org/pub/data_files/homo_sapiens/GRCh38/variation_genotype/ALL.chr11.phase3_shapeit2_mvncall_integrated_v3plus_nounphased.rsID.genotypes.GRCh38_dbSNP.vcf.gz

rs334的坐标：

GRCh37，chr11:5248232；GRCh38，chr11:5227002

tabix（bash命令）切割上述FTP上的VCF文件

GRCh37：

conda install tabix
tabix -h ftp://ftp.1000genomes.ebi.ac.uk/vol1/ftp/release/20130502/ALL.chr11.phase3_shapeit2_mvncall_integrated_v5b.20130502.genotypes.vcf.gz 11:5248232-5248232

GRCh38：

tabix -h ftp://ftp.ensembl.org/pub/data_files/homo_sapiens/GRCh38/variation_genotype/ALL.chr11.phase3_shapeit2_mvncall_integrated_v3plus_nounphased.rsID.genotypes.GRCh38_dbSNP.vcf.gz  11:5227002-5227002

wget（bash命令）下载上述http上的VCF文件

# GRCh37:
wget -c http://ftp.1000genomes.ebi.ac.uk/vol1/ftp/release/20130502/ALL.chr11.phase3_shapeit2_mvncall_integrated_v5b.20130502.genotypes.vcf.gz
# GRCh38:
wget -c http://ftp.ensembl.org/pub/data_files/homo_sapiens/GRCh38/variation_genotype/ALL.chr11.phase3_shapeit2_mvncall_integrated_v3plus_nounphased.rsID.genotypes.GRCh38_dbSNP.vcf.gz

统计样本数量：

zcat ALL.chr11.phase3_shapeit2_mvncall_integrated_v5b.20130502.genotypes.vcf.gz | head -n 3000 |   grep CHROM  | tr "\t" "\n" | tail -n +10 | sort -u | wc -l
# 2504
zcat ALL.chr11.phase3_shapeit2_mvncall_integrated_v3plus_nounphased.rsID.genotypes.GRCh38_dbSNP.vcf.gz | head -n 3000 |   grep CHROM  | tr "\t" "\n" | tail -n +10 | sort -u |wc -l
# 2504

至此，我们获得了1000G数据库phase3的11号染色体上、2504个样本的VCF文件，其中包含人群频率和各个样本的具体基因型等信息。

若想从中单独提取人群频率数据，或只提取个别样本的基因型，可通过一些工具或bash命令实现。

下载1000G的fastq文件

先查看上一篇推送中，下载的样本列表中的某个样本编号：

head -n 2 igsr_samples.tsv

下载当前1000G数据库的FTP站点上所有文件的下载链接“树”（tree）。其实就是个文本文件，告诉你哪些文件放在了哪里。

wget -c http://ftp.1000genomes.ebi.ac.uk/vol1/ftp/current.tree
grep -w HG00404  current.tree | cut -f 1,2,3 | grep fastq | sort 

grep -w HG00404  current.tree | cut -f 1,2,3 | grep fastq | sort | wc -l
# 63
# 共21组

这21个文件不是拆分了不同的染色体，而是一个样本分别在多个Run中的测序结果。详见FAQ中关于fastq文件：

很多样本都有多个fastq文件，是因为许多人都是用测序仪进行了不止一次的测序。

例如命名为ERR001268_1.filt.fastq.gz、ERR001268_2.filt.fastq.gz和ERR001268.filt.fastq.gz的每一组文件代表了一个测序Run的所有序列。

ERR001268为“a sequencing run accessions ”。“ _1”和“_2”代表双端测序文件（paired-end files），其它为单端测序（项目早期生成，或来自双端测序未匹配成功的reads）。

使用wget命令逐个下载一个Run的测序文件：

wget -c ftp://ftp.1000genomes.ebi.ac.uk/vol1/ftp/phase3/data/HG00404/sequence_read/SRR100420_1.filt.fastq.gz
wget -c ftp://ftp.1000genomes.ebi.ac.uk/vol1/ftp/phase3/data/HG00404/sequence_read/SRR100420_2.filt.fastq.gz
wget -c ftp://ftp.1000genomes.ebi.ac.uk/vol1/ftp/phase3/data/HG00404/sequence_read/SRR100420.filt.fastq.gz

关于测序数据如何批量下载，如何合并，以及1000G/IGSR推荐的下载方式，可查看其FAQ。把握好正确的数据来源是关键。

撰写：宋红卫

校对：宋红卫

更多人类遗传学知识、文献和分析技术

请关注和星标聊生信