R数据可视化: PCA和PCoA图, 2D和3D

主成分分析(Principal Components Analysis，PCA) ，也称主分量分析或主成分回归分析法，是一种无监督的数据降维方法。PCA通过线性变换将原始数据变换为一组各维度线性无关的表示，可用于提取数据的主要特征分量，常用于高维数据的 降维 。这种降维的思想首先减少数据集的维数，同时还保持数据集的对方差贡献最大的特征，最终使数据直观呈现在二维坐标系。

直观上，第一主成分轴 优于 第二主成分轴，具有最大可分性。

主坐标分析（Principal Coordinates Analysis，PCoA），即经典多维标度（Classical multidimensional scaling），用于研究数据间的相似性。

主成分分析（Principal components analysis，PCA）是一种统计分析、简化数据集的方法。它利用正交变换来对一系列可能相关的变量的观测值进行线性变换，从而投影为一系列线性不相关变量的值，这些不相关变量称为主成分（Principal Components）。具体地，主成分可以看做一个线性方程，其包含一系列线性系数来指示投影方向（如图）。PCA对原始数据的正则化或预处理敏感（相对缩放）。PCA是最简单的以特征量分析多元统计分布的方法。通常情况下，这种运算可以被看作是揭露数据的内部结构，从而更好的解释数据的变量的方法。

主坐标分析（Principal Coordinates Analysis，PCoA），即经典多维标度（Classical multidimensional scaling），用于研究数据间的相似性。PCoA与PCA都是降低数据维度的方法，但是差异在在于PCA是基于原始矩阵，而PCoA是基于通过原始矩阵计算出的距离矩阵。因此，PCA是尽力保留数据中的变异让点的位置不改动，而PCoA是尽力保证原本的距离关系不发生改变，也就是使得原始数据间点的距离与投影中即结果中各点之间的距离尽可能相关（如图）。

R中有很多包都提供了PCA和PCoA,比如常用的ade4包。本文将基于该包进行PCA和PCoA的分析，数据是自带的deug，该数据提供了104个学生9门课程的成绩（见截图）和综合评定。综合评定有以下几个等级：A+,A,B,B-,C-,D。

让我们通过PCA和PCoA来看一看这样的综合评定是否合理，是否确实依据这9门课把这104个学生合理分配到不同组（每个等级一个组）。

前文已经介绍了PCA是基于原始数据，所以直接进行PCA分析即可。相信大家都比较熟悉散点图的绘制方法，这里不再细讲，PCA分析完毕后我们直接作图展示结果。

整体看起来还不错，就是B-和C-的学生似乎难以区分。

有时候PCA和PCoA的结果差不多，有时候某种方法能够把样本有效分开而另一种可能效果不佳，这些都要看样本数据的特性。

除转录组研究以外，在16S微生物的研究中我们会根据物种丰度的文件对数据进行PCA或者PCoA分析，也是我们所说的β多样性分析。根据PCA或者PCoA的结果看感染组和对照组能否分开，以了解微生物组的总体变化情况。

β多样性分析的概念

Beta多样性指的是样本间多样性。在肠道菌群分析中，Beta多样性是衡量个体间微生物组成相似性的一个指标。通过计算样本间距离可以获得β多样性计算矩阵，后续一般会利用PCoA、进化树聚类等分析对此数值关系进行图形展示。主要基于OTU的群落比较方法，有欧式距离、bray curtis距离、Jaccard 距离，这些方法优势在于算法简单，考虑物种丰度（有无）和均度（相对丰度），但其没有考虑OTUs之间的进化关系，认为OTU之间不存在进化上的联系，每个OTU间的关系平等。另一种算法Unifrac距离法，是根据系统发生树进行比较，并根据16s的序列信息对OTU进行进化树分类， 一般有加权和非加权分析。

QIIME2中重要的Beta多样性指数：

Jaccard距离：群落差异的定性度量，即只考虑种类，不考虑丰度。

Bray-Curtis距离：群落差异的定量度量，较常用。

Unweighted UniFrac距离：包含特征之间的系统发育关系的群落差异定性度量。

Weighted UniFrac距离：包含特征之间的系统发育关系的群落差异定量度量。

解压缩通过qiime2输出的 .qza文件，获得绘图的matrix和pcoa结果文件

将pcoa结果整理成下表，保存为 ***_site.txt

注意没有legend，需要AI加入。

后期需要继续摸索，其实可以加legend的，只是目前自己的技术做不到。。。

PCA思想解析：

https://www.jianshu.com/p/09bae5cbdc53

第1步，DGGE图谱数据的整理

简单明了些，将Quantityone或者类似分辨匹配后得出的DGGE条带的相对百分亮度数据整理到excel表格中。格式如下图，即纵向从上往下为样品编号，条带编号则在第一行从左往右，表格中填入对应条带的相对百分亮度值。比如，条带3在泳道1中的亮度为0..此表是本人某论文的原始数据的一部分。

第2步，DGGE图谱数据的导入。

使用Canoco进行相关分析。Canoco在不少论坛，比如普兰塔等，都有破解版。大家运行的时候只要调整下电脑系统时间到2002年以前即可（最新破解版似乎不需要这么麻烦了）。Canoco不能直接识别excel格式的数据。必须将DGGE图谱的条带数据转化为能识别的.dta后缀的文件才可以进行后续的PCA等相关分析。具体 *** 作如下：在excel表格中，将条带数据连同样品编号和条带编号一起鼠标选择，Ctrl + C 拷贝。然后“开始”-“程序”下找到Canoco的子程序WcanoImp，这个子程序是专门导入数据的。程序会d出如下对话框：

如果你的条带矩阵数据有很多0值，可以将最后的“Save in Condensed Format”勾上，至于为什么，没有太深究。不过不勾上，似乎不影响分析结果。

点击“Save”，d出如下对话框：

这个时候，要自己选择条带数据矩阵的保存路径，并给该文件命名。如，我将此文件保存到“PCA教程”文件夹下，命名为spe.dta，注意这个“spe.dta”是我自己打出来的。点击“保存”，d出如下对话框：

点击“OK”，系统会提示你成功生成文件。这个文件就是我们后续PCA分析所需要的Canoco可识别的条带组成数据。这个spe.dta文件就保存在PCA教程”文件夹下。当然，我们还可以对理化因子数据进行同样的 *** 作。如果对物种生存环境对生物的影响方面的研究感兴趣，可以进行RDA，CCA等分析，这个就需要物种和环境两个数据矩阵了。

第3步，Canoco进行PCA数据分析。

接下来我们就可以使用Canoco对我们的spe.dta文件记性PCA分析了。

在“开始”-“程序”下找到Canoco for Windows 4.5 (当然你我版本可能不同)。

Close Canocod出的小广告。就进入Canoco的分析界面：

点击最左边的符号“”，d出如下对话框：

我们只有条带（即物种species）数据，所以只能选择最上面的“Only species dataavailable”。如果你有物种数据、物种所在环境的理化参数数据或者其他辅助数据，可以选择对应的选项。不同的选项，后续的 *** 作略有不同，可以进行不同的分析。比如，“Species and environmental data available”，我们可以进行计算不同环境因子对物种组成影响的RDA，CCA等分析。感兴趣，可以查阅下文献。

选择“下一步”：

“Browse”，浏览选择上步生成的“spe.dta”数据。如下：

然后在“Canoco solution file name”一栏中，将“spe.dta”数据的保存路径复制过来。举例，我这次就将下图的蓝色选择部分复制，然后拷贝到“Canoco solution file name”下面，自己再后面键盘敲进去将来要生成的solution文件的名称，比如PCA.sol，注意后缀一定是.sol。

这一步比较重要。

选择“下一步”：

选择“Linear”线性的分析模式下的PCA，“下一步”：

默认，下一步：

如果你的不同条带数据差异巨大，比如，最大的100，最小的0.1，你可以选择开平方处理下。默认，下一步：

默认，下一步：

默认，下一步：

选择“完成”，OK，Canoco分析完成了。d出如下对话框：

在“文件名”一栏中，敲进去“PCA.con”。如下图。当然，名字你随便起，后缀最好是.con。

下面就要看分析结果和绘图了。

第4步，Canoco PCA数据分析结果及绘图。

上一步在选择“保存”后，d出如下对话框：

选择“Analyze”对数据进行分析

我们可以看到PCA的大体分析结果，如下图，我们可以看到第一轴解释了31.3%的物种组成变化，第二轴解释了46.9-31.3=15.6%的物种组成。

这个时候该绘图了！选择“CanoDraw”按钮，d出如下对话框，命名PCA.cdw

然后“保存”。然后就可以输出图了，先关掉如下的Canoco小广告。

我们可以看到这样的 *** 作界面：

在“Create”菜单下选择“Scatter plots”—“samples”或者“species”，就可以得到样点或者物种组成的PCA图了。

Sub MergeSort(Array() As Integer, First As Integer, Last As Integer)

Dim mid As Integer = 0

If first<last Then

mid = (first+last)\ 2

MergeSort(Array, first, mid)

MergeSort(Array, mid+1, last)

Merge(Array, first, mid, last)

End If

End Sub

/*

以下示例代码实现了归并 *** 作。array[]是元素序列，其中从索引p开始到q位置，按照升序排列，同时，从(q+1)到r也已经按照升序排列，merge()函数将把这两个已经排序好的子序列合并成一个排序序列。结果放到array中。

*/

/**

* 0 <= p <= q <r, subarray array[p..q] and array[q+1..r] are already sorted.

* the merge() function merges the two sub-arrays into one sorted array.

*/

void Merge(int array[], int p, int q, int r)

{

int i,k

int begin1,end1,begin2,end2

int* temp = (int*)malloc((r-p+1)*sizeof(int))

begin1 = p

end1 = q

begin2 = q+1

end2 = r

k = 0

while((begin1 <= end1)&&( begin2 <= end2))

{

if(array[begin1] <= array[begin2]){

temp[k] = array[begin1]

begin1++

}

else

{

temp[k] = array[begin2]

begin2++

}

k++

}

while(begin1<=end1 || begin2<=end2)

{

if(begin1<=end1)

{

temp[k++] = array[begin1++]

}

if(begin2<=end2)

{

temp[k++] = array[begin2++]

}

}

for (i = 0i <=(r - p)i++)

array[p+i] = temp[i]

free(temp)

}

---