Why are all my p-values so low in logistic regression model?

Question

I'm trying to make a regression model to explain a dependent variable that follows a binomial distribution - I have data on the number of successes and the number of trials for each observation. The proportion of successes to trials ranges from about 0.2 to 0.97. I have ~30 independent variables I am choosing from for my model, some of which are continuous, some are discrete, some are percentages.

In researching how best to model this situation, I thought it made sense to use logistic regression. I am working in MATLAB and have implemented this using a generalized linear model for a binomial distribution with a logit link function. From what I can tell the fact that my independent variables vary in data type is okay.

I have tried doing this for many different combinations of independent variables, and my results keep showing that the p-values for all of the independent variables are 0 (out to many decimal places). Note also that SSE is quite high.

Generalized Linear regression model:
    logit(FC) ~ 1 + DDP + SCHabs + Train + Exp_hab + Exp_ratio + Sch_tot
    Distribution = Binomial

Estimated Coefficients:
                Estimate         SE        tStat           pValue   
               __________    __________    _______    ___________

(Intercept)       -0.4451     0.0060914    -73.071              0
DDP             0.0036999    8.0929e-05     45.718              0
SCHabs             1.7664      0.025142     70.259              0
Train          2.0887e-05    9.8255e-08     212.58              0
Exp_hab           -6.8407       0.25137    -27.214    4.4104e-163
Exp_ratio         -0.1891      0.016817    -11.245     2.4572e-29
Sch_tot        2.9644e-06    1.5103e-08     196.28              0

30 observations, 23 error degrees of freedom
Dispersion: 1
Chi^2-statistic vs. constant model: 3.05e+05, p-value = 0
SSE = 2.74 e+09

EDIT: Here is the same output when I run the regression in R:

Call:
glm(formula = cbind(FC, NotFC) ~ DDP + SCHabs + Train + Exp_hab + 
    Exp_ratio + Sch_tot, family = binomial, data = regdata)

Deviance Residuals: 
     Min        1Q    Median        3Q       Max  
-185.628   -28.537    -1.076    45.945   244.305  

Coefficients:
              Estimate Std. Error z value Pr(>|z|)    
(Intercept) -4.453e-01  6.089e-03  -73.13   <2e-16 ***
DDP          3.586e-03  8.056e-05   44.52   <2e-16 ***
SCHabs       1.741e+00  2.470e-02   70.49   <2e-16 ***
Train        2.091e-05  9.833e-08  212.65   <2e-16 ***
Exp_hab     -6.916e+00  2.511e-01  -27.54   <2e-16 ***
Exp_ratio   -1.761e-01  1.672e-02  -10.53   <2e-16 ***
Sch_tot      2.960e-06  1.508e-08  196.30   <2e-16 *** 
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

(Dispersion parameter for binomial family taken to be 1)  

Null deviance: 489866  on 29  degrees of freedom
Residual deviance: 184863  on 23  degrees of freedom
AIC: 185168

Number of Fisher Scoring iterations: 5

At first I thought maybe it was a multicollinearity problem, but the standard errors are low, none of the values in the correlation matrix are above 0.65, and a collintest did not flag any issues.

Next I thought maybe I have too many independent variables in my model, so I tried using just a single independent variable. I tried several of them, no improvement. I then tried using stepwise regression to minimize the SSE. This added in all of the variables I included in the data table, as well as many of the cross terms. Only 1 p value was insignificant, which still seems absurd.

Generalized Linear regression model:
FC ~ [Linear formula with 29 terms in 10 predictors]
Distribution = Binomial

Estimated Coefficients:
                          Estimate        SE           tStat        pValue   
                         ___________    __________    ________    ___________

(Intercept)               6.7475       0.94871      7.1123     1.1413e-12
DDP                      0.77076       0.11393      6.7654     1.3296e-11
LWE                      -1.7078       0.05188     -32.919    1.1664e-237
MinBlocks              -0.053431      0.005388     -9.9167     3.5224e-23
STHabs                   -7.3881        1.0538     -7.0107     2.3707e-12
SCHabs                    5.4708        0.7251      7.5449     4.5264e-14
Train                -0.00074337    0.00012089     -6.1494       7.78e-10
Exp_hab                  -705.19         85.84     -8.2152     2.1178e-16
Exp_ratio                -4.3504       0.88336     -4.9249     8.4421e-07
Sch_tot              -8.5885e-07    1.4783e-06    -0.58097        0.56126
Nonfuncp                  122.31        2.6776      45.678              0
DDP:LWE                  0.13659     0.0097325      14.034     9.6492e-45
DDP:MinBlocks          0.0074521     0.0005404       13.79     2.9306e-43
DDP:STHabs              0.080665      0.030503      2.6445      0.0081805
DDP:SCHabs              -0.57123       0.11228     -5.0873     3.6318e-07
DDP:Train            -2.1825e-06    2.8467e-07     -7.6667     1.7648e-14
DDP:Exp_hab              -67.035        10.265     -6.5307     6.5451e-11
LWE:MinBlocks           0.037295      0.001902      19.608     1.3131e-85
LWE:STHabs               0.60362      0.085094      7.0936      1.307e-12
LWE:SCHabs                -11.87       0.77467     -15.322     5.4435e-53
LWE:Train             3.0939e-05     4.983e-06       6.209      5.333e-10
LWE:Exp_hab              -112.45        22.495     -4.9988     5.7686e-07
LWE:Sch_tot           4.8253e-06    1.1147e-07      43.289              0
MinBlocks:Exp_hab         1.1291        0.2784      4.0556     5.0011e-05
STHabs:Train          0.00060167    6.8549e-05      8.7772     1.6763e-18
STHabs:Exp_hab            724.38        83.805      8.6437      5.442e-18
STHabs:Nonfuncp          -357.31        8.4226     -42.423              0
Train:Exp_hab           0.051368     0.0083675      6.1389     8.3074e-10
Exp_hab:Exp_ratio          769.8        79.383      9.6974      3.094e-22


30 observations, 1 error degrees of freedom
Dispersion: 1
Chi^2-statistic vs. constant model: 4.9e+05, p-value = 0
SSE = 637  

Finally, I decided to simulate a meaningless independent variable by drawing from the normal distribution and adding it to the regression to see if it was also significant. I did this many times, and every time this meaningless variable also had a p-value of zero, which obviously can't be right.

What's going on here? What am I doing wrong?

EDIT: Here are scatter plots of the key variables from the regression above. A few of the discrete variables, like DDP and LWE, have lots of zeros in them. Which are meaningful, real zeros. I thought about maybe coding them as binary variables instead. I also include a correlation matrix plot.

EDIT AGAIN: Here is a public dropbox link to a .csv file with my data. Note the dependent variable I am modeling is "FCperHab" which is calculated by dividing "FC" by "Hab". For context, these are data on states in India. I am trying to explain the percentage of habitations in each state that have full access to drinking water. https://www.dropbox.com/s/w5gzdznhuad59gp/regdata_CV.csv?dl=0

Answer 1

As pointed out in comments, the more complicated model has too many predictors to be taken seriously. I focus here on models with 6 predictors.

I have used MATLAB very occasionally in the past but not for any related purpose and am emphatically no expert. But on the face of it your MATLAB call does not feed a fractional response (outcome, dependent variable) to the function either directly or indirectly and it is amazing that it produces any output at all. If you are using the same naming conventions across different software, then FC is a count, not a binary response, and not a fit present for a logit function.

Here is the result of a calculation in Stata 14 of a logit model in which the response is treated as a continuous proportion:

. fracreg logit fcperhab ddp schabs train exp_hab  exp_ratio sch_tot

Iteration 0:   log pseudolikelihood = -19.110279  
Iteration 1:   log pseudolikelihood =  -18.06008  
Iteration 2:   log pseudolikelihood = -17.979059  
Iteration 3:   log pseudolikelihood = -17.976617  
Iteration 4:   log pseudolikelihood = -17.976616  

Fractional logistic regression                  Number of obs     =         30
                                                Wald chi2(6)      =      21.06
                                                Prob > chi2       =     0.0018
Log pseudolikelihood = -17.976616               Pseudo R2         =     0.0865

------------------------------------------------------------------------------
             |               Robust
    fcperhab |      Coef.   Std. Err.      z    P>|z|     [95% Conf. Interval]
-------------+----------------------------------------------------------------
         ddp |   .0033372   .0044366     0.75   0.452    -.0053583    .0120327
      schabs |   .0195366   .0162846     1.20   0.230    -.0123806    .0514538
       train |   .0000243   8.40e-06     2.90   0.004     7.85e-06    .0000408
     exp_hab |  -4.832934   10.20002    -0.47   0.636     -24.8246    15.15873
   exp_ratio |  -.0016535    .904611    -0.00   0.999    -1.774658    1.771351
     sch_tot |   1.15e-06   9.91e-07     1.16   0.248    -7.98e-07    3.09e-06
       _cons |  -.1102337   .3944965    -0.28   0.780    -.8834326    .6629653
------------------------------------------------------------------------------

For all that the results may seem disappointing, there is nothing pathological about the output. Non-Stata users can probably guess that _cons means the intercept. As context here are some summary statistics:

. su fcperhab ddp schabs train exp_hab  exp_ratio sch_tot

    Variable |        Obs        Mean    Std. Dev.       Min        Max
-------------+---------------------------------------------------------
    fcperhab |         30    .6355612    .2472515   .2056616   .9978048
         ddp |         30    7.766667    19.39016          0         85
      schabs |         30        8.83    10.14478          0       40.2
       train |         30    18741.83    26769.19          0     108643
     exp_hab |         30    .0152065    .0195379          0   .1046035
-------------+---------------------------------------------------------
   exp_ratio |         30    .3794279    .2442145          0   .7318189
     sch_tot |         30    203229.1    410059.5        251    2198181

Stata users who have not yet upgraded to 14 [released 7 April 2015] may note that glm fcperhab ddp schabs train exp_hab exp_ratio sch_tot, link(logit) f(binomial) vce(robust) gives the same calculation. But it is important here to spell out to Stata, in whatever version is used, that binomial is at best a convenient fiction here. In fracreg that is automatic; otherwise, vce(robust) is that signal and if we omit it results are quite different, but still free of pathological very high or very low P-values:

. glm fcperhab ddp schabs train exp_hab  exp_ratio sch_tot, link(logit) f(binomial)
note: fcperhab has noninteger values

Iteration 0:   log likelihood = -13.058365  
Iteration 1:   log likelihood = -12.952102  
Iteration 2:   log likelihood = -12.950226  
Iteration 3:   log likelihood = -12.950225  
Iteration 4:   log likelihood = -12.950225  

Generalized linear models                         No. of obs      =         30
Optimization     : ML                             Residual df     =         23
                                                  Scale parameter =          1
Deviance         =  5.476451096                   (1/df) Deviance =   .2381066
Pearson          =   4.86766495                   (1/df) Pearson  =   .2116376

Variance function: V(u) = u*(1-u/1)               [Binomial]
Link function    : g(u) = ln(u/(1-u))             [Logit]

                                                  AIC             =   1.330015
Log likelihood   = -12.95022455                   BIC             =  -72.75109

------------------------------------------------------------------------------
             |                 OIM
    fcperhab |      Coef.   Std. Err.      z    P>|z|     [95% Conf. Interval]
-------------+----------------------------------------------------------------
         ddp |   .0033372    .023125     0.14   0.885    -.0419869    .0486613
      schabs |   .0195366   .0430982     0.45   0.650    -.0649343    .1040075
       train |   .0000243   .0000275     0.88   0.377    -.0000297    .0000783
     exp_hab |  -4.832937   27.53167    -0.18   0.861    -58.79402    49.12815
   exp_ratio |  -.0016538   1.961583    -0.00   0.999    -3.846287    3.842979
     sch_tot |   1.15e-06   2.69e-06     0.43   0.671    -4.13e-06    6.43e-06
       _cons |  -.1102337    .785408    -0.14   0.888    -1.649605    1.429138
------------------------------------------------------------------------------

There are many other scientific and statistical issues not clear without further discussion and analysis, and I will raise just two:

1. Some predictors appear to be absolute counts or amounts, and it's not all clear why absolute values are natural here.

2. Some of the predictors may need or benefit from transformation.