Sinusoidal unit-specific time trends

Question

Suppose I have a panel dataset with monthly observations over 10 years. I have a simple dummy intervention, where some policy is put in place around the Spring in every year and only affects some people or groups. In other words, the units over time experience multiple shocks. Examples include something like emergency room re-entry or a crime policy affecting districts during warm weather. Note: the dummy is turning 'on' and 'off' over time and there is a seasonal component observed in the raw trends over time.

Let's assume I'm observing aggregate hospital admissions over time. To include hospital-specific linear time trends, I would multiply each hospital dummy by a continuous linear time trend variable, where $t = 1, 2, 3, ... , T$. To assess the robustness of the findings, it is common to estimate the following equation:

$$ y_{it} = \alpha_{0i} + \beta_i (t \times \alpha_{1i}) + \lambda_t + \delta D_{it} + u_{it}, $$

where $D_{it}$ is a treatment dummy equal to 1 if the facility was treated and is in a post-treatment time period, 0 otherwise. The parameters $\alpha_{0,i}$ and $\lambda_t$ denote hospital and month fixed effects, respectively. The interaction of the continuous time trend variable $t$ with each $\alpha_{1,i}$ gives each $i$ it's own unique linear time trend. This works well and is quite common in the literature with one intervention. But suppose the policy is removed and then reinstated in subsequent years. In other words, the intervention repeats every year around Spring time.

Here is what I am observing in practice. The toy dataset shows visits to a facility over 24 months. Every 12 months I see this inverted U-shape. The shaded regions represent the shocks observed over time.

# Loading the required packages

library(tidyr)
library(dplyr)

# Creating some fake data

set.seed(1987)

y <- c(8.4, 10.0, 11.8, 12.2, 13.1, 13.3, 12.0, 12.4, 12.0, 10.3, 10.0, 10.5, 9.3, 8.0, 8.1, 10.1 , 11.5, 12.1, 12.5, 12.1, 10.7, 8.8, 8.7, 7.0)
shocks <- c(5, 6, 18, 19)  # on/off shocks

df <- 
  tibble(unit = rep(1:10, each = 24),
         month = rep(1:24, 10),
         y = rep(y, 10) + rnorm(240, 1, 5)
         ) %>%
  mutate(group = ifelse(unit > 5, "Exposed", "Unexposed"),
         intervention = ifelse(unit > 5 & month %in% shocks, 1, 0),
         time = rep(1:10, 24),       # linear time trend
         sine = sin(2*pi*time/12),   # sin()
         cosine = cos(2*pi*time/12)  # cos()
         )

# Producing a fitted trend line for each group
# Shaded regions show transient shocks over time
# Some units become exposed to a policy and others do not

df %>%
  ggplot(., aes(x = factor(month), y = y, color = factor(group))) +
  geom_smooth(aes(group = group)) +
  labs(x = "Month (24 Periods)", y = "Mean Outcome") +
  annotate("rect", xmin = 5, xmax = 7, ymin = -Inf, ymax = Inf, fill = "maroon", alpha = .2) +
  annotate("rect", xmin = 18, xmax = 20, ymin = -Inf, ymax = Inf, fill = "maroon", alpha = .2) +
  scale_color_manual(name = "Group:", values = c("Exposed" = "maroon",
                                "Unexposed" = "royalblue")) +
  theme_classic() +
  theme(
    legend.position = "bottom"
    )

Questions:

(1) Is it valid to allow each unit to have their own sinusoidal trend, which is allowed to flexibly model the cyclical pattern over time between the two groups?

(2) Is there a better approach? I don't see how adding a quadratic time trend would help given that the intervention repeats. I would argue that including a higher order polynomial term would be too demanding. Or maybe not?

To be clear, this is more of a robustness check—not the main specification. I don't think it has been done before, in part because most quasi-experimental evaluations involve units with one treatment history. In applied work, the linear time trend is a way to adjust for the possibility that the treatment group and the control group were on somewhat different growth trajectories before the shock. But with multiple interventions, I wanted to adjust for any differential "cyclical" patterns observed over time, as they seem to repeat every year with the reintroduction of the policy.

Please review my R code below. I am somewhat new to including sinusoidal trends. I welcome any insight or criticism of this approach.

# I multiplied each unit by sin() and cos() separately

> summary(mod <- lm(y ~ as.factor(unit)*sine + as.factor(unit)*cosine + as.factor(month) + intervention, data = df))

Call:
lm(formula = y ~ as.factor(unit) * sine + as.factor(unit) * cosine + 
    as.factor(month) + intervention, data = df)

Residuals:
     Min       1Q   Median       3Q      Max 
-11.4567  -3.1386  -0.0271   3.0985  10.8461 

Coefficients:
                         Estimate Std. Error t value Pr(>|t|)    
(Intercept)               10.8550     2.8408   3.821 0.000181 ***
as.factor(unit)2          -0.6000     1.9757  -0.304 0.761700    
as.factor(unit)3           0.8304     1.7146   0.484 0.628729    
as.factor(unit)4          -1.5308     1.7592  -0.870 0.385329    
as.factor(unit)5           1.4748     1.9535   0.755 0.451217    
as.factor(unit)6           1.9856     1.6350   1.214 0.226108    
as.factor(unit)7           2.5895     2.0160   1.284 0.200584    
as.factor(unit)8           1.2193     1.7436   0.699 0.485217    
as.factor(unit)9           0.5477     1.8026   0.304 0.761602    
as.factor(unit)10          3.3756     1.9880   1.698 0.091192 .  
sine                       2.6086     2.8070   0.929 0.353920    
cosine                    -4.1625     3.3955  -1.226 0.221802    
as.factor(month)2          0.2706     2.7212   0.099 0.920883    
as.factor(month)3          3.8197     3.8195   1.000 0.318591    
as.factor(month)4         -3.7247     4.6494  -0.801 0.424083    
as.factor(month)5         -1.4893     5.5054  -0.271 0.787066    
as.factor(month)6         -1.3394     5.7128  -0.234 0.814891    
as.factor(month)7          2.3822     5.3802   0.443 0.658440    
as.factor(month)8          2.3526     4.8477   0.485 0.628043    
as.factor(month)9          3.2434     3.7984   0.854 0.394272    
as.factor(month)10         3.3325     2.9592   1.126 0.261546    
as.factor(month)11         2.8721     2.3530   1.221 0.223786    
as.factor(month)12         0.3874     2.7212   0.142 0.886934    
as.factor(month)13        -0.7625     3.8195  -0.200 0.841997    
as.factor(month)14        -7.9720     4.6494  -1.715 0.088077 .  
as.factor(month)15        -7.7604     5.4126  -1.434 0.153318    
as.factor(month)16        -2.0334     5.6234  -0.362 0.718065    
as.factor(month)17        -1.1144     5.3802  -0.207 0.836129    
as.factor(month)18         2.9657     4.9512   0.599 0.549912    
as.factor(month)19         1.2505     3.9296   0.318 0.750660    
as.factor(month)20         4.6338     2.9592   1.566 0.119072    
as.factor(month)21         1.9352     2.3530   0.822 0.411883    
as.factor(month)22        -3.5087     2.7212  -1.289 0.198861    
as.factor(month)23        -4.7884     3.8195  -1.254 0.211542    
as.factor(month)24        -6.3957     4.6494  -1.376 0.170601    
intervention              -0.7227     2.0136  -0.359 0.720054    
as.factor(unit)2:sine     -2.5287     5.0387  -0.502 0.616358    
as.factor(unit)3:sine     -0.1849     3.4525  -0.054 0.957344    
as.factor(unit)4:sine     -0.7943     3.5224  -0.226 0.821831    
as.factor(unit)5:sine     -2.1947     4.9591  -0.443 0.658593    
as.factor(unit)6:sine     -1.7032     2.0679  -0.824 0.411209    
as.factor(unit)7:sine     -6.4227     5.0401  -1.274 0.204138    
as.factor(unit)8:sine     -0.3358     3.4440  -0.098 0.922423    
as.factor(unit)9:sine     -4.5655     3.5362  -1.291 0.198271    
as.factor(unit)10:sine    -2.6118     4.9475  -0.528 0.598196    
as.factor(unit)2:cosine    1.4540     6.0715   0.239 0.810997    
as.factor(unit)3:cosine    6.8385     4.2671   1.603 0.110720    
as.factor(unit)4:cosine   -1.1452     4.1782  -0.274 0.784312    
as.factor(unit)5:cosine    7.7818     6.1802   1.259 0.209551    
as.factor(unit)6:cosine    1.2134     2.4830   0.489 0.625636    
as.factor(unit)7:cosine    5.3385     6.0575   0.881 0.379291    
as.factor(unit)8:cosine    3.0046     4.2863   0.701 0.484185    
as.factor(unit)9:cosine    0.2422     4.1670   0.058 0.953707    
as.factor(unit)10:cosine   6.4106     6.1764   1.038 0.300653    
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

Residual standard error: 5.262 on 186 degrees of freedom
Multiple R-squared:  0.2956,    Adjusted R-squared:  0.0949 
F-statistic: 1.473 on 53 and 186 DF,  p-value: 0.03175

All that matters is the estimate on intervention, at least for reporting purposes. I'm curious to see if effects hold even after the inclusion of such trends. I understand this is a very demanding procedure. To offer some perspective, I'm working with about 60 units and each unit is observed over 120 months.

Answer 1

As far as I understand the issue, your question boils down to how to appropriately model a non-linear time trend as the repeating U-shaped pattern may confound your treatment effect estimate. This is a valid concern and I agree that flexibly modeling the time trend as a robustness check is a good idea. Here are a few thoughts on a number of possible specifications and their pros and cons:

In general, I would not suggest to use a high-order polynomial time trend. While they are very flexible, high-order polynomials are known to be problematic when estimating treatment effects, see for instance this paper by Andrew Gelman and Gudio Imbens.

You could instead go for simple polynomials, i.e., quadratic terms. However, instead of using one quadratic term that is unit-specific (which does not make sense for repeated U-shaped patterns), you could try a quadratic trend that is unit-year-specific. This accounts for U-shapes where the characteristics of the U-shape can also change over time. However, depending on the specifics of your data set, this might heavily overfit and could well be too much to ask from your data.

You could also stay in your proposed framework where you use sine and cosine based time trends. However, be aware that this comes with its own set of assumptions, most importantly that the intensity of the sine and cosine does not change over time. This means that the characteristics of the U-shape are constant across years. Assuming that this specification will be used in an economics paper, you will need to lay out in detail why you think that the trend in your y-variable behaves like this. This may or may not make sense, depending on the specific context you are looking at.

In case you want to fully dive into advanced time series econometrics, you could look into state space models. Here, you could fit a local level model which basically allows you to fit a smooth time trend function to your data. This can capture much more flexible patterns in your data (e.g., the intensity of the U-shape can also change over time), but will require more tuning. There may also be the "problem" that there could be no R-packages available for estimating state space models with unit-specific local means. All in all, this may be a bit of an overkill for a robustness check.

Long story short, in case the assumption of a sinusoidal trend in your data makes sense in your context, go for it and make sure to argue why you think that this is a sensible assumption.

Answer 2

I would use $x_{ij}$ and $\beta_j$ and ensure that all the hospital dummies are the model, the time dummies, and all their interactions. You are trying to frontload everything into one model, however, and it may work out that a "divide and conquer" approach to break up a large problem in order to solve easier smaller ones may be a better approach. Quite often, lumping everything together when there are so many coefficients, ends up throwing the model off. Thus, the goodness of fit and coefficient significance can be drastically different for individual hospital-specific models -- but you won't know until you run a few. I would also look at use of GEE for the longitudinal regression model, since it doesn't require uniformly distributed times.