Traditionally, both researchers and practitioners rely on standard Erlang queueing models to analyze call center operations. In those models, service times are assumed to be independent and identically distributed exponential random variables with a constant mean. Going beyond such unrealistic assumptions has strong implications, as is evidenced by theoretical advances in the recent literature. However, there is very little empirical research, analyzing the statistical properties of service times in practice, to support that body of theoretical work. In this paper, we carry out a large-scale data-based investigation of service times in a call center with many heterogeneous agents and multiple call types. We observe that, for a given call type: (a) the service-time distribution depends strongly on the individual agent, (b) that it changes with time, and (c) that average service times are correlated across successive days or weeks. We develop stochastic models that account for these facts. In our proposed models, the service-time distribution is assumed to be lognormal with a mean that obeys a linear mixed-effects model with a weekly Gaussian random effect, and these successive weekly effects obey an autoregressive process of order one. We compare our models to simpler ones, e.g., where the mean service time depends only on the agent and call type, or only on the call type, and we find that our proposed models have a better goodness-of-fit, both in-sample and out-of-sample. We also perform simulation experiments to show that the choice of model can have a significant impact on the estimates of common measures of quality of service in the call center. Our study provides empirical support to the theoretical research that goes beyond standard modelling assumptions in service systems.