1""" Testing models module 2""" 3 4import numpy as np 5import pytest 6 7from numpy.testing import assert_array_almost_equal 8from nilearn.glm import OLSModel 9 10 11N = 10 12X = np.c_[np.linspace(- 1, 1, N), np.ones((N,))] 13Y = np.r_[range(5), range(1, 6)] 14MODEL = OLSModel(X) 15RESULTS = MODEL.fit(Y) 16 17""" R script 18 19:: 20 21 X = cbind(0:9 * 2/9 -1, 1) 22 Y = as.matrix(c(0:4, 1:5)) 23 results = lm(Y ~ X-1) 24 print(results) 25 print(summary(results)) 26 27gives:: 28 29 Call: 30 lm(formula = Y ~ X - 1) 31 32 Coefficients: 33 X1 X2 34 1.773 2.500 35 36 Residuals: 37 Min 1Q Median 3Q Max 38 -1.6970 -0.6667 0.0000 0.6667 1.6970 39 40 Coefficients: 41 Estimate Std. Error t value Pr(>|t|) 42 X1 1.7727 0.5455 3.250 0.0117 * 43 X2 2.5000 0.3482 7.181 9.42e-05 *** 44 --- 45 46 Residual standard error: 1.101 on 8 degrees of freedom 47 Multiple R-squared: 0.8859, Adjusted R-squared: 0.8574 48 F-statistic: 31.06 on 2 and 8 DF, p-value: 0.0001694 49""" 50 51 52def test_model(): 53 # Check basics about the model fit 54 # Check we fit the mean 55 assert_array_almost_equal(RESULTS.theta[1], np.mean(Y)) 56 # Check we get the same as R 57 assert_array_almost_equal(RESULTS.theta, [1.773, 2.5], 3) 58 percentile = np.percentile 59 pcts = percentile(RESULTS.resid, [0, 25, 50, 75, 100]) 60 assert_array_almost_equal(pcts, [-1.6970, -0.6667, 0, 0.6667, 1.6970], 4) 61 62 63def test_t_contrast(): 64 # Test individual t against R 65 assert_array_almost_equal(RESULTS.t(0), 3.25) 66 assert_array_almost_equal(RESULTS.t(1), 7.181, 3) 67 # And contrast 68 assert_array_almost_equal(RESULTS.Tcontrast([1, 0]).t, 3.25) 69 assert_array_almost_equal(RESULTS.Tcontrast([0, 1]).t, 7.181, 3) 70 # Input matrix checked for size 71 with pytest.raises(ValueError): 72 RESULTS.Tcontrast([1]) 73 with pytest.raises(ValueError): 74 RESULTS.Tcontrast([1, 0, 0]) 75 # And shape 76 with pytest.raises(ValueError): 77 RESULTS.Tcontrast(np.array([1, 0])[:, None]) 78 79 80def test_t_output(): 81 # Check we get required outputs 82 exp_t = RESULTS.t(0) 83 exp_effect = RESULTS.theta[0] 84 exp_sd = exp_effect / exp_t 85 res = RESULTS.Tcontrast([1, 0]) 86 assert_array_almost_equal(res.t, exp_t) 87 assert_array_almost_equal(res.effect, exp_effect) 88 assert_array_almost_equal(res.sd, exp_sd) 89 res = RESULTS.Tcontrast([1, 0], store=('effect',)) 90 assert res.t is None 91 assert_array_almost_equal(res.effect, exp_effect) 92 assert res.sd is None 93 res = RESULTS.Tcontrast([1, 0], store=('t',)) 94 assert_array_almost_equal(res.t, exp_t) 95 assert res.effect is None 96 assert res.sd is None 97 res = RESULTS.Tcontrast([1, 0], store=('sd',)) 98 assert res.t is None 99 assert res.effect is None 100 assert_array_almost_equal(res.sd, exp_sd) 101 res = RESULTS.Tcontrast([1, 0], store=('effect', 'sd')) 102 assert res.t is None 103 assert_array_almost_equal(res.effect, exp_effect) 104 assert_array_almost_equal(res.sd, exp_sd) 105 106 107def test_f_output(): 108 # Test f_output 109 res = RESULTS.Fcontrast([1, 0]) 110 exp_f = RESULTS.t(0) ** 2 111 assert_array_almost_equal(exp_f, res.F) 112 # Test arrays work as well as lists 113 res = RESULTS.Fcontrast(np.array([1, 0])) 114 assert_array_almost_equal(exp_f, res.F) 115 # Test with matrix against R 116 res = RESULTS.Fcontrast(np.eye(2)) 117 assert_array_almost_equal(31.06, res.F, 2) 118 # Input matrix checked for size 119 with pytest.raises(ValueError): 120 RESULTS.Fcontrast([1]) 121 with pytest.raises(ValueError): 122 RESULTS.Fcontrast([1, 0, 0]) 123 # And shape 124 with pytest.raises(ValueError): 125 RESULTS.Fcontrast(np.array([1, 0])[:, None]) 126 127 128def test_f_output_new_api(): 129 res = RESULTS.Fcontrast([1, 0]) 130 assert_array_almost_equal(res.effect, RESULTS.theta[0]) 131 assert_array_almost_equal(res.covariance, RESULTS.vcov()[0][0]) 132 133 134def test_conf_int(): 135 lower_, upper_ = RESULTS.conf_int() 136 assert (lower_ < upper_).all() 137 assert (lower_ > upper_ - 10).all() 138 lower_, upper_ = RESULTS.conf_int(cols=[1]).T 139 assert lower_ < upper_ 140 assert lower_ > upper_ - 10 141