R help - Apr 2020 - Changes to stats::glm function between R versions 3.4.0 and 3.5.1

Hi all,

Does anyone know whether there was a change to the algorithm of the glm function
between versions 3.4.0 and 3.5.1 of the stats package? I noticed the
introduction of the 'singular.ok' option, but I'm seeing more
fundamental differences in the output of Generalised Linear Models between the
two versions, particularly when the models don't converge.

In the later version, I'm seeing more variables 'blowing up' and
giving large or NA standard error values when a model doesn't converge, but
I'm using the same value of 'maxit' for both versions.

The numerical precision seems to be the same in versions 3.4.0 and 3.5.1 of R,
as far as I can tell, but perhaps there is some difference that is indirectly
affecting glm? Alternatively, there is a C function named Cdqrls that is called
by glm, and I wondered if this had changed?

I have rather limited control over the version of R that I use, so I'm
hoping I can produce results with 3.5.1 that are as similar as possible to those
of 3.4.0.

Many thanks,

Mark Purver
Statistician, UK Ministry of Justice
________________________________
 This e-mail and any attachments is intended only for the attention of the
addressee(s). Its unauthorised use, disclosure, storage or copying is not
permitted. If you are not the intended recipient, please destroy all copies and
inform the sender by return e-mail. Internet e-mail is not a secure medium. Any
reply to this message could be intercepted and read by someone else. Please bear
that in mind when deciding whether to send material in response to this message
by e-mail. This e-mail (whether you are the sender or the recipient) may be
monitored, recorded and retained by the Ministry of Justice. Monitoring /
blocking software may be used, and e-mail content may be read at any time. You
have a responsibility to ensure laws are not broken when composing or forwarding
e-mails and their contents.

>>>>> Purver, Mark 
>>>>>     on Thu, 16 Apr 2020 10:02:09 +0000 writes:
    > Hi all,
    > Does anyone know whether there was a change to the algorithm of the glm
function between versions 3.4.0 and 3.5.1 of the stats package? I noticed the
introduction of the 'singular.ok' option, but I'm seeing more
fundamental differences in the output of Generalised Linear Models between the
two versions, particularly when the models don't converge.

    > In the later version, I'm seeing more variables 'blowing
up' and giving large or NA standard error values when a model doesn't
converge, but I'm using the same value of 'maxit' for both versions.

    > The numerical precision seems to be the same in versions 3.4.0 and
3.5.1 of R, as far as I can tell, but perhaps there is some difference that is
indirectly affecting glm? Alternatively, there is a C function named Cdqrls that
is called by glm, and I wondered if this had changed?

    > I have rather limited control over the version of R that I use, so
I'm hoping I can produce results with 3.5.1 that are as similar as possible
to those of 3.4.0.

    > Many thanks,

    > Mark Purver
    > Statistician, UK Ministry of Justice
    > ________________________________

Dear Mark,
is there any chance you can add a reproducible example to what
you claim above?

In general: Yes, there are changes between R versions, but of
course they should be improvements (aka "bug fixes").  So, for
this case, we really need a "repr.ex.", ideally as small as
possible.

OTOH, both 3.4.x  and  3.5.y  are "very old" in our eyes,
and at the moment, the smallest R version we could possibly
change is 4.0.0, but almost surely your problem would not be
both grave and simple to fix enough, and so a change would be
for R 4.0.1

Best regards,
Martin Maechler
ETH Zurich  and  R Core team

Dear Martin,

Thanks so much for this response. I have attempted a reproducible example below.
The 'real' issue involves a much larger data set, but the very small
example below produces very small differences in the standard errors between
different R versions, which may have the same root cause. I hope to understand
whether such differences indicate any change to the 'glm' function
between version 3.4.0 of the 'stats' package and later versions. These
differences are only apparent when the model does not converge.

Thanks,
Mark


dat <- data.frame(
  a = c(1, 3, 8, 9, 2, 3, 4, 5),
  b = c(1, 3, 8, 8, 2, 3, 4, 5),
  x = c(1, 0, 0, 1, 0, 0, 0, 0)
)

m <- glm(x~a+b, family=binomial(), data=dat)

print(summary(m)))

# Standard errors of variables 'a' and 'b' are 580356.74 and
638232.37 in R 3.4.0, but 580356.25 and 638231.93 in R 3.4.2
________________________________
 This e-mail and any attachments is intended only for the attention of the
addressee(s). Its unauthorised use, disclosure, storage or copying is not
permitted. If you are not the intended recipient, please destroy all copies and
inform the sender by return e-mail. Internet e-mail is not a secure medium. Any
reply to this message could be intercepted and read by someone else. Please bear
that in mind when deciding whether to send material in response to this message
by e-mail. This e-mail (whether you are the sender or the recipient) may be
monitored, recorded and retained by the Ministry of Justice. Monitoring /
blocking software may be used, and e-mail content may be read at any time. You
have a responsibility to ensure laws are not broken when composing or forwarding
e-mails and their contents.