github.com/jingcheng-WU/gonum@v0.9.1-0.20210323123734-f1a2a11a8f7b/stat/cca_test.go (about)

     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  
     5  package stat_test
     6  
     7  import (
     8  	"testing"
     9  
    10  	"github.com/jingcheng-WU/gonum/floats"
    11  	"github.com/jingcheng-WU/gonum/mat"
    12  	"github.com/jingcheng-WU/gonum/stat"
    13  )
    14  
    15  func TestCanonicalCorrelations(t *testing.T) {
    16  tests:
    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.2574391924541896, -0.015847751662118038, -0.21221699346310258, -0.09457338038947205,
   121  				-0.48365944300184865, -0.3837101908138455, -0.14744483174159395, 0.6597324886718278,
   122  				-0.08007763658732961, -0.34935567428092285, -0.3287336458109394, -0.2862040444334662,
   123  				0.127758636038638, 0.7337427663667616, -0.4851134819036985, 0.22479648659701942,
   124  				-0.6969432006136685, 0.43417487760028844, 0.360287288763638, 0.029066160862628414,
   125  				-0.0990903250057202, -0.05034112154538474, -0.6384330631742202, 0.10223671362182897,
   126  				0.42604599637650303, -0.032333435130815824, 0.22895275160308087, 0.6419232947608798,
   127  			}),
   128  			wantqVecs: mat.NewDense(4, 4, []float64{
   129  				0.018166050236326788, 0.1583489460479047, 0.006672357764289544, -0.9871935400650647,
   130  				-0.23476990459861324, -0.9483314614936598, 0.14624205056313114, -0.1554470767919039,
   131  				-0.9700704038477144, 0.24060717410000537, 0.025183898422704167, 0.020913407435834964,
   132  				0.05930006823184807, 0.13304600030976868, 0.9889057151969495, 0.029116149472076858,
   133  			}),
   134  			wantphiVs: mat.NewDense(7, 4, []float64{
   135  				-0.002746223410819314, -0.009344451350088911, -0.04896439327142919, -0.015496718980582016,
   136  				-0.042856445527953785, 0.024170870211944927, -0.036072347209397136, 0.18389832305881182,
   137  				-1.2248435648802678, -5.603092136472504, -5.809414458379886, -4.792681219042103,
   138  				-0.00436848250946508, 0.34241011649776265, -0.4469961215717922, 0.11501618143536857,
   139  				-0.07415340695219563, 0.11931357949236807, 0.1115518305471455, 0.002163875832307984,
   140  				-0.023327032310162924, -0.1046330818178401, -0.38530459750774165, -0.016092787010290065,
   141  				0.00012930513878583387, -0.0004540746921447011, 0.0030296315865439264, 0.008189547797465318,
   142  			}),
   143  			wantpsiVs: mat.NewDense(4, 4, []float64{
   144  				0.030159336201738367, 0.3002219289647159, -0.08782173775936601, -1.9583226531517122,
   145  				-0.00654831040738931, -0.03922120867162458, 0.011757077620998818, -0.006111306448187141,
   146  				-0.0052075523350125505, 0.004577020045295936, 0.0022762313289591976, 0.0008441873006823151,
   147  				0.0020111735096325924, -0.0037352799829939247, 0.12925780716217938, 0.10377090563297825,
   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  }