Some more information is needed to figure out the best solution here, so I'm simply answering a number of scenarios with example
Modeling the outcome
If the outcome is binary, use
family = binomial(). If it is count data, use
family = poisson() since it's a fixed time interval. You could also consider aggregating binary data to counts and then use Poisson.
I'll just assume binomial from here and on.
Modeling a fixed hospital:time interaction
If you have few hospitals (say, less than five), there will be very little information to infer the random effects hyperparameters. In that case, they may just be modeled as fixed:
fit_full = glm(status ~ time * hospital, data = df, family = binomial())
Then the resulting inference on the
time:hospital interaction could be of interest.
Modeling random hospital:time interaction
One of the only practical implications of modeling a term as random is that it applies shrinkage. That is, data points far from the model's prediction are regarded as partially random fluctuations with the true value being closer to the mean. Read more here.
Using the same model as above, but allowing for random slopes for each hospital:
fit_full = lme4::glmer(status ~ time + (1 + time|hospital), family = binomial())
To test the random slope, you can do a Likelihood Ratio Test (LRT) by comparing to a (nested) model that does not contain this term:
fit_null = lme4::glmer(status ~ time + (1|hospital), family = binomial())
Personally, I have a preference for Bayesian inference and you could use the
brms package which is much like
glmer. As a quick fix, you can also compute a BIC-based Bayes Factor (cf. Wagenmakers et al. (2007):
exp((BIC(fit_full) - BIC(fit_null))/2)
Modeling time series
I know of no other packages that can do the above and model some autocorrelation than
brms (and perhaps
brms may take hours to fit, though. For AR(1), it would be something like:
fit = brms::brm(status ~ time + (1 + time|hospital), data = df, family = bernoulli(), autocor = cor_ar(~1, p = 1))
If your dates only come in integers (2013, 2014, 2015, 2016, 2017), then there is may be too little information to estimate autoregressive coefficient(s) and you may consider leaving it out. You do have a lot of data, so this may now be necessary. Your time variable would need a finer resolution for an autoregressive model to be identifiable.