This is an R Markdown Notebook. When you execute code within the notebook, the results appear beneath the code.

Try executing this chunk by clicking the Run button within the chunk or by placing your cursor inside it and pressing Cmd+Shift+Enter.

library(fpp2)
## Loading required package: ggplot2
## Loading required package: forecast
## Loading required package: fma
## Loading required package: expsmooth
### 6.1 Time series components

### additive decomposition
### yt = St + Tt + Rt



### multiplicative decomposition
### yt = St * Tt * Rt
### log yt = log(St) + log(Tt) + log(Rt)

### 說明性的例子

#在下面的時間序列中，除了兩年 - 1983年和1998年 - 之外，每年的價值增加20％，每個價值下降 20％。
#在原始數據的情節中，1998年似乎比1983年更不尋常。事實上，很容易錯過1983年是不尋常的一年。
#由於我們認為1983年和1998年同樣不同尋常，因為它們與前一年的比例變化相同，因此建議採用乘法模型。
#將滑塊拖動到Log scale。對數據應用對數變換。觀察兩個不尋常的年份在對數刻度上顯示相同。另請注意，轉換會使系列中的趨勢線性化，從而使模型更容易。
### 6.2 Moving Average

### yt = St * Tt * Rt

### Tt 變化

autoplot(elecsales, lwd=1) + xlab("Year") + ylab("GWh") +         ### 原始數據
  ggtitle("Annual electricity sales: South Australia")

ma(elecsales, 5)                                                  ### k=2 , m = 2k+1
## Time Series:
## Start = 1989 
## End = 2008 
## Frequency = 1 
##  [1]       NA       NA 2381.530 2424.556 2463.758 2552.598 2627.700
##  [8] 2750.622 2858.348 3014.704 3077.300 3144.520 3188.700 3202.320
## [15] 3216.940 3307.296 3398.754 3485.434       NA       NA
autoplot(elecsales, series="Data", lwd=2) +
  autolayer(ma(elecsales,5), series="5-MA", lwd=1) +
  xlab("Year") + ylab("GWh") +
  ggtitle("Annual electricity sales: South Australia") +  
  scale_colour_manual(values=c("Data"="grey","5-MA"="red"),
                      breaks=c("Data","5-MA"))
## Warning: Removed 4 rows containing missing values (geom_path).

# Moving averages of moving averages
beer2 <- window(ausbeer,start=1992)
ma4 <- ma(beer2, order=4, centre=FALSE)
ma2x4 <- ma(beer2, order=4, centre=TRUE)

autoplot(elecequip, series="Data") +
  autolayer(ma(elecequip, 12), series="12-MA") +
  xlab("Year") + ylab("New orders index") +
  ggtitle("Electrical equipment manufacturing (Euro area)") +
  scale_colour_manual(values=c("Data"="grey","12-MA"="red"),
                      breaks=c("Data","12-MA"))
## Warning: Removed 12 rows containing missing values (geom_path).

6.3 classical decomposition

# 1.3分解方法
# 經典分解分為加法分解和乘法分解
# 假設所用數據為季節性週期為m的時間序列,m也成為季節性指數
# 其步驟如下：
# 1.如果季節性週期是偶數m，使用2*m-MA來估計trend-cyclic T̂ tT^t如果季節性週期是奇數m，使用m-MA來估計trend-cyclic T̂ tT^t
# 2.計算去趨勢序列y − T̂ ty−T^t
# 3.計算季節趨勢中每一個元素的平均值，如以一年中12月為周期，那麼計算所有的1月的平均值作為1月的值，計算12個月的作為周期的Ŝ tS^t
# 4.計算Ê t= yt− T̂ t− Ŝ tE^t=yt−T^t−S^t
# 乘法未解與上面的模式類似
# 只不過是將減法變為除法，不再贅述
# 
# 1.3.2 經典分解的缺點
# 使用moving average時會使數據的前面和後面的數值缺失，變為null
# 季節性的趨勢可能會隨時間變化，而經典中使用了平均值，忽略了這一變化
# 其求平均值時，容易受到異常值的影響，如有一個月特別不正常，直接計算平均值容易得到不正確的結果
# ---------------------


### 6.3 Classical Decomposition
# multiplicative decomposition
elecequip %>% decompose(type="multiplicative") %>%
  autoplot() + xlab("Year") +
  ggtitle("Classical multiplicative decomposition
          of electrical equipment index")

elecequip %>% decompose(type="additive") %>%
  autoplot() + xlab("Year") +
  ggtitle("Classical additive decomposition
          of electrical equipment index")

6.4 X11 Decomposition

### X11 Decomposition

library(seasonal)
elecequip %>% seas(x11="") -> fit
autoplot(fit) +
  ggtitle("X11 decomposition of electrical equipment index")

autoplot(elecequip, series="Data", lwd=1) +
  autolayer(trendcycle(fit), series="Trend") +
  autolayer(seasadj(fit), series="Seasonally Adjusted") +
  xlab("Year") + ylab("New orders index") +
  ggtitle("Electrical equipment manufacturing (Euro area)") +
  scale_colour_manual(values=c("gray","blue","red"),
                      breaks=c("Data","Seasonally Adjusted","Trend")) +
  theme(legend.position="bottom")

# seasonal sub-series
fit %>% seasonal() %>% ggsubseriesplot() + ylab("Seasonal")

### 6.5 SEATS decomposition
library(seasonal)
elecequip %>% seas() %>%
  autoplot() +
  ggtitle("SEATS decomposition of electrical equipment index")

### 適用於每季於每月的切割方法
### 6.6 STL Decomposition
elecequip %>%
  stl(t.window=13, s.window="periodic", robust=TRUE) %>%
  autoplot()

###LOESS叫做局部迴歸，迴歸的平滑值是由二次權重最小平方迴歸，簡單來說就是對局部資料作一次迴歸擬和，但針對選定的x資料其附近的資料點作權重加乘處理之後再放入迴歸方程之中。
# 优点 
# 1. 可以适用于任何季节周期的数据
# 2. 季节性成分随时间的变化率可以被控制 
# 3. trend-cycle 的平滑可以被控制 
# 缺点 
# 1. 不能自动处理trading day or calendar variation 
# 2. 只有加法模式 
# 可以使用Box-Cox 来将乘法变为加法，将λλ设为0为乘法，1为加法。

6.7 Measuring strength of trend and seasonality

```

### 6.8 Forecasting with decomposition

fit <- stl(elecequip, t.window=13, s.window="periodic",
           robust=TRUE)
fit %>% seasadj() %>% naive() %>%
  autoplot() + ylab("New orders index") +
  ggtitle("Naive forecasts of seasonally adjusted data")

fit %>% forecast(method="naive") %>%
  autoplot(flwd=1) + ylab("New orders index")

stlf(elecequip, method='naive') %>% 
  autoplot(flwd=1) + ylab("New orders index")

Add a new chunk by clicking the Insert Chunk button on the toolbar or by pressing Cmd+Option+I.

When you save the notebook, an HTML file containing the code and output will be saved alongside it (click the Preview button or press Cmd+Shift+K to preview the HTML file).

The preview shows you a rendered HTML copy of the contents of the editor. Consequently, unlike Knit, Preview does not run any R code chunks. Instead, the output of the chunk when it was last run in the editor is displayed.