ggplot2中两大精髓的函数分别为qplot(快速作图quickplot)，类似于R中基本的plot和ggplot(更符合ggplot2图层式绘图理念，一层层添加修改)。
测试数据是ggplot2包中自带的diamond数据，每一行为一种钻石，每一列为钻石不同的属性，如carat(克拉),cut(切工),color(色泽),clarity(透明度)等。
绘制散点图，对x,y值取log
颜色、大小、性状和其他属性的设置
几何对象
qplot()函数配合不同的几何对象便可绘制出不同的图形:
绘制其他常见图形
以上为ggplot2包中常见图形的快速绘制(quickplot)即qplot函数的应用。qplot函数还有很多其他的参数,对xlim,ylim设置x,y轴的区间例如xlim=c(0,20);对轴取log值，log="x",对x轴取对数，log="xy"表示对x轴和y轴取对数;main：图形的主题main=“qplottitle”;#xlab,ylab:设置轴标签文字
本篇使用R的ggplot2包实现从原始数据读入到热图输出的过程，并在教程结束后提供一份封装好的命令行绘图工具，只需要提供矩阵，即可一键绘图。
上一篇讲述了Rstudio的使用作为R写作和编译环境的入门，后面的命令都可以拷贝到Rstudio中运行，或写成一个R脚本，使用Rscriptheatmap.r运行。我们还提供了Bash的封装，在不修改R脚本的情况下，改变参数绘制出不同的图形。
绘图首先需要数据。通过生成几组向量，转换为矩阵，得到想要的数据。
还有另外2种获取数值矩阵的方式。
melt：把正常矩阵转换为长表格模式的函数。工作原理是把全部的非id列的数值列转为1列(列名默认为value)；所有字符列转为一列，列名默认为variable。
热图出来了，但有点不对劲，横轴重叠一起了。一个办法是调整图像的宽度，另一个是旋转横轴标记。
设置想要的颜色。
调整legend的位置,legend.position,可以接受的值有top,bottom,left,right,和一个坐标c(0.05,0.8)(左上角，坐标是相对于图的左下角(即原点)计算的）
调整背景和背景格线以及X轴、Y轴的标题。
为了使横轴旋转45度，需要把这句话theme(axis.text.x=element_text(angle=45,hjust=1,vjust=1))放在theme_bw()的后面。
合并以上命令，就得到了下面这个看似复杂的绘图命令。
也可以只用Point
图形出来了，就得考虑存储了，一般输出为PDF格式，方便后期的修改。
至此，完成了简单的heatmap的绘图。但实际绘制时，经常会碰到由于数值变化很大，导致颜色过于集中，使得图的可读性下降很多。因此需要对数据进行一些处理，具体的下次再说。
上面的测试数据，数值的分布比较均一，相差不是太大，但是Gene_4和Gene_5由于整体的值低于其它的基因，从颜色上看，不仔细看，看不出差别。
图中只有右上角可以看到红色，其他地方就没了颜色的差异。这通常不是我们想要的。为了更好的可视化效果，需要对数据做些预处理，主要有对数转换，Z-score转换，抹去异常值，非线性颜色等方式。
假设下面的数据是基因表达数据，4个基因(a,b,c,d)和5个样品(Grp_1,Grp_2,Grp_3,Grp_4)，矩阵中的值代表基因表达FPKM值。
对数转换后的数据，看起来就清晰的多了。而且对数转换后，数据还保留着之前的变化趋势，不只是基因在不同样品之间的表达可比(同一行的不同列)，不同基因在同一样品的值也可比(同一列的不同行)(不同基因之间比较表达值存在理论上的问题，即便是按照长度标准化之后的FPKM也不代表基因之间是完全可比的)。
Z-score又称为标准分数，是一组数中的每个数减去这一组数的平均值再除以这一组数的标准差，代表的是原始分数距离原始平均值的距离，以标准差为单位。可以对不同分布的各原始分数进行比较，用来反映数据的相对变化趋势，而非绝对变化量。
Z-score转换后，颜色分布也相对均一了，每个基因在不同样品之间的表达的高低一目了然。但是不同基因之间就完全不可比了。
粗暴一点，假设检测饱和度为100，大于100的值都视为100对待。
虽然损失了一部分信息，但整体模式还是出来了。但是在选择异常值标准时需要根据实际确认。
正常来讲，颜色的赋予在最小值到最大值之间是均匀分布的。如果最小值到最大值之间用100个颜色区分，则其中每一个bin，不论其大小、有没有值都会赋予一个颜色。非线性颜色则是对数据比较小但密集的地方赋予更多颜色，数据大但分布散的地方赋予更少颜色，这样既能加大区分度，又最小的影响原始数值。通常可以根据数据模式，手动设置颜色区间。为了方便自动化处理，也可选择用四分位数的方式设置颜色区间。
如果想保持图中每一行的顺序与输入的数据框一致，需要设置因子的水平。这也是ggplot2中调整图例或横纵轴字符顺序的常用方式。
基于ggplot2的heatmap绘制到现在就差不多了，但总是这么画下去也会觉得有点累，有没有办法更简化呢？。
绘制热图除了使用ggplot2，还可以有其它的包或函数，比如pheatmap::pheatmap(pheatmap包中的pheatmap函数)、gplots::heatmap.2等。
相比于ggplot2作heatmap,pheatmap会更为简单一些，一个函数设置不同的参数，可以完成行列聚类、行列注释、Z-score计算、颜色自定义等。那我们来看看效果怎样。
虽然有点丑，但一步就出来了。
在heatmap美化篇提到的数据前期处理方式，都可以用于pheatmap的画图。此外Z-score计算在pheatmap中只要一个参数就可以实现。
有时可能不需要行或列的聚类，原始展示就可以了。
给矩阵(data)中行和列不同的分组注释。假如有两个文件，第一个文件为行注释，其第一列与矩阵中的第一列内容相同(顺序没有关系)，其它列为第一列的不同的标记，如下面示例中(假设行为基因，列为样品)的2,3列对应基因的不同类型(TForenzyme)和不同分组。第二个文件为列注释，其第一列与矩阵中第一行内容相同，其它列则为样品的注释。
自定义下颜色吧。
自己做个hclust传进去，顺序跟pheatmap默认是一样的
可以按任意指标排序，基因表达是一个例子。
样本量多时的自动较忧排序
箱线图是能同时反映数据统计量和整体分布，又很漂亮的展示图。在2014年的NatureMethod上有2篇Correspondence论述了使用箱线图的好处和一个在线绘制箱线图的工具。就这样都可以发两篇Naturemethod，没天理，但也说明了箱线图的重要意义。
下面这张图展示了Barplot、Boxplot、Volinplot和Beanplot对数据分布的反应。从Barplot上只能看到数据标准差或标准误不同；Boxplot可以看到数据分布的集中性不同；Violinplot和Beanplot展示的是数据真正的分布，尤其是对Biomodal数据的展示。
假设有这么一个基因表达矩阵，第一列为基因名字，第一行为样品名字，想绘制样品中基因表达的整体分布。
箱线图出来了，看上去还可以，再加点色彩(fill)。
再看看Violinplot
还有Jitterplot(这里使用的是ggbeeswarm包)
也可以用geom_jitter(aes(colour=factor(variable)))代替geom_quasirandom(aes(colour=factor(variable)))绘制抖动图，但个人认为geom_quasirandom给出的结果更有特色。
为了更好的展示效果，下面的矩阵增加了样品数量和样品的分组信息。
常规矩阵绘制箱线图要求必须是个方正的矩阵输入，而有时想比较的几个组里面检测的值数目不同。比如有三个组，GrpA组检测了6个病人，GrpB组检测了10个病人，GrpC组是12个正常人的检测数据。这时就很难形成一个行为检测值，列为样品的矩阵，就需要长表格模式解决这一问题。
长表格形式自身就是常规矩阵melt后的格式，这种用来绘制箱线图就很简单了，就不举例子了。
线图是反映趋势变化的一种方式，其输入数据一般也是一个矩阵。
假设有这么一个矩阵，第一列为转录起始位点及其上下游5kb的区域，第二列为H3K27ac修饰在这些区域的丰度，想绘制一张线图展示。
满心期待一个倒钟形曲线，结果，什么也没有。
仔细看，出来一段提示
图出来了，一条线，看一眼没问题；再仔细看，不对了，怎么还不是倒钟形，原来横坐标错位了。
检查下数据格式
图终于出来了，调了下legend的位置，看上去有点意思了。
有点难看，如果平滑下，会不会好一些，stat_smooth可以对绘制的线进行局部拟合。在不影响变化趋势的情况下，可以使用(但慎用)。
从图中看，趋势还是一致的，线条更优美了。另外一个方式是增加区间的数量，线也会好些，而且更真实。
stat_smooth和geom_line各绘制了一条线，只保留一条就好。
好了，终于完成了单条线图的绘制。
那么再来一个多线图的例子吧，只要给之前的数据矩阵多加几列就好了。
如果横轴是文本，又该怎么调整顺序呢？还记得之前热图旁的行或列的顺序调整吗？重新设置变量的factor水平就可以控制其顺序。
比较下位置信息做为数字(前面的线图)和位置信息横轴的差别。当为数值时，ggplot2会选择合适的几个刻度做标记，当为文本时，会全部标记。另外文本横轴，smooth效果不明显。
绘制散点图的输入一般都是规规矩矩的矩阵，可以让不同的列分别代表X轴、Y轴、点的大小、颜色、形状、名称等。
火山图用于展示基因表达差异的分布，横轴为Log2FoldChange，越偏离中心差异倍数越大；纵轴为(-1)*Log10P_adjust，值越大差异越显著。一般横轴越偏离中心的点其纵轴值也会比较大，因此呈现火山喷发的形状。
火山图需要的数据格式如下
说好的火山图的例子，但怎么也看不出喷发的态势。
对数据坐下预处理，差异大的基因padj小，先对其求取负对数，所谓负负得正，差异大的基因就会处于图的上方了。
有点意思了，数据太少不明显，下一步加上颜色看看。
利用现有的数据，基本上就是这个样子了。虽然还不太像，原理都已经都点到了。
六个点少了只剩下了3个，重叠在一起了，而且也不知道哪个点代表什么样品。这时需要把点抖动下，用到一个包ggbeeswarm，抖动图的神器。
假设有这么一个富集分析结果矩阵(文件名为GOenrichement.xls)存储了EHBIO样品和Baodian样品中各自上调的基因富集的通路。
通过这张图解释下，富集分析的结果怎么解读。富集分析实际是查找哪些通路里面包含的差异基因占总差异基因的比例显著高于通路中总基因占所有已经注释的基因的比例。这一显著性通常用多重假设检验矫正过的pvalue(又称qvalue,FDR或p.adjust)来表示。在图中体现为点的颜色。从绿到红富集显著性逐渐增高。点的大小表示对应通路中包含的差异基因的数目。点的形状代表了不同类型的基因，如EHBIO中上调的基因和Baodian中上调的基因。横轴表示对应通路包含的差异基因占总的差异基因的比例,本图中最高不过5%,这个值越大说明通路被影响的越多。
维恩图是用来反映不同集合之间的交集和并集情况的展示图。一般用于展示2-5个集合之间的交并关系。集合数目更多时，将会比较难分辨，更多集合的展示方式一般使用upSetView。
假设有这么一个矩阵，第一列为不同集合中的ID，第二列为集合的名字，无标题行，存储为venn.txt。
对于集合比较多的时候，包括上面提到的5个集合的交并集情况，如果只是为了展示个炫图，还可以，但如果想解释结果，就会比较头疼，难判断区域的归属。
因此对于这种多集合情况，推荐使用UpSetView展示，看效果如下。
测试数据，存储为upsetview.txt（第一行为集合的名，每个集合一列；每一行为一个ID，如果对应ID出现在这个集合则标记1，否则标记0）：
单个点表示特有，连起来的点表示共有，相当于venn图中重叠的部分。
垂直的柱子代表的是Venn图中的数字，看连接的点判断归属。
水平的柱子代表对应样品中Item的总数。
有如下4个基因在5组样品中的表达值
柱子有点多，也可以利用mean±SD的形式展现
首先计算平均值和标准差，使用group_by按gene分组，对每组做summarize
设置误差线的宽度和位置
每个基因的原始表达值堆积柱状图(只需要修改positon=stack)
堆积柱状图显示没问题，但文本标记错位了
指定下分组信息，位置计算就正确了
比较每组各个基因的相对表达(position=fill)
纵轴的显示改为百分比
在柱子中标记百分比值
首先计算百分比，同样是group_by(按照给定的变量分组，然后按组操作)和mutate两个函数(在当前数据表增加新变量)
每组里面各个基因的相对表达,纵轴的显示改为百分比
facet后，显示正常，不需要做特别的修改
在柱子中标记百分比值(计算百分比值需要注意了,文本显示位置还是跟之前一致)
文本显示位置没有问题，但柱子的位置有些奇怪，使得两组之间不可比。
先对数据做下排序，然后再标记文本(这部分可以看下视频，或者自己输出下数据看看为什么要排序。底层原因可能是这种新的数据结构的问题。)
这样两种条件下的比较更容易了。
ggplot2中修改图形字体。
#修改坐标轴和legend、标题的字体theme(text=element_text(family="Arial"))#或者theme_bw(base_family="Arial")#修改geom_text的字体geom_text(family="Arial")3.12.2ggplot2支持中文字体输出PDFshowtext包可给定字体文件，加载到R环境中，生成新的字体家族名字，后期调用这个名字设定字体，并且支持中文写入pdf不乱码
library(showtext)showtext.auto(enable=TRUE)font_path="FZSTK.TTF"font_name=tools::file_path_sans_ext(basename(font_path))font.add(font_name,font_path)#修改坐标轴和legend、标题的字体theme(text=element_text(family=font_name))#修改geom_text的字体geom_text(family=font_name)3.12.3系统可用字体Linux字体一般在/usr/share/fonts下，也可以使用fc-list列出所以加载的字体。
Windows字体在C:\Windows\Fonts\下，直接可以看到，也可以拷贝到Linux下使用。
通常情况下，作图的字体都是英文，ggplot2默认的或按需求加载一种字体就可以了。但如果中英文混合出现时，单个字体只能支持一种文字，最好的方式是合并两种字体，类似于Word中设置中英文分别使用不同的字体。
假如需要合并FZSTK.TTF(windows下获取)和Schoolbell-Regular.ttf(谷歌下载)，这两个都是手写字体。按如下，把字体文件和程序脚本mergefont.pe放在同一目录下，运行fontforge-scriptmergefont.pe即可获得合并后的字体FZ_School.ttf。
ct@ehbio$lsFZSTK.TTFmergefont.peSchoolbell-Regular.ttfct@ehbio$catmergefont.peOpen("FZSTK.TTF")SelectAll()ScaleToEm(1024)Generate("temp.ttf","",0x14)Close()#OpenEnglishfontandmergetotheChinesefontOpen("Schoolbell-Regular.ttf")SelectAll()ScaleToEm(1024)MergeFonts("temp.ttf")SetFontNames("FZ_School","FZST","Schoolbel","Regular","")Generate("FZ_School.ttf","",0x14)Close()ct@ehbio$fontforge-scriptmergefont.pect@ehbio$lsFZ_School.ttfFZSTK.TTFmergefont.peSchoolbell-Regular.ttf然后安装前面的介绍使用showtext导入即可使用。
字体文件自己从Windows获取，Schoolbell从Googlefonts获取。
library(showtext)##Addfontsthatareavailableoncurrentpath#方正字体+schoolebell(中英混合)font.add("FZ_School","font/FZ_School.ttf")#黑体font.add("simhei","font/simhei.ttf")font.add("Arial","font/arial.ttf")#黑体和Arial的合体font.add("HeiArial","font/HeiArial.ttf")showtext.auto()##automaticallyuseshowtextfornewdeviceslibrary(ggplot2)p=ggplot(NULL,aes(x=1:10,y=2^(1:10),group=1))+geom_line()+theme(axis.title.y=element_text(family="Arial"), axis.title.x=element_text(family="HeiArial"),plot.title=element_text(family="simhei"))+xlab("Daysspenton生信宝典")+ylab("Thingsyouhavelearned")+ggtitle("生信宝典，换个角度学生信")+annotate("text",7,300,family="FZ_School",size=8,label="收获曲线(Harvestcurve)",angle=15)#annotate指定的是文字的中间部分的位置#ggsave(p,filename="example-SXBD.pdf",width=7,height=4)##PDFdevicep
主成分是原有变量的线性组合，其数目不多于原始变量。组合之后，相当于我们获得了一批新的观测数据，这些数据的含义不同于原有数据，但包含了之前数据的大部分特征，并且有着较低的维度，便于进一步的分析。
在空间上，PCA可以理解为把原始数据投射到一个新的坐标系统，第一主成分为第一坐标轴，它的含义代表了原始数据中多个变量经过某种变换得到的新变量的变化区间；第二成分为第二坐标轴，代表了原始数据中多个变量经过某种变换得到的第二个新变量的变化区间。这样我们把利用原始数据解释样品的差异转变为利用新变量解释样品的差异。
这种投射方式会有很多，为了最大限度保留对原始数据的解释，一般会用最大方差理论或最小损失理论，使得第一主成分有着最大的方差或变异数(就是说其能尽量多的解释原始数据的差异)；随后的每一个主成分都与前面的主成分正交，且有着仅次于前一主成分的最大方差(正交简单的理解就是两个主成分空间夹角为90°，两者之间无线性关联，从而完成去冗余操作)。
简化运算。
比如说我们在样品的制备过程中，由于不完全一致的操作，导致样品的状态有细微的改变，从而造成一些持家基因也发生了相应的变化，但变化幅度远小于核心基因(一般认为噪音的方差小于信息的方差）。而PCA在降维的过程中滤去了这些变化幅度较小的噪音变化，增大了数据的信噪比。利用散点图实现多维数据可视化。
假设有一套数据集，包含100个样品中某一基因的表达量。如下所示，每一行为一个样品，每一列为基因的表达值。这也是做PCA分析的基本数据组织方式，每一行代表一个样品，每一列代表一组观察数据即一个变量。
那么如果有2个基因呢？
如果有3个基因呢？
当我们向由Gene1和Gene2构成的X-Y平面做垂线时，可以很明显的看出，Gene2所在的轴对样品的分类没有贡献。因为投射到X-Y屏幕上的点在Y轴几乎处于同一位置。
我们把坐标轴做一个转换，可以看到在由Gene1和Gene3构成的X-Y平面上，样品被分为了4类。Gene2对样品的分类几乎没有贡献，因为几乎所有样品在Gene2维度上的值都一样。
在上述例子中，我们可以很容易的区分出Gene1和Gene3可以作为分类的主成分，而Gene2则对分类没有帮助，可以在计算中去除。
但是如果我们测序了几万个基因的表达时，就很难通过肉眼去看，或者作出一个图供我们筛选哪些基因对样本分类贡献大。这时我们应该怎么做呢？
其中有一个方法是，在这个基因表达矩阵中选出3个变化最大的基因做为3个主成分对样品进行分类。我们试验下效果怎么样。
我们看到图中的样品并没有按照预先设定的标签完全分开。当然我们也可以通过其他方法筛选变异最大的三个基因，最终的分类效果不会相差很大。因为不管怎么筛选，我们都只用到了3个基因的表达量。
假如我们把这个数据用PCA来分类，结果是怎样的呢？
采用3个主成分获得的分类效果优于2个主成分，因为这样保留的原始信息更多。
方差和协方差的计算公式如下：
\[方差Var(X)=\frac{\sum_{i=1}^{n}(X_i-\barX)^2}{n-1}\]
\[协方差cov(X,Y)=\frac{\sum_{i=1}^{n}(X_i-\barX)(Y_i-\barY)}{n-1}\]
如果数据的均值为0，这个公式可以进一步简化。简化后的公式把计算协方差转变为了矩阵乘法运算。这也是为什么PCA需要中心化数据。
\[方差Var(X)=\frac{\sum_{i=1}^{n}{X_i}^2}{n-1}\]
\[协方差cov(X,Y)=\frac{\sum_{i=1}^{n}X_iY_i}{n-1}\]
\[协方差矩阵cov(X,Y)=\frac{X_{n,m}^TY_{n,m}}{n-1}\]
假如我们有一个矩阵如下，
mat中5组基因的表达值的方差计算如下：
\[\mathbf{cov(mat)}=\frac{1}{3}\left[\begin{array}{rrrr}0.16&0.83&-0.72&-0.26\\-0.41&-0.04&-0.52&0.97\\-0.22&-0.02&-0.82&1.05\\-0.35&-0.39&-1.56&2.3\\0.16&2.24&-0.43&-1.96\\\end{array}\right]\left[\begin{array}{rrrrr}0.16&-0.41&-0.22&-0.35&0.16\\0.83&-0.04&-0.02&-0.39&2.24\\-0.72&-0.52&-0.82&-1.56&-0.43\\-0.26&0.97&1.05&2.3&-1.96\\\end{array}\right]=\left[\begin{array}{rrrrr}0.43&0.01&0.09&0.04&0.9\\0.01&0.46&0.51&1.07&-0.61\\0.09&0.51&0.61&1.26&-0.6\\0.04&1.07&1.26&2.67&-1.59\\0.9&-0.61&-0.6&-1.59&3.02\\\end{array}\right]\]
而做PCA分析就是产生一组新的变量，使得新变量的协方差矩阵为对角阵，满足上面的要求。从而达到去冗余的目的。然后再选取方差大的变量，实现降维和去噪。
如果正向推导，这种组合可能会有很多种，一一计算会比较麻烦。那反过来看呢？我们不去寻找这种组合，而是计算如何使原变量的协方差矩阵变为对角阵。
数学推导中谨记的两个概念：
前面提到，新变量(\(Y_{m,k}\))是原始变量(\(X_{m,n}\))(原始变量的协方差矩阵为(\(C_{n,n}\)))的线性组合，那么假设我们找到了这么一个线性组合(命名为特征矩阵(\(P_{n,k}\)))，得到一组新变量\(Y_{m,k}=X_{m,n}P_{n,k}\)，并且新变量的协方差矩阵(\(D_{k,k}\))为对角阵。那么这个特征矩阵(\(P_{n,k}\))需要符合什么条件呢？
从矩阵运算可以看出，最终的特征矩阵(\(P_{n,k}\))需要把原变量协方差矩阵(\(C_{n,n}\))转换为对角阵(因为新变量的协方差矩阵(\(D_{k,k}\))为对角阵)，并且对角元素从大到小排列（保证每个主成分的贡献度依次降低）。
现在就把求解新变量的任务转变为了求解原变量协方差矩阵的对角化问题了。在线性代数中，矩阵对角化的问题就是求解矩阵的特征值和特征向量的问题。
我们举一个例子讲述怎么求解特征值和特征向量。
假设\(A_{n,n}\)为n阶对称阵，如存在\(\lambda\)和非零向量\(x\)，使得\(Ax=\lambdax\)，则称\(\lambda\)为矩阵\(A_{n,n}\)的特征值，非零向量\(x\)为为矩阵\(A_{n,n}\)对应于特征值\(\lambda\)的特征向量。
根据这个定义可以得出\((A_{n,n}-\lambdaE)x=0\)，由于\(x\)为非零向量，所以行列式\(|A-\lambdaE|=0\)。
由此求解出n个根\(\lambda_{1},\lambda_{2},...,\lambda_{3}\)就是矩阵\(A\)的特征值。
回顾下行列式的计算：
假如我们有一个矩阵\(\mathbf{A}=\left[\begin{array}{cc}3&-1\\-1&3\end{array}\right]\)，如何计算它的特征值和特征向量呢？
则\(\lambda\)的值为2或4。
对\(\lambda_{1}=2\)时，求解\((A-2E)x=\left|\begin{array}{cc}1&-1\\-1&1\end{array}\right|x=0\)，得\(x=k\left|\begin{array}{c}1\\1\end{array}\right|\)，则对应于\(\lambda_{1}=2\)时的特征向量\(p_{1}=\left|\begin{array}{c}1\\1\end{array}\right|\)
对\(\lambda_{2}=4\)时，求解\((A-2E)x=\left|\begin{array}{cc}-1&-1\\-1&-1\end{array}\right|x=0\)，得\(x=k\left|\begin{array}{c}1\\-1\end{array}\right|\)，则对应于\(\lambda_{2}=4\)时的特征向量\(p_{2}=\left|\begin{array}{c}1\\-1\end{array}\right|\)
以上就完成了PCA的数学推导。
我们使用前面用到的数据data3来演示下如何用R函数实现PCA的计算，并与R中自带的prcomp做个比较。
利用prcomp进行主成分分析
自定义PCA计算函数
prcomp函数会返回主成分的标准差、特征向量和主成分构成的新矩阵。接下来，探索下不同主成分对数据差异的贡献和主成分与原始变量的关系。
这两个信息可以判断主成分分析的质量：
指导选择主成分的数目：
鉴定核心变量和变量间的隐性关系:
在测试数据中，scale后，三个主成分对数据差异的贡献度大都在30%左右，而未scale的数据，三个主成分对数据差异的贡献度相差很大。这是因为三个基因由于自身表达量级所引起的方差的差异导致它们各自对数据的权重差异，从而使主成分偏向于数值大的变量。
前面用到一组比较大的测试数据集，并做了PCA分析，现在测试不同的处理对结果的影响。
首先回顾下我们用到的数据。
如果首先提取500个变化最大的基因，再执行PCA分析会怎样呢？
一般说来，在PCA之前原始数据需要中心化（centering，数值减去平均值）。中心化的方法很多，除了平均值中心化（mean-centering）外，还包括其它更稳健的方法，比如中位数中心化等。
除了中心化以外，定标(Scale,数值除以标准差)也是数据前处理中需要考虑的一点。如果数据没有定标，则原始数据中方差大的变量对主成分的贡献会很大。数据的方差与其量级成指数关系，比如一组数据(1,2,3,4)的方差是1.67，而(10,20,30,40)的方差就是167,数据变大10倍，方差放大了100倍。
但是定标(scale)可能会有一些负面效果，因为定标后变量之间的权重就是变得相同。如果我们的变量中有噪音的话，我们就在无形中把噪音和信息的权重变得相同，但PCA本身无法区分信号和噪音。在这样的情形下，我们就不必做定标。
一般而言，对于度量单位不同的指标或是取值范围彼此差异非常大的指标不直接由其协方差矩阵出发进行主成分分析，而应该考虑对数据的标准化。比如度量单位不同，有万人、万吨、万元、亿元，而数据间的差异性也非常大，小则几十大则几万，因此在用协方差矩阵求解主成分时存在协方差矩阵中数据的差异性很大。在后面提取主成分时发现，只提取了一个主成分，而此时并不能将所有的变量都解释到，这就没有真正起到降维的作用。此时就需要对数据进行定标(scale)，这样提取的主成分可以覆盖更多的变量，这就实现主成分分析的最终目的。但是对原始数据进行标准化后更倾向于使得各个指标的作用在主成分分析构成中相等。对于数据取值范围不大或是度量单位相同的指标进行标准化处理后，其主成分分析的结果与仍由协方差矩阵出发求得的结果有较大区别。这是因为对数据标准化的过程实际上就是抹杀原有变量离散程度差异的过程，标准化后方差均为1，而实际上方差是对数据信息的重要概括形式，也就是说，对原始数据进行标准化后抹杀了一部分重要信息，因此才使得标准化后各变量在主成分构成中的作用趋于相等。因此，对同度量或是取值范围在同量级的数据还是直接使用非定标数据求解主成分为宜。
中心化和定标都会受数据中离群值（outliers）或者数据不均匀（比如数据被分为若干个小组）的影响，应该用更稳健的中心化和定标方法。
生存分析一个常用的方法是寿命表法。
读入数据
PAM50是通过50个基因的表达量把乳腺癌分为四种类型(LuminalA,LuminalB,HER2-enriched,andBasal-like)作为预后的标志。根据PAM50属性对病人进行分组，评估比较两组之间生存率的差别。
简化Stage信息，先只查看大的阶段
目前这部分已不再更新，具体见后面提到的在线绘图网站（在绘制每个图的同时，也提供了这个图的绘制代码）
一步作图相比于直接写R代码有什么好处呢？
模块化好。数据处理和可视化分开；一步作图，只是作图，不做无关的处理，更随意。
易用性强。敲代码，总不如给改数来得快。
重用性强。假如我要做十几个分开的基因集合呢？一段段复制代码？改错了或忘记改某个地方了怎么办？
快速出图。先快速出个原型，再接着调整。
s-plot可以绘制的图的类型还有一些，列举如下；
为了简化绘图、维持脚本的一致，我用bash对R做了一个封装，然后就可以通过修改命令好参数绘制不同的图了。
先看一看怎么使用
首先把测试数据存储到文件中方便调用。数据矩阵存储在heatmap_data.xls文件中；行注释存储在heatmap_row_anno.xls文件中；列注释存储在heatmap_col_anno.xls文件中。
字有点小，是因为图太大了，把图的宽和高缩小下试试。
横轴的标记水平放置
使用正常矩阵默认参数绘制箱线图
#-f:指定输入的矩阵文件，第一列为行名字，第一行为header 列数不限，列名字不限；行数不限，行名字默认为文本sp_boxplot.sh-fboxplot.normal.data箱线图出来了，但有点小乱。
#-f:指定输入的矩阵文件，第一列为行名字，第一行为header 列数不限，列名字不限；行数不限，行名字默认为文本#-P:none,去掉legend(uppercaseP)#-b:X-axis旋转45度#-V:TRUE绘制小提琴图sp_boxplot.sh-fboxplot.normal.data-Pnone-b45-VTRUE
绘制单个基因的小提琴图加抖动图
使用melted矩阵默认参数绘箱线图
#-f:指定输入文件#-mTRUE:指定输入的矩阵为meltedformat#-dExpr：指定表达值所在的列#-FRep:指定子类所在列，也就是legend#-aGroup：指定X轴分组信息#-jTRUE:jitterplotsp_boxplot.sh-fboxplot.melt.data-mTRUE-dExpr-FRep-aGroup-jTRUE
#如果没有子类，则-a和-F指定为同一值#-RTRUE:旋转boxplotsp_boxplot.sh-fboxplot.melt.data-mTRUE-dExpr-aGroup-FGroup-JTRUE-RTRUE
参数中最需要注意的是引号的使用：
#-f:指定输入的矩阵文件，第一列为行名字，第一行为header 列数不限，列名字不限；行数不限，行名字默认为文本#-AFALSE:指定行名为数字sp_lines.sh-fline_data.xls-AFALSE
使用melted矩阵默认参数绘制个线图(除需要改变文件格式，指定-mTRUE,-axvariable外其它与正常矩阵一样)
完整的图
完整参数列表如下：
Samp Gene1 Gene2 Color Size GC_quality Base_qualitya 1 1 grp1 10 PASS PASSb 2 2 grp1 10 PASS PASSc 1 3 grp1 10 WARN PASSd 3 1 grp2 15 WARN WARNe 2 2 grp2 15 PASS WARNf 3 3 grp3 5 PASS PASSg 2 1 grp3 5 WARN PASS想绘制样品在这两个Gene为轴的空间的分布，并标记样品的属性，只需要运行如下命令
如果横纵轴为字符串，且有重复,则需指定参数-JTRUE以错开重叠的点，具体如下
只有想不到，没有做不到，sp_scatterplot2.sh还可以完成更多你想做的散点图，而且只需调参数，无需改代码，简单可重用。