vg
vg是一个基于图论的泛基因组分析软件,可用于存储、处理和分析大规模基因组序列数据。它支持多种数据格式和算法,并提供了丰富的功能接口,方便用户进行高效的基因组数据分析和可视化。
vg主要有以下几个方面的功能:
1.数据导入和转换:可以将多种基因组和变异数据格式导入到vg中,并通过内置的工具和API进行转换和处理。vg支持FASTA,FASTQ,GFA,VCF,SAM/BAM等常见数据格式。
2.图构建和索引:vg可以快速构建泛基因组的组装图,并对其进行索引以支持快速查询和搜索。
3.基因组比对和变异检测:vg可以根据构建的图结构对多个基因组进行比对,并检测其间的变异信息。支持多种比对算法和变异检测方法,如BWA-MEM,minimap2,Giraffe等。
4.可视化和交互式分析:vg提供了多种快速可视化和交互式分析工具,包括Tview、dotplot、GBrowse等,可以帮助用户深入理解基因组数据和变异信息。
代码库:https://github.com/vgteam/vg
官方Wiki:https://github.com/vgteam/vg/wiki
泛基因组分析中用到的方法:Practical Graphical Pangenomics
vg更新迭代比较快,经常有各种问题,运行遇到问题多在github issue中搜一下,https://github.com/vgteam/vg/issues
安装¶
直接下载
$ wget https://github.com/vgteam/vg/releases/download/v1.53.0/vg          
$ chmod +x vg
$ module load vg/1.53.0
$ vg 
vg: variation graph tool, version v1.46.0 "Altamura"
usage: vg <command> [options]
main mapping and calling pipeline:
  -- autoindex     mapping tool-oriented index construction from interchange formats
  -- construct     graph construction
  -- rna           construct splicing graphs and pantranscriptomes
  -- index         index graphs or alignments for random access or mapping
  -- map           MEM-based read alignment
  -- giraffe       fast haplotype-aware short read alignment
  -- mpmap         splice-aware multipath alignment of short reads
  -- augment       augment a graph from an alignment
  -- pack          convert alignments to a compact coverage index
  -- call          call or genotype VCF variants
  -- help          show all subcommands
For more commands, type `vg help`.
For technical support, please visit: https://www.biostars.org/tag/vg/
构建泛基因组¶
从vcf构建¶
使用变异数据(vcf)和参考基因组构建泛基因组,vcf需要使用tabix构建索引(tabix -p vcf x.vcf.gz)。
$ vg construct -r small/x.fa -v small/x.vcf.gz >x.vg
从基因组构建¶
也使用组装好的基因组来构建泛基因组,只需提供不同个体的基因组fasta文件列表即可,使用的软件为Minigraph-Cactus,具体可参考 The Minigraph-Cactus Pangenome Pipeline
cactus-pangenome ./jobstore genome.txt \
  --outDir pangenome  --workDir `pwd`/tmp \
  --outName pangenome --chop --haplo --permissiveContigFilter --gbz full clip  filter \
  --gfa full clip  filter \
  --vcf full clip  filter \
  --giraffe  full clip  filter \
  --chrom-vg  full clip  filter \
  --viz full clip   \
  --odgi full clip  filter \
  --chrom-og full clip  filter  \
  --binariesMode local  \
  --filter 2 \
  --mgCores 50 --mapCores  8   --consCores 32  --indexCores 32 \
  --refContigs chr1 chr2 chr3 chr4 chr5 --reference Ath  --logFile cactus-pangenome.log   #--restart
变异检测¶
详细文档参考 SV Genotyping and variant calling
创建索引¶
# 从 fasta 文件创建
vg autoindex -p pangenome -w giraffe -t 10 -r ../pangenome.fasta
pangenome.distpangenome.giraffe.gbzpangenome.min比对¶
https://github.com/vgteam/vg/wiki/Mapping-short-reads-with-Giraffe
https://github.com/vgteam/vg/wiki/Giraffe-best-practices
vg 提供了3种比对工具
- vg giraffe:用于快速将二代度短序列reads比对到带有单倍型信息的图基因组
- vg map:通用的比对工具
- vg mpmap:用于转录组比对
一般使用二代数据检测变异使用 vg giraffe ,-t设置线程数。详细用法见 Mapping short reads with Giraffe
$ vg giraffe -m pangenome.min \
    -d pangenome.dist \
    -Z pangenome.gbz \
    -f $sample.r1.fq.gz \
    -f $sample.r2.fq.gz  \
    -N $sample  -t 32  >$sample.gam
比对结果统计¶
$ vg stats -a $sample.gam > $sample.stats
变异检测¶
# -Q 5: ignore mapping and base qualitiy < 5
$ vg pack -x  pangenome.gbz -g $sample.sorted.gam -Q 5 -o $sample.pack -t 32
# Compute the snarls (using fewer threads with `-t` can reduce memory at the cost of increased runtime)
$ vg snarls pangenome.gbz > pangenome.snarls
# Genotype the graph 
$ vg call pangenome.gbz -r pangenome.snarls -k  $sample.pack -a -A -t 32 -s $sample > $sample.vcf
运行问题¶
- vg giraffe运行缓慢
集群上运行 vg giraffe 时即使设置了使用多个线程,cputime 和 runtime 相等,进程为 D 状态,运行速度非常缓慢。如下所示:
top - 10:43:24 up 11 days, 11:25,  1 user,  load average: 19.42, 12.95, 14.74
Tasks: 292 total,   2 running, 290 sleeping,   0 stopped,   0 zombie
%Cpu(s):  3.4 us,  2.9 sy,  0.0 ni, 93.6 id,  0.0 wa,  0.0 hi,  0.0 si,  0.0 st
KiB Mem : 13183338+total, 93854384 free, 13359480 used, 24619524 buff/cache
KiB Swap:  4194300 total,  4194300 free,        0 used. 11196654+avail Mem
  PID USER      PR  NI    VIRT    RES    SHR S  %CPU %MEM     TIME+ COMMAND
12957 user      20   0 8862844   4.9g   4204 D 101.7  3.9  10:38.59 vg
 1187 root      39  19       0      0      0 R   2.6  0.0 396:49.96 kipmi0
13114 root      20   0  166356   2560   1672 R   0.7  0.0   0:00.47 top
    1 root      20   0  192444   5444   2612 S   0.0  0.0   5:17.11 systemd
    2 root      20   0       0      0      0 S   0.0  0.0   0:00.21 kthreadd
    3 root      20   0       0      0      0 S   0.0  0.0   0:06.06 ksoftirqd/0
    5 root       0 -20       0      0      0 S   0.0  0.0   0:00.00 kworker/0:0H
    8 root      rt   0       0      0      0 S   0.0  0.0   0:03.98 migration/0
    9 root      20   0       0      0      0 S   0.0  0.0   0:00.00 rcu_bh
   10 root      20   0       0      0      0 S   0.0  0.0   9:49.19 rcu_sched
   11 root       0 -20       0      0      0 S   0.0  0.0   0:00.00 lru-add-drain
   12 root      rt   0       0      0      0 S   0.0  0.0   0:02.97 watchdog/0
   13 root      rt   0       0      0      0 S   0.0  0.0   0:03.63 watchdog/1
   14 root      rt   0       0      0      0 S   0.0  0.0   0:04.42 migration/1
dist 文件。解决办法为,每个样本都拷贝一份临时 dist 文件然后再运行 vg giraffe,运行完成后删除临时 dist 文件。 #BSUB -J vg
#BSUB -n 10
#BSUB -o %J.out
#BSUB -e %J.err
#BSUB -q normal
module load vg/1.53.0
# 生成随机字符串
tmp_dist=`mktemp -u dist_XXXXX`
# 拷贝一份临时 dist 文件
cp ref.dist ${tmp_dist}
vg giraffe -p -t 10 -Z ref.giraffe.gbz -d ${tmp_dist} -m ref.min -f sample_clean_1.fq.gz -f sample_clean_2.fq.gz > sample.gam
#清理临时 dist 文件
rm -r ${tmp_dist}
泛基因组¶
泛基因组文件¶
泛基因组文件(Pan-genome file)是指包含多个相关物种或同一物种的不同个体的基因组序列的综合文件。泛基因组文件汇集了多个基因组的信息,能够更全面地反映一个物种的遗传多样性和基因组变异。
泛基因组文件与基因组文件的区别在于其所包含的内容和结构:
1.内容:基因组文件通常只包含单一个体的基因组序列,而泛基因组文件则包含多个个体的基因组序列,用于研究它们之间的差异和共性。
2.结构:基因组文件一般以染色体为单位存储,每个染色体由连续的碱基序列组成,并伴随其他相关信息,如注释等,文件为线性单一结构。而泛基因组文件则以基因族(gene family)为单位组织数据,每个基因族中包含多个拷贝基因(gene copy),可以是来自不同个体的结构变异单倍型结构,整体的文件格式为图形结构存储。
泛基因组文件的文件格式¶
https://github.com/vgteam/vg/wiki/File-Types
- GFA文件(可读)
GFA文件(Graphical Fragment Assembly file):GFA文件是一种用于存储和表示基因组组装图的文件格式。它包含了由测序数据生成的所有片段(fragments)及其相互之间的关系(边)。GFA文件的数据结构主要包括以下几个部分:
- 线或线段是图的节点。它们由一个唯一的非空字符串名称和一个非空字符串标签组成。 
- 线或链接是图的边。它们连接两个片段访问。我们不允许在连接的部分之间有任何重叠。 
- 线或路径在图中称为路径。它们由一个唯一的非空字符串名称和一个非空的段访问序列组成。 
- -line或walk是具有适合表示组装基因组的元数据的路径。它们由样本名称(非空字符串)、单倍型标识符(整数)、序列名称(非空字符串)、命名序列中位置的可选间隔以及段访问的非空序列组成 
- GBWT文件(不可直接读)
GBWT(Generalized Burrows-Wheeler Transform file) 是一个索引文件,可以有效地存储类似的线程(轻量级路径)。在大多数应用中,它被用于存储(真实的或人工的)单体型。gbwt文件的数据结构主要包括以下几个部分:
- 分区(Partition):将基因组序列和变异信息切分成多个区域,每个区域包含一个或多个线性链(linear chain)。 
- 线性链(Linear Chain):由连续的变异位置和对应的基因组字符组成,用于表示一条基因组序列的一部分。 
- 元数据(Metadata):存储与gbwt文件相关的附加信息,如版本号、索引等。 
- xg文件(不可直接读) 
XG(XG Light Graph)文件是一种用于存储DNA序列的索引信息的文件格式。它是基于图形数据结构的索引,可以支持快速的基因组序列查询和路径遍历。xg文件的数据结构主要包括以下几个部分:
- 节点(Node):表示基因组序列的节点,包含唯一标识符和基因组字符等信息。 
- 边(Edge):表示节点之间的连接关系,例如前驱节点和后继节点等。 
- 序列(Sequence):存储与节点相关联的基因组序列信息。 
- 属性(Attribute):存储附加于节点或边上的其他信息,如注释等。 
参考¶
泛基因组call SNP比普通基因组call SNP更有准确
本文阅读量 次本站总访问量 次