# Sample time series data
library(TTR) # You must run this line first
library(DLL)
library(forecast)
## Registered S3 method overwritten by 'quantmod':
## method from
## as.zoo.data.frame zoo
library(fpp3)
## ── Attaching packages ────────────────────────────────────────────── fpp3 0.5 ──
## ✔ tibble 3.2.1 ✔ tsibble 1.1.4
## ✔ dplyr 1.1.4 ✔ tsibbledata 0.4.1
## ✔ tidyr 1.3.1 ✔ feasts 0.3.2
## ✔ lubridate 1.9.3 ✔ fable 0.3.4
## ✔ ggplot2 3.5.1 ✔ fabletools 0.4.2
## ── Conflicts ───────────────────────────────────────────────── fpp3_conflicts ──
## ✖ lubridate::date() masks base::date()
## ✖ dplyr::filter() masks stats::filter()
## ✖ tsibble::intersect() masks base::intersect()
## ✖ tsibble::interval() masks lubridate::interval()
## ✖ dplyr::lag() masks stats::lag()
## ✖ tsibble::setdiff() masks base::setdiff()
## ✖ tsibble::union() masks base::union()
library(readxl)
library(DataCombine)
library(dplyr)
library(Hmisc)
##
## Attaching package: 'Hmisc'
## The following objects are masked from 'package:dplyr':
##
## src, summarize
## The following objects are masked from 'package:base':
##
## format.pval, units
salesdata <- read.csv("Sales.csv")
print(salesdata)
## Year Quarter Period Sales
## 1 90 1 1 362
## 2 90 2 2 385
## 3 90 3 3 432
## 4 90 4 4 341
## 5 91 1 5 382
## 6 91 2 6 409
## 7 91 3 7 498
## 8 91 4 8 387
## 9 92 1 9 473
## 10 92 2 10 513
## 11 92 3 11 582
## 12 92 4 12 474
## 13 93 1 13 544
## 14 93 2 14 582
## 15 93 3 15 681
## 16 93 4 16 557
## 17 94 1 17 628
## 18 94 2 18 707
## 19 94 3 19 773
## 20 94 4 20 592
## 21 95 1 21 627
## 22 95 2 22 725
## 23 95 3 23 854
## 24 95 4 24 661
tsales<-ts(salesdata)
plot.ts(tsales)

# Deterministic model for Trend and Seasonality
fit <- lm(Sales ~ Period + Quarter, data=salesdata)
summary(fit)
##
## Call:
## lm(formula = Sales ~ Period + Quarter, data = salesdata)
##
## Residuals:
## Min 1Q Median 3Q Max
## -90.408 -41.546 -8.792 53.775 120.358
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 347.542 37.059 9.378 5.89e-09 ***
## Period 18.087 1.901 9.514 4.60e-09 ***
## Quarter -9.971 11.770 -0.847 0.406
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 63.62 on 21 degrees of freedom
## Multiple R-squared: 0.8125, Adjusted R-squared: 0.7947
## F-statistic: 45.51 on 2 and 21 DF, p-value: 2.323e-08
coefficients(fit)
## (Intercept) Period Quarter
## 347.541667 18.087500 -9.970833
fit <- lm(Sales ~ Period + factor(Quarter), data=salesdata)
summary(fit)
##
## Call:
## lm(formula = Sales ~ Period + factor(Quarter), data = salesdata)
##
## Residuals:
## Min 1Q Median 3Q Max
## -56.542 -18.500 0.417 17.312 44.975
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 303.7042 14.5356 20.894 1.44e-14 ***
## Period 18.0875 0.8287 21.827 6.47e-15 ***
## factor(Quarter)2 32.7458 16.0333 2.042 0.05524 .
## factor(Quarter)3 97.8250 16.0974 6.077 7.63e-06 ***
## factor(Quarter)4 -54.9292 16.2037 -3.390 0.00307 **
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 27.73 on 19 degrees of freedom
## Multiple R-squared: 0.9678, Adjusted R-squared: 0.961
## F-statistic: 142.6 on 4 and 19 DF, p-value: 6.87e-14
coefficients(fit)
## (Intercept) Period factor(Quarter)2 factor(Quarter)3
## 303.70417 18.08750 32.74583 97.82500
## factor(Quarter)4
## -54.92917
# Random Walk model (linear model of first difference)
# Seasonal Difference linear model
# Seasonal Means model
smeans <- aggregate(Sales~Quarter, data=salesdata, mean)
# Simple Moving Averages
kings <- scan("kings.dat",skip=3)
kingstimeseries <- ts(kings)
kingstimeseries
## Time Series:
## Start = 1
## End = 42
## Frequency = 1
## [1] 60 43 67 50 56 42 50 65 68 43 65 34 47 34 49 41 13 35 53 56 16 43 69 59 48
## [26] 59 86 55 68 51 33 49 67 77 81 67 71 81 68 70 77 56
plot.ts(kingstimeseries)

kingstimeseriesSMA3 <- SMA(kingstimeseries,n=3)
plot.ts(kingstimeseriesSMA3)

# Decomposing Time series
births <- scan("nybirths.dat",skip=3)
birthstimeseries <- ts(births, frequency=12, start=c(1946,1))
plot.ts(birthstimeseries)

birthstimeseriescomponents <- decompose(birthstimeseries)
plot(birthstimeseriescomponents)

# Simple Exponetial Smoothing
#If you believe there is no seasonality or trend just modelling error terms
#Simple exponential smoothing is predicting time series with no trend or seasonsality
#In simple exponential we are estimating Alpha only
#The coefficient is the level or baseline while the alphas is smoothing parameter
#YO=LO
#Y1=ALPHA*LO
rain <- scan("precip1.dat",skip=1)
rainseries <- ts(rain,start=c(1813))
plot.ts(rainseries)

rainseriesforecasts <- HoltWinters(rainseries, beta=FALSE, gamma=FALSE)
rainseriesforecasts
## Holt-Winters exponential smoothing without trend and without seasonal component.
##
## Call:
## HoltWinters(x = rainseries, beta = FALSE, gamma = FALSE)
##
## Smoothing parameters:
## alpha: 0.02412151
## beta : FALSE
## gamma: FALSE
##
## Coefficients:
## [,1]
## a 24.67819
rainseriesforecasts$fitted
## Time Series:
## Start = 1814
## End = 1912
## Frequency = 1
## xhat level
## 1814 23.56000 23.56000
## 1815 23.62054 23.62054
## 1816 23.57808 23.57808
## 1817 23.76290 23.76290
## 1818 23.76017 23.76017
## 1819 23.76306 23.76306
## 1820 23.82691 23.82691
## 1821 23.79900 23.79900
## 1822 23.98935 23.98935
## 1823 23.98623 23.98623
## 1824 23.98921 23.98921
## 1825 24.19282 24.19282
## 1826 24.17032 24.17032
## 1827 24.13171 24.13171
## 1828 24.10442 24.10442
## 1829 24.19549 24.19549
## 1830 24.22261 24.22261
## 1831 24.24329 24.24329
## 1832 24.32812 24.32812
## 1833 24.21938 24.21938
## 1834 24.23290 24.23290
## 1835 24.13369 24.13369
## 1836 24.13867 24.13867
## 1837 24.21782 24.21782
## 1838 24.10257 24.10257
## 1839 24.04293 24.04293
## 1840 24.12608 24.12608
## 1841 24.01280 24.01280
## 1842 24.18448 24.18448
## 1843 24.15808 24.15808
## 1844 24.19889 24.19889
## 1845 24.16153 24.16153
## 1846 24.12748 24.12748
## 1847 24.18133 24.18133
## 1848 24.02499 24.02499
## 1849 24.16454 24.16454
## 1850 24.13476 24.13476
## 1851 24.01621 24.01621
## 1852 23.93453 23.93453
## 1853 24.20964 24.20964
## 1854 24.25018 24.25018
## 1855 24.11509 24.11509
## 1856 24.08964 24.08964
## 1857 24.04430 24.04430
## 1858 23.99933 23.99933
## 1859 23.87319 23.87319
## 1860 23.97780 23.97780
## 1861 24.17710 24.17710
## 1862 24.13110 24.13110
## 1863 24.21405 24.21405
## 1864 24.15075 24.15075
## 1865 23.97658 23.97658
## 1866 24.10933 24.10933
## 1867 24.29001 24.29001
## 1868 24.33729 24.33729
## 1869 24.31468 24.31468
## 1870 24.34134 24.34134
## 1871 24.26847 24.26847
## 1872 24.28659 24.28659
## 1873 24.51752 24.51752
## 1874 24.47295 24.47295
## 1875 24.33660 24.33660
## 1876 24.43558 24.43558
## 1877 24.47717 24.47717
## 1878 24.56625 24.56625
## 1879 24.79573 24.79573
## 1880 25.01341 25.01341
## 1881 25.14045 25.14045
## 1882 25.20750 25.20750
## 1883 25.25411 25.25411
## 1884 25.23351 25.23351
## 1885 25.11571 25.11571
## 1886 25.15248 25.15248
## 1887 25.19729 25.19729
## 1888 25.05286 25.05286
## 1889 25.11768 25.11768
## 1890 25.08710 25.08710
## 1891 24.99407 24.99407
## 1892 25.07019 25.07019
## 1893 25.01085 25.01085
## 1894 24.88515 24.88515
## 1895 24.95884 24.95884
## 1896 24.87469 24.87469
## 1897 24.84201 24.84201
## 1898 24.79420 24.79420
## 1899 24.62284 24.62284
## 1900 24.57259 24.57259
## 1901 24.54141 24.54141
## 1902 24.48421 24.48421
## 1903 24.39631 24.39631
## 1904 24.72686 24.72686
## 1905 24.62852 24.62852
## 1906 24.58852 24.58852
## 1907 24.58059 24.58059
## 1908 24.54271 24.54271
## 1909 24.52166 24.52166
## 1910 24.57541 24.57541
## 1911 24.59433 24.59433
## 1912 24.59905 24.59905
plot(rainseriesforecasts)

rainseriesforecasts$SSE
## [1] 1828.855
#simple exponential
HoltWinters(rainseries, beta=FALSE, gamma=FALSE, l.start=23.56)
## Holt-Winters exponential smoothing without trend and without seasonal component.
##
## Call:
## HoltWinters(x = rainseries, beta = FALSE, gamma = FALSE, l.start = 23.56)
##
## Smoothing parameters:
## alpha: 0.02412151
## beta : FALSE
## gamma: FALSE
##
## Coefficients:
## [,1]
## a 24.67819
# predicting 8 years ahead
library(forecast)
rainseriesforecasts2 <- forecast(rainseriesforecasts, h=8)
rainseriesforecasts2
## Point Forecast Lo 80 Hi 80 Lo 95 Hi 95
## 1913 24.67819 19.17493 30.18145 16.26169 33.09470
## 1914 24.67819 19.17333 30.18305 16.25924 33.09715
## 1915 24.67819 19.17173 30.18465 16.25679 33.09960
## 1916 24.67819 19.17013 30.18625 16.25434 33.10204
## 1917 24.67819 19.16853 30.18785 16.25190 33.10449
## 1918 24.67819 19.16694 30.18945 16.24945 33.10694
## 1919 24.67819 19.16534 30.19105 16.24701 33.10938
## 1920 24.67819 19.16374 30.19265 16.24456 33.11182
rainseriesforecasts2$residuals
## Time Series:
## Start = 1813
## End = 1912
## Frequency = 1
## [1] NA 2.5100000 -1.7605450 7.6619220 -0.1128951 0.1198281
## [7] 2.6469377 -1.1569105 7.8909960 -0.1293468 0.1237733 8.4407877
## [13] -0.9328169 -1.6003159 -1.1317139 3.7755848 1.1245120 0.8573870
## [19] 3.5167056 -4.5081227 0.5606200 -4.1129030 0.2063065 3.2813300
## [25] -4.7778206 -2.4725723 3.4470698 -4.6960787 7.1171978 -1.0944797
## [31] 1.6919208 -1.5488909 -1.4115293 2.2325189 -6.4813328 5.7850067
## [37] -1.2345364 -4.9147575 -3.3862062 11.4054742 1.6803570 -5.6001757
## [43] -1.0550911 -1.8796407 -1.8643009 -5.2293312 4.3368082 8.2621979
## [49] -1.9070988 3.4389033 -2.6240483 -7.2207523 5.5034232 7.4906723
## [55] 1.9599860 -0.9372918 1.1053171 -3.0213449 0.7515345 9.5734064
## [61] -1.8475186 -5.6529537 4.1034041 1.7244238 3.6928281 9.5137515
## [67] 9.0242655 5.2665866 2.7795486 1.9325016 -0.8541132 -4.8835107
## [73] 1.5242869 1.8575188 -5.9872873 2.6871351 -1.2676827 -3.8571042
## [79] 3.1559349 -2.4601910 -5.2108475 3.0548460 -3.4888415 -1.3546853
## [85] -1.9820083 -7.1041993 -2.0828352 -1.2925941 -2.3714148 -3.6442127
## [91] 13.7036912 -4.0768625 -1.6585224 -0.3285164 -1.5705920 -0.8727070
## [97] 2.2283440 0.7845930 0.1956674 3.2809476
acf(rainseriesforecasts2$residuals, lag.max=20, na.action = na.pass)

Box.test(rainseriesforecasts2$residuals, lag=20, type="Ljung-Box")
##
## Box-Ljung test
##
## data: rainseriesforecasts2$residuals
## X-squared = 17.401, df = 20, p-value = 0.6268
plot.ts(rainseriesforecasts2$residuals)

#Holt's Exponential Smoothing (Double)
skirts <- scan("skirts.dat",skip=5)
skirtsseries <- ts(skirts,start=c(1866))
plot.ts(skirtsseries)

#double exponential smoothing
skirtsseriesforecasts <- HoltWinters(skirtsseries, gamma=FALSE)
skirtsseriesforecasts
## Holt-Winters exponential smoothing with trend and without seasonal component.
##
## Call:
## HoltWinters(x = skirtsseries, gamma = FALSE)
##
## Smoothing parameters:
## alpha: 0.8383481
## beta : 1
## gamma: FALSE
##
## Coefficients:
## [,1]
## a 529.308585
## b 5.690464
skirtsseriesforecasts$fitted
## Time Series:
## Start = 1868
## End = 1911
## Frequency = 1
## xhat level trend
## 1868 626.0000 617.0000 9.0000000
## 1869 633.3233 625.1617 8.1616519
## 1870 645.9730 635.5673 10.4056551
## 1871 674.8676 655.2175 19.6501514
## 1872 721.5669 688.3922 33.1747101
## 1873 765.5280 726.9601 38.5679053
## 1874 835.0679 781.0140 54.0538890
## 1875 857.1507 819.0824 38.0683919
## 1876 926.8236 872.9530 53.8706278
## 1877 1017.8773 945.4151 72.4621621
## 1878 1055.3346 1000.3749 54.9597120
## 1879 1062.7858 1031.5803 31.2054806
## 1880 1067.5233 1049.5518 17.9714705
## 1881 1054.4368 1051.9943 2.4425157
## 1882 995.7858 1023.8901 -28.1042388
## 1883 1009.9581 1016.9241 -6.9660004
## 1884 1006.4158 1011.6699 -5.2541636
## 1885 1020.5848 1016.1274 4.4574647
## 1886 975.4375 995.7824 -20.3449378
## 1887 932.5620 964.1722 -31.6102483
## 1888 881.5658 922.8690 -41.3031877
## 1889 822.5470 872.7080 -50.1610258
## 1890 783.2057 827.9569 -44.7511318
## 1891 803.5008 815.7288 -12.2280226
## 1892 818.9369 817.3329 1.6040432
## 1893 792.1429 804.7379 -12.5949893
## 1894 827.9325 816.3352 11.5973092
## 1895 824.5527 820.4440 4.1087620
## 1896 772.4038 796.4239 -24.0200954
## 1897 704.1126 750.2682 -46.1556418
## 1898 694.6555 722.4619 -27.8063625
## 1899 610.4191 666.4405 -56.0213818
## 1900 570.4620 618.4512 -47.9892465
## 1901 551.8787 585.1649 -33.2862822
## 1902 545.6230 565.3940 -19.7709906
## 1903 519.7774 542.5857 -22.8082962
## 1904 549.3199 545.9528 3.3671422
## 1905 538.7371 542.3449 -3.6078666
## 1906 565.7506 554.0478 11.7028276
## 1907 561.1046 557.5762 3.5284391
## 1908 519.1868 538.3815 -19.1946847
## 1909 514.7692 526.5753 -11.8061821
## 1910 526.8238 526.6996 0.1242218
## 1911 520.5367 523.6181 -3.0814354
plot(skirtsseriesforecasts)

HoltWinters(skirtsseries, gamma=FALSE, l.start=608, b.start=9)
## Holt-Winters exponential smoothing with trend and without seasonal component.
##
## Call:
## HoltWinters(x = skirtsseries, gamma = FALSE, l.start = 608, b.start = 9)
##
## Smoothing parameters:
## alpha: 0.8346775
## beta : 1
## gamma: FALSE
##
## Coefficients:
## [,1]
## a 529.278637
## b 5.670129
skirtsseriesforecasts2 <- forecast(skirtsseriesforecasts, h=5)
skirtsseriesforecasts2
## Point Forecast Lo 80 Hi 80 Lo 95 Hi 95
## 1912 534.9990 509.5521 560.4460 496.0813 573.9168
## 1913 540.6895 491.0105 590.3685 464.7120 616.6670
## 1914 546.3800 465.3613 627.3987 422.4726 670.2874
## 1915 552.0704 434.4020 669.7388 372.1122 732.0287
## 1916 557.7609 398.9412 716.5806 314.8671 800.6547
plot(skirtsseriesforecasts2)

acf(skirtsseriesforecasts2$residuals, lag.max=20, na.action = na.pass)

Box.test(skirtsseriesforecasts2$residuals, lag=20, type="Ljung-Box")
##
## Box-Ljung test
##
## data: skirtsseriesforecasts2$residuals
## X-squared = 19.731, df = 20, p-value = 0.4749
plot.ts(skirtsseriesforecasts2$residuals) # make a time plot

#Holt-Winters Exponential Smoothing (triple)
souvenir <- scan("fancy.dat")
souvenir
## [1] 1664.81 2397.53 2840.71 3547.29 3752.96 3714.74 4349.61
## [8] 3566.34 5021.82 6423.48 7600.60 19756.21 2499.81 5198.24
## [15] 7225.14 4806.03 5900.88 4951.34 6179.12 4752.15 5496.43
## [22] 5835.10 12600.08 28541.72 4717.02 5702.63 9957.58 5304.78
## [29] 6492.43 6630.80 7349.62 8176.62 8573.17 9690.50 15151.84
## [36] 34061.01 5921.10 5814.58 12421.25 6369.77 7609.12 7224.75
## [43] 8121.22 7979.25 8093.06 8476.70 17914.66 30114.41 4826.64
## [50] 6470.23 9638.77 8821.17 8722.37 10209.48 11276.55 12552.22
## [57] 11637.39 13606.89 21822.11 45060.69 7615.03 9849.69 14558.40
## [64] 11587.33 9332.56 13082.09 16732.78 19888.61 23933.38 25391.35
## [71] 36024.80 80721.71 10243.24 11266.88 21826.84 17357.33 15997.79
## [78] 18601.53 26155.15 28586.52 30505.41 30821.33 46634.38 104660.67