1// Copyright ©2016 The Gonum Authors. All rights reserved. 2// Use of this source code is governed by a BSD-style 3// license that can be found in the LICENSE file. 4 5package stat_test 6 7import ( 8 "testing" 9 10 "gonum.org/v1/gonum/floats" 11 "gonum.org/v1/gonum/mat" 12 "gonum.org/v1/gonum/stat" 13) 14 15func TestCanonicalCorrelations(t *testing.T) { 16tests: 17 for i, test := range []struct { 18 xdata mat.Matrix 19 ydata mat.Matrix 20 weights []float64 21 wantCorrs []float64 22 wantpVecs *mat.Dense 23 wantqVecs *mat.Dense 24 wantphiVs *mat.Dense 25 wantpsiVs *mat.Dense 26 epsilon float64 27 }{ 28 // Test results verified using R. 29 { // Truncated iris data, Sepal vs Petal measurements. 30 xdata: mat.NewDense(10, 2, []float64{ 31 5.1, 3.5, 32 4.9, 3.0, 33 4.7, 3.2, 34 4.6, 3.1, 35 5.0, 3.6, 36 5.4, 3.9, 37 4.6, 3.4, 38 5.0, 3.4, 39 4.4, 2.9, 40 4.9, 3.1, 41 }), 42 ydata: mat.NewDense(10, 2, []float64{ 43 1.4, 0.2, 44 1.4, 0.2, 45 1.3, 0.2, 46 1.5, 0.2, 47 1.4, 0.2, 48 1.7, 0.4, 49 1.4, 0.3, 50 1.5, 0.2, 51 1.4, 0.2, 52 1.5, 0.1, 53 }), 54 wantCorrs: []float64{0.7250624174504773, 0.5547679185730191}, 55 wantpVecs: mat.NewDense(2, 2, []float64{ 56 0.0765914610875867, 0.9970625597666721, 57 0.9970625597666721, -0.0765914610875868, 58 }), 59 wantqVecs: mat.NewDense(2, 2, []float64{ 60 0.3075184850910837, 0.9515421069649439, 61 0.9515421069649439, -0.3075184850910837, 62 }), 63 wantphiVs: mat.NewDense(2, 2, []float64{ 64 -1.9794877596804641, 5.2016325219025124, 65 4.5211829944066553, -2.7263663170835697, 66 }), 67 wantpsiVs: mat.NewDense(2, 2, []float64{ 68 -0.0613084818030103, 10.8514169865438941, 69 12.7209032660734298, -7.6793888180353775, 70 }), 71 epsilon: 1e-12, 72 }, 73 // Test results compared to those results presented in examples by 74 // Koch, Inge. Analysis of multivariate and high-dimensional data. 75 // Vol. 32. Cambridge University Press, 2013. ISBN: 9780521887939 76 { // ASA Car Exposition Data of Ramos and Donoho (1983) 77 // Displacement, Horsepower, Weight 78 xdata: carData.Slice(0, 392, 0, 3), 79 // Acceleration, MPG 80 ydata: carData.Slice(0, 392, 3, 5), 81 wantCorrs: []float64{0.8782187384352336, 0.6328187219216761}, 82 wantpVecs: mat.NewDense(3, 2, []float64{ 83 0.3218296374829181, 0.3947540257657075, 84 0.4162807660635797, 0.7573719053303306, 85 0.8503740401982725, -0.5201509936144236, 86 }), 87 wantqVecs: mat.NewDense(2, 2, []float64{ 88 -0.5161984172278830, -0.8564690269072364, 89 -0.8564690269072364, 0.5161984172278830, 90 }), 91 wantphiVs: mat.NewDense(3, 2, []float64{ 92 0.0025033152994308, 0.0047795464118615, 93 0.0201923608080173, 0.0409150208725958, 94 -0.0000247374128745, -0.0026766435161875, 95 }), 96 wantpsiVs: mat.NewDense(2, 2, []float64{ 97 -0.1666196759760772, -0.3637393866139658, 98 -0.0915512109649727, 0.1077863777929168, 99 }), 100 epsilon: 1e-12, 101 }, 102 // Test results compared to those results presented in examples by 103 // Koch, Inge. Analysis of multivariate and high-dimensional data. 104 // Vol. 32. Cambridge University Press, 2013. ISBN: 9780521887939 105 { // Boston Housing Data of Harrison and Rubinfeld (1978) 106 // Per capita crime rate by town, 107 // Proportion of non-retail business acres per town, 108 // Nitric oxide concentration (parts per 10 million), 109 // Weighted distances to Boston employment centres, 110 // Index of accessibility to radial highways, 111 // Pupil-teacher ratio by town, Proportion of blacks by town 112 xdata: bostonData.Slice(0, 506, 0, 7), 113 // Average number of rooms per dwelling, 114 // Proportion of owner-occupied units built prior to 1940, 115 // Full-value property-tax rate per $10000, 116 // Median value of owner-occupied homes in $1000s 117 ydata: bostonData.Slice(0, 506, 7, 11), 118 wantCorrs: []float64{0.9451239443886021, 0.6786622733370654, 0.5714338361583764, 0.2009739704710440}, 119 wantpVecs: mat.NewDense(7, 4, []float64{ 120 -0.2574391924541903, 0.0158477516621194, 0.2122169934631024, -0.0945733803894706, 121 -0.4836594430018478, 0.3837101908138468, 0.1474448317415911, 0.6597324886718275, 122 -0.0800776365873296, 0.3493556742809252, 0.3287336458109373, -0.2862040444334655, 123 0.1277586360386374, -0.7337427663667596, 0.4851134819037011, 0.2247964865970192, 124 -0.6969432006136684, -0.4341748776002893, -0.3602872887636357, 0.0290661608626292, 125 -0.0990903250057199, 0.0503411215453873, 0.6384330631742202, 0.1022367136218303, 126 0.4260459963765036, 0.0323334351308141, -0.2289527516030810, 0.6419232947608805, 127 }), 128 wantqVecs: mat.NewDense(4, 4, []float64{ 129 0.0181660502363264, -0.1583489460479038, -0.0066723577642883, -0.9871935400650649, 130 -0.2347699045986119, 0.9483314614936594, -0.1462420505631345, -0.1554470767919033, 131 -0.9700704038477141, -0.2406071741000039, -0.0251838984227037, 0.0209134074358349, 132 0.0593000682318482, -0.1330460003097728, -0.9889057151969489, 0.0291161494720761, 133 }), 134 wantphiVs: mat.NewDense(7, 4, []float64{ 135 -0.0027462234108197, 0.0093444513500898, 0.0489643932714296, -0.0154967189805819, 136 -0.0428564455279537, -0.0241708702119420, 0.0360723472093996, 0.1838983230588095, 137 -1.2248435648802380, 5.6030921364723980, 5.8094144583797025, -4.7926812190419676, 138 -0.0043684825094649, -0.3424101164977618, 0.4469961215717917, 0.1150161814353696, 139 -0.0741534069521954, -0.1193135794923700, -0.1115518305471460, 0.0021638758323088, 140 -0.0233270323101624, 0.1046330818178399, 0.3853045975077387, -0.0160927870102877, 141 0.0001293051387859, 0.0004540746921446, -0.0030296315865440, 0.0081895477974654, 142 }), 143 wantpsiVs: mat.NewDense(4, 4, []float64{ 144 0.0301593362017375, -0.3002219289647127, 0.0878217377593682, -1.9583226531517062, 145 -0.0065483104073892, 0.0392212086716247, -0.0117570776209991, -0.0061113064481860, 146 -0.0052075523350125, -0.0045770200452960, -0.0022762313289592, 0.0008441873006821, 147 0.0020111735096327, 0.0037352799829930, -0.1292578071621794, 0.1037709056329765, 148 }), 149 epsilon: 1e-12, 150 }, 151 } { 152 var cc stat.CC 153 var corrs []float64 154 var pVecs, qVecs mat.Dense 155 var phiVs, psiVs mat.Dense 156 for j := 0; j < 2; j++ { 157 err := cc.CanonicalCorrelations(test.xdata, test.ydata, test.weights) 158 if err != nil { 159 t.Errorf("%d use %d: unexpected error: %v", i, j, err) 160 continue tests 161 } 162 163 corrs = cc.CorrsTo(corrs) 164 cc.LeftTo(&pVecs, true) 165 cc.RightTo(&qVecs, true) 166 cc.LeftTo(&phiVs, false) 167 cc.RightTo(&psiVs, false) 168 169 if !floats.EqualApprox(corrs, test.wantCorrs, test.epsilon) { 170 t.Errorf("%d use %d: unexpected variance result got:%v, want:%v", 171 i, j, corrs, test.wantCorrs) 172 } 173 if !mat.EqualApprox(&pVecs, test.wantpVecs, test.epsilon) { 174 t.Errorf("%d use %d: unexpected CCA result got:\n%v\nwant:\n%v", 175 i, j, mat.Formatted(&pVecs), mat.Formatted(test.wantpVecs)) 176 } 177 if !mat.EqualApprox(&qVecs, test.wantqVecs, test.epsilon) { 178 t.Errorf("%d use %d: unexpected CCA result got:\n%v\nwant:\n%v", 179 i, j, mat.Formatted(&qVecs), mat.Formatted(test.wantqVecs)) 180 } 181 if !mat.EqualApprox(&phiVs, test.wantphiVs, test.epsilon) { 182 t.Errorf("%d use %d: unexpected CCA result got:\n%v\nwant:\n%v", 183 i, j, mat.Formatted(&phiVs), mat.Formatted(test.wantphiVs)) 184 } 185 if !mat.EqualApprox(&psiVs, test.wantpsiVs, test.epsilon) { 186 t.Errorf("%d use %d: unexpected CCA result got:\n%v\nwant:\n%v", 187 i, j, mat.Formatted(&psiVs), mat.Formatted(test.wantpsiVs)) 188 } 189 } 190 } 191} 192