This brief analysis looks at daily deaths between 1st Jan 1995 and 31st December 2014.The source of the data is the births and deaths extract held by PHE derived from ONS.

First read in the data and plot the number of deaths per day.

autoplot(deaths.ts) + ggtitle("Daily deaths 1995 - 2014") + ylab('No deaths')
Figure 1: Daily deaths

Figure 1: Daily deaths

There is clearly the well known seasonal pattern with winter peaks and summer troughs. The shape of the peaks varies and the size of the winter peak has reduced over time. In some years the the peaks are sharp e.g. 2000, and in others e.g. 2013, they are more ‘spread out’. Looking at the number of deaths per year suggests a downward linear trend.

ggplot(yearlydeaths, aes(year, deaths)) + geom_point() + geom_line() +ggtitle ('Yearly deaths in England. 1995-2014') + geom_smooth(method = 'lm')
Figure 2: Yearly deaths

Figure 2: Yearly deaths

Seasonal variation

We can decompose the time series into the ‘trend component’, seasonal component and noise

autoplot(deathscomp1)
Figure 3: Seasonal decomposition

Figure 3: Seasonal decomposition

This suggests and general downward trend in the number of deaths per day from 1500 to 1350, a consistent seasonal pattern, and larger random variation coinciding with the largest data peaks in the winters of 1996, 1997, 1999 and 2000. Removing the seasonal component reveals the underlying trend.

Simple seasonal adjustment

autoplot(deaths.tsadj, ylab = 'Seasonally adjusted no deaths')
Figure 4: Seasonally adjusted number of deaths

Figure 4: Seasonally adjusted number of deaths

Projecting ahead and detecting change

The next step is to project ahead. In this example forecasting ahead 2 years using HoltWinters modelling.

library('forecast')
deathforecast <- HoltWinters(deaths.ts)

##autoplot(deathforecast) + ggtitle('Holt-Winters model for daily deaths')
deathforecast$SSE
## [1] 28116574
deathforecast2 <- forecast.HoltWinters(deathforecast, t = 50)
##deathforecast2
autoplot(deathforecast2)
Figure 5: Forecast of daily deaths

Figure 5: Forecast of daily deaths

Changepoint analysis cna be used to detect shifts in mean daily deaths and variability in dail deaths. This suggests a shift in 2004 and again in 2014.

library(changepoint)

autoplot(cpt.mean(deaths.ts))
Figure 6: Breaks in mean no of daily deaths

Figure 6: Breaks in mean no of daily deaths 

library(changepoint)

set.seed(10)
cpts <- cpt.mean(deaths.ts, method = 'BinSeg')
## Warning in BINSEG(sumstat, pen = pen.value, cost_func = costfunc, minseglen
## = minseglen, : The number of changepoints identified is Q, it is advised to
## increase Q to make sure changepoints have not been missed.
plot(cpts, pch = 20, col ='grey', cpt.col = "black", type = "p")

ncpts(cpts)
## [1] 5
cpts(cpts)
## [1] 1552 1810 1841 1876 3370