1
$\begingroup$

I have a series of observations as y_actual and fit a glm model with binomial family using logit link function using glm() of R. The series y_fit gives the model estimate. Below are my data and fit

y_actual = c(0.136334692, 0.37840196, 0.390388878, 0.332306853, 0.015579242, 
    0.004312221, 0.197257627, 0.003849912, 0.102556055, 0.034236727, 
    0.044338145, 0.323566006, 0.002906722, 0.365485532, 0.497487763, 
    0.407797227, 0.125140868, 0.003000614, 0.010255274, 0.085995523, 
    0.250887837, 0.348495864, 0.490360423, 0.00639617, 0.306278475, 
    0.188225809, 0.084915355, 0.248019023, 0.393584319, 0.248291519, 
    0.102701397, 0.003875779, 0.019110771, 0.141684518, 0.161157541, 
    0.152268462, 0.232824583, 0.211531935, 0.49862286, 0.013384138, 
    0.005328948, 0.002979605, 0.103723205, 0.084474023, 0.012880473, 
    0.405629168, 0.001864907, 0.345550022, 0.43418713, 0.008840066, 
    0.990839128, 0.879488777, 0.874808353, 0.990273562, 0.978276244, 
    0.757710621, 0.653928755, 0.947022794, 0.724590969, 0.996377637, 
    0.735576744, 0.684193226, 1, 0.92141491, 0.857815327, 1, 1, 0.585688824, 
    0.995244609, 0.937017936, 0.721010869, 0.950997678, 1, 0.994234475, 
    0.654585458, 0.595042818, 0.847643574, 0.580531937, 1, 0.947160812, 
    0.978015774, 0.528595922, 0.942787953, 1, 0.833256332, 0.542124525, 
    0.662855368, 0.99478029, 0.929398164, 0.674319383, 0.719835478, 
    0.998391873, 0.629446018, 0.945572679, 0.542700111, 0.849182625, 
    0.571305647, 0.650508505, 0.69688673, 0.632782458)
    
y_fit = c(0.5108746, 0.510187, 0.5119261, 0.5099491, 0.5060829, 0.5114042, 
    0.5058282, 0.5099703, 0.5106101, 0.5115868, 0.513446, 0.5115868, 
    0.5120128, 0.5086123, 0.5142434, 0.5139871, 0.5086469, 0.5107573, 
    0.5120035, 0.5092959, 0.5101969, 0.5101844, 0.512246, 0.5102176, 
    0.5104995, 0.5132078, 0.507861, 0.5096159, 0.5058775, 0.5119135, 
    0.5083143, 0.5092259, 0.512577, 0.5070328, 0.510359, 0.5087454, 
    0.5135267, 0.5106125, 0.5146122, 0.5145157, 0.5083908, 0.5098692, 
    0.5073968, 0.5112421, 0.5122689, 0.5075709, 0.5130123, 0.512024, 
    0.5102805, 0.512763, 0.5135497, 0.5100924, 0.514118, 0.5067507, 
    0.5135497, 0.5100924, 0.514118, 0.5067507, 0.5137501, 0.513939, 
    0.5131861, 0.5127353, 0.5074683, 0.5111746, 0.5088329, 0.5140566, 
    0.5101618, 0.5131737, 0.5095153, 0.514056, 0.5110791, 0.5083455, 
    0.5057883, 0.5072669, 0.5134591, 0.506769, 0.5075491, 0.5094409, 
    0.5137069, 0.5111754, 0.5141377, 0.5128511, 0.5144231, 0.5063635, 
    0.5077987, 0.5071685, 0.5076486, 0.5083421, 0.512321, 0.5145854, 
    0.509603, 0.5074661, 0.5107244, 0.5100801, 0.5123697, 0.5088153, 
    0.5131287, 0.5138794, 0.5076785, 0.5104665)

Based on above fit, I estimated quantile residuals as below,

set.seed(1); qError = apply(cbind(y_actual, 
    y_fit), 1, function(Y) 
    qnorm(runif(1, pbinom(Y[1] - 1, 1, Y[2]), 
    pbinom(Y[1], 1, Y[2]))))

Below is the distribution of this error estimates,

enter image description here

This histogram is providing some interesting observations. This distribution is highly negatively skewed and left portion appears tapering off smoothly from the peak at around -0.10, however the right portion falls suddenly, overall it is not a well behaved bell shaped curve.

I also calculated the pearson residuals

pError = apply(cbind(y_actual, y_fit), 1, 
     function(Y) (Y[1] - Y[2]) / 
        sqrt(Y[2] * (1 - Y[2])))

The histogram of this residual is as below,

enter image description here

This histogram reveals completely different story about the model fit. Important difference is that the distribution of pearson residuals appears to be somewhat symmetric.

Also, the response residual i.e. actual - fit is also somewhat symmetrical. Therefore from response resid and pearson resid, I can infer that the model fit does not have any preference on under-estimation or over-estimation., But if I go with inference from quantile residual, I would infer that the model is fairly over-estimation as in this case, distribution of residuals is visibly negatively skewed.

Furthermore, I don't think this issue has anything to do with small sample size (100 observation), because I have taken out a subset of my original data with 5,000 observations, and overall inference remain same.

In this present case, the y_actual represents the success rate so, I basically fitted a logistic regression using glm for ratio data.

Could you please help me on how can I interpret this scenario with the distribution of quantile residuals? I know it is bad fit, but what more can be inferred about the quality of overall model looking at such asymmetrical distribution of quantile residuals?

Any insight and pointer will be very helpful.

$\endgroup$
2
  • $\begingroup$ It sounds like you fitted a binomial regression directly to the proportions y_actual. This doesn't make sense. [Did you get the warning: non-integer #successes in a binomial glm?] If you plot(y_actual, y_fit) you can see that there is nothing to learn from these results, other than what you already know -- this is the wrong model for the data, such as it is. $\endgroup$
    – dipetkov
    Apr 24, 2022 at 20:15
  • 1
    $\begingroup$ There is warning. But this is standard approach to fit logistic type regression when response data is fractional type. Some insights can be obtained from stackoverflow.com/questions/37584715/…. There you would see, myglm <- glm(prate ~ mrate + totemp + age + sole, data = df, family = quasibinomial('logit')) and myglm1 <- glm(prate ~ mrate + totemp + age + sole, data = df, family = binomial) would give same estimate (of-course an warning will pop-up in the second case I you stated) $\endgroup$ Apr 24, 2022 at 22:17

0

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge you have read our privacy policy.