R: regression coefficients and lubridate

Question

I'm getting some odd coefficients when I apply lm to dates that have been processed and rounded using the lubridate package. MWE:

library(ggplot2)
library(lubridate)
library(dplyr)

lakers$month <- ymd(lakers$date) %>% round_date(unit = 'month')
items_by_month <- lakers %>% group_by(month) %>% summarize(count = n()) %>%
    mutate(count = count / 1000)

ggplot(data = items_by_month, aes(x = month, y = count)) + 
    geom_line() +
    stat_smooth(method = 'lm', data = items_by_month)

model <- lm(data = items_by_month, count ~ month)
summary(model)
time <- max(items_by_month$month) - min(items_by_month$month)
coef(model)['month'] * as.numeric(time)

The plot indicates that ggplot, at least, understands what's going on with the regression model.

But in summary(model) the coefficient on month is on the order of 10^-7, which is about 5 orders of magnitude too small: the plot shows an increase of about 2.5 between the first and last dates, but the last line shows an increase of about 2.5 * 10^-5.

Note that I've divided the count column by 10^3, in order to get values that are easier to read (and closer to my actual use case). But that shouldn't effect either the plot or lm. Also, I know there are more sophisticated techniques than linear regression for analyzing time series data; but I'm just looking at gross trends over time, not factor out seasonal patterns, etc.

Answer 1

Possible solution

It would help if you were to report specific quantities in the question, but even so one can make reasonable guesses. My eye says the slope of the fitted line is around $2$ per 5 months, or $5$ per year. If your output is "on the order of" $10^{-7}$, that means it is around $5/10^{-7}$ times what you expected. That is close to the number of seconds per year (equal to $10^7\pi$ to a good approximation), suggesting that the internal numerical value of your "month" variable is in terms of seconds rather than years. All you need to do, therefore, is convert it back from a rate of change per second into a rate of change per year. The conversion factor is approximately

$$60\text{ seconds/minute}\times 60\text{ minutes/hour}\times 24\text{ hours/day}\times 365.2422\text{ days/year} = 3.1556926 \times 10^7\text{ seconds/year}.$$

General advice and comments

I have provided this answer, rather than migrating the question to StackOverflow, because such problems with dates are common: they occur on almost every computing platfrom, from Excel through R. Experience with a large number of platforms suggests some simple principles to follow:

1. Use a system's internal date datatype to store dates, perform date-specific manipulation (such as finding the day of the week, etc.), and to produce good labels on graphics.

2. For statistical analysis, circumvent the system's default by computing a numeric equivalent of the dates. It is often best to establish a project-specific date origin and to represent all dates in terms of days, months, or years from that particular origin. This achieves several important things:

   • You will not make mistakes concerning the units in which dates are represented.
   • Your statistical output, such as regression coefficients, will be readily interpretable.
   • Your calculations will tend to be more numerically stable, because they will involve numbers of reasonable size. (R's internal date values, which are in seconds since the end of 1969, are in the many billions: even in double-precision computation, the sums of squares involved in many statistical procedures cause catastrophic loss of precision. See https://stats.stackexchange.com/a/318516/919 for a discussion.)

3. Allow an exception to rule (2) when working with time series procedures that "understand" how to cope with varying-length months, how annual, monthly, and weekly seasons work, etc.