Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

Sign up for GitHub

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Jump to bottom

BUG: df in wald_test_terms for single column items #5475

Closed
josef-pkt opened this issue Jan 19, 2019 · 7 comments · Fixed by #5684
Closed

BUG: df in wald_test_terms for single column items #5475

josef-pkt opened this issue Jan 19, 2019 · 7 comments · Fixed by #5684
Labels
comp-base comp-stats prio-high type-bug-wrong serious bugs that silently return incorrect numbers
Milestone
0.10

Comments

@josef-pkt
Copy link
Member

I'm not sure what's going on with the df in wald_test_terms results
I'm working with a GAM example in a notebook, and the df in the results is equal to the number of coefficients and not 1 as it should be for a single constraint.

The constraint matrix (returned in temp) are 1-D for single column constraints and not 2-D as for joint hypothesis. So there might be a incorrect df = constraint.shape[0] or missing atleast_2d.

Why doesn't this show up in the unit tests?

res_bsf.wald_test_terms(skip_single=False, combine_terms=['fuel', 'drive', 'weight', 'hp'])
<class 'statsmodels.stats.contrast.WaldTestResults'>
                    chi2         P>chi2  df constraint
Intercept     677.786660  2.024193e-149             24
fuel[T.gas]    63.823523   1.360790e-15             24
drive[T.fwd]    2.887571   8.926542e-02             24
drive[T.rwd]    1.594575   2.066736e-01             24
weight_s0      13.741523   2.097653e-04             24
weight_s1      26.966268   2.070372e-07             24
weight_s2      50.869767   9.870223e-13             24
weight_s3      59.216023   1.412788e-14             24
weight_s4      64.379100   1.026412e-15             24
weight_s5      66.065909   4.360938e-16             24
weight_s6      66.597275   3.330462e-16             24
weight_s7      69.019581   9.749068e-17             24
weight_s8      66.816994   2.979188e-16             24
weight_s9      71.769472   2.418654e-17             24
weight_s10     72.534338   1.641508e-17             24
hp_s0           6.717964   9.544654e-03             24
hp_s1          11.429164   7.230019e-04             24
hp_s2          22.246778   2.397608e-06             24
hp_s3          28.663073   8.613056e-08             24
hp_s4          43.307451   4.677986e-11             24
hp_s5          42.858347   5.885046e-11             24
hp_s6          35.971865   2.001873e-09             24
hp_s7          18.913675   1.367687e-05             24
hp_s8          21.947246   2.802488e-06             24
fuel           63.823523   1.360790e-15              1
drive           2.982567   2.250835e-01              2
weight         97.733732   5.013118e-16             11
hp             74.837143   1.701728e-12              9
@josef-pkt josef-pkt added prio-high comp-base comp-stats type-bug-wrong serious bugs that silently return incorrect numbers labels Jan 19, 2019
@josef-pkt josef-pkt added this to the 0.10 milestone Jan 19, 2019
@josef-pkt
Copy link
Member Author

Needs investigating to see whether this is the case in general or whether there is something specific to my model. (There should not be anything specific to GAM and my GAM example for this.)

@josef-pkt
Copy link
Member Author

confirmed bug in results table for wald_test_terms

I'm using constraint.shape[0] for df_constraint instead of the value that wald_test returns.
The p-values are correct and the correct test is used because wald_test handles 1-D constraint correctly.
But the returned df_constraint is wrong if constraint is 1-D.

        for name, constraint in constraints + combined_constraints + extra_constraints:
            wt = result.wald_test(constraint)
            row = [wt.statistic.item(), wt.pvalue.item(), constraint.shape[0]]
            if use_t:
                row.append(wt.df_denom)
            res_wald.append(row)
            index.append(name)

unit tests in base test_generic_methods

compare_waldres checks the df and uses df = c.shape[0] if c.ndim == 2 else 1
However, AFAICS there is no test class with a single column variable, i.e. with single, simple constraint.

@josef-pkt
Copy link
Member Author

names for df
I don't like df_num, df_denom because I never remember what they are and which is which

it looks like chi2 has df_constraints in df_denom, and F has it in df_num

vars(res_bsf.wald_test(row, use_f=True))
{'df_denom': 189.13394496903189,
 'df_num': 1.0,
 'dist': <scipy.stats._continuous_distns.f_gen at 0xb73c908>,
 'dist_args': (1.0, 189.13394496903189),
 'distribution': 'F',
 'effect': None,
 'fvalue': array([[ 21.9472461]]),
 'pvalue': array(5.343056530859978e-06),
 'statistic': array([[ 21.9472461]])}

vars(res_bsf.wald_test(row, use_f=False))
{'df_denom': 1.0,
 'dist': <scipy.stats._continuous_distns.chi2_gen at 0xb737080>,
 'dist_args': (),
 'distribution': 'chi2',
 'effect': None,
 'pvalue': array(2.802487582963751e-06),
 'sd': None,
 'statistic': array([[ 21.9472461]]),
 'tvalue': array([[ 21.9472461]])}

wald_test_terms uses df_constraints and additionally df_denom in case of F

wtt3 = res_bsf.wald_test_terms(skip_single=True, combine_terms=['fuel', 'drive', 'weight', 'hp'])
vars(wtt3)
{'df_constraints': array([ 1,  2, 11,  9], dtype=int64),
 'dist_args': None,
 'distribution': 'chi2',
 'pvalues': array([  1.36079010e-15,   2.25083535e-01,   5.01311810e-16,
          1.70172799e-12]),
 'statistic': array([ 63.82352279,   2.98256735,  97.73373168,  74.83714306]),
 'table':         statistic        pvalue  df_constraint
 fuel    63.823523  1.360790e-15              1
 drive    2.982567  2.250835e-01              2
 weight  97.733732  5.013118e-16             11
 hp      74.837143  1.701728e-12              9,
 'temp': [('fuel',
   array([[ 0.,  1.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,
            0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.]])),
  ('drive',
   array([[ 0.,  0.,  1.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,
            0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.],
          [ 0.,  0.,  0.,  1.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,
            0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.]])),
  ('weight',
   array([[ 0.,  0.,  0.,  0.,  1.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,
            0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.],
          [ 0.,  0.,  0.,  0.,  0.,  1.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,
            0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.],
          [ 0.,  0.,  0.,  0.,  0.,  0.,  1.,  0.,  0.,  0.,  0.,  0.,  0.,
            0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.],
          [ 0.,  0.,  0.,  0.,  0.,  0.,  0.,  1.,  0.,  0.,  0.,  0.,  0.,
            0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.],
          [ 0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  1.,  0.,  0.,  0.,  0.,
            0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.],
          [ 0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  1.,  0.,  0.,  0.,
            0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.],
          [ 0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  1.,  0.,  0.,
            0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.],
          [ 0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  1.,  0.,
            0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.],
          [ 0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  1.,
            0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.],
          [ 0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,
            1.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.],
          [ 0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,
            0.,  1.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.]])),
  ('hp',
   array([[ 0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,
            0.,  0.,  1.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.],
          [ 0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,
            0.,  0.,  0.,  1.,  0.,  0.,  0.,  0.,  0.,  0.,  0.],
          [ 0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,
            0.,  0.,  0.,  0.,  1.,  0.,  0.,  0.,  0.,  0.,  0.],
          [ 0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,
            0.,  0.,  0.,  0.,  0.,  1.,  0.,  0.,  0.,  0.,  0.],
          [ 0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,
            0.,  0.,  0.,  0.,  0.,  0.,  1.,  0.,  0.,  0.,  0.],
          [ 0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,
            0.,  0.,  0.,  0.,  0.,  0.,  0.,  1.,  0.,  0.,  0.],
          [ 0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,
            0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  1.,  0.,  0.],
          [ 0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,
            0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  1.,  0.],
          [ 0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,
            0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  1.]]))]}

wald_test_terms does not have a use_f option.
so, cheating


res_bsf._results.use_t = True
wtt3 = res_bsf.wald_test_terms(skip_single=True, combine_terms=['fuel', 'drive', 'weight', 'hp'])
res_bsf._results.use_t = False
vars(wtt3)

{'df_constraints': array([ 1,  2, 11,  9], dtype=int64),
 'df_denom': array([ 189.13394497,  189.13394497,  189.13394497,  189.13394497]),
 'dist_args': None,
 'distribution': 'F',
 'pvalues': array([  1.29954437e-13,   2.27717973e-01,   1.13680464e-12,
          2.07948129e-10]),
 'statistic': array([ 63.82352279,   1.49128368,   8.8848847 ,   8.31523812]),
 'table':         statistic        pvalue  df_constraint    df_denom
 fuel    63.823523  1.299544e-13              1  189.133945
 drive    1.491284  2.277180e-01              2  189.133945
 weight   8.884885  1.136805e-12             11  189.133945
 hp       8.315238  2.079481e-10              9  189.133945,
 'temp': [('fuel',
   array([[ 0.,  1.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,
            0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.]])),
  ('drive',
   array([[ 0.,  0.,  1.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,
            0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.],
          [ 0.,  0.,  0.,  1.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,
            0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.]])),
  ('weight',
   array([[ 0.,  0.,  0.,  0.,  1.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,
            0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.],
          [ 0.,  0.,  0.,  0.,  0.,  1.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,
            0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.],
          [ 0.,  0.,  0.,  0.,  0.,  0.,  1.,  0.,  0.,  0.,  0.,  0.,  0.,
            0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.],
          [ 0.,  0.,  0.,  0.,  0.,  0.,  0.,  1.,  0.,  0.,  0.,  0.,  0.,
            0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.],
          [ 0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  1.,  0.,  0.,  0.,  0.,
            0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.],
          [ 0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  1.,  0.,  0.,  0.,
            0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.],
          [ 0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  1.,  0.,  0.,
            0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.],
          [ 0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  1.,  0.,
            0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.],
          [ 0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  1.,
            0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.],
          [ 0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,
            1.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.],
          [ 0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,
            0.,  1.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.]])),
  ('hp',
   array([[ 0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,
            0.,  0.,  1.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.],
          [ 0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,
            0.,  0.,  0.,  1.,  0.,  0.,  0.,  0.,  0.,  0.,  0.],
          [ 0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,
            0.,  0.,  0.,  0.,  1.,  0.,  0.,  0.,  0.,  0.,  0.],
          [ 0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,
            0.,  0.,  0.,  0.,  0.,  1.,  0.,  0.,  0.,  0.,  0.],
          [ 0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,
            0.,  0.,  0.,  0.,  0.,  0.,  1.,  0.,  0.,  0.,  0.],
          [ 0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,
            0.,  0.,  0.,  0.,  0.,  0.,  0.,  1.,  0.,  0.,  0.],
          [ 0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,
            0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  1.,  0.,  0.],
          [ 0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,
            0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  1.,  0.],
          [ 0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,
            0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  1.]]))]}

​```

@josef-pkt
Copy link
Member Author

quickfix for wald_test_terms: make simple constraint matrix into 2-D array (one 2-D row)

cleaner version:
unify and add/rename df_constraints in wald_test, reuse df_constraint for wald_test_terms return table.

@josef-pkt
Copy link
Member Author

I was using this uncommitted change in the GAM PR #5370 (I will not include it there. PR not recently rebased)
This worked fine in my examples, but there are no unit tests yet.

@@ -1814,11 +1814,11 @@ class LikelihoodModelResults(Results):
                         combined[cname].append(constraint_matrix)
 
                 if skip_single:
                     continue
 
-                constraints.append((name, constraint_matrix))
+                constraints.append((name, np.atleast_2d(constraint_matrix)))
 
             combined_constraints = []
             for cname in combine_terms:
                 combined_constraints.append((cname, np.vstack(combined[cname])))

@josef-pkt josef-pkt mentioned this issue May 9, 2019
Closed
@josef-pkt
Copy link
Member Author

simple example with problem using GLM with data from gam tests

wrong in model created with pandas data

res1a.wald_test_terms(combine_terms=['drive', 'fuel'])
​
<class 'statsmodels.stats.contrast.WaldTestResults'>
                    chi2         P>chi2  df constraint
Intercept     553.340987  2.361532e-122              6
fuel[T.gas]    72.626605   1.566529e-17              6
drive[T.fwd]    2.440157   1.182644e-01              6
drive[T.rwd]    1.063788   3.023527e-01              6
weight         83.153705   7.591496e-20              6
hp             23.361817   1.342164e-06              6
drive           2.718371   2.568699e-01              2
fuel           72.626605   1.566529e-17              1

model created with formulas shows correct df

modf = GLM.from_formula('city_mpg ~ fuel + drive + weight + hp', df_autos)
res1f = modf.fit()  # default is pirls with use_t=False
f
res1f.wald_test_terms(combine_terms=['drive', 'fuel'])
<class 'statsmodels.stats.contrast.WaldTestResults'>
                 chi2         P>chi2  df constraint
Intercept  553.340987  2.361532e-122              1
fuel        72.626605   1.566529e-17              1
drive        2.718371   2.568699e-01              2
weight      83.153705   7.591496e-20              1
hp          23.361817   1.342164e-06              1
drive        2.718371   2.568699e-01              2
fuel        72.626605   1.566529e-17              1

@josef-pkt
Copy link
Member Author

josef-pkt commented May 9, 2019 •
edited
Loading

In the formula path of wald_test_terms all constraint matrices are 2-dim,
so adding atleast_2d to the non-formula path makes this the same

from vars(wtf)

'temp': [('Intercept', array([[ 1.,  0.,  0.,  0.,  0.,  0.]])),
  ('fuel', array([[ 0.,  1.,  0.,  0.,  0.,  0.]])),
  ('drive', array([[ 0.,  0.,  1.,  0.,  0.,  0.],
          [ 0.,  0.,  0.,  1.,  0.,  0.]])),
  ('weight', array([[ 0.,  0.,  0.,  0.,  1.,  0.]])),
  ('hp', array([[ 0.,  0.,  0.,  0.,  0.,  1.]])),
  ('drive', array([[ 0.,  0.,  1.,  0.,  0.,  0.],
          [ 0.,  0.,  0.,  1.,  0.,  0.]])),
  ('fuel', array([[ 0.,  1.,  0.,  0.,  0.,  0.]]))]

josef-pkt added a commit to josef-pkt/statsmodels that referenced this issue May 10, 2019
@josef-pkt
BUG: fix df in summary for single constraint in wald_test_terms closes …
02b9eaa 
…statsmodels#5475
@josef-pkt josef-pkt mentioned this issue May 10, 2019
Merged
4 tasks
josef-pkt added a commit to josef-pkt/statsmodels that referenced this issue May 10, 2019
@josef-pkt
BUG: fix df in summary for single constraint in wald_test_terms closes …
b9c27d9 
…statsmodels#5475
@jbrockmendel jbrockmendel mentioned this issue May 12, 2019
Merged
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment
Assignees
No one assigned
Labels
comp-base comp-stats prio-high type-bug-wrong serious bugs that silently return incorrect numbers
Projects
None yet
Milestone
0.10
Development

Successfully merging a pull request may close this issue.

BUG: fix df in summary for single constraint in wald_test_terms josef-pkt/statsmodels
1 participant
@josef-pkt