> I noticed that the 3.2.1 release cycle is about to start. Is there any
> chance that this fix will make it into the next version of R?
>
> This bug is fairly serious: getting the wrong variance estimate leads to
> the wrong log-likelihood and the wrong AIC, BIC etc, which can and does
> lead to suboptimal model selection. If it's not fixed, this issue will
> affect every student taking our time series course in Fall 2015 (and
> probably lots of other students in other time series courses). When I
> taught time series in Spring 2015, I had to teach students how to work
> around the bug, which wasted class time and shook student confidence in R.
> It'd be great if we didn’t have to deal with this issue next semester.
>
> Again, the fix is trivial:
>
> --- a/src/library/stats/R/arima.R
> +++ b/src/library/stats/R/arima.R
> @@ -211,8 +211,10 @@ arima <- function(x, order = c(0L, 0L, 0L),
> if(fit$rank == 0L) {
> ## Degenerate model. Proceed anyway so as not to break old code
> fit <- lm(x ~ xreg - 1, na.action = na.omit)
> + n.used <- sum(!is.na(resid(fit))) - length(Delta)
> + } else {
> + n.used <- sum(!is.na(resid(fit)))
> }
> - n.used <- sum(!is.na(resid(fit))) - length(Delta)
> init0 <- c(init0, coef(fit))
> ses <- summary(fit)$coefficients[, 2L]
> parscale <- c(parscale, 10 * ses)
>
Yes, such a change *is* small in the source code.
But we have to be sure about its desirability.
In another post about this you mention "REML", and I think we
really are discussing if variance estimates should use a
denominator of 'n' or 'n - p' in this case.
> The patch that introduced the bug (
> https://github.com/wch/r-source/commit/32f633885a903bc422537dc426644f743cc645e0
> ) was designed to change the initialization for the optimization routine.
> The proposed fix leaves the deliberate part of the patch unchanged (it
> preserves the value of "init0").
I can confirm this... a change introduced in R 3.0.2.
I'm about to commit changes ... after also adding a proper
regression test.
Martin Maechler, ETH Zurich
______________________________________________
R-devel@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/r-devel
