第一种方案绘制的Heatmap需要借助于Corrplot包，我们求出dat特征的相关系数矩阵，进一步利用corrplot函数来画图，在该包的官方文档中，对于corrplot函数的参数描述可以说是非常多，这里我给出几种常用的参数：

method 表示热力图中每一块所展示的形状，可选值有： "circle", "square", "ellipse", "number", "shade", "color", "pie"；

type 表示相关系数矩阵展示的方式，比如只展示上三角或下三角或者全部展示，可选值有：“full”，“upper"，"lower”；

tl.pos 指定文本标签(变量名称)的位置，当type=full时，默认标签位置在左边和顶部(lt)，当type=lower时，默认标签在左边和对角线(ld)，当type=upper时，默认标签在顶部和对角线，d表示对角线，n表示不添加文本标签；

diag 表示对角线上取值，默认为FALSE

cl.pos 表示图例位置，当type=upper或full时，图例在右方，当type=lower时，图例在底部，不需要图例时，需指定该参数为n；

...

下面给出方法运用：

当然，利用corrplot函数画图可以实现图层的叠加，上面这张图就是分上下两部分完成的，其中默认的颜色样式个人觉得还是可以的，只不过对于相关系数值会根据高低颜色深浅会发生变化，对于一些相关性低的值颜色会非常浅，所以看得不是很明显。

关于这个包具体的的使用方法可以参考这位博主写的文章：

这里我们还是用方案一随机生成的矩阵，介绍用pheatmap包来绘制热力图。pheatmap包里关于绘制热力图的参数相对来说比较少，可以帮助我们快速的绘图，这里给出一些常用参数：

color 设置渐变的颜色，通常借助于colorRampPalette函数，比如说设置红黄蓝渐变，并在这之间分成50个等级，我们可以设置color=colorRampPalette(c("red","yellow","blue"))(50)

cluster_cols &cluster_rows 表示是否按行或列聚类，默认值为FALSE

clustering_method 表示聚类方法，默认是complete,此外还"ward.D",“single”,“average”,等；

display_numbers 表示是否在heatmap里面显示数值，默认是FALSE；

show_rownames &show_colnames 表示是否显示行名或列名；

file 设置图片保存位置

...

下面给出方法运用

去掉边框线可能会好看一点：

由于是随机生成的数据，就不显示聚类的效果（只需要把cluster_row和cluster_col删掉即可），总体来说用pheatmap绘制热图会相对简单一点，但是毫不逊色于其他包绘制的热图。此外，如果想对于行或列来显示一些注释信息（annotation），比如将特征分成2类，每一类是不同的颜色，这里就可以先生成一个行名是特征，列名是分类结果的数据框，然后利用annotation_row(或col)参数，将生成的数据框赋给它即可，具体可以参考这位博主的文章：

接下来介绍的ComplexHeatmap包就比较全面了，他可以兼容pheatmap函数的所有功能，可以说是pheatmap包的加强版，能够创建更加复杂的热力图，如果你会pheatmap包的应用，那么在ComplexHeatmap包里面，你只需要指明是该包下的pheatmap函数即可使用（ComplexHeatmap::pheatmap()）。接下来列举出一些常用参数：

name 、column_title、row_title设置图例、列标题与行标题的名字；

column_title_side &row_title_side 设置列标题与行标题的位置，之注意：列标题只能跟"top"或"buttom"参数，行标题只能跟"left"或"right"参数；

column_names_side &row_names_side 设置行名与列名的位置，后面跟的是位置参数，如"left"、"top"等；

column_names_rot &row_names_rot 设置行名与列名的倾斜角度，后面跟的是角度，如0、30、90等；

column_names_gp &row_names_gp 设置行名与列名的颜色，比如 column_names_gp =gpar(col=rep("red",5))

column_title_gp &row_title_gp 设置列与行标题的颜色，注意：这个需要和聚类分割的数量来决定，要指定row_split &column_split，颜色的设置才能生效；

col 设置渐变的颜色向量参数，这里推荐用RColorBrewer包中的颜色，比如 col = rev(brewer.pal(n = 7, name ="RdYlBu"))

cluster_rows &cluster_columns 表示是否对行列进行聚类，默认是TRUE

cluster_rows &cluster_columns 表示是否对行列进行聚类，默认是TRUE，如果是特定值，则表示对聚类树进行处理；

row_dend_reorder &column_dend_reorder 表示将行或列进行排序，默认是TRUE，所以我们在利用这个包绘制相关系数热力图时，会看到对角线不是1，那么我们就需要检查是否设置了这个参数；

show_column_dend &show_row_dend 表示是否展示行与列的聚类树；

...

下面利用上述随机生成的数据来绘制heatmap：

最值得一提的是，cluster_rows参数，可以结合hlust函数来使用，并通过color_branches函数来为不同类别设置颜色，使得整个heatmap看起来更加美观。如果我们要显示聚类后的数据分割并命名，我们可以这样：

如若想得到更加详细的说明，可以看ComplexHeatmap包的官方文档，或者参见这位博主的文章：

当然，画heatmap怎么能少的了ggplot2呢，我们在利用ggplot画图时，只需要设置scale_fill_gradient即可，例如：scale_fill_gradient(low = "yellow", high = "red") 表示颜色从黄色到红色渐变。注意要把数据处理成ggplot所需要的样式！下面来绘制heatmap:

如果要实现聚类树在heatmap上，我们需要利用ggtree函数，分别绘制聚类树与热力图，最后用aplot包进行拼接即可。

由于ComplexHeatmap包绘制的热力图是一个Heatmap对象，故他与其他图形不同，自身可以与其他Heatmap对象结合，我们只需要利用"+"号或者"%v%"连接符对多个Heatmap对象进行水平或垂直连接就可以了。

当我们需要将pheatmap包绘制的热力图与ggplot画的其他图贴在一起时，我们可以利用ggplotify包来实现，具体操作流程为：

我们用上回利用iris数据集画组合小提琴图的例子，进一步组合heatmap：

当然ggplot也可以画heatmap，这里不再阐述，对于上面几种绘图方案，我们只需选取一种最美观，最有效的方式来画heatmap即可。

此处的ld.x是一个x轴和y轴一样的矩阵，代码中选取了前20行和20列进行画图。这个赋值一个snp.position就能得到snp在染色体的位置。

这里新学到一个函数pa。可以看到LD.plot有一个graphical.par参数，是对图性进行修改的，但是我当时没理解这里具体怎么修改，修改哪些成份。后来发现R语言画图有par函数。直接输入par就可以看到其中有哪些参数了。这些参数和graphical.par是对应的。

要理解LD衰减图，我们就必须先理解连锁不平衡（Linkage

disequilibrium，LD）的概念。连锁不平衡是由两个名词构成，连锁+不平衡。前者，很容易让我们产生概念混淆；后者，让这个概念变得愈加晦涩。因此从一个类似的概念入手，大家可能更容易理解LD的概念，那就是基因的共表达。

基因的共表达，通常指的是两个基因的表达量呈现相关性。比较常见的例子就是：转录组因子和靶基因间的关系。因为转录因子对它的靶基因有正调控作用，所以转录因子的表达量提高会导致靶基因的表达量也上调，两者往往存在正相关关系。这个正相关关系，可以使用相关系数r^2来度量，这个数值在-1~1之间。总而言之，相关性可以理解为两个元素共同变化，步调一致。

类似的，连锁不平衡（LD）就是度量两个分子标记的基因型变化是否步调一致，存在相关性的指标。如果两个SNP标记位置相邻，那么在群体中也会呈现基因型步调一致的情况。比如有两个基因座，分别对应A/a和B/b两种等位基因。如果两个基因座是相关的，我们将会看到某些基因型往往共同遗传，即某些单倍型的频率会高于期望值。

例如在下图中，在群体中（A，a，B，b）各个基因型的频率已知的情况下，各种单倍型的期望频率（AB、Ab、aB、ab）都是可以计算出来。例如，AB的频率=（A的频率）X（B的频率）。但我们实际统计群体中各个单倍型的频率的时候，会观察到某些单倍型的频率会大于期望值，例如下图中的单倍型AB的理论频率是0.12，但观察到的实际频率是0.29。那么说明，基因型A更倾向于基因型B共同遗传。

这一般往往是由于在祖先的基因组中，A和B就是位于同一条染色体上，在传代过程中，这种共同遗传的关系被保留了下来。位点间的这种相关性，在杂交家系中一般被称为连锁（孟德尔老师豌豆实验中的发现），在自然群体中则一般被称为连锁不平衡。所以连锁不平衡中的“不平衡”，我认为可以理解为单倍型的频率分布偏离期望值，偏离了平衡。

这种不同基因组间的相关性，用一个数值来衡量就是D值（图中有计算公式）。类似相关系数是标准化后的协方差，LD系数（r^2）则是标准化后的D值（图中有计算公式），这个数值在0~1波动。r^2=0就是两个位点完全不相关，群体中单倍型分布是随机的（观测值=期望值）。r^2=1就是两个位点完全相关，某些基因型（A）只与特定的基因型（B）共同出现。

一般而言，两个位点在基因组上离得越近，相关性就越强，LD系数就越大。反之，LD系数越小。也就是说，随着位点间的距离不断增加，LD系数通常情况下会慢慢下降。这个规律，通常就会使用LD衰减图来呈现。

LD衰减图就是利用曲线图来呈现基因组上分子标记间的平均LD系数随着标记间距离增加而降低的过程。大概的计算原理就是先统计基因组上两两标记间的LD系数大小，再按照标记间的距离对LD系数进行分类，最终可以计算出一定距离的分子标记间的平均LD系数大小。第一个图是黄瓜重测序文章中统计各个亚群体的LD衰减速度的图形。横坐标是物理距离（kb），纵坐标是LD系数（r^2）。

从图中我们可以看出，西双版纳这个亚群体（紫色线）在基因组上50kb距离的平均LD系数大小约为0.4，但到了100kb的距离，对应的平均LD系数大小则降低到了不到0.3。而且，我们从图中也可以观察到LD系数的衰减速度在不同的亚群体快慢不同，衰减速度是

india >East Asian&Eurasian >

Xishuanbanna。那说明india群体的LD衰减距离最大，可能是india这个群体遗传多样性最高导致。

LD衰减距离

实际上，LD衰减的速度在不同物种间或同物种的不同亚群体间，往往差异非常巨大。所以，通常会使用1个标准——“LD衰减距离”来描述LD衰减速度的快慢。

LD衰减距离通常指的是：当平均LD系数衰减到一定大小的时候，对应的物理距离。

“一定大小”是这个定义的关键点，但没有特别统一的标准，在不同文章中标准不同。常见的标准包括：

a）LD系数降低到最大值的一半；

b）LD系数降低到0.5以下；

c）LD系数降低到0.1以下；

d）LD系数降低到基线水平（但注意，不同材料的基线值是不同的。比如图3黄瓜群体的基线大概是0.1）。

所以，下次你在文章中看到“LDdecay distance is

XXkb”的时候，别忘了看看作者使用的标准是什么。

LD衰减会受什么的影响？

如第一个图所示， LD系数衰退速度会受到不同因素的影响而有所不同。常见的因素包括：

1）物种类型

LD存在的本质是两个位点的连锁遗传导致的相关性。但这种相关性理论上会随着世代的增加、重组次数的增加而不断下降。所以，那些繁殖力强、时代间隔短的物种（例如，昆虫），其LD衰减的速度是非常快的。例如在家蚕和野蚕群体中，LD系数下降到最大值的1/2仅仅需要46bp和7bp的距离。

2）群体类型

相同物种的不同群体，由于其遗传背景不同，LD衰减速度也存在很大的差异。驯化选择，会导致群体遗传多样性下降，位点间的相关性（连锁程度）加强。所以，通常驯化程度越高，选择强度越大的群体，LD衰减速度是最慢的。例如，栽培稻比野生稻通常更大的LD衰减距离。类似的，自然选择、遗传漂变导致的群体遗传多样性下降，也会减慢LD衰减的速度。

3）在染色体的位置

染色体不同区域的LD衰减距离而是不同的。通常着丝粒区更易重组，所以LD衰减更慢。而基因组上那些受选择的区域相比普通的区域，LD衰减速度也是更慢的。

LD衰减距离的应用

LD衰减速度，在群体遗传分析中本身是对群体特性的评估，与群体类型的特性（自然群体还是驯化群体，选择强度大小）是相关的。但在其他研究中还有更多的应用价值。

基于分子标记（例如，SNP芯片，GBS测序）的GWAS分析，其实并没有检测到功能突变，本质就是利用标记和功能突变的相关性（LD关系），来检测与性状相关的功能突变的位置。一般而言，LD系数大于0.8就是强相关。如果LD系数小于0.1，则可以认为没有相关性。如果LD衰减到0.1这么大的区间内都没有标记覆盖的话，即使这个区间有一个效应很强的功能突变，也是检测不到关联信号的。所以，通常可以通过比较LD衰减（到0.1）距离和标记间的平均距离，来判断标记是否对全基因组有足够的覆盖度。

而如果GWAS检测到显著关联的区间后，则可以通过进一步绘制局部的LD单体型块图，来进一步判断显著相关的SNP和目标基因间是否存在强LD关系。这个图形我们下一篇文章会介绍。

再提一个应用的例子。在进行STRUCTURE分析的时候理论上必须输入不相关的位点。那么，就可以通过预估LD衰减到0.1的距离，来判断标记间的距离必须大于多少才能保证标记间不具相关性（LD<0.1）。

LD衰减图的绘制，实际上有两个步骤：

1）计算marker间两两的LD系数大小

这个可以使用haploview软件完成。计算的时候，只要设定一个关键的参数：区间大小。例如设定为5Mb，那么软件就会计算基因组上所有距离<5Mb的两两位点间的LD系数。实际上这个参数设定更大也没有意义，一般情况下位点间的相关性不会延伸到大于5Mb这么远的距离。

2）绘图

将LD系数按照对应的两个marker间的距离进行分类，例如：距离按照区间大0~5k,5k~10k,10k~15k…..分别分类。如果重测序的数据，SNP标记密度较大，这个分类区间可以设置小一些；如果是简化基因组数据，SNP标记较为稀疏，则分类区间可以适当加大。然后计算每种距离分类的LD系数的均值。最后在利用均值绘制曲线图就ok了。这一步的绘图，使用excel或R语言都可以轻松完成。

引用http://blog.sina.com.cn/s/blog_4ab0b3390102wfnm.html