探索数据集

本例使用了DeRisi等人1997年发表的酵母基因表达微阵列研究[1]中的数据。作者使用DNA微阵列研究了酵母中几乎所有基因的时间基因表达酿酒酵母在从发酵到呼吸的代谢转换过程中。在二次转换过程中的七个时间点测量表达水平。完整数据集可从基因表达综合网站下载，https://www.ncbi.nlm.nih.gov/geo/query/acc.cgi?acc=GSE28．

垫子锉yeastdata.mat包含表达式值（比率的log2）CH2DN_平均值和CH1DN_MEAN)，从实验的七个时间步骤，基因的名称，和一组表达水平测量的时间。

负载yeastdata.mat

来了解你可以使用的数据的大小元素个数(基因)显示数据集中包含多少个基因。

元素个数(基因)

ans = 6400

您可以通过索引变量来访问与实验相关的基因名称基因，表示基因名称的单元格数组。例如，中的第15个元素基因是YAL054C。这表示变量的第15行YeastValue包含YAL054C的表达式级别。

基因{15}

ans='YAL054C'

可以使用一个简单的绘图来显示此ORF的表达式配置文件。

情节(次yeastvalues(15日:))包含(‘时间（小时）’）;ylabel ('Log2相对表达水平'）;

您还可以绘制实际的表达式比率，而不是经过log2转换的值。

绘图（时间，2.^yeastvalues（15，：））xlabel(‘时间（小时）’）;ylabel (的相对表达水平的）;

与这个ORF相关的基因，ACS1，似乎在双auxic转移中被强烈上调。你可以通过在同一个图上画多条线来比较这个基因和其他基因的表达。

持有在…上绘图（时间，2.^yeastvalues（16:26，：）'）xlabel(‘时间（小时）’）;ylabel (的相对表达水平的)；头衔(“配置文件表达式级别”）;

筛选基因

通常，基因表达数据集包含与实验期间未显示任何有趣变化的基因相对应的信息。为了更容易找到有趣的基因，您可以将数据集的大小减少到仅包含最重要基因的子集。

如果你仔细查看基因列表，你会看到几个标记为“空的”这些都是数组中的空点，尽管它们可能有与它们相关的数据，但是为了这个例子，您可以认为这些点是噪声。字符串比较函数函数，并使用索引命令从数据集中删除。

emptySpots = strcmp (“空的”,基因);yeastvalues (emptySpots:) = [];基因(emptySpots) = [];元素个数(基因)

ans=6314

还可以看到数据集中的一些地方将表达式级别标记为楠．这表明在特定的时间步长没有收集该点的数据。处理这些缺失值的一种方法是使用特定基因随时间变化的数据的平均值或中值来推断它们。这个例子使用了一种不那么严格的方法，即简单地丢弃任何没有测量一个或多个表达水平的基因的数据。这个函数isnan用于识别缺失数据的基因，索引命令用于删除缺失数据的基因。

NanIndexes=any（isnan（yeastvalues），2）；yeastvalues（NanIndexes，：）=[]；genes（NanIndexes）=[]；numel（genes）

ans = 6276

如果你要绘制所有剩余基因型的表达谱，你会发现大多数基因型是平坦的，与其他基因型没有显著差异。这个平坦的数据显然是有用的，因为它表明与这些基因型相关的基因不会受到二价移位的显著影响；但是，在这个例子中，你是感兴趣的基因在表达上有很大的变化伴随着二价键的转移。你可以使用生物信息学工具箱中的过滤功能™ 删除不提供有关受代谢变化影响的基因的有用信息的各种类型的基因。

你可以使用genevarfilter功能是过滤出随时间变化较小的基因。该函数返回与变量大小相同的逻辑数组（即掩码）基因1对应的行YeastValue方差大于第10个百分位且0对应于阈值以下的值。您可以使用掩码索引到值中并删除筛选的基因。

掩码= genevarfilter (yeastvalues);: yeastvalues = yeastvalues(面具);基因=基因(面具);元素个数(基因)

ans = 5648

这个函数genelowvalfilter删除绝对表达值非常低的基因。请注意，这些过滤函数还可以自动计算过滤后的数据和名称，因此无需使用掩码为原始数据编制索引。

(面具,yeastvalues,基因)= genelowvalfilter (yeastvalues,基因,“absval”log2 (3));元素个数(基因)

ans=822

最后，您可以使用该函数基因过滤器删除信息熵较低的基因，例如数据中第15百分位的信息熵水平。

[mask，yeastvalues，genes]=geneentropyfilter（yeastvalues，genes，“prctile”15);元素个数(基因)

ans=614

聚类分析

现在您已经有了一个可管理的基因列表，您可以使用统计和机器学习工具箱中的一些不同的聚类技术来查找配置文件之间的关系™. 对于层次聚类，函数pdist计算纵断面和纵断面之间的成对距离链接创建层次集群树。

corrDist=pList（年值，“相关系数”）;clusterTree =连杆(corrDist,“平均”）;

的集群函数根据截止距离或最大簇数计算簇最大值选项用于标识16个不同的集群。

集群=集群(clusterTree“maxclust”,16);

这些基因簇中的基因图谱可以用简单的循环和子地块命令。

图形对于c=1:16子图（4,4，c）；图（时间，年值（（聚类==c），：）'）；轴牢固的终止sgtitle (“配置文件的层次聚类”）;

统计和机器学习工具箱也有一个K-means聚类功能。同样，发现了16个聚类，但由于算法不同，这些聚类不一定与通过分层聚类发现的聚类相同。

初始化随机数生成器的状态，以确保这些命令生成的数字与本示例HTML版本中的数字匹配。

rng(“默认”)；[cidx，ctrs]=kmeans（yeastvalues，16，“dist”，“相关系数”，“代表”5.“disp”，“决赛”)；数字对于C = 1:16 subplot(4,4, C);情节(次,yeastvalues ((cidx = = c):) ');轴牢固的终止sgtitle (“K-均值配置文件聚类”）;

重复1次，21次迭代，总距离之和= 23.4699。重复2次，22次迭代，总距离之和= 23.5615。重复3次，10次迭代，总距离之和= 24.823。重复4次，28次迭代，总距离之和= 23.4501。重复5次，19次迭代，总距离之和= 23.5109。距离的最佳总和= 23.4501

您可以只打印质心，而不是打印所有轮廓。

图形对于c=1:16子地块（4,4，c）；地块（时间，CTR（c，：）'）；轴牢固的轴关终止sgtitle (“K-均值配置文件聚类”）;

你可以使用clustergram函数创建表达式级别的热图和层次聚类输出的树状图。

cgObj=聚类图（yeastvalues（：，2:end），“行标签”的基因,“ColumnLabels”,乘以(2:结束));

主成分分析

主成分分析(PCA)是一种有用的技术，可用于降低大数据集的维数，例如来自微阵列的数据集。PCA也可以用于在有噪声的数据中寻找信号。这个函数mapcaplot计算数据集的主要组成部分并创建结果的散点图。您可以从其中一个图中以交互方式选择数据点，这些点将在另一个图中自动高亮显示。这使您可以同时可视化多个维度。

h = mapcaplot (yeastvalues,基因);

请注意，前两个主成分得分的散点图显示有两个不同的区域。这并不意外，因为过滤过程去除了许多低方差或低信息的基因。这些基因可能出现在散点图的中间。

如果你想看看主分量的值主成分分析函数用于计算数据集的主成分。

[pc, zscores, pcvars] = pca(酵母值);

第一个输出,个人计算机的主成分矩阵YeastValue数据。矩阵的第一列是第一个主成分，第二列是第二个主成分，依此类推。第二个输出，zscores，由主成分得分组成，即在主成分空间中酵母值的表示。第三个输出,pcvars，包含主成分方差，它给出了每个主成分占数据方差多少的度量。

很明显，第一个主成分占了模型中方差的大部分。如下所示，您可以计算每个组件的方差的确切百分比。

pcvars. /笔(pcvars) * 100

ans=79.8316 9.5858 4.0781 2.6486 2.1723 0.9747 0.7089

这意味着前两个主成分占了几乎90%的方差累加命令查看方差的累积和。

累计金额（pcvars./sum（pcvars）*100

ans=79.8316 89.4174 93.4955 96.1441 98.3164 99.2911 100.0000

如果您想对主分量的绘制有更多的控制，您可以使用分散作用

图形散点（zscores（：，1），zscores（：，2））；xlabel(“第一主成分”）;ylabel (第二主成分的)；头衔(“主成分散点图”）;

创建散点图的另一种方法是使用函数gscatter来自统计和机器学习工具箱。gscatter创建分组散点图，其中每组的点具有不同的颜色或标记。您可以使用聚类数据，或任何其他群集功能，以对点进行分组。

图pcclusters=clusterdata（zscores（：，1:2），“maxclust”8“链接”，“av”)；gscatter（zscores（：，1），zscores（：，2），pcclusters）xlabel(“第一主成分”）;ylabel (第二主成分的)；头衔(“带有彩色簇的主成分散点图”）;

自组织映射

如果你有深度学习工具箱™, 可以使用自组织映射（SOM）对数据进行聚类。

%检查是否安装了深度学习工具箱如果~exist（哪个(“自组织地图”),“文件”)disp(“本例中的自组织地图部分需要使用深度学习工具箱。”）回来终止

的自组织映射函数创建一个新的SOM网络对象。此示例将使用前两个主组件生成一个SOM。

P = zscores(: 1:2)”;Net = selforgmap([4 4]);

使用默认参数训练网络。

网=火车(净,P);

使用plotsom在数据的散点图上显示网络。注意，SOM算法使用随机的起始点，因此每次运行的结果都会不同。

曲线图（P（1，：），P（2，：），“.g”，“markersize”, 20)在…上net.layers plotsom (net.iw {1}, {1} .distances)关

可以通过查找数据集中每个点最近的节点，使用SOM分配簇。

距离=dist（P'，net.IW{1}'）；[d，cndx]=min（距离，[]，2）；%cndx包含群集索引图gscatter（P（1，：），P（2，：），cndx；图例关持有在…上net.layers plotsom (net.iw {1}, {1} .distances);持有关

关闭所有数字。

关闭(“全部”)；删除（cgObj）；删除（h）；

工具书类

[1] DeRisi，J.L.，Iyer，V.R.和Brown，P.O.，“探索基因组尺度上基因表达的代谢和遗传控制”，《科学》，278（5338）：680-61997。