#--------------------------------- # Derby Data #--------------------------------- mydata = read.csv("http://faculty.chicagobooth.edu/nicholas.polson/teaching/41000/Kentucky_Derby_2018.csv",header=T) head(mydata) tail(mydata) attach(mydata) ts1 <- ts(speedmph, start=1875, end=2016, frequency=1) plot(ts1,xlab="year",ylab="speed mph",main="Kentucky Derby",col=2, lty=2, lwd=2) # we can see obvious time trend in the time series plot # here is a simple way to extract the time trend # after removing the time trend, we can talk about the ARMA model model = lm(speedmph~year_num) speed_tt = model\$coef[1] + model\$res ts2 <- ts(speed_tt, start=1875, end=2016, frequency=1) plot(ts2,xlab="year",ylab="speed mph",main="time trend removed",col=2, lty=2, lwd=2) # to find the outlier out = which(abs(model\$res)>1.645) # label those outliers plot(ts2,xlab="year",ylab="speed mph",main="time trend removed",col=2, lty=2, lwd=2) text(x=out+1874,y=ts2[out],labels=winner[out]) # remove the outliers ts3 = ts2[-out] winner2 = winner[-out] years = c(1:142)[-out]+1874 # plot the acf and pacf plots # an AR model might not be appropriate # but we can consider an MA model # acf(ts3) # pacf(ts3) model3 = arima(ts3,c(0,0,1)) out2 = which(abs(model3\$res)>1) winner2[out2] # MA(1) prediction + the time trend ts4 = ts3-model3\$res + c(1:142)[-out]*model\$coef[2] plot(years,ts4,type="l",ylim=c(32.5,37.5),main="MA(1) predictions, outliers removed", xlab="year",ylab="speed mph",col=2, lty=2, lwd=2) lines(years,ts1[-out],type="l",col=3, lty=3, lwd=3) text(x=years[out2],y=ts1[-out][out2],labels=winner2[out2])