R01 STA1341

1. Package

Menginstall dan mengaktivasi package yang akan digunakan:

##Cara Install
#install.packages("TSA")

#Cara Aktifkan
library(forecast)
library(graphics)
library(TTR)
library(TSA)

2. Import File ke R dari Berbagai Ekstensi

Format xlsx

library(readxl)
data1 <- read_excel("D:/MATERI KULIAH S2 IPB/ASPRAK/RESPONSI 1/Data_P1.xlsx", sheet = "Sheet1")
head(data1)
## # A tibble: 6 x 1
##      Yt
##   <dbl>
## 1  48.7
## 2  45.8
## 3  46.4
## 4  46.2
## 5  44  
## 6  53.8

Format .csv

data2 <- read.csv("D:/MATERI KULIAH S2 IPB/ASPRAK/RESPONSI 1/Data P1.csv", sep = ";", header = TRUE)
head(data2)
##     Yt
## 1 48.7
## 2 45.8
## 3 46.4
## 4 46.2
## 5 44.0
## 6 53.8

Format .sav

library(foreign)
data3 <- read.spss("D:/MATERI KULIAH S2 IPB/ASPRAK/RESPONSI 1/Data P1.sav", to.data.frame = TRUE)
head(data3)
##     Yt
## 1 48.7
## 2 45.8
## 3 46.4
## 4 46.2
## 5 44.0
## 6 53.8

Format .txt

data4 <- read.table("D:/MATERI KULIAH S2 IPB/ASPRAK/RESPONSI 1/data1.txt", header = F)
head(data4)
##     V1
## 1 48.7
## 2 45.8
## 3 46.4
## 4 46.2
## 5 44.0
## 6 53.8

3. Mendefinisikan Data Deret Waktu

#Melihat Data
data1
## # A tibble: 50 x 1
##       Yt
##    <dbl>
##  1  48.7
##  2  45.8
##  3  46.4
##  4  46.2
##  5  44  
##  6  53.8
##  7  47.6
##  8  47  
##  9  47.6
## 10  51.1
## # ... with 40 more rows
## # i Use `print(n = ...)` to see more rows
#Struktur Data
str(data1)
## tibble [50 x 1] (S3: tbl_df/tbl/data.frame)
##  $ Yt: num [1:50] 48.7 45.8 46.4 46.2 44 53.8 47.6 47 47.6 51.1 ...
#Dimensi Data
dim(data1)
## [1] 50  1
#Menjadikan data terbaca sebagai data deret waktu
data1.ts <- ts(data1)
head(data1.ts)
## Time Series:
## Start = 1 
## End = 6 
## Frequency = 1 
##        Yt
## [1,] 48.7
## [2,] 45.8
## [3,] 46.4
## [4,] 46.2
## [5,] 44.0
## [6,] 53.8
#Menambah keterangan waktu (misal: bulan)
data2.ts <- ts(data1, start = c(2020,1), frequency = 12)
head(data2.ts)
##       Jan  Feb  Mar  Apr  May  Jun
## 2020 48.7 45.8 46.4 46.2 44.0 53.8
#ringkasan data
summary(data1.ts)
##        Yt       
##  Min.   :43.30  
##  1st Qu.:46.40  
##  Median :48.65  
##  Mean   :48.92  
##  3rd Qu.:51.55  
##  Max.   :54.80

4. Eksplorasi Data Deret Waktu

#time series plot 
ts.plot(data1.ts, xlab="Waktu", ylab="Yt", main="Time Series Plot")
points(data1.ts)

#time series plot 
ts.plot(data2.ts, xlab="Waktu", ylab="Yt", main="Time Series Plot")
points(data2.ts)

Contoh lain:

data1_1 <- read_excel("D:/MATERI KULIAH S2 IPB/ASPRAK/RESPONSI 1/Data P1_1.xlsx", sheet = "Sheet1")
head(data1_1)
## # A tibble: 6 x 4
##    week aktual pemulusan peramalan
##   <dbl>  <dbl>     <dbl>     <dbl>
## 1     1 10618.       NA        NA 
## 2     2 10538.       NA        NA 
## 3     3 10209.    10455.       NA 
## 4     4 10553     10433.    10455.
## 5     5  9935.    10232.    10433.
## 6     6 10534.    10341.    10232.
summary(data1_1)
##       week            aktual        pemulusan       peramalan    
##  Min.   :  1.00   Min.   : 9815   Min.   :10078   Min.   :10078  
##  1st Qu.: 30.75   1st Qu.:10210   1st Qu.:10267   1st Qu.:10268  
##  Median : 60.50   Median :10392   Median :10398   Median :10399  
##  Mean   : 60.50   Mean   :10379   Mean   :10379   Mean   :10380  
##  3rd Qu.: 90.25   3rd Qu.:10535   3rd Qu.:10477   3rd Qu.:10477  
##  Max.   :120.00   Max.   :10827   Max.   :10663   Max.   :10663  
##                                   NA's   :2       NA's   :3
#mengubah data ke dalam time series
data1_1.ts <- ts(data1_1$aktual)
data1_2.ts <- ts(data1_1$pemulusan)
data1_3.ts <- ts(data1_1$peramalan)
ts.plot(data1_1.ts, xlab="Waktu", ylab="Yt", main= "Time Series Plot", ylim=c(9000,11000))
points(data1_1.ts)
lines(data1_2.ts,col="green",lwd=2)
lines(data1_3.ts,col="red",lwd=2)
legend("bottomleft",c("data aktual","data pemulusan","data peramalan"), lty=8,
col=c("black","green","red"), cex=0.8)

5. Single Moving Average (SMA)

Data

data1 <- read_excel("D:/MATERI KULIAH S2 IPB/ASPRAK/RESPONSI 1/Data_P1.xlsx", sheet = "Sheet1")
head(data1)
## # A tibble: 6 x 1
##      Yt
##   <dbl>
## 1  48.7
## 2  45.8
## 3  46.4
## 4  46.2
## 5  44  
## 6  53.8
data1.ts <- ts(data1)
head(data1.ts)
## Time Series:
## Start = 1 
## End = 6 
## Frequency = 1 
##        Yt
## [1,] 48.7
## [2,] 45.8
## [3,] 46.4
## [4,] 46.2
## [5,] 44.0
## [6,] 53.8

SMA N=3

#Single Moving  Average dengan N=3
df_ts_sma <- TTR::SMA(data1.ts,  n=3)
ramal_sma <- c(NA,df_ts_sma)

df_sma <- cbind(df_aktual=c(data1.ts,NA), pemulusan=c(df_ts_sma,NA), ramal_sma)
head(df_sma)
##      df_aktual pemulusan ramal_sma
## [1,]      48.7        NA        NA
## [2,]      45.8        NA        NA
## [3,]      46.4  46.96667        NA
## [4,]      46.2  46.13333  46.96667
## [5,]      44.0  45.53333  46.13333
## [6,]      53.8  48.00000  45.53333
tail(df_sma)
##       df_aktual pemulusan ramal_sma
## [46,]      52.0  52.76667  50.10000
## [47,]      50.6  51.70000  52.76667
## [48,]      48.7  50.43333  51.70000
## [49,]      51.4  50.23333  50.43333
## [50,]      47.7  49.26667  50.23333
## [51,]        NA        NA  49.26667

Plot SMA N=3

#Plot 
ts.plot(data1.ts, xlab="periode  waktu", ylab="Yt", col="blue", lty=3, ylim=c(40,65))
points(data1.ts)
lines(df_ts_sma, col="red", lwd=2)
lines(ramal_sma, col="black", lwd= 2)
title("Rataan Bergerak Sederhana N=3", cex.main=1, font.main=4 ,col.main="black")
legend("topleft", c("Data aktual","Pemulusan SMA","Ramalan SMA"),lty=1:3,col=c ("blue","red","black"))

Ukuran Keakuratan Ramalan

# Ukuran Keakuratan
error.sma <- df_sma[, 1] - df_sma[, 3]

## SSE (Sum Square Error)
SSE.sma <- sum(error.sma^2, na.rm = T)

## MSE (Mean Squared Error)
MSE.sma <- mean(error.sma^2, na.rm = T)

## RMSE (Root Mean Square Error)
RMSE.sma <- sqrt(mean(error.sma^2, na.rm = T))

## MAD (Mean Absolute Deviation)
MAD.sma <- mean(abs(error.sma), na.rm = T)

## MAPE (Mean Absolute Percentage Error)
r.error.sma <- (error.sma/df_sma[, 1])*100 # Relative Error
MAPE.sma <- mean(abs(r.error.sma), na.rm = T)


akurasi <- data.frame(
          "Ukuran Keakuratan" = c("SSE", "MSE", "MAPE", "RMSE", "MAD"), 
          "Simple Moving Average N=2" = c(SSE.sma, MSE.sma, MAPE.sma, RMSE.sma, MAD.sma))
akurasi
##   Ukuran.Keakuratan Simple.Moving.Average.N.2
## 1               SSE                545.863333
## 2               MSE                 11.614113
## 3              MAPE                  5.562152
## 4              RMSE                  3.407949
## 5               MAD                  2.739716

Single Moving Average dengan fungsi sma

library(smooth)
data1ts_sma<-sma(data1.ts, order=3, ic=c("AIC","BIC"))
summary(data1ts_sma)
## Time elapsed: 0.01 seconds
## Model estimated: SMA(3)
## Initial values were produced using backcasting.
## 
## Loss function type: MSE; Loss function value: 11.0842
## Error standard deviation: 3.398
## Sample size: 50
## Number of estimated parameters: 2
## Number of degrees of freedom: 48
## Information criteria:
##      AIC     AICc      BIC     BICc 
## 266.1701 266.4254 269.9941 270.4935
forecast(data1ts_sma)
## Time Series:
## Start = 51 
## End = 60 
## Frequency = 1 
##    Point forecast Lower bound (2.5%) Upper bound (97.5%)
## 51       49.26667           42.43462            56.09871
## 52       49.45556           42.25395            56.65716
## 53       48.80741           40.99183            56.62298
## 54       49.17654           40.37458            57.97851
## 55       49.14650           39.80766            58.48535
## 56       49.04348           39.10430            58.98267
## 57       49.12218           38.57682            59.66753
## 58       49.10405           38.03973            60.16838
## 59       49.08990           37.50763            60.67218
## 60       49.10538           37.02518            61.18558
plot(forecast(data1ts_sma))

6. Double Moving Average (DMA)

DMA N=3

#Single Moving  Average dengan N=3
df_ts_sma <- TTR::SMA(data1.ts,  n=3)

#Double Moving  Average dengan N=3
df_ts_dma <- TTR::SMA(df_ts_sma, n=3)
At <- 2*df_ts_sma-df_ts_dma
Bt <- df_ts_sma-df_ts_dma
pemulusan_dma <- At+Bt
ramal_dma <- c(NA, pemulusan_dma)
df_dma <- cbind(df_aktual=c(data1.ts,NA), pemulusan_dma=c(pemulusan_dma,NA), ramal_dma)
head(df_dma)
##      df_aktual pemulusan_dma ramal_dma
## [1,]      48.7            NA        NA
## [2,]      45.8            NA        NA
## [3,]      46.4            NA        NA
## [4,]      46.2            NA        NA
## [5,]      44.0      44.17778        NA
## [6,]      53.8      50.88889  44.17778
tail(df_dma)
##       df_aktual pemulusan_dma ramal_dma
## [46,]      52.0      56.98889  52.74444
## [47,]      50.6      52.05556  56.98889
## [48,]      48.7      48.03333  52.05556
## [49,]      51.4      49.12222  48.03333
## [50,]      47.7      47.84444  49.12222
## [51,]        NA            NA  47.84444

Plot DMA

#Plot
ts.plot(data1.ts, xlab="periode waktu", ylab="Yt", col="blue", lty=3, ylim=c(40,65))
points(data1.ts)
lines(pemulusan_dma,col="red",lwd=2)
lines(ramal_dma,col="black",lwd= 2)
title("Rataan Bergerak  Berganda N=3", cex.main=1, font.main=4 ,col.main="black")
legend("topleft", c("Data Aktual","Pemulusan","Ramalan"),lty=1:3,col=c ("blue","red","black"))

## Ukuran Keakuratan Ramalan DMA 3

# Ukuran Keakuratan
error.dma <- df_dma[, 1] - df_dma[, 3]

## SSE (Sum Square Error)
SSE.dma <- sum(error.dma^2, na.rm = T)

## MSE (Mean Squared Error)
MSE.dma <- mean(error.dma^2, na.rm = T)

## RMSE (Root Mean Square Error)
RMSE.dma <- sqrt(mean(error.dma^2, na.rm = T))

## MAD (Mean Absolute Deviation)
MAD.dma <- mean(abs(error.dma), na.rm = T)

## MAPE (Mean Absolute Percentage Error)
r.error.dma <- (error.dma/df_dma[, 1])*100 # Relative Error
MAPE.dma <- mean(abs(r.error.dma), na.rm = T)


akurasi <- data.frame(
          "Ukuran Keakuratan" = c("SSE", "MSE", "MAPE", "RMSE", "MAD"), 
          "Double Moving Average N=2" = c(SSE.dma, MSE.dma, MAPE.dma, RMSE.dma, MAD.dma))
akurasi
##   Ukuran.Keakuratan Double.Moving.Average.N.2
## 1               SSE                874.330370
## 2               MSE                 19.429564
## 3              MAPE                  7.126068
## 4              RMSE                  4.407898
## 5               MAD                  3.511111

Perbandingan SMA dan DMA

Plot SMA dan DMA N=3

#perbandingan SMA dan DMA
ts.plot(data1.ts, xlab="periode waktu", ylab="Yt", col="blue", lty=3, ylim=c(40,65))
points(data1.ts)
lines(pemulusan_dma,col="red",lwd=2)
lines(df_ts_sma,col="black",lwd= 2)
title("Perbandingan SMA dan DMA",cex.main=1,font.main=4 ,col.main="black")
legend("topleft", c("Data aktual","Pemulusan SMA","Pemulusan DMA"),lty=1:3,col=c ("blue","black","red"))

Ukuran Keakuratan Ramalan SMA dan DMA

akurasi <- data.frame(
          "Ukuran Keakuratan" = c("SSE", "MSE", "MAPE", "RMSE", "MAD"), 
          "Simple Moving Average N=3" = c(SSE.sma, MSE.sma, MAPE.sma, RMSE.sma, MAD.sma),
          "Double Moving Average N=3" = c(SSE.dma, MSE.dma, MAPE.dma, RMSE.dma, MAD.dma))
akurasi
##   Ukuran.Keakuratan Simple.Moving.Average.N.3 Double.Moving.Average.N.3
## 1               SSE                545.863333                874.330370
## 2               MSE                 11.614113                 19.429564
## 3              MAPE                  5.562152                  7.126068
## 4              RMSE                  3.407949                  4.407898
## 5               MAD                  2.739716                  3.511111