10
$\begingroup$

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.
Plot with monthly totals and regression line

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.

$\endgroup$

1 Answer 1

9
$\begingroup$

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.

$\endgroup$
2
  • 1
    $\begingroup$ The lubridate package is designed to facilitate working with dates. See the vignette here, and note the x-axis labels in the plot. However, your observation that months might be internally stored as seconds led me to look more carefully at the underlying class, which is the R base POSIXct, and which does indeed use seconds from 1970-01-01 00:00:00 UTC. So the multiplier on the coefficient did the job. $\endgroup$
    – Dan Hicks
    Mar 17, 2016 at 17:55
  • 2
    $\begingroup$ That's right: but "working with dates" is not the same thing as representing dates numerically for statistical calculations. Those aims are so profoundly different that they often require totally different ways of representing dates. For instance, I have seen well-known R packages for maximum likelihood regression models utterly fail when dates are regressors (because the software does not internally standardize the values). For dates past the year 2000 or so this can create a particularly toxic error in which most, but not all, precision is lost--and so goes undetected. $\endgroup$
    – whuber
    Mar 17, 2016 at 18:19

Your Answer

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

Not the answer you're looking for? Browse other questions tagged or ask your own question.