Post-GWAS-LD图与基因，关联结果一起可视化

image.png

IntAssoplot包可以将GWAS结果和LD，基因结构图整合在一起画出。
介绍来自：https://github.com/whweve/IntAssoPlot

1下载IntAssoplot包

#### 前期准备包
install.packages(c("devtools","remotes"))

#install depended packages, including ggplot2, SNPRelate, ggrepel, gdsfmt and reshape2 #ggplot2, ggrepel, and reshape2 are installed from CRAN

install.packages(c("ggplot2","ggrepel","reshape2"))

#SNPRelate and gdsfmt are installed from Bioconductor

if (!requireNamespace("BiocManager", quietly = TRUE)) install.packages("BiocManager")

BiocManager::install(c("SNPRelate","gdsfmt"))
##### 下载IntAssoPlot from Github
library(remotes) # version 2.1.0

#download, build, and install IntAssoPlot without creating vignette

install_github("whweve/IntAssoPlot")
或者
#download, build, and install IntAssoPlot with creating vignette

install_github("whweve/IntAssoPlot",build=TRUE,build_vignettes = TRUE)

2 数据准备

2.1 GWAS的结果文件

包含关联分析文件的数据集。目前，我们针对相应标记的物理位置绘制了 -log10 转换后的 p 值，该值可以从基因组扫描中得出。因此，需要一个包含标记名称及其 p 值的数据框。需要的变量如下：Marker（分子标记名称）、Locus（标记的染色体）、Site（标记的位置）、p（标记的p值）这里是关联数据格式的快速演示：

head(association)
#>     Marker  Locus     Site         p
#> 1 SNP37049     9 90042244 0.2073684
#> 2 SNP37050     9 90042408 0.5028760
#> 3 SNP37052     9 90210602 0.1561271
#> 4 SNP37053     9 90255652 0.2814653
#> 5 SNP37056     9 90290939 0.6392511
#> 6 SNP44597     9 90367887 0.6582597
#the attribute of each column could be viewed as:
str(association$Marker)
#>  chr [1:2316] "SNP37049" "SNP37050" "SNP37052" "SNP37053" "SNP37056" ...
str(association$Locus)
#>  int [1:2316] 9 9 9 9 9 9 9 9 9 9 ...
str(association$Site)
#>  int [1:2316] 90042244 90042408 90210602 90255652 90290939 90367887 90377686 90385837 90560965 90759605 ...
str(association$p)
#>  num [1:2316] 0.207 0.503 0.156 0.281 0.639 ..

2.2 gene结构注释的文件格式

包含注释文件的数据集，通常是 gtf 文件，没有列名。当对一个物种的基因组进行测序时，可以在基因组注释网站上找到 gtf 文件。对于已发布的基因组，请参阅 www.ensembl.org。由于大多数 gtf 文件缺少列名，因此可以使用 read.table(?where is your file?,header=FALSE) 读取注释文件。这是 gtf 数据格式的快速演示：

head(gtf)
#>   V1             V2          V3       V4       V5 V6 V7 V8
#> 1  9 protein_coding        exon 90030791 90030847  .  -  .
#> 2  9 protein_coding        exon 90030635 90030701  .  -  .
#> 3  9 protein_coding         CDS 90030635 90030681  .  -  0
#> 4  9 protein_coding start_codon 90030679 90030681  .  -  0
#> 5  9 protein_coding        exon 90029368 90029471  .  -  .
#> 6  9 protein_coding         CDS 90029368 90029471  .  -  1
#>                                                                                                             V9
#> 1                gene_id "GRMZM2G416644"; transcript_id "GRMZM2G416644_T01"; exon_number "1"; seqedit "false";
#> 2                gene_id "GRMZM2G416644"; transcript_id "GRMZM2G416644_T01"; exon_number "2"; seqedit "false";
#> 3 gene_id "GRMZM2G416644"; transcript_id "GRMZM2G416644_T01"; exon_number "2"; protein_id "GRMZM2G416644_P01";
#> 4                                 gene_id "GRMZM2G416644"; transcript_id "GRMZM2G416644_T01"; exon_number "2";
#> 5                gene_id "GRMZM2G416644"; transcript_id "GRMZM2G416644_T01"; exon_number "3"; seqedit "false";
#> 6 gene_id "GRMZM2G416644"; transcript_id "GRMZM2G416644_T01"; exon_number "3"; protein_id "GRMZM2G416644_P01";
#the attribute of each column could be viewed as:
str(gtf$V1)
#>  chr [1:4665] "9" "9" "9" "9" "9" "9" "9" "9" "9" "9" "9" "9" "9" "9" "9" ...
str(gtf$V2)
#>  chr [1:4665] "protein_coding" "protein_coding" "protein_coding" ...
str(gtf$V3)
#>  chr [1:4665] "exon" "exon" "CDS" "start_codon" "exon" "CDS" "exon" "CDS" ...
str(gtf$V4)
#>  int [1:4665] 90030791 90030635 90030635 90030679 90029368 90029368 90028806 90029227 90029224 90030697 ...
str(gtf$V5)
#>  int [1:4665] 90030847 90030701 90030681 90030681 90029471 90029471 90029282 90029282 90029226 90031013 ...
str(gtf$V6)
#>  chr [1:4665] "." "." "." "." "." "." "." "." "." "." "." "." "." "." "." ...
str(gtf$V7)
#>  chr [1:4665] "-" "-" "-" "-" "-" "-" "-" "-" "-" "-" "-" "-" "-" "-" "-" ...
str(gtf$V8)
#>  chr [1:4665] "." "." "0" "0" "." "1" "." "2" "0" "." "0" "0" "." "0" "0" ...
str(gtf$V9)
#>  chr [1:4665] "gene_id \"GRMZM2G416644\"; transcript_id \"GRMZM2G416644_T01\"; exon_number \"1\"; seqedit \"false\";" ...

2.3 基因型数据的格式

包含基因型文件（通常是 hapmap 文件）和列名的数据集。 Hapmap 是一种常用的存储基因信息的格式。可以对重要材料进行高通量测序，将序列与参考基因组对齐，提取 SNP/InDels。另一方面，人们可以对一个特定的基因进行重新测序。为此，设计重叠扩增引物、扩增基因组或转录组片段、多重比对序列、提取 SNP/InDels。这是关联数据格式的快速演示：

#only 20 column of the genotype markers are shown.
head(zmvpp1_hapmap[,1:20])
#>          rs       allele   chrom  pos strand assembly center protLSID assayLSID
#> 1 INDEL-665    -/+     9 94183611      +       NA     NA       NA        NA
#> 2   snp-663    A/T     9 94183597      +       NA     NA       NA        NA
#> 3   snp-649    A/G     9 94183594      +       NA     NA       NA        NA
#> 4   snp-646    C/T     9 94183585      +       NA     NA       NA        NA
#> 5   snp-637    G/T     9 94183584      +       NA     NA       NA        NA
#> 6   snp-636    A/G     9 94183575      +       NA     NA       NA        NA
#>   panel QCcode B73 CIMBL55 CIMBL32 CML298 CML170 CIMBL12 CIMBL95 CIMBL1 CIMBL56
#> 1    NA     NA  --      ++      ++     ++     ++      ++      ++     --      --
#> 2    NA     NA  AA      TT      TT     TT     TT      TT      TT     NN      NN
#> 3    NA     NA  AA      GG      GG     GG     GG      GG      GG     NN      NN
#> 4    NA     NA  CC      TT      TT     TT     TT      TT      TT     NN      NN
#> 5    NA     NA  GG      TT      TT     TT     TT      TT      TT     NN      NN
#> 6    NA     NA  AA      AA      AA     AA     AA      AA      AA     NN      NN

3 绘制与注释和 LD 矩阵的关联

IntAssoPlot 是将关联分析、基因结构和 LD 矩阵的结果自动整合到一个视图中。事实上，这项任务很困难，因为散点图 - log10P 与 SNP 的物理位置相对，而 LD 矩阵没有物理位置。该包的工作其实就是为了解决上述问题。在这里，我们使用先前发布的数据（Wang 等人，2016 年）展示了一个示例来展示 IntAssoPlot 的用法。

3.1 具有一组基因型标记的区域综合图

在跨越 400 kbp 区域的区域绘制关联结果，并使用与关联映射相同的 SNP 标记绘制 LD 矩阵。

IntRegionalPlot(chr=9,left=94178074-200000,right=94178074+200000,gtf=gtf,association=association,hapmap=hapmap_am368,hapmap_ld=hapmap_am368,threshold=5,leadsnp_size=2)
#> Linkage Disequilibrium (LD) estimation on genotypes:
#>     # of samples: 368
#>     # of SNPs: 270
#>     using 1 thread
#>     method: correlation
#> LD matrix:    the sum of all selected genotypes (0,1,2) = 143276
#> Linkage Disequilibrium (LD) estimation on genotypes:
#>     # of samples: 368
#>     # of SNPs: 270
#>     using 1 thread
#>     method: correlation
#> LD matrix:    the sum of all selected genotypes (0,1,2) = 143276

结果图：

image.png

3.2 用从浅灰色到深灰色的颜色绘制 LD 值

IntRegionalPlot(chr=9,left=94178074-200000,right=94178074+200000,gtf=gtf,association=association,hapmap=hapmap_am368,hapmap_ld=hapmap_am368,threshold=5,leadsnp_size=2,colour02 = "gray1",colour04 = "gray21",colour06 = "gray41",colour08 = "gray61",colour10 = "gray81",)
#> Linkage Disequilibrium (LD) estimation on genotypes:
#>     # of samples: 368
#>     # of SNPs: 270
#>     using 1 thread
#>     method: correlation
#> LD matrix:    the sum of all selected genotypes (0,1,2) = 143276
#> Linkage Disequilibrium (LD) estimation on genotypes:
#>     # of samples: 368
#>     # of SNPs: 270
#>     using 1 thread
#>     method: correlation
#> LD matrix:    the sum of all selected genotypes (0,1,2) = 143276

image.png

3.3 用从白色到红色的颜色绘制 LD 值

#get five colors ranging from white to red
pal <- colorRampPalette(c("white", "red"))
IntRegionalPlot(chr=9,left=94178074-200000,right=94178074+200000,gtf=gtf,association=association,hapmap=hapmap_am368,hapmap_ld=hapmap_am368,threshold=5,leadsnp_size=2,colour02 = pal(5)[1],colour04 = pal(5)[2],colour06 = pal(5)[3],colour08 = pal(5)[4],colour10 = pal(5)[5])
#> Linkage Disequilibrium (LD) estimation on genotypes:
#>     # of samples: 368
#>     # of SNPs: 270
#>     using 1 thread
#>     method: correlation
#> LD matrix:    the sum of all selected genotypes (0,1,2) = 143276
#> Linkage Disequilibrium (LD) estimation on genotypes:
#>     # of samples: 368
#>     # of SNPs: 270
#>     using 1 thread
#>     method: correlation
#> LD matrix:    the sum of all selected genotypes (0,1,2) = 143276

image.png

3.4 用从白色到红色的颜色绘制 LD 值并标记基因名称

#get five colors ranging from white to red
pal <- colorRampPalette(c("white", "red"))
IntRegionalPlot(chr=9,left=94178074-200000,right=94178074+200000,gtf=gtf,association=association,hapmap=hapmap_am368,hapmap_ld=hapmap_am368,threshold=5,leadsnp_size=2,colour02 = pal(5)[1],colour04 = pal(5)[2],colour06 = pal(5)[3],colour08 = pal(5)[4],colour10 = pal(5)[5],label_gene_name = TRUE)
#> Linkage Disequilibrium (LD) estimation on genotypes:
#>     # of samples: 368
#>     # of SNPs: 270
#>     using 1 thread
#>     method: correlation
#> LD matrix:    the sum of all selected genotypes (0,1,2) = 143276
#> Linkage Disequilibrium (LD) estimation on genotypes:
#>     # of samples: 368
#>     # of SNPs: 270
#>     using 1 thread
#>     method: correlation
#> LD matrix:    the sum of all selected genotypes (0,1,2) = 143276

image.png

3.5 具有两组基因型标记的区域综合图

在跨越 200 kbp 区域的区域绘制关联结果，并使用与关联映射不同的 SNP 标记器绘制 LD 矩阵。此功能允许研究人员在更广泛的标记范围内研究 LD 结构。

IntRegionalPlot(chr=9,left=94178074-100000,right=94178074+100000,gtf=gtf,association=association,hapmap=hapmap_am368,hapmap_ld=hapmap2,threshold=5,leadsnp_size=2)
#> Linkage Disequilibrium (LD) estimation on genotypes:
#>     # of samples: 368
#>     # of SNPs: 108
#>     using 1 thread
#>     method: correlation
#> LD matrix:    the sum of all selected genotypes (0,1,2) = 56746
#> Linkage Disequilibrium (LD) estimation on genotypes:
#>     # of samples: 104
#>     # of SNPs: 1,081
#>     using 1 thread
#>     method: correlation
#> LD matrix:    the sum of all selected genotypes (0,1,2) = 114672

image.png

3.6 具有一组基因型标记的相对较小的区域综合图

IntRegionalPlot(chr=9,left=94178074-2000,right=94178074+5000,gtf=gtf,association=association,hapmap=hapmap_am368,hapmap_ld=hapmap_am368,threshold=5,leadsnp_size=2)
#> Linkage Disequilibrium (LD) estimation on genotypes:
#>     # of samples: 368
#>     # of SNPs: 41
#>     using 1 thread
#>     method: correlation
#> LD matrix:    the sum of all selected genotypes (0,1,2) = 21412
#> Linkage Disequilibrium (LD) estimation on genotypes:
#>     # of samples: 368
#>     # of SNPs: 41
#>     using 1 thread
#>     method: correlation
#> LD matrix:    the sum of all selected genotypes (0,1,2) = 21412

image.png

3.7 单基因水平图

绘制给定基因的关联结果，并使用与关联映射相同的 SNP 标记物绘制 LD 矩阵。还指定的标记通过各种形状和颜色突出显示。

3.7.1 基本版

IntGenicPlot('GRMZM2G170927_T01',gtf,association=zmvpp1_association,hapmap=zmvpp1_hapmap,hapmap_ld = zmvpp1_hapmap,threshold=8,leadsnpLD = FALSE)
#> Linkage Disequilibrium (LD) estimation on genotypes:
#>     # of samples: 141
#>     # of SNPs: 53
#>     using 1 thread
#>     method: correlation
#> LD matrix:    the sum of all selected genotypes (0,1,2) = 9992
#> Linkage Disequilibrium (LD) estimation on genotypes:
#>     # of samples: 141
#>     # of SNPs: 53
#>     using 1 thread
#>     method: correlation
#> LD matrix:    the sum of all selected genotypes (0,1,2) = 9992

image.png

3.7.2从基因的上游/下游扩展区域

IntGenicPlot('GRMZM2G170927_T01',gtf,association=zmvpp1_association,hapmap=zmvpp1_hapmap,hapmap_ld = zmvpp1_hapmap,threshold=8,up=500,down=600,leadsnpLD = FALSE)
#> Linkage Disequilibrium (LD) estimation on genotypes:
#>     # of samples: 141
#>     # of SNPs: 124
#>     using 1 thread
#>     method: correlation
#> LD matrix:    the sum of all selected genotypes (0,1,2) = 23488
#> Linkage Disequilibrium (LD) estimation on genotypes:
#>     # of samples: 141
#>     # of SNPs: 124
#>     using 1 thread
#>     method: correlation
#> LD matrix:    the sum of all selected genotypes (0,1,2) = 23488

image.png

3.7.3 突出显示选定的标记，在数据框中指定颜色和形状：标记 2 突出显示

IntGenicPlot('GRMZM2G170927_T01',gtf,association=zmvpp1_association,hapmap=zmvpp1_hapmap,hapmap_ld = zmvpp1_hapmap,threshold=8,up=500,down=600,leadsnpLD = FALSE,marker2highlight=marker2highlight)
#> Linkage Disequilibrium (LD) estimation on genotypes:
#>     # of samples: 141
#>     # of SNPs: 124
#>     using 1 thread
#>     method: correlation
#> LD matrix:    the sum of all selected genotypes (0,1,2) = 23488
#> Linkage Disequilibrium (LD) estimation on genotypes:
#>     # of samples: 141
#>     # of SNPs: 124
#>     using 1 thread
#>     method: correlation
#> LD matrix:    the sum of all selected genotypes (0,1,2) = 23488

image.png

3.7.4 添加连接线

IntGenicPlot('GRMZM2G170927_T01',gtf,association=zmvpp1_association,hapmap=zmvpp1_hapmap,hapmap_ld = zmvpp1_hapmap,threshold=8,up=500,down=600,leadsnpLD = FALSE,marker2highlight=marker2highlight,link2gene=marker2link,link2LD=marker2link)
#> Linkage Disequilibrium (LD) estimation on genotypes:
#>     # of samples: 141
#>     # of SNPs: 124
#>     using 1 thread
#>     method: correlation
#> LD matrix:    the sum of all selected genotypes (0,1,2) = 23488
#> Linkage Disequilibrium (LD) estimation on genotypes:
#>     # of samples: 141
#>     # of SNPs: 124
#>     using 1 thread
#>     method: correlation
#> LD matrix:    the sum of all selected genotypes (0,1,2) = 23488

image.png

3.7.5 为高亮标记添加名称

IntGenicPlot('GRMZM2G170927_T01',gtf,association=zmvpp1_association,hapmap=zmvpp1_hapmap,hapmap_ld = zmvpp1_hapmap,threshold=8,up=500,down=600,leadsnpLD = FALSE,marker2highlight=marker2highlight,link2gene=marker2link,link2LD=marker2link,marker2label=marker2link,marker2label_angle=60,marker2label_size=2)
#> Linkage Disequilibrium (LD) estimation on genotypes:
#>     # of samples: 141
#>     # of SNPs: 124
#>     using 1 thread
#>     method: correlation
#> LD matrix:    the sum of all selected genotypes (0,1,2) = 23488
#> Linkage Disequilibrium (LD) estimation on genotypes:
#>     # of samples: 141
#>     # of SNPs: 124
#>     using 1 thread
#>     method: correlation
#> LD matrix:    the sum of all selected genotypes (0,1,2) = 23488

image.png

最后编辑于：2021.10.30 20:51:16

人面猴
序言：七十年代末，一起剥皮案震惊了整个滨河市，随后出现的几起案子，更是在滨河造成了极大的恐慌，老刑警刘岩，带你破解...
沈念sama阅读 206,378评论 6赞 481
死咒
序言：滨河连续发生了三起死亡事件，死亡现场离奇诡异，居然都是意外死亡，警方通过查阅死者的电脑和手机，发现死者居然都...
沈念sama阅读 88,356评论 2赞 382
救了他两次的神仙让他今天三更去死
文/潘晓璐我一进店门，熙熙楼的掌柜王于贵愁眉苦脸地迎上来，“玉大人，你说我怎么就摊上这事。” “怎么了？”我有些...
开封第一讲书人阅读 152,702评论 0赞 342
道士缉凶录：失踪的卖姜人
文/不坏的土叔我叫张陵，是天一观的道长。经常有香客问我，道长，这世上最难降的妖魔是什么？我笑而不...
开封第一讲书人阅读 55,259评论 1赞 279
港岛之恋（遗憾婚礼）
正文为了忘掉前任，我火速办了婚礼，结果婚礼上，老公的妹妹穿的比我还像新娘。我一直安慰自己，他们只是感情好，可当我...
茶点故事阅读 64,263评论 5赞 371
恶毒庶女顶嫁案：这布局不是一般人想出来的
文/花漫我一把揭开白布。她就那样静静地躺着，像睡着了一般。火红的嫁衣衬着肌肤如雪。梳的纹丝不乱的头发上，一...
开封第一讲书人阅读 49,036评论 1赞 285
城市分裂传说
那天，我揣着相机与录音，去河边找鬼。笑死，一个胖子当着我的面吹牛，可吹牛的内容都是我干的。我是一名探鬼主播，决...
沈念sama阅读 38,349评论 3赞 400
双鸳鸯连环套：你想象不到人心有多黑
文/苍兰香墨我猛地睁开眼，长吁一口气：“原来是场噩梦啊……” “哼！你这毒妇竟也来了？” 一声冷哼从身侧响起，我...
开封第一讲书人阅读 36,979评论 0赞 259
万荣杀人案实录
序言：老挝万荣一对情侣失踪，失踪者是张志新（化名）和其女友刘颖，没想到半个月后，有当地人在树林里发现了一具尸体，经...
沈念sama阅读 43,469评论 1赞 300
护林员之死
正文独居荒郊野岭守林人离奇死亡，尸身上长有42处带血的脓包…… 初始之章·张勋以下内容为张勋视角年9月15日...
茶点故事阅读 35,938评论 2赞 323
白月光启示录
正文我和宋清朗相恋三年，在试婚纱的时候发现自己被绿了。大学时的朋友给我发了我未婚夫和他白月光在一起吃饭的照片。...
茶点故事阅读 38,059评论 1赞 333
活死人
序言：一个原本活蹦乱跳的男人离奇死亡，死状恐怖，灵堂内的尸体忽然破棺而出，到底是诈尸还是另有隐情，我是刑警宁泽，带...
沈念sama阅读 33,703评论 4赞 323
日本核电站爆炸内幕
正文年R本政府宣布，位于F岛的核电站，受9级特大地震影响，放射性物质发生泄漏。R本人自食恶果不足惜，却给世界环境...
茶点故事阅读 39,257评论 3赞 307
男人毒药：我在死后第九天来索命
文/蒙蒙一、第九天我趴在偏房一处隐蔽的房顶上张望。院中可真热闹，春花似锦、人声如沸。这庄子的主人今日做“春日...
开封第一讲书人阅读 30,262评论 0赞 19
一桩弑父案，背后竟有这般阴谋
文/苍兰香墨我抬头看了看天上的太阳。三九已至，却和暖如春，着一层夹袄步出监牢的瞬间，已是汗流浃背。一阵脚步声响...
开封第一讲书人阅读 31,485评论 1赞 262
情欲美人皮
我被黑心中介骗来泰国打工，没想到刚下飞机就差点儿被人妖公主榨干…… 1. 我叫王不留，地道东北人。一个月前我还...
沈念sama阅读 45,501评论 2赞 354
代替公主和亲
正文我出身青楼，却偏偏与公主长得像，于是被迫代替她去往敌国和亲。传闻我的和亲对象是个残疾皇子，可洞房花烛夜当晚...
茶点故事阅读 42,792评论 2赞 345

Post-GWAS-LD图与基因，关联结果一起可视化

1下载IntAssoplot包

2 数据准备

2.1 GWAS的结果文件

2.2 gene结构注释的文件格式

2.3 基因型数据的格式

3 绘制与注释和 LD 矩阵的关联

3.1 具有一组基因型标记的区域综合图

3.2 用从浅灰色到深灰色的颜色绘制 LD 值

3.3 用从白色到红色的颜色绘制 LD 值

3.4 用从白色到红色的颜色绘制 LD 值并标记基因名称

3.5 具有两组基因型标记的区域综合图

3.6 具有一组基因型标记的相对较小的区域综合图

3.7 单基因水平图

3.7.1 基本版

3.7.2从基因的上游/下游扩展区域

3.7.3 突出显示选定的标记，在数据框中指定颜色和形状：标记 2 突出显示

3.7.4 添加连接线

3.7.5 为高亮标记添加名称

推荐阅读更多精彩内容