stat/mds: reduce work in TorgersonScaling and fix k return
The fix for k papers over inaccuracies in near zero Eigenvalues. So this
may want further work there.
MODIFY
stat/mds/mds.go
stat/mds/mds.go
@@ -15,6 +15,7 @@ import (
// Euclidean coordinates. TorgersonScaling places the coordinates in dst and
// returns it and the number of positive Eigenvalues if successful.
// If the scaling is not successful, dst is returned, but will not be a valid scaling.
+// When the scaling is successful, mds will be resized to k columns wide.
// Eigenvalues will be copied into eigdst and returned as eig if it is provided.
//
// If dst is nil, a new mat.Dense is allocated. If dst is not a zero matrix,
@@ -62,6 +63,7 @@ func TorgersonScaling(dst *mat.Dense, eigdst []float64, dis mat.Symmetric) (k in
vals := ed.Values(nil)
reverse(vals, dst.RawMatrix())
copy(eigdst, vals)
+
for i, v := range vals {
if v < 0 {
vals[i] = 0
@@ -71,9 +73,10 @@ func TorgersonScaling(dst *mat.Dense, eigdst []float64, dis mat.Symmetric) (k in
vals[i] = math.Sqrt(v)
}
- // TODO(kortschak): Make use of the knowledge
- // of k to avoid doing unnecessary work.
- dst.Mul(dst, mat.NewDiagDense(len(vals), vals))
+ var tmp mat.Dense
+ tmp.Mul(dst, mat.NewDiagonalRect(n, k, vals[:k]))
+ dst = dst.Slice(0, n, 0, k).(*mat.Dense)
+ dst.Copy(&tmp)
return k, dst, eigdst
}

MODIFY
stat/mds/mds_test.go
stat/mds/mds_test.go
@@ -71,29 +71,32 @@ var torgersonScalingTests = []struct {
3927, 2868, 1616, 1786, 2196, 1403, 650, 2068, 3886, 949, 1500, 3231, 2108, 3188, 2428, 2187, 1754, 1827, 2707, 0, 2105,
1991, 1802, 1175, 1381, 1588, 937, 1455, 1019, 2974, 1155, 1205, 2937, 1157, 2409, 1363, 898, 428, 1249, 1209, 2105, 0,
}),
+ // Note that k here is 12 despite the result from R's cmdscale.
+ // This is due to disparity between the ED performed by R and Gonum.
+ // See https://github.com/gonum/gonum/issues/768
wantK: 12,
- want: mat.NewDense(21, 11, []float64{
- 2290.274680, 1798.80293, 53.79314, -103.826958, -156.955115, 54.755434, -47.6768205, 1.241284, -14.893196, -6.366664, 4.818373,
- -825.382790, 546.81148, -113.85842, 84.585831, 291.440759, -33.046236, -74.5267190, 3.766233, 225.620420, -21.270973, 22.050484,
- 59.183341, -367.08135, 177.55291, 38.797514, -95.620447, 40.058268, 2.3212184, 34.351715, -2.262151, -129.298820, 75.568686,
- -82.845973, -429.91466, 300.19274, 106.353695, -180.446140, 31.336985, 88.6540176, 5.102633, 87.657362, 86.867202, 156.097411,
- -352.499435, -290.90843, 457.35294, 111.449150, -417.496682, -138.972847, 46.7759317, 23.806769, 108.307203, 99.578699, -103.623270,
- 293.689633, -405.31194, 360.09323, -636.202379, 159.392662, -9.560652, 33.1089500, 9.193680, 19.417287, 7.925230, -1.272483,
- 681.931545, -1108.64478, 26.09257, 151.693056, 254.114892, 197.254208, -176.0359152, -6.508596, -59.829991, 64.904841, -26.466935,
- -9.423364, 240.40600, -344.20659, 656.121110, -138.082259, -19.891841, -0.3128244, -20.167733, -18.250209, -11.175662, 2.160076,
- -2048.449113, 642.45854, 167.86631, 78.621423, 239.881353, -51.009028, 60.7527717, 125.492010, -144.940590, 35.754867, 13.633563,
- 561.108970, -773.36929, 80.91722, 48.548472, -129.512374, 20.426493, 207.7920058, 23.646986, -166.184382, -64.027105, -5.187863,
- 164.921799, -549.36704, 270.82327, 116.886334, 62.711250, 117.489888, -252.4906676, -47.193881, -62.053055, 75.466644, 4.907670,
- -1935.040811, 49.12514, -483.02056, -315.241752, -277.696589, -92.444183, -76.2549266, -28.529861, -26.279840, -2.503630, 7.050184,
- -226.423236, 187.08779, -358.43234, -257.737009, -190.609088, -110.666015, -119.7638108, -43.118326, -76.066370, 89.289877, 18.820576,
- -1423.353697, 305.87513, 253.26763, 2.478812, 5.366334, 205.060613, -45.8020725, -75.244576, -4.303057, -150.778068, -23.118399,
- -299.498710, 388.80726, -109.17417, 12.651217, 231.405239, 154.749482, 237.4507512, -227.372584, 25.449825, 71.404170, -35.167507,
- 260.878046, 416.67381, -171.52428, 20.926369, 194.345825, -19.740272, 121.3353507, 184.415363, -5.563374, 90.507909, -22.562758,
- 587.675679, 81.18224, -75.88485, 13.080496, 114.688714, -185.462926, -47.4707495, 209.170525, 39.452744, -65.352092, -23.127035,
- -156.836257, -211.13911, 131.30852, 27.089432, -100.800743, 6.794791, -25.6913210, -7.904762, 39.731234, -79.091817, -75.116800,
- 709.413282, 1109.36665, -179.83052, -109.895049, -90.243055, 316.665218, 9.9070412, 9.976867, 6.316941, -13.285987, 1.496742,
- 839.445911, -1836.79055, -541.35188, -108.755016, 25.155430, 41.399927, 71.8793196, 9.438407, 73.554358, -35.415426, -6.021569,
- 911.230500, 205.93020, 98.02313, 62.375253, 198.960032, -525.197308, -13.9515306, -183.562155, -44.881159, -43.133195, 15.060856,
+ want: mat.NewDense(21, 12, []float64{
+ 2290.274680, 1798.80293, 53.79314, -103.826958, -156.955115, 54.755434, -47.6768205, 1.241284, -14.893196, -6.366664, 4.818373, 0,
+ -825.382790, 546.81148, -113.85842, 84.585831, 291.440759, -33.046236, -74.5267190, 3.766233, 225.620420, -21.270973, 22.050484, 0,
+ 59.183341, -367.08135, 177.55291, 38.797514, -95.620447, 40.058268, 2.3212184, 34.351715, -2.262151, -129.298820, 75.568686, 0,
+ -82.845973, -429.91466, 300.19274, 106.353695, -180.446140, 31.336985, 88.6540176, 5.102633, 87.657362, 86.867202, 156.097411, 0,
+ -352.499435, -290.90843, 457.35294, 111.449150, -417.496682, -138.972847, 46.7759317, 23.806769, 108.307203, 99.578699, -103.623270, 0,
+ 293.689633, -405.31194, 360.09323, -636.202379, 159.392662, -9.560652, 33.1089500, 9.193680, 19.417287, 7.925230, -1.272483, 0,
+ 681.931545, -1108.64478, 26.09257, 151.693056, 254.114892, 197.254208, -176.0359152, -6.508596, -59.829991, 64.904841, -26.466935, 0,
+ -9.423364, 240.40600, -344.20659, 656.121110, -138.082259, -19.891841, -0.3128244, -20.167733, -18.250209, -11.175662, 2.160076, 0,
+ -2048.449113, 642.45854, 167.86631, 78.621423, 239.881353, -51.009028, 60.7527717, 125.492010, -144.940590, 35.754867, 13.633563, 0,
+ 561.108970, -773.36929, 80.91722, 48.548472, -129.512374, 20.426493, 207.7920058, 23.646986, -166.184382, -64.027105, -5.187863, 0,
+ 164.921799, -549.36704, 270.82327, 116.886334, 62.711250, 117.489888, -252.4906676, -47.193881, -62.053055, 75.466644, 4.907670, 0,
+ -1935.040811, 49.12514, -483.02056, -315.241752, -277.696589, -92.444183, -76.2549266, -28.529861, -26.279840, -2.503630, 7.050184, 0,
+ -226.423236, 187.08779, -358.43234, -257.737009, -190.609088, -110.666015, -119.7638108, -43.118326, -76.066370, 89.289877, 18.820576, 0,
+ -1423.353697, 305.87513, 253.26763, 2.478812, 5.366334, 205.060613, -45.8020725, -75.244576, -4.303057, -150.778068, -23.118399, 0,
+ -299.498710, 388.80726, -109.17417, 12.651217, 231.405239, 154.749482, 237.4507512, -227.372584, 25.449825, 71.404170, -35.167507, 0,
+ 260.878046, 416.67381, -171.52428, 20.926369, 194.345825, -19.740272, 121.3353507, 184.415363, -5.563374, 90.507909, -22.562758, 0,
+ 587.675679, 81.18224, -75.88485, 13.080496, 114.688714, -185.462926, -47.4707495, 209.170525, 39.452744, -65.352092, -23.127035, 0,
+ -156.836257, -211.13911, 131.30852, 27.089432, -100.800743, 6.794791, -25.6913210, -7.904762, 39.731234, -79.091817, -75.116800, 0,
+ 709.413282, 1109.36665, -179.83052, -109.895049, -90.243055, 316.665218, 9.9070412, 9.976867, 6.316941, -13.285987, 1.496742, 0,
+ 839.445911, -1836.79055, -541.35188, -108.755016, 25.155430, 41.399927, 71.8793196, 9.438407, 73.554358, -35.415426, -6.021569, 0,
+ 911.230500, 205.93020, 98.02313, 62.375253, 198.960032, -525.197308, -13.9515306, -183.562155, -44.881159, -43.133195, 15.060856, 0,
}),
wantVals: []float64{
1.953837708954e+07, 1.185655533400e+07, 1.528844467987e+06, 1.118741950509e+06,
@@ -114,14 +117,12 @@ func TestTorgersonScaling(t *testing.T) {
t.Error("unexpected scaling failure")
continue
}
- _, wc := test.want.Dims()
- n := test.dis.Symmetric()
if gotK != test.wantK {
t.Errorf("unexpected k for test %d: got:%d want:%d", i, gotK, test.wantK)
}
- if !mat.EqualApprox(colAbs{got.Slice(0, n, 0, wc)}, colAbs{test.want}, 1e-6) {
- t.Errorf("unexpected result for test %d:\ngot:\n%.4f\nwant:\n%.4f",
- i, mat.Formatted(got.Slice(0, n, 0, wc)), mat.Formatted(test.want))
+ if !mat.EqualApprox(colAbs{got}, colAbs{test.want}, 1e-5) {
+ t.Errorf("unexpected result for test %d:\ngot:\n%.5f\nwant:\n%.5f",
+ i, mat.Formatted(got), mat.Formatted(test.want))
}
if !floats.EqualApprox(gotVals, test.wantVals, 1e-10) {
t.Errorf("unexpected Eigenvalues for test %d:\ngot: %.12e\nwant:%.12e", i, gotVals, test.wantVals)

stat/mds: reduce work in TorgersonScaling and fix k return

The fix for k papers over inaccuracies in near zero Eigenvalues. So this
may want further work there.