0.1 Zadanie 1.

Zestaw danych plastiki składa się z miesięcznej sprzedaży (w tysiącach) produktu A dla producenta tworzyw sztucznych przez pięć lat.

  1. Wykreśl szereg czasowy sprzedaży produktu A.
# wykres:
  1. Czy możesz zidentyfikować sezonowe wahania i/lub trend-cykl?
  1. Oblicz i wykreśl dane dostosowane sezonowo.

0.2 Zadanie 2.

Rozważmy miesięczne dane dotyczące sprzedaży i reklamy dla firmy produkującej części samochodowe (zestaw danych advert).

  1. Wykreśl dane używając autoplot. Dlaczego warto ustawić facets=TRUE?
  1. Ponieważ dane są w różnych skalach, nie powinny być wykreślane w tym samym panelu. Ustawienie facets=TRUE umieszcza je w różnych panelach.
  1. Dopasuj standardowy model regresji \(y_t = a + b*x_t + \a_t\), gdzie \(y_t\) oznacza sprzedaż, a \(x_t\) reklamę, używając funkcji tslm().
  1. Wykazać, że reszty mają znaczną autokorelację.

0.3 Zadanie 3.

W ćwiczeniu wykorzystano dane cangas (miesięczna produkcja gazu w Kanadzie w miliardach metrów sześciennych, styczeń 1960 - luty 2005).

  1. Wykreśl dane używając autoplot, ggsubseriesplot i ggseasonplot, aby spojrzeć na efekt zmieniającej się sezonowości w czasie. Jak myślisz, co powoduje, że zmienia się ona tak bardzo?
  1. Wykonaj dekompozycję STL danych. Musisz wybrać s.window, aby wychwycić zmienny kształt składnika sezonowego.
  1. Porównaj wyniki z tymi uzyskanymi przy użyciu SEATS i X11. Czym się różnią?

0.4 Grafika

0.4.1 Quandl

Pobieramy dane z Quandl i jednocześnie je dekomponujemy miesięcznie oraz wykreślamy:

plot(stl(Quandl("WIKI/GOOG",type="ts",collapse="monthly")[,11],s.window="per"))

NSE <- Quandl("NSE/OIL", type="xts")
candleChart(NSE)

getSymbols("P9O.F", src='yahoo', from = "2019-01-01", to = "2023-03-29")
## [1] "P9O.F"
tail(P9O.F)
##            P9O.F.Open P9O.F.High P9O.F.Low P9O.F.Close P9O.F.Volume
## 2023-03-21      5.664      5.664     5.664       5.664            0
## 2023-03-22      5.928      5.928     5.928       5.928            0
## 2023-03-23      5.828      5.828     5.810       5.810          173
## 2023-03-24      5.792      5.792     5.686       5.686          800
## 2023-03-27      5.668      5.668     5.668       5.668            0
## 2023-03-28      5.748      5.748     5.748       5.748            0
##            P9O.F.Adjusted
## 2023-03-21          5.664
## 2023-03-22          5.928
## 2023-03-23          5.810
## 2023-03-24          5.686
## 2023-03-27          5.668
## 2023-03-28          5.748
plot(P9O.F$P9O.F.Close)

p=100
plot(rollmean(P9O.F$P9O.F.Close,p),main='Prosta vs Trójkątna Średnia Ruchoma')

lines(rollmean(P9O.F$P9O.F.Close,10),col='red')

lines(rollmean(rollmean(P9O.F$P9O.F.Close,5),5),col='blue')

Aby tworzyć prognozy przy użyciu prostego wygładzania wykładniczego w R, możemy dopasować prosty model predykcyjny wygładzania wykładniczego przy użyciu funkcji “HoltWinters()” w R. Aby użyć funkcji HoltWinters() do prostego wygładzania wykładniczego, musimy ustawić parametry beta=FALSE i gamma=FALSE w funkcji HoltWinters() (parametry beta i gamma są używane do wygładzania wykładniczego Holta, lub wygładzania wykładniczego Holta-Wintersa, jak opisano poniżej).

Funkcja HoltWinters() zwraca obiekt typu lista, która zawiera kilka nazwanych elementów.

wynik <- HoltWinters(P9O.F$P9O.F.Close, beta=FALSE, gamma=FALSE)
wynik
## Holt-Winters exponential smoothing without trend and without seasonal component.
## 
## Call:
## HoltWinters(x = P9O.F$P9O.F.Close, beta = FALSE, gamma = FALSE)
## 
## Smoothing parameters:
##  alpha: 0.80824
##  beta : FALSE
##  gamma: FALSE
## 
## Coefficients:
##     [,1]
## a 5.7343

Wynik działania funkcji HoltWinters() mówi nam, że szacowana wartość parametru alfa wynosi około 0.8. Jest to wartość bardzo bliska 1, co mówi nam, że prognozy są oparte zarówno na ostatnich, jak i mniej aktualnych obserwacjach.

Możemy wykreślić oryginalny szereg czasowy w stosunku do prognoz wpisując:

plot(wynik)

Wykres pokazuje oryginalny szereg czasowy w kolorze czarnym, a prognozy jako czerwoną linię. Szereg czasowy prognoz jest znacznie gładszy niż szereg czasowy oryginalnych danych.

Jako miarę dokładności prognoz możemy obliczyć sumę kwadratów błędów dla błędów prognozy w próbie, czyli błędów prognozy dla okresu czasu objętego naszym oryginalnym szeregiem czasowym.

Suma błędów kwadratowych jest przechowywana w nazwanym elemencie zmiennej listowej “wynik” o nazwie “SSE”, więc możemy uzyskać jej wartość wpisując:

wynik$SSE
## [1] 44.84

Prognozy dla kolejnych punktów czasowych możemy wykonać za pomocą funkcji “forecast.HoltWinters()” z pakietu R “forecast”.

wynik2 <- forecast(wynik, h=20)
plot(wynik2)

Na podstawie: https://a-little-book-of-r-for-time-series.readthedocs.io/en/latest/src/timeseries.html

LS0tCnRpdGxlOiAiUmFwb3J0IDMuIgphdXRob3I6ICJUd29qZSBpbWnEmSBpIG5hendpc2tvIgpvdXRwdXQ6CiAgaHRtbF9kb2N1bWVudDoKICAgIHRoZW1lOiBjZXJ1bGVhbgogICAgaGlnaGxpZ2h0OiB0ZXh0bWF0ZQogICAgZm9udHNpemU6IDhwdAogICAgdG9jOiB0cnVlCiAgICBudW1iZXJfc2VjdGlvbnM6IHRydWUKICAgIGNvZGVfZG93bmxvYWQ6IHRydWUKICAgIHRvY19mbG9hdDoKICAgICAgY29sbGFwc2VkOiBmYWxzZQplZGl0b3Jfb3B0aW9uczogCiAgbWFya2Rvd246IAogICAgd3JhcDogNzIKLS0tCgpgYGB7ciBzZXR1cCwgaW5jbHVkZT1GQUxTRX0Ka25pdHI6Om9wdHNfY2h1bmskc2V0KGVjaG8gPSBUUlVFLCBjYWNoZT1UUlVFKQpsaWJyYXJ5KGZwcDIpCmxpYnJhcnkoUXVhbmRsKQpvcHRpb25zKGRpZ2l0cz01KQpsaWJyYXJ5KGdncGxvdDIpCnRoZW1lX3NldCh0aGVtZV9taW5pbWFsKCkpCmxpYnJhcnkocXVhbnRtb2QpCmxpYnJhcnkoImZvcmVjYXN0IikKYGBgCgojIyBaYWRhbmllIDEuCgpaZXN0YXcgZGFueWNoIGBwbGFzdGlraWAgc2vFgmFkYSBzacSZIHogbWllc2nEmWN6bmVqIHNwcnplZGHFvHkgKHcgdHlzacSFY2FjaCkgcHJvZHVrdHUgQSBkbGEgcHJvZHVjZW50YSB0d29yenl3IHN6dHVjem55Y2ggcHJ6ZXogcGnEmcSHIGxhdC4KCj4gYS4gV3lrcmXFm2wgc3plcmVnIGN6YXNvd3kgc3ByemVkYcW8eSBwcm9kdWt0dSBBLiAKCmBgYHtyIHd5a3Jlc30KIyB3eWtyZXM6CgpgYGAKCj4gYi4gQ3p5IG1vxbxlc3ogemlkZW50eWZpa293YcSHIHNlem9ub3dlIHdhaGFuaWEgaS9sdWIgdHJlbmQtY3lrbD8KCmBgYHtyIHRyZW5kX2N5a2x9CgoKYGBgCgo+IGMuIE9ibGljeiBpIHd5a3JlxZtsIGRhbmUgZG9zdG9zb3dhbmUgc2V6b25vd28uCgpgYGB7cn0KCgpgYGAKCiMjIFphZGFuaWUgMi4gCgpSb3p3YcW8bXkgbWllc2nEmWN6bmUgZGFuZSBkb3R5Y3rEhWNlIHNwcnplZGHFvHkgaSByZWtsYW15IGRsYSBmaXJteSBwcm9kdWt1asSFY2VqIGN6xJnFm2NpIHNhbW9jaG9kb3dlICh6ZXN0YXcgZGFueWNoIGBhZHZlcnRgKS4KCj4gYS4gV3lrcmXFm2wgZGFuZSB1xbx5d2FqxIVjIGBhdXRvcGxvdGAuIERsYWN6ZWdvIHdhcnRvIHVzdGF3acSHIGBmYWNldHM9VFJVRWA/CgpgYGB7cn0KCmBgYAoKPiBiLiBQb25pZXdhxbwgZGFuZSBzxIUgdyByw7PFvG55Y2ggc2thbGFjaCwgbmllIHBvd2lubnkgYnnEhyB3eWtyZcWbbGFuZSB3IHR5bSBzYW15bSBwYW5lbHUuIFVzdGF3aWVuaWUgYGZhY2V0cz1UUlVFYCB1bWllc3pjemEgamUgdyByw7PFvG55Y2ggcGFuZWxhY2guCgpgYGB7cn0KCmBgYAoKPiBjLiBEb3Bhc3VqIHN0YW5kYXJkb3d5IG1vZGVsIHJlZ3Jlc2ppICR5X3QgPSBhICsgYip4X3QgKyBcYV90JCwgZ2R6aWUgJHlfdCQgb3puYWN6YSBzcHJ6ZWRhxbwsIGEgJHhfdCQgcmVrbGFtxJksIHXFvHl3YWrEhWMgZnVua2NqaSBgdHNsbSgpYC4KCmBgYHtyfQoKYGBgCgo+IGQuIFd5a2F6YcSHLCDFvGUgcmVzenR5IG1hasSFIHpuYWN6bsSFIGF1dG9rb3JlbGFjasSZLgoKYGBge3J9CgpgYGAKCiMjIFphZGFuaWUgMy4KClcgxId3aWN6ZW5pdSB3eWtvcnp5c3Rhbm8gZGFuZSBgY2FuZ2FzYCAobWllc2nEmWN6bmEgcHJvZHVrY2phIGdhenUgdyBLYW5hZHppZSB3IG1pbGlhcmRhY2ggbWV0csOzdyBzemXFm2NpZW5ueWNoLCBzdHljemXFhCAxOTYwIC0gbHV0eSAyMDA1KS4KCj4gYS4gV3lrcmXFm2wgZGFuZSB1xbx5d2FqxIVjIGBhdXRvcGxvdGAsIGBnZ3N1YnNlcmllc3Bsb3RgIGkgYGdnc2Vhc29ucGxvdGAsIGFieSBzcG9qcnplxIcgbmEgZWZla3Qgem1pZW5pYWrEhWNlaiBzacSZIHNlem9ub3dvxZtjaSB3IGN6YXNpZS4gSmFrIG15xZtsaXN6LCBjbyBwb3dvZHVqZSwgxbxlIHptaWVuaWEgc2nEmSBvbmEgdGFrIGJhcmR6bz8KCmBgYHtyfQoKYGBgCgo+IGIuIFd5a29uYWogZGVrb21wb3p5Y2rEmSBTVEwgZGFueWNoLiBNdXNpc3ogd3licmHEhyBgcy53aW5kb3dgLCBhYnkgd3ljaHd5Y2nEhyB6bWllbm55IGtzenRhxYJ0IHNrxYJhZG5pa2Egc2V6b25vd2Vnby4KCmBgYHtyfQoKYGBgCgo+IGMuIFBvcsOzd25haiB3eW5pa2kgeiB0eW1pIHV6eXNrYW55bWkgcHJ6eSB1xbx5Y2l1IFNFQVRTIGkgWDExLiBDenltIHNpxJkgcsOzxbxuacSFPwoKYGBge3J9CgoKYGBgCgojIyBHcmFmaWthCgojIyMgUXVhbmRsCgpQb2JpZXJhbXkgZGFuZSB6IFF1YW5kbCBpIGplZG5vY3plxZtuaWUgamUgZGVrb21wb251amVteSBtaWVzacSZY3puaWUgb3JheiB3eWtyZcWbbGFteToKCmBgYHtyfQpwbG90KHN0bChRdWFuZGwoIldJS0kvR09PRyIsdHlwZT0idHMiLGNvbGxhcHNlPSJtb250aGx5IilbLDExXSxzLndpbmRvdz0icGVyIikpCmBgYAoKYGBge3J9Ck5TRSA8LSBRdWFuZGwoIk5TRS9PSUwiLCB0eXBlPSJ4dHMiKQpjYW5kbGVDaGFydChOU0UpCmBgYAoKCgpgYGB7cn0KZ2V0U3ltYm9scygiUDlPLkYiLCBzcmM9J3lhaG9vJywgZnJvbSA9ICIyMDE5LTAxLTAxIiwgdG8gPSAiMjAyMy0wMy0yOSIpCnRhaWwoUDlPLkYpCmBgYAoKCmBgYHtyfQpwbG90KFA5Ty5GJFA5Ty5GLkNsb3NlKQpgYGAKCmBgYHtyfQpwPTEwMApwbG90KHJvbGxtZWFuKFA5Ty5GJFA5Ty5GLkNsb3NlLHApLG1haW49J1Byb3N0YSB2cyBUcsOzamvEhXRuYSDFmnJlZG5pYSBSdWNob21hJykKbGluZXMocm9sbG1lYW4oUDlPLkYkUDlPLkYuQ2xvc2UsMTApLGNvbD0ncmVkJykKbGluZXMocm9sbG1lYW4ocm9sbG1lYW4oUDlPLkYkUDlPLkYuQ2xvc2UsNSksNSksY29sPSdibHVlJykKYGBgCgoKQWJ5IHR3b3J6ecSHIHByb2dub3p5IHByenkgdcW8eWNpdSBwcm9zdGVnbyB3eWfFgmFkemFuaWEgd3lrxYJhZG5pY3plZ28gdyBSLCBtb8W8ZW15IGRvcGFzb3dhxIcgcHJvc3R5IG1vZGVsIHByZWR5a2N5am55IHd5Z8WCYWR6YW5pYSB3eWvFgmFkbmljemVnbyBwcnp5IHXFvHljaXUgZnVua2NqaSAiSG9sdFdpbnRlcnMoKSIgdyBSLiBBYnkgdcW8ecSHIGZ1bmtjamkgSG9sdFdpbnRlcnMoKSBkbyBwcm9zdGVnbyB3eWfFgmFkemFuaWEgd3lrxYJhZG5pY3plZ28sIG11c2lteSB1c3Rhd2nEhyBwYXJhbWV0cnkgYmV0YT1GQUxTRSBpIGdhbW1hPUZBTFNFIHcgZnVua2NqaSBIb2x0V2ludGVycygpIChwYXJhbWV0cnkgYmV0YSBpIGdhbW1hIHPEhSB1xbx5d2FuZSBkbyB3eWfFgmFkemFuaWEgd3lrxYJhZG5pY3plZ28gSG9sdGEsIGx1YiB3eWfFgmFkemFuaWEgd3lrxYJhZG5pY3plZ28gSG9sdGEtV2ludGVyc2EsIGphayBvcGlzYW5vIHBvbmnFvGVqKS4KCkZ1bmtjamEgSG9sdFdpbnRlcnMoKSB6d3JhY2Egb2JpZWt0IHR5cHUgbGlzdGEsIGt0w7NyYSB6YXdpZXJhIGtpbGthIG5hendhbnljaCBlbGVtZW50w7N3LgoKYGBge3J9Cnd5bmlrIDwtIEhvbHRXaW50ZXJzKFA5Ty5GJFA5Ty5GLkNsb3NlLCBiZXRhPUZBTFNFLCBnYW1tYT1GQUxTRSkKd3luaWsKYGBgCgpXeW5payBkemlhxYJhbmlhIGZ1bmtjamkgSG9sdFdpbnRlcnMoKSBtw7N3aSBuYW0sIMW8ZSBzemFjb3dhbmEgd2FydG/Fm8SHIHBhcmFtZXRydSBhbGZhIHd5bm9zaSBva2/Fgm8gMC44LiBKZXN0IHRvIHdhcnRvxZvEhyBiYXJkem8gYmxpc2thIDEsIGNvIG3Ds3dpIG5hbSwgxbxlIHByb2dub3p5IHPEhSBvcGFydGUgemFyw7N3bm8gbmEgb3N0YXRuaWNoLCBqYWsgaSBtbmllaiBha3R1YWxueWNoIG9ic2Vyd2FjamFjaC4KCk1vxbxlbXkgd3lrcmXFm2xpxIcgb3J5Z2luYWxueSBzemVyZWcgY3phc293eSB3IHN0b3N1bmt1IGRvIHByb2dub3ogd3Bpc3VqxIVjOgoKYGBge3J9CnBsb3Qod3luaWspCmBgYAoKCld5a3JlcyBwb2thenVqZSBvcnlnaW5hbG55IHN6ZXJlZyBjemFzb3d5IHcga29sb3J6ZSBjemFybnltLCBhIHByb2dub3p5IGpha28gY3plcndvbsSFIGxpbmnEmS4gU3plcmVnIGN6YXNvd3kgcHJvZ25veiBqZXN0IHpuYWN6bmllIGfFgmFkc3p5IG5pxbwgc3plcmVnIGN6YXNvd3kgb3J5Z2luYWxueWNoIGRhbnljaC4KCkpha28gbWlhcsSZIGRva8WCYWRub8WbY2kgcHJvZ25veiBtb8W8ZW15IG9ibGljennEhyBzdW3EmSBrd2FkcmF0w7N3IGLFgsSZZMOzdyBkbGEgYsWCxJlkw7N3IHByb2dub3p5IHcgcHLDs2JpZSwgY3p5bGkgYsWCxJlkw7N3IHByb2dub3p5IGRsYSBva3Jlc3UgY3phc3Ugb2JqxJl0ZWdvIG5hc3p5bSBvcnlnaW5hbG55bSBzemVyZWdpZW0gY3phc293eW0uIAoKU3VtYSBixYLEmWTDs3cga3dhZHJhdG93eWNoIGplc3QgcHJ6ZWNob3d5d2FuYSB3IG5hendhbnltIGVsZW1lbmNpZSB6bWllbm5laiBsaXN0b3dlaiAid3luaWsiIG8gbmF6d2llICJTU0UiLCB3acSZYyBtb8W8ZW15IHV6eXNrYcSHIGplaiB3YXJ0b8WbxIcgd3Bpc3VqxIVjOgoKYGBge3J9Cnd5bmlrJFNTRQpgYGAKClByb2dub3p5IGRsYSBrb2xlam55Y2ggcHVua3TDs3cgY3phc293eWNoIG1vxbxlbXkgd3lrb25hxIcgemEgcG9tb2PEhSBmdW5rY2ppICJmb3JlY2FzdC5Ib2x0V2ludGVycygpIiB6IHBha2lldHUgUiAiZm9yZWNhc3QiLgoKYGBge3J9Cnd5bmlrMiA8LSBmb3JlY2FzdCh3eW5paywgaD0yMCkKcGxvdCh3eW5pazIpCmBgYAoKTmEgcG9kc3Rhd2llOiBodHRwczovL2EtbGl0dGxlLWJvb2stb2Ytci1mb3ItdGltZS1zZXJpZXMucmVhZHRoZWRvY3MuaW8vZW4vbGF0ZXN0L3NyYy90aW1lc2VyaWVzLmh0bWwKCgoKCg==