MDS

library(biplotEZ)

This vignette deals with biplots based on Multidimensional scaling (MDS).

1 What is MDS

In general, multidimensional scaling deals with constructing a low dimensional map of \(n\) samples such that the Euclidean distances between the samples match a given set of dissimilarities \(\mathbf{\Delta}:n \times n\) as closely as possible. Here the focus is on Principal Coordinate Analysis (PCO), an approach based on the singular value decomposition of a matrix. However, the regression biplot provides a general structure for fitting any 2D map of samples with biplot axes.

2 Regression biplot

The function regress accepts as arguments an object of class biplot. The call to the function biplot should contain the data set that will be used to construct the biplot axes. The second argument Z contains the coordinates of the samples in the MDS map. By default linear regression axes will be fitted to the plot. If \(\mathbf{X}\) denote the data matrix, the directions of the biplot axes are found by solving the regression equation

\[ \mathbf{X} = \mathbf{ZH}' + \mathbf{E}. \] The matrix \(\mathbf{H}'=(\mathbf{Z'Z})^{-1}\mathbf{Z'X}\) and calibration of the axes proceed analogous to equation (2) in the biplotEZ vignette. The coordinates of the marker \(/mu\) on biplot axis \(j\) is given as follows

\[ p_{\mu} = \frac{\mu}{\mathbf{h}_{(j)}'\mathbf{h}_{(j)}}\mathbf{h}_{(j)}. \]

biplot(rock) |> 
  regress(Z = MASS::sammon(dist(scale(rock), method="manhattan"))$points) |> 
  plot()
#> Initial stress        : 0.02554
#> stress after  10 iters: 0.01094, magic = 0.500
#> stress after  20 iters: 0.01080, magic = 0.500
#> stress after  30 iters: 0.01078, magic = 0.500

More flexibility is provided by approximating the biplot axes using B-splines. A calibrated trajectory represented in a matrix \(\mathbf{H}:m \times 2\) snakes through the samples points \(\mathbf{Z}\) such that the marker on the trajectory which is closest to a particular sample gives the predicted value of that sample for the particular variable. To retain smoothness of the trajectory \(\mathbf{H}\), B-splines are used.

biplot(rock) |> 
  regress(Z = MASS::sammon(dist(scale(rock), method="manhattan"))$points,
          axes = "splines") |> 
  plot()
#> Initial stress        : 0.02554
#> stress after  10 iters: 0.01094, magic = 0.500
#> stress after  20 iters: 0.01080, magic = 0.500
#> stress after  30 iters: 0.01078, magic = 0.500

#> Calculating spline axis for variable 1 
#>     area     peri     shape   perm
#> 1   4990 2791.900 0.0903296    6.3
#> 2   7002 3892.600 0.1486220    6.3
#> 3   7558 3930.660 0.1833120    6.3
#> 4   7352 3869.320 0.1170630    6.3
#> 5   7943 3948.540 0.1224170   17.1
#> 6   7979 4010.150 0.1670450   17.1
#> 7   9333 4345.750 0.1896510   17.1
#> 8   8209 4344.750 0.1641270   17.1
#> 9   8393 3682.040 0.2036540  119.0
#> 10  6425 3098.650 0.1623940  119.0
#> 11  9364 4480.050 0.1509440  119.0
#> 12  8624 3986.240 0.1481410  119.0
#> 13 10651 4036.540 0.2285950   82.4
#> 14  8868 3518.040 0.2316230   82.4
#> 15  9417 3999.370 0.1725670   82.4
#> 16  8874 3629.070 0.1534810   82.4
#> 17 10962 4608.660 0.2043140   58.6
#> 18 10743 4787.620 0.2627270   58.6
#> 19 11878 4864.220 0.2000710   58.6
#> 20  9867 4479.410 0.1448100   58.6
#> 21  7838 3428.740 0.1138520  142.0
#> 22 11876 4353.140 0.2910290  142.0
#> 23 12212 4697.650 0.2400770  142.0
#> 24  8233 3518.440 0.1618650  142.0
#> 25  6360 1977.390 0.2808870  740.0
#> 26  4193 1379.350 0.1794550  740.0
#> 27  7416 1916.240 0.1918020  740.0
#> 28  5246 1585.420 0.1330830  740.0
#> 29  6509 1851.210 0.2252140  890.0
#> 30  4895 1239.660 0.3412730  890.0
#> 31  6775 1728.140 0.3116460  890.0
#> 32  7894 1461.060 0.2760160  890.0
#> 33  5980 1426.760 0.1976530  950.0
#> 34  5318  990.388 0.3266350  950.0
#> 35  7392 1350.760 0.1541920  950.0
#> 36  7894 1461.060 0.2760160  950.0
#> 37  3469 1376.700 0.1769690  100.0
#> 38  1468  476.322 0.4387120  100.0
#> 39  3524 1189.460 0.1635860  100.0
#> 40  5267 1644.960 0.2538320  100.0
#> 41  5048  941.543 0.3286410 1300.0
#> 42  1016  308.642 0.2300810 1300.0
#> 43  5605 1145.690 0.4641250 1300.0
#> 44  8793 2280.490 0.4204770 1300.0
#> 45  3475 1174.110 0.2007440  580.0
#> 46  1651  597.808 0.2626510  580.0
#> 47  5514 1455.880 0.1824530  580.0
#> 48  9718 1485.580 0.2004470  580.0
#> attr(,"scaled:center")
#>          area          peri         shape          perm 
#> -7187.7291667 -2682.2119375    -0.2181104  -415.4500000 
#> Calculating spline axis for variable 2 
#>     area     peri     shape   perm
#> 1   4990 2791.900 0.0903296    6.3
#> 2   7002 3892.600 0.1486220    6.3
#> 3   7558 3930.660 0.1833120    6.3
#> 4   7352 3869.320 0.1170630    6.3
#> 5   7943 3948.540 0.1224170   17.1
#> 6   7979 4010.150 0.1670450   17.1
#> 7   9333 4345.750 0.1896510   17.1
#> 8   8209 4344.750 0.1641270   17.1
#> 9   8393 3682.040 0.2036540  119.0
#> 10  6425 3098.650 0.1623940  119.0
#> 11  9364 4480.050 0.1509440  119.0
#> 12  8624 3986.240 0.1481410  119.0
#> 13 10651 4036.540 0.2285950   82.4
#> 14  8868 3518.040 0.2316230   82.4
#> 15  9417 3999.370 0.1725670   82.4
#> 16  8874 3629.070 0.1534810   82.4
#> 17 10962 4608.660 0.2043140   58.6
#> 18 10743 4787.620 0.2627270   58.6
#> 19 11878 4864.220 0.2000710   58.6
#> 20  9867 4479.410 0.1448100   58.6
#> 21  7838 3428.740 0.1138520  142.0
#> 22 11876 4353.140 0.2910290  142.0
#> 23 12212 4697.650 0.2400770  142.0
#> 24  8233 3518.440 0.1618650  142.0
#> 25  6360 1977.390 0.2808870  740.0
#> 26  4193 1379.350 0.1794550  740.0
#> 27  7416 1916.240 0.1918020  740.0
#> 28  5246 1585.420 0.1330830  740.0
#> 29  6509 1851.210 0.2252140  890.0
#> 30  4895 1239.660 0.3412730  890.0
#> 31  6775 1728.140 0.3116460  890.0
#> 32  7894 1461.060 0.2760160  890.0
#> 33  5980 1426.760 0.1976530  950.0
#> 34  5318  990.388 0.3266350  950.0
#> 35  7392 1350.760 0.1541920  950.0
#> 36  7894 1461.060 0.2760160  950.0
#> 37  3469 1376.700 0.1769690  100.0
#> 38  1468  476.322 0.4387120  100.0
#> 39  3524 1189.460 0.1635860  100.0
#> 40  5267 1644.960 0.2538320  100.0
#> 41  5048  941.543 0.3286410 1300.0
#> 42  1016  308.642 0.2300810 1300.0
#> 43  5605 1145.690 0.4641250 1300.0
#> 44  8793 2280.490 0.4204770 1300.0
#> 45  3475 1174.110 0.2007440  580.0
#> 46  1651  597.808 0.2626510  580.0
#> 47  5514 1455.880 0.1824530  580.0
#> 48  9718 1485.580 0.2004470  580.0
#> attr(,"scaled:center")
#>          area          peri         shape          perm 
#> -7187.7291667 -2682.2119375    -0.2181104  -415.4500000 
#> Calculating spline axis for variable 3 
#>     area     peri     shape   perm
#> 1   4990 2791.900 0.0903296    6.3
#> 2   7002 3892.600 0.1486220    6.3
#> 3   7558 3930.660 0.1833120    6.3
#> 4   7352 3869.320 0.1170630    6.3
#> 5   7943 3948.540 0.1224170   17.1
#> 6   7979 4010.150 0.1670450   17.1
#> 7   9333 4345.750 0.1896510   17.1
#> 8   8209 4344.750 0.1641270   17.1
#> 9   8393 3682.040 0.2036540  119.0
#> 10  6425 3098.650 0.1623940  119.0
#> 11  9364 4480.050 0.1509440  119.0
#> 12  8624 3986.240 0.1481410  119.0
#> 13 10651 4036.540 0.2285950   82.4
#> 14  8868 3518.040 0.2316230   82.4
#> 15  9417 3999.370 0.1725670   82.4
#> 16  8874 3629.070 0.1534810   82.4
#> 17 10962 4608.660 0.2043140   58.6
#> 18 10743 4787.620 0.2627270   58.6
#> 19 11878 4864.220 0.2000710   58.6
#> 20  9867 4479.410 0.1448100   58.6
#> 21  7838 3428.740 0.1138520  142.0
#> 22 11876 4353.140 0.2910290  142.0
#> 23 12212 4697.650 0.2400770  142.0
#> 24  8233 3518.440 0.1618650  142.0
#> 25  6360 1977.390 0.2808870  740.0
#> 26  4193 1379.350 0.1794550  740.0
#> 27  7416 1916.240 0.1918020  740.0
#> 28  5246 1585.420 0.1330830  740.0
#> 29  6509 1851.210 0.2252140  890.0
#> 30  4895 1239.660 0.3412730  890.0
#> 31  6775 1728.140 0.3116460  890.0
#> 32  7894 1461.060 0.2760160  890.0
#> 33  5980 1426.760 0.1976530  950.0
#> 34  5318  990.388 0.3266350  950.0
#> 35  7392 1350.760 0.1541920  950.0
#> 36  7894 1461.060 0.2760160  950.0
#> 37  3469 1376.700 0.1769690  100.0
#> 38  1468  476.322 0.4387120  100.0
#> 39  3524 1189.460 0.1635860  100.0
#> 40  5267 1644.960 0.2538320  100.0
#> 41  5048  941.543 0.3286410 1300.0
#> 42  1016  308.642 0.2300810 1300.0
#> 43  5605 1145.690 0.4641250 1300.0
#> 44  8793 2280.490 0.4204770 1300.0
#> 45  3475 1174.110 0.2007440  580.0
#> 46  1651  597.808 0.2626510  580.0
#> 47  5514 1455.880 0.1824530  580.0
#> 48  9718 1485.580 0.2004470  580.0
#> attr(,"scaled:center")
#>          area          peri         shape          perm 
#> -7187.7291667 -2682.2119375    -0.2181104  -415.4500000 
#> Calculating spline axis for variable 4 
#>     area     peri     shape   perm
#> 1   4990 2791.900 0.0903296    6.3
#> 2   7002 3892.600 0.1486220    6.3
#> 3   7558 3930.660 0.1833120    6.3
#> 4   7352 3869.320 0.1170630    6.3
#> 5   7943 3948.540 0.1224170   17.1
#> 6   7979 4010.150 0.1670450   17.1
#> 7   9333 4345.750 0.1896510   17.1
#> 8   8209 4344.750 0.1641270   17.1
#> 9   8393 3682.040 0.2036540  119.0
#> 10  6425 3098.650 0.1623940  119.0
#> 11  9364 4480.050 0.1509440  119.0
#> 12  8624 3986.240 0.1481410  119.0
#> 13 10651 4036.540 0.2285950   82.4
#> 14  8868 3518.040 0.2316230   82.4
#> 15  9417 3999.370 0.1725670   82.4
#> 16  8874 3629.070 0.1534810   82.4
#> 17 10962 4608.660 0.2043140   58.6
#> 18 10743 4787.620 0.2627270   58.6
#> 19 11878 4864.220 0.2000710   58.6
#> 20  9867 4479.410 0.1448100   58.6
#> 21  7838 3428.740 0.1138520  142.0
#> 22 11876 4353.140 0.2910290  142.0
#> 23 12212 4697.650 0.2400770  142.0
#> 24  8233 3518.440 0.1618650  142.0
#> 25  6360 1977.390 0.2808870  740.0
#> 26  4193 1379.350 0.1794550  740.0
#> 27  7416 1916.240 0.1918020  740.0
#> 28  5246 1585.420 0.1330830  740.0
#> 29  6509 1851.210 0.2252140  890.0
#> 30  4895 1239.660 0.3412730  890.0
#> 31  6775 1728.140 0.3116460  890.0
#> 32  7894 1461.060 0.2760160  890.0
#> 33  5980 1426.760 0.1976530  950.0
#> 34  5318  990.388 0.3266350  950.0
#> 35  7392 1350.760 0.1541920  950.0
#> 36  7894 1461.060 0.2760160  950.0
#> 37  3469 1376.700 0.1769690  100.0
#> 38  1468  476.322 0.4387120  100.0
#> 39  3524 1189.460 0.1635860  100.0
#> 40  5267 1644.960 0.2538320  100.0
#> 41  5048  941.543 0.3286410 1300.0
#> 42  1016  308.642 0.2300810 1300.0
#> 43  5605 1145.690 0.4641250 1300.0
#> 44  8793 2280.490 0.4204770 1300.0
#> 45  3475 1174.110 0.2007440  580.0
#> 46  1651  597.808 0.2626510  580.0
#> 47  5514 1455.880 0.1824530  580.0
#> 48  9718 1485.580 0.2004470  580.0
#> attr(,"scaled:center")
#>          area          peri         shape          perm 
#> -7187.7291667 -2682.2119375    -0.2181104  -415.4500000

3 Principal Coordinate Analysis (PCO) biplots

Let \(\mathbf{\tilde\Delta}\) be the \(n \times n\) matrix containing the values \(-\frac{1}{2}\delta_{ij}^2\) where \(\delta_{ij}\) represent the dissimilarity between objects \(i\) and \(j\). If it is possible to find a set of coordinates \(\mathbf{Y}:n\times m\) where typically \(m = n-1\) such that the Euclidean distances between the rows of \(\mathbf{Y}\) exactly match the dissimilarities \(\delta_{ij}\), the dissimilarity is known as a Euclidean embeddable distance metric. J. C. Gower (1982) shows that if the distances \(\delta_{ij}\) are Euclidean embeddable, then

\[ (\mathbf{I}-\mathbf{1s}')\mathbf{\tilde\Delta} (\mathbf{I}-\mathbf{s1}') \]

is positive semi-definite.

The Euclidean representation of the samples is obtained as \(\mathbf{Y=V\Lambda}^{\frac{1}{2}}\) where \((\mathbf{I}-\mathbf{1s}')\mathbf{\tilde\Delta} (\mathbf{I}-\mathbf{s1}') = \mathbf{V \Lambda V'}\). Since the coordinates in \(\mathbf{Y}\) is already referred to their principal axes with \(\mathbf{Y'Y=\Lambda}\), the representation of the samples in a 2D biplot is obtained from the first two columns of \(\mathbf{Y}\).

In addition to the exact biplot axis representations discussed in J. C. Gower, Lubbe, and Roux (2011) approximate axes can be obtained. Linear axes are fitted with the regression method. The variables in the data matrix \(\mathbf{X}:n \times p\) can be represented as biplot axes in the PCO biplot with the sample points \(\mathbf{Z=V\Lambda}^{\frac{1}{2}}\mathbf{J}_2\) according to the regression method discussed in section 2 above.

biplot(rock, scale = TRUE) |> PCO() |> plot()

Using B-splines instead of linear regression provides the user with more flexibility. This is achieved by setting the argument axes = "splines".

biplot(rock, scale = TRUE) |> PCO(axes = "splines") |> plot()

#> Calculating spline axis for variable 1 
#>     area     peri     shape   perm
#> 1   4990 2791.900 0.0903296    6.3
#> 2   7002 3892.600 0.1486220    6.3
#> 3   7558 3930.660 0.1833120    6.3
#> 4   7352 3869.320 0.1170630    6.3
#> 5   7943 3948.540 0.1224170   17.1
#> 6   7979 4010.150 0.1670450   17.1
#> 7   9333 4345.750 0.1896510   17.1
#> 8   8209 4344.750 0.1641270   17.1
#> 9   8393 3682.040 0.2036540  119.0
#> 10  6425 3098.650 0.1623940  119.0
#> 11  9364 4480.050 0.1509440  119.0
#> 12  8624 3986.240 0.1481410  119.0
#> 13 10651 4036.540 0.2285950   82.4
#> 14  8868 3518.040 0.2316230   82.4
#> 15  9417 3999.370 0.1725670   82.4
#> 16  8874 3629.070 0.1534810   82.4
#> 17 10962 4608.660 0.2043140   58.6
#> 18 10743 4787.620 0.2627270   58.6
#> 19 11878 4864.220 0.2000710   58.6
#> 20  9867 4479.410 0.1448100   58.6
#> 21  7838 3428.740 0.1138520  142.0
#> 22 11876 4353.140 0.2910290  142.0
#> 23 12212 4697.650 0.2400770  142.0
#> 24  8233 3518.440 0.1618650  142.0
#> 25  6360 1977.390 0.2808870  740.0
#> 26  4193 1379.350 0.1794550  740.0
#> 27  7416 1916.240 0.1918020  740.0
#> 28  5246 1585.420 0.1330830  740.0
#> 29  6509 1851.210 0.2252140  890.0
#> 30  4895 1239.660 0.3412730  890.0
#> 31  6775 1728.140 0.3116460  890.0
#> 32  7894 1461.060 0.2760160  890.0
#> 33  5980 1426.760 0.1976530  950.0
#> 34  5318  990.388 0.3266350  950.0
#> 35  7392 1350.760 0.1541920  950.0
#> 36  7894 1461.060 0.2760160  950.0
#> 37  3469 1376.700 0.1769690  100.0
#> 38  1468  476.322 0.4387120  100.0
#> 39  3524 1189.460 0.1635860  100.0
#> 40  5267 1644.960 0.2538320  100.0
#> 41  5048  941.543 0.3286410 1300.0
#> 42  1016  308.642 0.2300810 1300.0
#> 43  5605 1145.690 0.4641250 1300.0
#> 44  8793 2280.490 0.4204770 1300.0
#> 45  3475 1174.110 0.2007440  580.0
#> 46  1651  597.808 0.2626510  580.0
#> 47  5514 1455.880 0.1824530  580.0
#> 48  9718 1485.580 0.2004470  580.0
#> attr(,"scaled:center")
#>          area          peri         shape          perm 
#> -7187.7291667 -2682.2119375    -0.2181104  -415.4500000 
#> attr(,"scaled:scale")
#>         area         peri        shape         perm 
#> 3.725992e-04 6.984893e-04 1.197656e+01 2.284053e-03 
#> Calculating spline axis for variable 2 
#>     area     peri     shape   perm
#> 1   4990 2791.900 0.0903296    6.3
#> 2   7002 3892.600 0.1486220    6.3
#> 3   7558 3930.660 0.1833120    6.3
#> 4   7352 3869.320 0.1170630    6.3
#> 5   7943 3948.540 0.1224170   17.1
#> 6   7979 4010.150 0.1670450   17.1
#> 7   9333 4345.750 0.1896510   17.1
#> 8   8209 4344.750 0.1641270   17.1
#> 9   8393 3682.040 0.2036540  119.0
#> 10  6425 3098.650 0.1623940  119.0
#> 11  9364 4480.050 0.1509440  119.0
#> 12  8624 3986.240 0.1481410  119.0
#> 13 10651 4036.540 0.2285950   82.4
#> 14  8868 3518.040 0.2316230   82.4
#> 15  9417 3999.370 0.1725670   82.4
#> 16  8874 3629.070 0.1534810   82.4
#> 17 10962 4608.660 0.2043140   58.6
#> 18 10743 4787.620 0.2627270   58.6
#> 19 11878 4864.220 0.2000710   58.6
#> 20  9867 4479.410 0.1448100   58.6
#> 21  7838 3428.740 0.1138520  142.0
#> 22 11876 4353.140 0.2910290  142.0
#> 23 12212 4697.650 0.2400770  142.0
#> 24  8233 3518.440 0.1618650  142.0
#> 25  6360 1977.390 0.2808870  740.0
#> 26  4193 1379.350 0.1794550  740.0
#> 27  7416 1916.240 0.1918020  740.0
#> 28  5246 1585.420 0.1330830  740.0
#> 29  6509 1851.210 0.2252140  890.0
#> 30  4895 1239.660 0.3412730  890.0
#> 31  6775 1728.140 0.3116460  890.0
#> 32  7894 1461.060 0.2760160  890.0
#> 33  5980 1426.760 0.1976530  950.0
#> 34  5318  990.388 0.3266350  950.0
#> 35  7392 1350.760 0.1541920  950.0
#> 36  7894 1461.060 0.2760160  950.0
#> 37  3469 1376.700 0.1769690  100.0
#> 38  1468  476.322 0.4387120  100.0
#> 39  3524 1189.460 0.1635860  100.0
#> 40  5267 1644.960 0.2538320  100.0
#> 41  5048  941.543 0.3286410 1300.0
#> 42  1016  308.642 0.2300810 1300.0
#> 43  5605 1145.690 0.4641250 1300.0
#> 44  8793 2280.490 0.4204770 1300.0
#> 45  3475 1174.110 0.2007440  580.0
#> 46  1651  597.808 0.2626510  580.0
#> 47  5514 1455.880 0.1824530  580.0
#> 48  9718 1485.580 0.2004470  580.0
#> attr(,"scaled:center")
#>          area          peri         shape          perm 
#> -7187.7291667 -2682.2119375    -0.2181104  -415.4500000 
#> attr(,"scaled:scale")
#>         area         peri        shape         perm 
#> 3.725992e-04 6.984893e-04 1.197656e+01 2.284053e-03 
#> Calculating spline axis for variable 3 
#>     area     peri     shape   perm
#> 1   4990 2791.900 0.0903296    6.3
#> 2   7002 3892.600 0.1486220    6.3
#> 3   7558 3930.660 0.1833120    6.3
#> 4   7352 3869.320 0.1170630    6.3
#> 5   7943 3948.540 0.1224170   17.1
#> 6   7979 4010.150 0.1670450   17.1
#> 7   9333 4345.750 0.1896510   17.1
#> 8   8209 4344.750 0.1641270   17.1
#> 9   8393 3682.040 0.2036540  119.0
#> 10  6425 3098.650 0.1623940  119.0
#> 11  9364 4480.050 0.1509440  119.0
#> 12  8624 3986.240 0.1481410  119.0
#> 13 10651 4036.540 0.2285950   82.4
#> 14  8868 3518.040 0.2316230   82.4
#> 15  9417 3999.370 0.1725670   82.4
#> 16  8874 3629.070 0.1534810   82.4
#> 17 10962 4608.660 0.2043140   58.6
#> 18 10743 4787.620 0.2627270   58.6
#> 19 11878 4864.220 0.2000710   58.6
#> 20  9867 4479.410 0.1448100   58.6
#> 21  7838 3428.740 0.1138520  142.0
#> 22 11876 4353.140 0.2910290  142.0
#> 23 12212 4697.650 0.2400770  142.0
#> 24  8233 3518.440 0.1618650  142.0
#> 25  6360 1977.390 0.2808870  740.0
#> 26  4193 1379.350 0.1794550  740.0
#> 27  7416 1916.240 0.1918020  740.0
#> 28  5246 1585.420 0.1330830  740.0
#> 29  6509 1851.210 0.2252140  890.0
#> 30  4895 1239.660 0.3412730  890.0
#> 31  6775 1728.140 0.3116460  890.0
#> 32  7894 1461.060 0.2760160  890.0
#> 33  5980 1426.760 0.1976530  950.0
#> 34  5318  990.388 0.3266350  950.0
#> 35  7392 1350.760 0.1541920  950.0
#> 36  7894 1461.060 0.2760160  950.0
#> 37  3469 1376.700 0.1769690  100.0
#> 38  1468  476.322 0.4387120  100.0
#> 39  3524 1189.460 0.1635860  100.0
#> 40  5267 1644.960 0.2538320  100.0
#> 41  5048  941.543 0.3286410 1300.0
#> 42  1016  308.642 0.2300810 1300.0
#> 43  5605 1145.690 0.4641250 1300.0
#> 44  8793 2280.490 0.4204770 1300.0
#> 45  3475 1174.110 0.2007440  580.0
#> 46  1651  597.808 0.2626510  580.0
#> 47  5514 1455.880 0.1824530  580.0
#> 48  9718 1485.580 0.2004470  580.0
#> attr(,"scaled:center")
#>          area          peri         shape          perm 
#> -7187.7291667 -2682.2119375    -0.2181104  -415.4500000 
#> attr(,"scaled:scale")
#>         area         peri        shape         perm 
#> 3.725992e-04 6.984893e-04 1.197656e+01 2.284053e-03 
#> Calculating spline axis for variable 4 
#>     area     peri     shape   perm
#> 1   4990 2791.900 0.0903296    6.3
#> 2   7002 3892.600 0.1486220    6.3
#> 3   7558 3930.660 0.1833120    6.3
#> 4   7352 3869.320 0.1170630    6.3
#> 5   7943 3948.540 0.1224170   17.1
#> 6   7979 4010.150 0.1670450   17.1
#> 7   9333 4345.750 0.1896510   17.1
#> 8   8209 4344.750 0.1641270   17.1
#> 9   8393 3682.040 0.2036540  119.0
#> 10  6425 3098.650 0.1623940  119.0
#> 11  9364 4480.050 0.1509440  119.0
#> 12  8624 3986.240 0.1481410  119.0
#> 13 10651 4036.540 0.2285950   82.4
#> 14  8868 3518.040 0.2316230   82.4
#> 15  9417 3999.370 0.1725670   82.4
#> 16  8874 3629.070 0.1534810   82.4
#> 17 10962 4608.660 0.2043140   58.6
#> 18 10743 4787.620 0.2627270   58.6
#> 19 11878 4864.220 0.2000710   58.6
#> 20  9867 4479.410 0.1448100   58.6
#> 21  7838 3428.740 0.1138520  142.0
#> 22 11876 4353.140 0.2910290  142.0
#> 23 12212 4697.650 0.2400770  142.0
#> 24  8233 3518.440 0.1618650  142.0
#> 25  6360 1977.390 0.2808870  740.0
#> 26  4193 1379.350 0.1794550  740.0
#> 27  7416 1916.240 0.1918020  740.0
#> 28  5246 1585.420 0.1330830  740.0
#> 29  6509 1851.210 0.2252140  890.0
#> 30  4895 1239.660 0.3412730  890.0
#> 31  6775 1728.140 0.3116460  890.0
#> 32  7894 1461.060 0.2760160  890.0
#> 33  5980 1426.760 0.1976530  950.0
#> 34  5318  990.388 0.3266350  950.0
#> 35  7392 1350.760 0.1541920  950.0
#> 36  7894 1461.060 0.2760160  950.0
#> 37  3469 1376.700 0.1769690  100.0
#> 38  1468  476.322 0.4387120  100.0
#> 39  3524 1189.460 0.1635860  100.0
#> 40  5267 1644.960 0.2538320  100.0
#> 41  5048  941.543 0.3286410 1300.0
#> 42  1016  308.642 0.2300810 1300.0
#> 43  5605 1145.690 0.4641250 1300.0
#> 44  8793 2280.490 0.4204770 1300.0
#> 45  3475 1174.110 0.2007440  580.0
#> 46  1651  597.808 0.2626510  580.0
#> 47  5514 1455.880 0.1824530  580.0
#> 48  9718 1485.580 0.2004470  580.0
#> attr(,"scaled:center")
#>          area          peri         shape          perm 
#> -7187.7291667 -2682.2119375    -0.2181104  -415.4500000 
#> attr(,"scaled:scale")
#>         area         peri        shape         perm 
#> 3.725992e-04 6.984893e-04 1.197656e+01 2.284053e-03

By default the distance metric used for the analysis is Euclidean distance. The user can also specify the distance matrix Dmat as an \(n \times n\) matrix of a dist object. As illustration of a metric that is Euclidean embeddable, the Clark’s distance

\[ \delta_{ij}^2 = \sum_{k=1}^{p}\left(\frac{x_{ik}-x_{jk}}{x_{ik}+x_{jk}}\right)^2 \]

as defined by John C. Gower and Ngouenet (2005) is calculated below and used for constructing a PCO biplot. Note that the metric is only defined for strictly positive values, so that the data is scaled to values between \(1\) and \(2\).

Clark.dist <- function(X)
{
  n <- nrow(X)
  p <- ncol(X)
  Dmat <- matrix (0, nrow=n, ncol=n)
  for (i in 1:(n-1))
    for (j in (i+1):n)
      Dmat[i,j] <- sum(((X[i,] - X[j,])/(X[i,] + X[j,]))^2)
  sqrt(Dmat + t(Dmat))
}
my.data <- scale(rock, center=apply(rock,2,min), scale=diff(apply(rock,2,range)))+1
biplot(rock) |> PCO(Dmat = Clark.dist (rock), axes = "splines") |> plot()

#> Calculating spline axis for variable 1 
#>     area     peri     shape   perm
#> 1   4990 2791.900 0.0903296    6.3
#> 2   7002 3892.600 0.1486220    6.3
#> 3   7558 3930.660 0.1833120    6.3
#> 4   7352 3869.320 0.1170630    6.3
#> 5   7943 3948.540 0.1224170   17.1
#> 6   7979 4010.150 0.1670450   17.1
#> 7   9333 4345.750 0.1896510   17.1
#> 8   8209 4344.750 0.1641270   17.1
#> 9   8393 3682.040 0.2036540  119.0
#> 10  6425 3098.650 0.1623940  119.0
#> 11  9364 4480.050 0.1509440  119.0
#> 12  8624 3986.240 0.1481410  119.0
#> 13 10651 4036.540 0.2285950   82.4
#> 14  8868 3518.040 0.2316230   82.4
#> 15  9417 3999.370 0.1725670   82.4
#> 16  8874 3629.070 0.1534810   82.4
#> 17 10962 4608.660 0.2043140   58.6
#> 18 10743 4787.620 0.2627270   58.6
#> 19 11878 4864.220 0.2000710   58.6
#> 20  9867 4479.410 0.1448100   58.6
#> 21  7838 3428.740 0.1138520  142.0
#> 22 11876 4353.140 0.2910290  142.0
#> 23 12212 4697.650 0.2400770  142.0
#> 24  8233 3518.440 0.1618650  142.0
#> 25  6360 1977.390 0.2808870  740.0
#> 26  4193 1379.350 0.1794550  740.0
#> 27  7416 1916.240 0.1918020  740.0
#> 28  5246 1585.420 0.1330830  740.0
#> 29  6509 1851.210 0.2252140  890.0
#> 30  4895 1239.660 0.3412730  890.0
#> 31  6775 1728.140 0.3116460  890.0
#> 32  7894 1461.060 0.2760160  890.0
#> 33  5980 1426.760 0.1976530  950.0
#> 34  5318  990.388 0.3266350  950.0
#> 35  7392 1350.760 0.1541920  950.0
#> 36  7894 1461.060 0.2760160  950.0
#> 37  3469 1376.700 0.1769690  100.0
#> 38  1468  476.322 0.4387120  100.0
#> 39  3524 1189.460 0.1635860  100.0
#> 40  5267 1644.960 0.2538320  100.0
#> 41  5048  941.543 0.3286410 1300.0
#> 42  1016  308.642 0.2300810 1300.0
#> 43  5605 1145.690 0.4641250 1300.0
#> 44  8793 2280.490 0.4204770 1300.0
#> 45  3475 1174.110 0.2007440  580.0
#> 46  1651  597.808 0.2626510  580.0
#> 47  5514 1455.880 0.1824530  580.0
#> 48  9718 1485.580 0.2004470  580.0
#> attr(,"scaled:center")
#>          area          peri         shape          perm 
#> -7187.7291667 -2682.2119375    -0.2181104  -415.4500000 
#> Calculating spline axis for variable 2 
#>     area     peri     shape   perm
#> 1   4990 2791.900 0.0903296    6.3
#> 2   7002 3892.600 0.1486220    6.3
#> 3   7558 3930.660 0.1833120    6.3
#> 4   7352 3869.320 0.1170630    6.3
#> 5   7943 3948.540 0.1224170   17.1
#> 6   7979 4010.150 0.1670450   17.1
#> 7   9333 4345.750 0.1896510   17.1
#> 8   8209 4344.750 0.1641270   17.1
#> 9   8393 3682.040 0.2036540  119.0
#> 10  6425 3098.650 0.1623940  119.0
#> 11  9364 4480.050 0.1509440  119.0
#> 12  8624 3986.240 0.1481410  119.0
#> 13 10651 4036.540 0.2285950   82.4
#> 14  8868 3518.040 0.2316230   82.4
#> 15  9417 3999.370 0.1725670   82.4
#> 16  8874 3629.070 0.1534810   82.4
#> 17 10962 4608.660 0.2043140   58.6
#> 18 10743 4787.620 0.2627270   58.6
#> 19 11878 4864.220 0.2000710   58.6
#> 20  9867 4479.410 0.1448100   58.6
#> 21  7838 3428.740 0.1138520  142.0
#> 22 11876 4353.140 0.2910290  142.0
#> 23 12212 4697.650 0.2400770  142.0
#> 24  8233 3518.440 0.1618650  142.0
#> 25  6360 1977.390 0.2808870  740.0
#> 26  4193 1379.350 0.1794550  740.0
#> 27  7416 1916.240 0.1918020  740.0
#> 28  5246 1585.420 0.1330830  740.0
#> 29  6509 1851.210 0.2252140  890.0
#> 30  4895 1239.660 0.3412730  890.0
#> 31  6775 1728.140 0.3116460  890.0
#> 32  7894 1461.060 0.2760160  890.0
#> 33  5980 1426.760 0.1976530  950.0
#> 34  5318  990.388 0.3266350  950.0
#> 35  7392 1350.760 0.1541920  950.0
#> 36  7894 1461.060 0.2760160  950.0
#> 37  3469 1376.700 0.1769690  100.0
#> 38  1468  476.322 0.4387120  100.0
#> 39  3524 1189.460 0.1635860  100.0
#> 40  5267 1644.960 0.2538320  100.0
#> 41  5048  941.543 0.3286410 1300.0
#> 42  1016  308.642 0.2300810 1300.0
#> 43  5605 1145.690 0.4641250 1300.0
#> 44  8793 2280.490 0.4204770 1300.0
#> 45  3475 1174.110 0.2007440  580.0
#> 46  1651  597.808 0.2626510  580.0
#> 47  5514 1455.880 0.1824530  580.0
#> 48  9718 1485.580 0.2004470  580.0
#> attr(,"scaled:center")
#>          area          peri         shape          perm 
#> -7187.7291667 -2682.2119375    -0.2181104  -415.4500000 
#> Calculating spline axis for variable 3 
#>     area     peri     shape   perm
#> 1   4990 2791.900 0.0903296    6.3
#> 2   7002 3892.600 0.1486220    6.3
#> 3   7558 3930.660 0.1833120    6.3
#> 4   7352 3869.320 0.1170630    6.3
#> 5   7943 3948.540 0.1224170   17.1
#> 6   7979 4010.150 0.1670450   17.1
#> 7   9333 4345.750 0.1896510   17.1
#> 8   8209 4344.750 0.1641270   17.1
#> 9   8393 3682.040 0.2036540  119.0
#> 10  6425 3098.650 0.1623940  119.0
#> 11  9364 4480.050 0.1509440  119.0
#> 12  8624 3986.240 0.1481410  119.0
#> 13 10651 4036.540 0.2285950   82.4
#> 14  8868 3518.040 0.2316230   82.4
#> 15  9417 3999.370 0.1725670   82.4
#> 16  8874 3629.070 0.1534810   82.4
#> 17 10962 4608.660 0.2043140   58.6
#> 18 10743 4787.620 0.2627270   58.6
#> 19 11878 4864.220 0.2000710   58.6
#> 20  9867 4479.410 0.1448100   58.6
#> 21  7838 3428.740 0.1138520  142.0
#> 22 11876 4353.140 0.2910290  142.0
#> 23 12212 4697.650 0.2400770  142.0
#> 24  8233 3518.440 0.1618650  142.0
#> 25  6360 1977.390 0.2808870  740.0
#> 26  4193 1379.350 0.1794550  740.0
#> 27  7416 1916.240 0.1918020  740.0
#> 28  5246 1585.420 0.1330830  740.0
#> 29  6509 1851.210 0.2252140  890.0
#> 30  4895 1239.660 0.3412730  890.0
#> 31  6775 1728.140 0.3116460  890.0
#> 32  7894 1461.060 0.2760160  890.0
#> 33  5980 1426.760 0.1976530  950.0
#> 34  5318  990.388 0.3266350  950.0
#> 35  7392 1350.760 0.1541920  950.0
#> 36  7894 1461.060 0.2760160  950.0
#> 37  3469 1376.700 0.1769690  100.0
#> 38  1468  476.322 0.4387120  100.0
#> 39  3524 1189.460 0.1635860  100.0
#> 40  5267 1644.960 0.2538320  100.0
#> 41  5048  941.543 0.3286410 1300.0
#> 42  1016  308.642 0.2300810 1300.0
#> 43  5605 1145.690 0.4641250 1300.0
#> 44  8793 2280.490 0.4204770 1300.0
#> 45  3475 1174.110 0.2007440  580.0
#> 46  1651  597.808 0.2626510  580.0
#> 47  5514 1455.880 0.1824530  580.0
#> 48  9718 1485.580 0.2004470  580.0
#> attr(,"scaled:center")
#>          area          peri         shape          perm 
#> -7187.7291667 -2682.2119375    -0.2181104  -415.4500000 
#> Calculating spline axis for variable 4 
#>     area     peri     shape   perm
#> 1   4990 2791.900 0.0903296    6.3
#> 2   7002 3892.600 0.1486220    6.3
#> 3   7558 3930.660 0.1833120    6.3
#> 4   7352 3869.320 0.1170630    6.3
#> 5   7943 3948.540 0.1224170   17.1
#> 6   7979 4010.150 0.1670450   17.1
#> 7   9333 4345.750 0.1896510   17.1
#> 8   8209 4344.750 0.1641270   17.1
#> 9   8393 3682.040 0.2036540  119.0
#> 10  6425 3098.650 0.1623940  119.0
#> 11  9364 4480.050 0.1509440  119.0
#> 12  8624 3986.240 0.1481410  119.0
#> 13 10651 4036.540 0.2285950   82.4
#> 14  8868 3518.040 0.2316230   82.4
#> 15  9417 3999.370 0.1725670   82.4
#> 16  8874 3629.070 0.1534810   82.4
#> 17 10962 4608.660 0.2043140   58.6
#> 18 10743 4787.620 0.2627270   58.6
#> 19 11878 4864.220 0.2000710   58.6
#> 20  9867 4479.410 0.1448100   58.6
#> 21  7838 3428.740 0.1138520  142.0
#> 22 11876 4353.140 0.2910290  142.0
#> 23 12212 4697.650 0.2400770  142.0
#> 24  8233 3518.440 0.1618650  142.0
#> 25  6360 1977.390 0.2808870  740.0
#> 26  4193 1379.350 0.1794550  740.0
#> 27  7416 1916.240 0.1918020  740.0
#> 28  5246 1585.420 0.1330830  740.0
#> 29  6509 1851.210 0.2252140  890.0
#> 30  4895 1239.660 0.3412730  890.0
#> 31  6775 1728.140 0.3116460  890.0
#> 32  7894 1461.060 0.2760160  890.0
#> 33  5980 1426.760 0.1976530  950.0
#> 34  5318  990.388 0.3266350  950.0
#> 35  7392 1350.760 0.1541920  950.0
#> 36  7894 1461.060 0.2760160  950.0
#> 37  3469 1376.700 0.1769690  100.0
#> 38  1468  476.322 0.4387120  100.0
#> 39  3524 1189.460 0.1635860  100.0
#> 40  5267 1644.960 0.2538320  100.0
#> 41  5048  941.543 0.3286410 1300.0
#> 42  1016  308.642 0.2300810 1300.0
#> 43  5605 1145.690 0.4641250 1300.0
#> 44  8793 2280.490 0.4204770 1300.0
#> 45  3475 1174.110 0.2007440  580.0
#> 46  1651  597.808 0.2626510  580.0
#> 47  5514 1455.880 0.1824530  580.0
#> 48  9718 1485.580 0.2004470  580.0
#> attr(,"scaled:center")
#>          area          peri         shape          perm 
#> -7187.7291667 -2682.2119375    -0.2181104  -415.4500000

Alternatively, the user can specify a function which computes a distance matrix or dist object. The Manhattan distance is not Euclidean embeddable, but the square root of the distance is. The function sqrtManhattan is included as an example of a function computing a Euclidean embeddable dist object.

sqrtManhattan
#> function (X) 
#> {
#>     sqrt(stats::dist(X, method = "manhattan"))
#> }
#> <bytecode: 0x000002271155f128>
#> <environment: namespace:biplotEZ>
biplot(rock, scaled = TRUE) |> PCO(dist.func = sqrtManhattan) |> plot()

References

Gower, J. C. 1982. “Euclidean Distance Geometry.” Mathematical Scientist, 1–14.
Gower, J. C., S. Lubbe, and N. J. le Roux. 2011. Understanding Biplots. Wiley.
Gower, John C, and Roger F Ngouenet. 2005. “Nonlinearity Effects in Multidimensional Scaling.” Journal of Multivariate Analysis 94 (2): 344–65.