# 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