32
$\begingroup$

The following code evaluates the similarity between two time series:

set.seed(10) RandData <- rnorm(8760*2) America <- rep(c('NewYork','Miami'),each=8760) Date = seq(from=as.POSIXct("1991-01-01 00:00"), to=as.POSIXct("1991-12-31 23:00"), length=8760) DatNew <- data.frame(Loc = America, Doy = as.numeric(format(Date,format = "%j")), Tod = as.numeric(format(Date,format = "%H")), Temp = RandData, DecTime = rep(seq(1, length(RandData)/2) / (length(RandData)/2), 2)) require(mgcv) mod1 <- gam(Temp ~ Loc + s(Doy) + s(Doy,by = Loc) + s(Tod) + s(Tod,by = Loc),data = DatNew, method = "ML")

Here, gam is used to evaluate how the temperature at New York and Miami vary from the mean temperature (of both locations) at different times of the day. The problem that I now have is that I need to include an interaction term which shows how the temperature of each location varies throughout at the day for different days of the year. I eventually hope to display all of this information on one graph (for each location). So, for Miami I hope to have one graph that shows how the temperature varies from the mean during different times of the day and different times of the year (3d plot?)

Improve this question
$\endgroup$
1

2 Answers 2

Reset to default
29
$\begingroup$

The "a" in "gam" stands for "additive" which means no interactions, so if you fit interactions you are really not fitting a gam model any more.

That said, there are ways to get some interaction like terms within the additive terms in a gam, you are already using one of those by using the by argument to s. You could try extending this to having the argument by be a matrix with a function (sin, cos) of doy or tod. You could also just fit smoothing splines in a regular linear model that allows interactions (this does not give the backfitting that gam does, but could still be useful).

You might also look at projection pursuit regression as another fitting tool. Loess or more parametric models (with sin and/or cos) might also be useful.

Part of decision on what tool(s) to use is what question you are trying to answer. Are you just trying to find a model to predict future dates and times? are you trying to test to see if particular predictors are significant in the model? are you trying to understand the shape of the relationship between a predictor and the outcome? Something else?

Improve this answer
$\endgroup$
5
  • 6
    $\begingroup$ Suppose you have two predictors $x_1,x_2$ - isn't $$y = f_1(x_1) + f_2(x_2) + f_3(x_1x_2) + \varepsilon$$ still an gam? This could be thought of as an 'interaction' in some sense. Also, I think the gam package allows you to fit models like $$y = f_1(x_1) + f_2(x_2) + f_3(x_1)x_2 + f_4(x_2)x_1 + \varepsilon$$ (which I think is what the by argument is doing). Is that still a gam? $\endgroup$ Commented Jul 21, 2012 at 17:43
  • 1
    $\begingroup$ @Macro, It depends on if you want to split hairs or not (well technically what you wrote would maybe be am's since you are not really using the g part). $\endgroup$ Commented Jul 21, 2012 at 18:02
  • 2
    $\begingroup$ @Macro, I'm sure more has been said on the subject, but see page 4 on GAMs of this Venable's article, Exegeses on Linear Models. It isn't clear exactly how non-additive main effects and interaction effects are simultaneously identified. $\endgroup$ Commented Jul 22, 2012 at 2:20
  • $\begingroup$ many thanks for your comments. My main aim here is to not predict future values but to simply see how each time series alters from the mean of both series. For example, the outcome of mod1 shows how the time series at each location varies from the mean firstly at the time scale of day of year (Doy) and then time of day (Tod). From this I would like to see how each series varies at a function of Doy and Tod which I expect the time series to differ greatly during the summer period. $\endgroup$ Commented Jul 24, 2012 at 9:40
  • 7
    $\begingroup$ @AndyW It is important to note that GAM's have come a long way since Venable's was commenting on them - developments in the last decade or so with the penalised splines sensu Simon Wood (as implemented in mgcv) and the fully bayesian treatments address issues such as smoothness selection, interactions & how to fit them (tensor products of marginal bases being one approach) in the additive model framework. I'm reasonably confident Venable's & Cox's objections to GAMs as outlined in the former's Exegeses have to large degree been addressed by these recent developments in GAM theory. $\endgroup$ Commented Apr 11, 2017 at 23:24
38
$\begingroup$

For two continuous variables then you can do what you want (whether this is an interaction or not I'll leave others to discuss as per comments to @Greg's Answer) using:

mod1 <- gam(Temp ~ Loc + s(Doy, bs = "cc", k = 5) + s(Doy, bs = "cc", by = Loc, k = 5, m = 1) + s(Tod, bs = "cc", k = 5) + s(Tod, bs = "cc", by = Loc, k = 5, m = 1) + te(Tod, Doy, by = Loc, bs = rep("cc",2)), data = DatNew, method = "ML")

The simpler model should then be nested within the more complex model above. That simpler model is:

mod0 <- gam(Temp ~ Loc + s(Doy, bs = "cc", k = 5) + s(Doy, bs = "cc", by = Loc, k = 5, m = 1) + s(Tod, bs = "cc", k = 5) + s(Tod, bs = "cc", by = Loc, k = 5, m = 1), data = DatNew, method = "ML")

Note two things here:

  1. The basis-type for each smoother is stated. In this case we would expect that there are no discontinuities in Temp between 23:59 and 00:00 for Tod nor between Doy == 1 and Doy == 365.25. Hence cyclic cubic splines are appropriate, indicated here via bs = "cc".
  2. The basis dimension is stated explicitly (k = 5). This matches the default basis dimension for each smooth in a te() term.

Together these features ensure that the simpler model really is nested within the more complex model.

For more see ?gam.models in mgcv.

Improve this answer
$\endgroup$
9
  • $\begingroup$ Related to your 2nd point - in addition to the specification of k, should one also fix the number of knots (e.g. fx=TRUE). If not, the resulting model shows varying edf for each term. $\endgroup$ Commented Aug 3, 2015 at 14:55
  • $\begingroup$ I should update this answer a little given some new functionality in the mgcv package in terms of splines for marginal bases. That said, I don't agree that you need to fix the degrees of freedom for the spline. The key is to ensure that the bases for the models are appropriately nested. Then differences between models are possible by setting some of the coefficients for the basis functions to zero, just as would happen in a linear model with non-spline terms. $\endgroup$ Commented Aug 3, 2015 at 15:02
  • 4
    $\begingroup$ Hoping someone are still watching this thread and may answer. In these models, how come you must specify both s(Doy...) and s(Doy, by =Loc...)? I thought the first one would be nested into the latter and thus be unneccessary to specify? $\endgroup$ Commented Jan 25, 2016 at 13:48
  • 6
    $\begingroup$ No, the first smooth is the global function and the by smooth represent site-specific differences between it and the global smooth. The by smooths really need m = 1 added to them to put the penalty on the first derivative for the difference smooths though. $\endgroup$ Commented Jan 25, 2016 at 14:05
  • 2
    $\begingroup$ @JoshuaRosenberg If you mean te(), it depends what you are including in the tensor product? The interactions described here are factor-smooth interactions, but te() would imply two or more continuous variables. If you want global terms and subject-specific deviations, then yes, te(DoY, Year, by = Loc, m = 1) could be use alongside te(DoY, Year), although there are other ways to achieve similar things using random effect-like factor-smooth interaction and te() terms containing a random effect spline. $\endgroup$ Commented Jan 11, 2017 at 2:05

Start asking to get answers

Find the answer to your question by asking.

Ask question

Explore related questions

See similar questions with these tags.