33import numpy as np
44import pandas as pd
55import pytest
6- from scipy .stats import norm
76
87import doubleml as dml
9- from doubleml .utils ._estimation import _aggregate_coefs_and_ses
108
119from ._utils_blp_manual import blp_confint , fit_blp
1210
@@ -31,43 +29,64 @@ def use_t(request):
3129 return request .param
3230
3331
32+ @pytest .fixture (scope = "module" , params = [1 , 3 ])
33+ def n_rep (request ):
34+ return request .param
35+
36+
3437@pytest .fixture (scope = "module" )
35- def dml_blp_fixture (ci_joint , ci_level , cov_type , use_t ):
38+ def dml_blp_fixture (ci_joint , ci_level , cov_type , use_t , n_rep ):
3639 n = 50
3740 kwargs = {"cov_type" : cov_type , "use_t" : use_t }
3841
3942 np .random .seed (42 )
4043 random_basis = pd .DataFrame (np .random .normal (0 , 1 , size = (n , 3 )))
41- random_signal = np .random .normal (0 , 1 , size = (n ,))
44+ random_signal = np .random .normal (0 , 1 , size = (n , n_rep ))
4245
4346 blp = dml .DoubleMLBLP (random_signal , random_basis )
4447
4548 blp_obj = copy .copy (blp )
4649 blp .fit (** kwargs )
47- blp_manual = fit_blp (random_signal , random_basis , ** kwargs )
50+ blp_manual = []
51+ for i in range (n_rep ):
52+ blp_manual .append (fit_blp (random_signal [:, i ], random_basis , ** kwargs ))
4853
4954 np .random .seed (42 )
5055 ci_1 = blp .confint (random_basis , joint = ci_joint , level = ci_level , n_rep_boot = 1000 )
5156 np .random .seed (42 )
5257 ci_2 = blp .confint (joint = ci_joint , level = ci_level , n_rep_boot = 1000 )
53- expected_ci_2 = np .vstack (
54- (
55- blp .blp_model .conf_int (alpha = (1 - ci_level ))[0 ],
56- blp .blp_model .params ,
57- blp .blp_model .conf_int (alpha = (1 - ci_level ))[1 ],
58+ expected_ci_2 = []
59+ for i in range (n_rep ):
60+ expected_ci_2 .append (
61+ np .vstack (
62+ (
63+ blp .blp_model [i ].conf_int (alpha = (1 - ci_level ))[0 ],
64+ blp .blp_model [i ].params ,
65+ blp .blp_model [i ].conf_int (alpha = (1 - ci_level ))[1 ],
66+ )
67+ ).T
5868 )
59- ). T
69+ expected_ci_2 = np . median ( np . array ( expected_ci_2 ), axis = 0 )
6070
6171 np .random .seed (42 )
62- ci_manual = blp_confint (blp_manual , random_basis , joint = ci_joint , level = ci_level , n_rep_boot = 1000 )
72+ ci_manual = []
73+ for i in range (n_rep ):
74+ ci_manual .append (blp_confint (blp_manual [i ], random_basis , joint = ci_joint , level = ci_level , n_rep_boot = 1000 ))
75+ ci_manual = np .median (np .array ([ci_manual [i ].to_numpy () for i in range (n_rep )]), axis = 0 )
76+
77+ coef_manual = np .median (np .array ([blp_manual [i ].params for i in range (n_rep )]), axis = 0 )
78+ omega_manual = np .transpose (np .array ([blp_manual [i ].cov_params ().to_numpy () for i in range (n_rep )]), (1 , 2 , 0 ))
79+
80+ fittedvalues_manual = np .array ([blp_manual [i ].fittedvalues for i in range (n_rep )])
81+ fittedvalues = np .array ([blp .blp_model [i ].fittedvalues for i in range (n_rep )])
6382
6483 res_dict = {
65- "coef" : blp .blp_model . params ,
66- "coef_manual" : blp_manual . params ,
67- "values" : blp . blp_model . fittedvalues ,
68- "values_manual" : blp_manual . fittedvalues ,
84+ "coef" : blp .coef ,
85+ "coef_manual" : coef_manual ,
86+ "values" : fittedvalues ,
87+ "values_manual" : fittedvalues_manual ,
6988 "omega" : blp .blp_omega ,
70- "omega_manual" : blp_manual . cov_params (). to_numpy () ,
89+ "omega_manual" : omega_manual ,
7190 "basis" : blp .basis ,
7291 "signal" : blp .orth_signal ,
7392 "ci_1" : ci_1 ,
@@ -123,49 +142,11 @@ def test_dml_blp_defaults():
123142 blp = dml .DoubleMLBLP (random_signal , random_basis )
124143 blp .fit ()
125144
126- assert np .allclose (blp .blp_omega , blp .blp_model .cov_HC0 , rtol = 1e-9 , atol = 1e-4 )
145+ assert np .allclose (blp .blp_omega [:, :, 0 ], blp .blp_model [ 0 ] .cov_HC0 , rtol = 1e-9 , atol = 1e-4 )
127146
128147 assert blp ._is_gate is False
129148
130149
131- @pytest .mark .ci
132- def test_dml_blp_multi_rep ():
133- n = 50
134- n_rep = 3
135- level = 0.9
136- np .random .seed (42 )
137- random_basis = pd .DataFrame (np .random .normal (0 , 1 , size = (n , 3 )))
138- random_signal = np .random .normal (0 , 1 , size = (n , n_rep ))
139-
140- blp = dml .DoubleMLBLP (random_signal , random_basis )
141- blp .fit ()
142-
143- expected_coef , expected_se = _aggregate_coefs_and_ses (blp .all_coef , blp .all_se )
144-
145- assert blp .n_rep == n_rep
146- assert blp .all_coef .shape == (random_basis .shape [1 ], n_rep )
147- assert blp .all_se .shape == (random_basis .shape [1 ], n_rep )
148- assert blp .coef .shape == (random_basis .shape [1 ],)
149- assert blp .se .shape == (random_basis .shape [1 ],)
150- assert blp .blp_omega .shape == (random_basis .shape [1 ], random_basis .shape [1 ], n_rep )
151-
152- assert np .allclose (blp .coef , expected_coef , rtol = 1e-9 , atol = 1e-4 )
153- assert np .allclose (blp .se , expected_se , rtol = 1e-9 , atol = 1e-4 )
154-
155- ci = blp .confint (random_basis )
156- assert isinstance (ci , pd .DataFrame )
157- assert ci .shape [0 ] == n
158-
159- ci_coef = blp .confint (level = level )
160- critical_value = norm .ppf (1 - (1 - level ) / 2 )
161- expected_ci_lower = blp .coef - critical_value * blp .se
162- expected_ci_upper = blp .coef + critical_value * blp .se
163- expected_ci_coef = np .vstack ((expected_ci_lower , blp .coef , expected_ci_upper )).T
164- assert np .allclose (ci_coef .to_numpy (), expected_ci_coef , rtol = 1e-9 , atol = 1e-4 )
165-
166- assert isinstance (blp .summary , pd .DataFrame )
167-
168-
169150@pytest .mark .ci
170151def test_doubleml_exception_blp ():
171152 random_basis = pd .DataFrame (np .random .normal (0 , 1 , size = (2 , 3 )))
0 commit comments