← Back to branch summary

~scipystats/+junk/skipper-revjp

« back to all changes in this revision

Viewing changes to nipy/fixes/scipy/stats/models/tests/test_regression.py

  • Committer: joep
  • Date: 2009-08-17 03:43:53 UTC
  • mfrom: (1799.1.26 skipper-working)
  • Revision ID: josef.pktd@gmail.com-20090817034353-ktg588ju3qhmruvq
merge from Skipper, trying to clear up move, keeping history

Show diffs side-by-side

added added

removed removed

Lines of Context:
nipy/fixes/scipy/stats/models/tests/test_regression.py
8
8
from models.tools import add_constant
9
9
from models.regression import OLS, AR, WLS, GLS, yule_walker
10
10
import models
 
11
from models import tools
11
12
from check_for_rpy import skip_rpy
12
13
from nose import SkipTest
 
14
from scipy.stats import t
13
15
 
14
16
W = standard_normal
15
17
DECIMAL = 4
20
22
skipR = skip_rpy()
21
23
if not skipR:
22
24
    from rpy import r
 
25
    from rmodelwrap import RModel
23
26
 
24
27
 
25
28
class check_regression_results(object):
68
71
        assert_almost_equal(self.res1.llf, self.res2.llf, DECIMAL)
69
72
 
70
73
    def test_AIC(self):
71
 
        aic = self.res1.information_criteria()['aic']
72
 
        assert_almost_equal(aic, self.res2.AIC, DECIMAL)
 
74
        assert_almost_equal(self.res1.aic, self.res2.AIC, DECIMAL)
73
75
 
74
76
    def test_BIC(self):
75
 
        bic = self.res1.information_criteria()['bic']
76
 
        assert_almost_equal(bic, self.res2.BIC, DECIMAL)
77
 
    
78
 
    @dec.skipif(True, "Results not included yet")
79
 
    def test_resids(self):
80
 
        pass
 
77
        assert_almost_equal(self.res1.bic, self.res2.BIC, DECIMAL)
 
78
 
 
79
#TODO: add residuals from another package or R
 
80
#    def test_resids(self):
 
81
#        assert_almost_equal(self.res1.resid, self.res2.resid, DECIMAL)
81
82
 
82
83
class test_ols(check_regression_results):
83
84
    def __init__(self):
88
89
        results = OLS(data.endog, data.exog).fit()
89
90
        self.res1 = results
90
91
        self.res2 = longley()
91
 
 
 
92
       
92
93
    def check_confidenceintervals(self, conf1, conf2):
93
94
        for i in range(len(conf1)):
94
95
            assert_approx_equal(conf1[i][0], conf2[i][0], 6)
98
99
    def check_params(self, params1, params2):
99
100
        assert_almost_equal(params1, params2, DECIMAL)
100
101
 
 
102
class TestFtest(object):
 
103
    def __init__(self):
 
104
        from models.datasets.longley.data import load
 
105
        data = load()
 
106
        self.data = data
 
107
        self.data.exog = add_constant(self.data.exog)
 
108
        self.res1 = OLS(data.endog, data.exog).fit()
 
109
        self.R = np.identity(7)[:-1,:]
 
110
        self.Ftest = self.res1.Ftest(self.R)
 
111
        
 
112
    def test_F(self):
 
113
        assert_almost_equal(self.Ftest.F, self.res1.F, DECIMAL)
 
114
 
 
115
    def test_p(self):
 
116
        assert_almost_equal(self.Ftest.p_val, self.res1.F_p, DECIMAL)
 
117
        
 
118
    def test_Df_denom(self):
 
119
        assert_equal(self.Ftest.df_denom, self.res1.df_resid)
 
120
 
 
121
    def test_Df_num(self):
 
122
        assert_equal(self.Ftest.df_num, tools.rank(self.R))
 
123
 
 
124
class TestFTest2(TestFtest):
 
125
    '''
 
126
    The first set of tests is a test of the fully specified model.
 
127
    The second set of tests is a joint test that the coefficient on 
 
128
    GNP = the coefficient on UNEMP  and that the coefficient on 
 
129
    POP = the coefficient on YEAR for the Longley dataset.
 
130
    '''
 
131
 
 
132
    def setup(self):
 
133
        if skipR:
 
134
            raise SkipTest, "Rpy not installed"
 
135
        try:
 
136
            r.library('car')
 
137
            self.R2 = [[0,1,-1,0,0,0,0],[0, 0, 0, 0, 1, -1, 0]]
 
138
            self.Ftest2 = self.res1.Ftest(self.R2)
 
139
            self.R_Results = RModel(self.data.endog, self.data.exog, r.lm).robj
 
140
            self.F = r.linear_hypothesis(self.R_Results,
 
141
                    r.c('x.2 = x.3', 'x.5 = x.6'))
 
142
        except:
 
143
            raise SkipTest, "car library not installed for R"
 
144
    
 
145
    def test_F(self):
 
146
        assert_almost_equal(self.Ftest2.F, self.F['F'][1], DECIMAL)
 
147
 
 
148
    def test_p(self):
 
149
        assert_almost_equal(self.Ftest2.p_val, self.F['Pr(>F)'][1], DECIMAL)
 
150
        
 
151
    def test_Df_denom(self):
 
152
        assert_equal(self.Ftest2.df_denom, self.F['Res.Df'][0])
 
153
 
 
154
    def test_Df_num(self):
 
155
        self.F['Res.Df'].reverse()
 
156
        assert_equal(self.Ftest2.df_num, np.subtract.reduce(self.F['Res.Df']))
 
157
        
 
158
class TestTtest(object):
 
159
    '''
 
160
    Test individual t-tests.  Ie., are the coefficients significantly
 
161
    different than zero.
 
162
    '''
 
163
    def __init__(self):
 
164
        from models.datasets.longley.data import load
 
165
        data = load()
 
166
        data.exog = add_constant(data.exog)
 
167
        self.res1 = OLS(data.endog, data.exog).fit()
 
168
 
 
169
    def setup(self):
 
170
        if skipR:
 
171
            raise SkipTest, "Rpy not installed"
 
172
        else:
 
173
            self.R = np.identity(len(self.res1.params))
 
174
            self.Ttest = self.res1.Ttest(self.R)
 
175
            
 
176
    def test_T(self):
 
177
        assert_almost_equal(np.diag(self.Ttest.t), self.res1.t(), DECIMAL)
 
178
 
 
179
    def test_sd(self):
 
180
        assert_almost_equal(np.diag(self.Ttest.sd), self.res1.bse, DECIMAL)
 
181
 
 
182
    def test_p_val(self):
 
183
        assert_almost_equal(np.diag(self.Ttest.p_val), 
 
184
                t.sf(np.abs(self.res1.t()),self.res1.df_resid), DECIMAL)
 
185
        
 
186
    def test_Df_denom(self):
 
187
        assert_equal(self.Ttest.df_denom, self.res1.df_resid)
 
188
 
 
189
    def test_effect(self):
 
190
        assert_almost_equal(self.Ttest.effect, self.res1.params)
 
191
 
 
192
 
 
193
class TestTtest2(TestTtest):
 
194
    '''
 
195
    Tests the hypothesis that the coefficients on POP and YEAR
 
196
    are equal.
 
197
    '''
 
198
    def setup(self):
 
199
        if skipR:
 
200
            raise SkipTest, "Rpy not installed"
 
201
        try:
 
202
            r.library('car')
 
203
            R = np.zeros(len(self.res1.params))
 
204
            R[4:6] = [1,-1]
 
205
            self.R = R
 
206
            self.Ttest1 = self.res1.Ttest(self.R)
 
207
            self.Ttest2 = r.linear_hypothesis(self.R_Results['F'], 'x.5 = x.6')
 
208
            t = np.inner(self.R, self.res1.params)/\
 
209
                (np.sign(np.inner(self.R, self.res1.params))*\
 
210
                np.sqrt(self.Ttest2['F'][1]))
 
211
            self.t = t
 
212
            effect = np.inner(self.R, self.res1.params)
 
213
            self.effect = effect
 
214
        except:
 
215
            raise SkipTest, "car library not installed for R"
 
216
 
 
217
    def test_T(self):
 
218
        assert_equal(self.Ttest1.t, self.t, DECIMAL)
 
219
 
 
220
    def test_sd(self):
 
221
        assert_almost_equal(self.Ttest1.sd, effect/t, DECIMAL)
 
222
 
 
223
    def test_p_val(self):
 
224
        assert_almost_equal(self.Ftest2.p_val, t.sf(self.t, self.F['Res.Df'][0]),
 
225
            DECIMAL)
 
226
        
 
227
    def test_Df_denom(self):
 
228
        assert_equal(self.Ftest2.df_denom, self.F['Res.Df'][0])
 
229
 
 
230
    def test_effect(self):
 
231
        assert_equal(self.Ttest1.effect, self.effect, DECIMAL)
 
232
 
 
233
 
101
234
class test_gls(object):
102
235
    '''
103
236
    These test results were obtained by replication with R.
139
272
        self.res1 = gls_res
140
273
        self.res2 = ols_res
141
274
        self.res2.conf_int = self.res2.conf_int()
142
 
        self.res2.BIC = self.res2.information_criteria()['bic']
143
 
        self.res2.AIC = self.res2.information_criteria()['aic']
 
275
        self.res2.BIC = self.res2.bic
 
276
        self.res2.AIC = self.res2.aic
144
277
 
145
278
    def check_confidenceintervals(self, conf1, conf2):
146
279
        assert_almost_equal(conf1, conf2, DECIMAL)