多维尺度变换(multidimensional scaling, MDS)

　　多维尺度变换(multidimensional scaling, MDS)是在低维空间去展示高维多元数据的一种可视化方法。该方法看起来类似于利用主成分得分作图，或者对两个线性判别量的得分作图。与上述方法不同的是，多维尺度变换的基本目标是将原始数据“拟合”到一个低维坐标系中，使得由降维所引起的任何变形最小。多维尺度变换涉及到的问题可以描述为：当n个项目中各对项目之间的相似性(或距离)给定时，求这些项目在低维空间的表示，并使得项目间的接近程度与原先的相似性(或距离) “大体匹配”。
　　根据样本是否可以计量，多维尺度变换可分为度量型多维尺度变换(metric MDS)和非度量型多维尺度变换(non-metric MDS)。在可计量的多维尺度变换中，依据度量的距离标准是否是欧氏距离，分为经典多维尺度变换(classical MDS)和非经典多维尺度变换(non-classical MDS)。

度量型多维尺度变换(metric MDS)

经典多维尺度变换(classical MDS)

　　多维尺度变换就是从经典多维尺度变换这个问题开始的。这个问题起源于当我们仅能获取到物体之间的距离的时候，如何由此去重构它的欧几里得坐标(Torgerson, W.S. 1958. Theory & Methods of Scaling. New York: Wiley)。譬如，在以前，对一个国家的许多城市而言，并不能确定各个城市的经纬度信息，却知道所有城市之间的两两距离。任意两个城市的距离可表示为

$$\mathbf{D}=\left[ \begin{matrix} {{d}_{11}} & {{d}_{12}} & \cdots & {{d}_{1n}} \\ {{d}_{21}} & {{d}_{21}} & \cdots & {{d}_{2n}} \\ \vdots & \vdots & \ddots & \vdots \\ {{d}_{n1}} & {{d}_{n1}} & \cdots & {{d}_{nn}} \\ \end{matrix} \right] (1)$$
经典多维尺度变换方法便是利用这些距离数据计算出各个城市在二维空间上的地理分布，即重构出城市的位置坐标 $\mathbf{\tilde{x}}$(原始真实的位置坐标记为 x)，使得第 i 个和第 j 个城市的距离 $\left\| {{{\mathbf{\tilde{x}}}}_{i}}-{{{\mathbf{\tilde{x}}}}_{j}} \right\|\approx {{d}_{ij}}$。这样多维尺度变换问题可以看成一个优化问题，优化的目标函数为
$$\underset{{{{\mathbf{\tilde{x}}}}_{1}},\cdots ,{{{\mathbf{\tilde{x}}}}_{n}}}{\mathop{\min }}\,{{\left( {{\sum\limits_{i<j}{\left( \left\| {{{\mathbf{\tilde{x}}}}_{i}}-{{{\mathbf{\tilde{x}}}}_{j}} \right\|-{{d}_{ij}} \right)}}^{2}} \right)}^{{1}/{2}\;}} (2)$$
这里 $\left\| . \right\|$是向量的范数。在经典多维尺度变换中，该范数是欧氏距离，但广义的讲，这个范数可以是任意函数。

Matlab实现

%原始数据
cities = ...
{'Atl','Chi','Den','Hou','LA','Mia','NYC','SF','Sea','WDC'};
D = [    0  587 1212  701 1936  604  748 2139 2182   543;
       587    0  920  940 1745 1188  713 1858 1737   597;
      1212  920    0  879  831 1726 1631  949 1021  1494;
       701  940  879    0 1374  968 1420 1645 1891  1220;
      1936 1745  831 1374    0 2339 2451  347  959  2300;
       604 1188 1726  968 2339    0 1092 2594 2734   923;
       748  713 1631 1420 2451 1092    0 2571 2408   205;
      2139 1858  949 1645  347 2594 2571    0  678  2442;
      2182 1737 1021 1891  959 2734 2408  678    0  2329;
       543  597 1494 1220 2300  923  205 2442 2329     0];
[Y,eigvals] = cmdscale(D);
%绘图
plot(Y(:,1),Y(:,2),'.')
text(Y(:,1)+25,Y(:,2),cities)
xlabel('Miles')
ylabel('Miles')

　　经典多维尺度变换是基于欧氏距离假设的，所以不能直接应用在一般相似性矩阵情况上(wiki: Multidimensional scaling)。虽然Matlab官方文档中说D矩阵不必要是欧氏距离矩阵(Euclidean distance matrix)，但应用非欧氏矩阵或者更一般的非相似性矩阵（dissimilarity matrix）时，特征值向量eigvals中会出现负值。需要注意的是，当eigvals中负值较大的情况下计算出来的低维几何表示就没有意义了。

　　在前面的问题中，各个城市的坐标信息是未知的。然而即使样本的坐标信息是已知的，但样本若是处于高维空间中，依然需要使用多维尺度变换，将其转化到低维可视化空间上展示。
　　给定n 个q 维(q 值较大)的样本数据：

$$\mathbf{X}=\left[ \begin{matrix} {{\mathbf{x}}_{1}} \\ {{\mathbf{x}}_{2}} \\ \vdots \\ {{\mathbf{x}}_{n}} \\ \end{matrix} \right]=\left[ \begin{matrix} {{x}_{11}} & {{x}_{12}} & \cdots & {{x}_{1q}} \\ {{x}_{21}} & {{x}_{21}} & \cdots & {{x}_{2q}} \\ \vdots & \vdots & \ddots & \vdots \\ {{x}_{n1}} & {{x}_{n1}} & \cdots & {{x}_{nq}} \\ \end{matrix} \right] (3)$$
那么， n个样本间的欧氏距离 D可以表示为公式(1)的形式，其中 ${{d}_{ij}}=\left\| {{\mathbf{x}}_{i}}-{{\mathbf{x}}_{j}} \right\|=\sqrt{\sum\limits_{k=1}^{q}{{{\left( {{x}_{ik}}-{{x}_{jk}} \right)}^{2}}}}$。
　　在经典多维尺度变换问题中，假定矩阵 D是欧氏距离(如城市间的距离)，在计算的过程中，不需要采用数值计算的方法去优化目标函数(2)。可以利用主成分分析的思想，可以简化目标函数，得到准确的解析解。或者说可以根据欧氏距离矩阵准确的重构原坐标的位置矩阵。这是基于一个事实，原坐标矩阵可以由中心化后的内积矩阵进行特征值分解得到(参考 wiki, Multidimensional scaling)。经典解析解的详细推导参考《多元统计分析》 1 何晓群编著P311-314。

经典多维尺度变换解析算法步骤

1. 由距离矩阵$\mathbf{D}=\left[ {{d}_{ij}} \right]\in {{\Re }^{n\times n}}$，构造矩阵$\mathbf{A}=\left[ {{a}_{ij}} \right]=\left[ -\frac{1}{2}d_{ij}^{2} \right]$。
2. 计算内积阵$\mathbf{B}=\left[ {{b}_{ij}} \right]=\left[ {{a}_{ij}}-{{{\bar{a}}}_{i.}}-{{{\bar{a}}}_{.}}_{j}+{{{\bar{a}}}_{..}} \right]$。
3. 求B的特征值${{\lambda }_{1}}\ge {{\lambda }_{2}}\ge \cdots \ge {{\lambda }_{n}}$和特征向量。这里的距离矩阵时欧氏矩阵，所以特征值都为非负，若出现负的特征值，说明矩阵D一定不是欧氏型。令
   $${{s}_{k}}=\frac{\sum\limits_{i=1}^{k}{{{\lambda }_{i}}}}{\sum\limits_{i=1}^{n}{\left| {{\lambda }_{i}} \right|}} (4)$$
   其中，s_k 相当于主成分分析中的累计贡献率，当然我们希望k不要取太大，而使累计贡献率比较大。在前面提出的预测城市坐标的是示例中，由于x_i 是二维坐标数据，所以k值取2时，累计贡献率就几乎等于100%。但对于一般的数据来说，可能保留二维数据累计贡献率并不能达到100%。需要根据情况作出判断。
4. 确定一个k值后，重构$\mathbf{\tilde{X}}={{\mathbf{E}}_{k}}\mathbf{\Lambda }_{k}^{{1}/{2}\;}$, 其中${{\mathbf{E}}_{k}}$是矩阵B保留下来的前k个特征向量组成的矩阵，${{\mathbf{\Lambda }}_{k}}$是k个特征值组成的对角矩阵。
   从上面的推导中可以看出，经典多维尺度变换和主成分分析本质上是一样的，区别在于多维尺度变换是以样本为分析对象，主成分分析是以变量为分析对象。而且可以证明经过经典多维尺度变换后的$\mathbf{\tilde{X}}$的k维主坐标正好是X中心化后的利用主成分分析得到的前k个主成分的值（参考<多元统计分析>何晓群编著P317）。

非经典多维尺度变换(non-classical MDS)

　　当矩阵D不是欧氏距离矩阵，或者只是一般的相似性矩阵时，不满足经典多维度变化假设，此时不存在解析解，需要采用优化算法求解。为了能够直观的衡量重构出的低维空间的几何表示偏离原始数据的程度，目标函数(2)可改为百分比形式的偏离度度量

$$Stress({{\mathbf{\tilde{x}}}_{1}},{{\mathbf{\tilde{x}}}_{2}},\cdots ,{{\mathbf{\tilde{x}}}_{n}})={{\left( {{{\sum\limits_{i<j}{\left( \left\| {{{\mathbf{\tilde{x}}}}_{i}}-{{{\mathbf{\tilde{x}}}}_{j}} \right\|-{{d}_{ij}} \right)}}^{2}}}/{\sum\limits_{i<j}{d_{ij}^{2}}}\; \right)}^{{1}/{2}\;}} (5)$$
除了这种偏离度度量，还有另一种偏离度的度量，记为Squared stress
$$SStress({{\mathbf{\tilde{x}}}_{1}},{{\mathbf{\tilde{x}}}_{2}},\cdots ,{{\mathbf{\tilde{x}}}_{n}})={{\left( {{{\sum\limits_{i<j}{\left( {{\left\| {{{\mathbf{\tilde{x}}}}_{i}}-{{{\mathbf{\tilde{x}}}}_{j}} \right\|}^{2}}-d_{ij}^{2} \right)}}^{2}}}/{\sum\limits_{i<j}{d_{ij}^{4}}}\; \right)}^{{1}/{2}\;}} (6)$$
关于该目标函数有很多数值优化方法，参考 MDS (multidimensional scaling)多维尺度分析。

Matlab实现

[Y, stress, disparities] = mdscale(D, p, 'criterion', 'metricstress')

　　与经典多维尺度变换cmdscale不同的是，这里需要指定希望的维度和重构的方法。多维尺度变换目的是根据距离矩阵寻找高维样本数据X(q 维)在低维空间中的表示，也就是说多维尺度变换试图找到高维空间${{\Re }^{q}}$的一个子空间。所以优化目标函数(5)求解低维空间的坐标之前，需给定子空间的维度p。如果给定p 维子空间中，经过优化后得到的stress依然较大(一般2.5%称为非常好)，需要增加子空间的维度重新进行优化计算，直至得到满意的stress值(更多的最优维度确定方法请参考《实用多元统计分析》2陆璇 译 P553)。当偏离度度量选择公式(5)时，重构的方法’criterion’选择’metricstress’，而当偏离度度量选择公式(6)，’criterion’选择’metricsstress’。

非度量型多维尺度变换(non-metric MDS)

　　上面的度量多维尺度变换根据差异(或距离)矩阵创建了一系列点，使得这些点的距离矩阵与原始距离矩阵近似，这个要求有时是比较严格的，因为有些差异是不易量化的。非度量多维尺度变换就是为了放松这一要求。与前面近似距离矩阵不同，非度量多维尺度变换是近似一个非线性但单调的变换，只用距离矩阵D的顺序信息而不用其数值。它只是试图保存差异性的次序。因此，在不同的尺度可能有距离的压缩或展开。
类似于度量型多维尺度变换，非度量型尺度变换也有两种偏离度的度量标准

$$Stress({{\mathbf{\tilde{x}}}_{1}},{{\mathbf{\tilde{x}}}_{2}},\cdots ,{{\mathbf{\tilde{x}}}_{n}})={{\left( {{{\sum\limits_{i<j}{\left( \left\| {{{\mathbf{\tilde{x}}}}_{i}}-{{{\mathbf{\tilde{x}}}}_{j}} \right\|-{{{\hat{d}}}_{ij}} \right)}}^{2}}}/{\sum\limits_{i<j}{{{\left\| {{{\mathbf{\tilde{x}}}}_{i}}-{{{\mathbf{\tilde{x}}}}_{j}} \right\|}^{2}}}}\; \right)}^{{1}/{2}\;}} (7)$$
Squared stress
$$SStress({{\mathbf{\tilde{x}}}_{1}},{{\mathbf{\tilde{x}}}_{2}},\cdots ,{{\mathbf{\tilde{x}}}_{n}})={{\left( {{{\sum\limits_{i<j}{\left( {{\left\| {{{\mathbf{\tilde{x}}}}_{i}}-{{{\mathbf{\tilde{x}}}}_{j}} \right\|}^{2}}-\hat{d}_{ij}^{2} \right)}}^{2}}}/{\sum\limits_{i<j}{{{\left\| {{{\mathbf{\tilde{x}}}}_{i}}-{{{\mathbf{\tilde{x}}}}_{j}} \right\|}^{4}}}}\; \right)}^{{1}/{2}\;}} (8)$$

需要注意的是，这里的${{\hat{d}}_{ij}}$与度量型多维尺度变换方法中的${{d}_{ij}}$是不一样的。在度量型尺度变换方法中，${{d}_{ij}}$表示距离矩阵D对应的值。而非度量型尺度变换不考虑距离矩阵的元素的值的大小，只关心其排序。将距离矩阵D所有上对角元素按照严格上升的顺序排成一列

$${{d}_{{{i}_{1}}{{j}_{1}}}}\le {{d}_{{{i}_{2}}{{j}_{2}}}}\le \cdots \le {{d}_{{{i}_{m}}{{j}_{m}}}},\ m={n\left( n-1 \right)}/{2}\; (9)$$
这里 ${{d}_{{{i}_{1}}{{j}_{1}}}}$是所有样本间距离最小值。 ${{\hat{d}}_{ij}}$满足公式(9)，即它们与 ${{d}_{ij}}$有一样的排序的单调关系，其具体值与 ${{d}_{ij}}$无关。因此 ${{\hat{d}}_{ij}}$并不是真正的距离，仅仅是用来评价与真实距离 $\left\| {{{\mathbf{\tilde{x}}}}_{i}}-{{{\mathbf{\tilde{x}}}}_{j}} \right\|$的非单调性的一些参考值。
　　非度量多维尺度变换的算法确定 ${{\mathbf{\tilde{x}}}_{i}}$和 ${{\hat{d}}_{ij}}$是一个双重优化过程，首先根据预先给定的点寻找最优的单调变换，然后再反过来优化点的位置结构，具体步骤为(参考 wiki: Multidimensional scaling)：

1. 找到一组随机点结构${{\mathbf{\tilde{x}}}_{i}}$
2. 计算这些点之间的距离，得到低维空间重构点之间的距离矩阵${{\tilde{d}}_{ij}}$
3. 根据原距离矩阵${{d}_{ij}}$的顺序，对新重构的距离矩阵${{\tilde{d}}_{ij}}$进行保序回归(isotonic regression)，得到${{\hat{d}}_{ij}}$
4. 将得到${{\hat{d}}_{ij}}$带入公式(7)或者(8)，通过优化计算，找到新的点的位置${{\mathbf{\tilde{x}}}_{i}}$
5. 比较stress值是否满足一些指标，如果不满足将回到步骤2

Matlab实现

[Y, stress, disparities] = mdscale(D, p, 'criterion', 'stress')

两种多维尺度变换的比较

　　对于Metric MDS，这个方法以样本间相似度作为实际输入，需要样本是等距（interval）比例（ratio）尺度，优点是精确，可以根据多个准则评估样本间差异，缺点是计算成本高，耗时。对于很多应用问题，样本不费可计量，需要使用NonMetric MDS，这种方法接受样本的顺序尺度作为输入，并以此自动计算相似值。样本尺度要求是顺序的（ordinal），较简便，直观，从非计量的样本导出计量的分析结果，应用范围更广，但没法知道评估准则，效果较差。参考MDS（multidimensional scaling）多维尺度分析

---

参考资料

---

1. 《多元统计分析》何晓群 编著 中国人民大学出版社 第三版 ↩
2. 《实用多元统计分析》陆璇 译 清华大学出版社 第六版 ↩