Data Tidak Stasioner

Angga Fathan Rofiqy

02 October, 2023

Kode di Hide dalam default, untuk menampilkan kode, klik Code .

#                      -=( Install & Load Package Function )=-
install_load <- function (package1, ...)  {   

   # convert arguments to vector
   packages <- c(package1, ...)

   # start loop to determine if each package is installed
   for(package in packages){

       # if package is installed locally, load
       if(package %in% rownames(installed.packages()))
          do.call('library', list(package))

       # if package is not installed locally, download, then load
       else {
          install.packages(package)
          do.call("library", list(package))
       }
   } 
}
path <- function(){
  gsub  ( "\\\\",  "/",  readClipboard ()  )
}
#Copy path, Panggil function di console
#Copy r path, paste ke var yang diinginkan
#Export chart
export.chart <- "C:/Users/Fathan/Documents/Obsidian Vault/2. Kuliah/Smt 5/6. Metode Peramalan Deret Waktu/@Proj/STA1341-MPDW/Pertemuan 5/Chart"

1 A. Tentang Data

Dataset ini dibangun dari dataset awal yang terdiri dari \(600.000\) data transaksi yang dikumpulkan dalam waktu 6 tahun (periode 2014-2019), yang mencatat tanggal dan waktu penjualan, nama merek obat farmasi, dan jumlah yang terjual. Data ini diekspor dari sistem Point-of-Sale di apotek individu. Kelompok obat yang dipilih dari dataset (57 obat) diklasifikasikan ke dalam kategori berikut dalam Sistem Klasifikasi Anatomi Terapeutik Kimia (ATC):

  1. M01AB - Produk antiinflamasi dan antirheumatik, non-steroid, Turunan asam asetat dan substansi terkait
  2. M01AE - Produk antiinflamasi dan antirheumatik, non-steroid, Turunan asam propionat
  3. N02BA - Analgesik dan antipiretik lainnya, Asam salisilat dan turunannya
  4. N02BE/B - Analgesik dan antipiretik lainnya, Pirazolona dan Anilida
  5. N05B - Obat psikoleptik, Obat ansiolitik
  6. N05C - Obat psikoleptik, Obat hipnotik dan sedatif
  7. R03 - Obat untuk penyakit saluran napas obstruktif
  8. R06 - Antihistamin untuk penggunaan sistemik

Data penjualan diambil ulang dalam periode per jam, per hari, per minggu, dan per bulan. Data sudah diproses sebelumnya, termasuk deteksi dan penanganan outlier serta pengisian data yang hilang. Data yang akan saya gunakan adalah data dalam periode per bulan.

Beberapa waktu yang lalu saya merasa sulit untuk tidur. Kebetulan ada obat N05C yang merupakan obat tidur dalam dataset saya. N05C merupakan obat penenang, sifat tenang ini yang mungkin membantu seseorang untuk tidur. Ternyata ada obat tipe penenang lainnya, yakni N05B.

Kedua obat ini memiliki kesamaan yakni termasuk dalam klasifikasi ATC N (Nervous System), yang berarti keduanya berhubungan dengan sistem saraf. Kedua obat ini termasuk dalam kelompok “Psycholeptics drugs,” yang berarti keduanya dapat memengaruhi sistem saraf dan suasana hati pasien.

Namun ada beberapa perbedaan, diantaranya: N05B lebih fokus pada obat penenang dan pengobatan gangguan kecemasan. Sementara N05C lebih fokus pada obat tidur dan pengobatan gangguan tidur. Contoh obat-obatan yang termasuk dalam kategori ini berbeda, dengan N05B memiliki contoh obat-obatan penenang seperti benzodiazepin, sementara N05C memiliki obat tidur seperti zolpidem.

Tujuan Praktikum kali ini adalah membuat Model Regresi dengan Peubah Lag terbaik. Lag secara bahasa sendiri memiliki arti keterlambatan. Bisa dikatakan Lag mengacu pada nilai suatu variabel pada waktu sebelumnya dalam rangkaian data waktu.

Kali ini saya akan melihat “Apakah penjualan obat penenang N05B pada periode waktu saat inidipengaruhi oleh penjualan obat tidur N05C pada periode waktu sebelumnya?” Atau “Apakah penjualan obat tidur N05Cperiodesebelumnya memiliki pengaruh signifikan terhadap penjualan obat penenang N05B pada periode waktu saat ini ?.

1.1 Data Preparation

1.1.1 Import Data

install_load('rio')
raw.data <- import("https://raw.githubusercontent.com/Zen-Rofiqy/STA1341-MPDW/main/Data/MAANG%20Stock%20Prices.csv")

1.1.2 Data Checking

Cek Tipe data.

str(raw.data)
## 'data.frame':    8812 obs. of  8 variables:
##  $ Name     : chr  "AMZN" "AMZN" "AMZN" "AMZN" ...
##  $ Date     : chr  "7/30/18" "7/31/18" "8/1/18" "8/2/18" ...
##  $ Open     : chr  "91.366501" "89.324501" "89.199997" "89.438499" ...
##  $ High     : chr  "91.474998" "90.091499" "89.921997" "91.828003" ...
##  $ Low      : chr  "88.301003" "86.966003" "88.801003" "89.300003" ...
##  $ Close    : chr  "88.960999" "88.872002" "89.858498" "91.716499" ...
##  $ Adj Close: chr  "88.960999" "88.872002" "89.858498" "91.716499" ...
##  $ Volume   : chr  "131246000" "114774000" "83062000" "87094000" ...

Semua data Karakter, harus diubah.

Cek Data kosong.

sum(is.na(raw.data))
## [1] 0

Tidak ada data kosong.

1.1.3 Penyesuaian Tipe Data

Semua tipe data masih berupa character. Harus diubah menjadi tipe data yang sesuai.

install_load('dplyr')
data <- raw.data %>%  
  mutate(
    Date = as.Date(raw.data[, 2], format = "%m/%d/%y"), #Mengubah menjadi Date 
    across(3:ncol(raw.data), as.numeric)                #Mengubah menjadi Numerik
  )
str(data)
## 'data.frame':    8812 obs. of  8 variables:
##  $ Name     : chr  "AMZN" "AMZN" "AMZN" "AMZN" ...
##  $ Date     : Date, format: "2018-07-30" "2018-07-31" ...
##  $ Open     : num  91.4 89.3 89.2 89.4 91.9 ...
##  $ High     : num  91.5 90.1 89.9 91.8 92.1 ...
##  $ Low      : num  88.3 87 88.8 89.3 91.1 ...
##  $ Close    : num  89 88.9 89.9 91.7 91.2 ...
##  $ Adj Close: num  89 88.9 89.9 91.7 91.2 ...
##  $ Volume   : num  1.31e+08 1.15e+08 8.31e+07 8.71e+07 6.92e+07 ...

1.1.4 Rechecking Data 

Cek kembali data kosong.

cat('Banyaknya Data Kosong', sum(is.na(data)))
## Banyaknya Data Kosong 42

Melihat baris, kolom mana data yang kosong.

# Mencari indeks baris dan kolom yang mengandung NA
na.idx <- which(is.na(data), arr.ind = TRUE)

# Menampilkan data raw dengan baris dan kolom yang mengandung NA
install_load('kableExtra','DT')
## 
## Attaching package: 'kableExtra'
## The following object is masked from 'package:dplyr':
## 
##     group_rows
## Warning: package 'DT' was built under R version 4.2.3
datatable(raw.data[                        # Subsetting
                unique(na.idx[, 1]),   # Vektor indeks baris yang mengandung NA
                unique(na.idx[, 2])  ] # Vektor indeks kolom yang mengandung NA
          )

Ternyata pada baris tersebut ada data karakter text yang merupakan label dari tiap kolomnya. Sehingga ketika diubah ke numerik akan menjadi NA. Maka saya akan menghapus baris tersebut.

1.1.5 Data Cleaned

data <- data %>%
  filter(!row_number() %in% unique(na.idx[, 1]))
datatable(data, filter = 'top', 
          options = list(pageLength = 5))

1.2 Pemilihan Variabel Saham

1.2.1 Time Series MAANG

Melihat keseluruhan Time Series data saham.

install_load('ggplot2','extrafont')
# font_import(); loadfonts() #Run ini sekali aja
theme1 <- list(
  guides(fill="none"), #No Legends
  theme(
  text = element_text(size = 33),
  axis.title = element_text(size=15),
  axis.title.x = element_text(size=33),
  axis.title.y = element_text(size=33),
  axis.text.y = element_text(vjust = .5, face = "bold"),
  plot.title = element_text(hjust = 0.5, size=40),
  panel.background = element_rect(fill = 'transparent'),
  plot.background = element_rect(fill='transparent', color=NA),
  panel.grid.major = element_line(colour = "grey90"),
  axis.line = element_line(linewidth = 2, colour = "grey90"))
)
install_load('viridis','ggrepel')
#Plot
cts.maang <-
ggplot(data, aes(x=Date, y=`Adj Close`)) + #Data
  geom_line(aes(color=Name), linewidth=1) + #Timeseries
  #Color
    scale_color_viridis(alpha = 0.75, #Opacity
                     begin = 0, #Color pallte scale begins
                     end = 0.9, #Color pallte scale ends
                     direction = -1, #Flip color scale
                     discrete = T, #Discrete Value
                     option = "D") + #Color Palette
  theme1 + #THeme
  labs(x = "\nTahun", y = "Harga Saham (USD)\n") + #Label X & Y
  ggtitle("Time Series MAANG") +
  # Label / legend
  geom_text_repel(
    data=data[data$Date == max(data$Date),], #Posisi di ujung data
    aes(color = Name, label = Name), #Warna garis & label saham
    size = 8, #Ukuran text
    nudge_x = 80, #Posisi Text (kanan 50)
    hjust = 0, #Ujung
    segment.size = 1,               #Ukuran garis
    segment.alpha = .75,             #transparasi garis
    segment.linetype = "dotted",    #Time garis
    box.padding = .4, #Biar label saham nggak dempetan
    segment.curvature = -0.1, #biar garis mulus
    segment.ncp = 8, 
    segment.angle = 60 
  ) +
  #Axis
    coord_cartesian(clip = "off"
  ) +
    scale_x_date( #Sumbu x
    date_breaks = "1 year",  # Menampilkan label setiap tahun
    date_labels = "%Y",  # Format label tahun
    limits = c(as.Date("2018-07-30"), as.Date("2023-12-28"))
    #Tampilin lebih dari 20023-07-28 agar label saham bisa masuk
  ) +
    scale_y_continuous( #Sumbu y
    labels = scales::dollar_format(prefix = "$") #tambahin dolar
  ) +
    annotate( #Buat nandain batas data
    "text", x = as.Date("2023-7-28"), y = 50, 
    label = "28 Juli", size=6
  ) +
  geom_vline( #Buat garis batas data
    xintercept = as.numeric(as.Date("2023-07-28")), 
             linetype = "dotted", color = "red")
cts.maang

#Export Chart
ggsave("01_Time Series MAANG.png", cts.maang, path = export.chart,
        dpi = 300, height = 9, width = 16)

Jika dilihat dari tahun 2019-2022, semua saham cenderung memiliki pola trend naik. Lalu dari 2021-2023 polanya cenderung trend turun. Untuk tugas praktikum kali ini, saya hanya akan menggunakan rentang tahun 2022-2023 dengan tren cenderung turun. Agar pengerjaannya tidak terlalu sulit, karena masih tahap awal pembelajaran.

Ada yang menarik perhatian saya. Kenapa dulu Pendiri Amazon, Jeff Bezos yang pernah menjadi orang terkaya di dunia pada tahun 2017 lalu, harga saham sekarang tidak setinggi yang saya kira?. Oleh karena itu saya memutuskan untuk menggunakan data saham Amazon untuk praktikum kali ini.

amzn <- data %>%
  select(1, 2, 7) %>%  # Memilih kolom 1, 2, dan 7
  filter(Name == "AMZN", Date >= as.Date("2022-01-01"))  # Filter data saham Amazon tahun 2022 ke atas

rownames(amzn) <- NULL
str(amzn)
## 'data.frame':    394 obs. of  3 variables:
##  $ Name     : chr  "AMZN" "AMZN" "AMZN" "AMZN" ...
##  $ Date     : Date, format: "2022-01-03" "2022-01-04" ...
##  $ Adj Close: num  170 168 164 163 163 ...
datatable(amzn)

Mengubah Ajd Close Menjadi Time series.

amzn.ts <- ts(amzn[,3])

Ringkasan Data Ajd CLose.

summary(amzn.ts)
##    Min. 1st Qu.  Median    Mean 3rd Qu.    Max. 
##   81.82  102.01  115.77  119.81  134.88  170.40

1.2.2 Time Series Amazon

ts.plot(amzn.ts, xlab="Time Period", ylab="Harga Saham", 
        main = "Time Series Amazon", col='orange', lwd=2)
points(amzn.ts, col='orange', lwd=1.5)

library(ggplot2)
library(tsibble)
## Warning: package 'tsibble' was built under R version 4.2.3
## 
## Attaching package: 'tsibble'
## The following objects are masked from 'package:base':
## 
##     intersect, setdiff, union
library(tseries)
## Warning: package 'tseries' was built under R version 4.2.3
## Registered S3 method overwritten by 'quantmod':
##   method            from
##   as.zoo.data.frame zoo

2 Stasioner dalam Rataan dan Ragam

2.1 Plot Time Series

amzn <- cbind(index = 1:nrow(amzn), amzn)
install_load('ggplot2')
plot_stas <- amzn |> 
  ggplot(aes(x = index, y = `Adj Close`)) + geom_line() + theme_bw() +
  xlab("Bulan") + ylab("Penjualan Obat N02BA")
plot_stas

mean(amzn.ts)
## [1] 119.8065

Plot deret waktu di atas menunjukkan bahwa data tidak stasioner dalam rataan, ditandai dengan data yang tidak menyebar di sekitar nilai tengahnya (119.8065) dan tidak stasioner dalam ragam, ditandai dengan lebar pita yang cenderung tidak sama.

2.2 Plot ACF

acf(amzn.ts)

Berdasarkan plot ACF, terlihat bahwa plot ACF pada data tersebut cenderung tails off dan tidak membentuk gelombang sinus.

2.3 Uji ADF

tseries::adf.test(amzn.ts)
## 
##  Augmented Dickey-Fuller Test
## 
## data:  amzn.ts
## Dickey-Fuller = -1.7347, Lag order = 7, p-value = 0.6893
## alternative hypothesis: stationary

\(H_0\) : Data tidak stasioner dalam rataan

\(H_1\) : Data stasioner dalam rataan

Berdasarkan uji ADF tersebut, didapat p-value sebesar 0.6893 yang lebih besar dari taraf nyata 5% sehingga tak tolak \(H_0\) dan menandakan bahwa data tidak stasioner dalam rataan. Hal ini sesuai dengan hasil eksplorasi menggunakan plot time series dan plot ACF.

2.4 Plot Box-Cox

install_load('MASS')
## 
## Attaching package: 'MASS'
## The following object is masked from 'package:dplyr':
## 
##     select
index <- seq(1:nrow(amzn))
bc = boxcox(amzn.ts~index, lambda = seq(0,4,by=0.01))

#Nilai Rounded Lambda
lambda <- bc$x[which.max(bc$y)]
lambda
## [1] 0.21
#SK
bc$x[bc$y > max(bc$y) - 1/2 * qchisq(.95,1)]
##  [1] 0.00 0.01 0.02 0.03 0.04 0.05 0.06 0.07 0.08 0.09 0.10 0.11 0.12 0.13 0.14
## [16] 0.15 0.16 0.17 0.18 0.19 0.20 0.21 0.22 0.23 0.24 0.25 0.26 0.27 0.28 0.29
## [31] 0.30 0.31 0.32 0.33 0.34 0.35 0.36 0.37 0.38 0.39 0.40 0.41 0.42 0.43 0.44
## [46] 0.45 0.46 0.47 0.48 0.49 0.50 0.51 0.52 0.53 0.54 0.55 0.56 0.57 0.58 0.59
## [61] 0.60 0.61 0.62 0.63 0.64 0.65 0.66 0.67 0.68 0.69 0.70 0.71 0.72 0.73 0.74

Gambar di atas menunjukkan nilai rounded value (\(\lambda\)) optimum sebesar 0.21 dan pada selang kepercayaan 95% nilai memiliki batas bawah 0 dan batas atas 0.74. Selang tersebut tidak memuat nilai satu sehingga dapat dikatakan bahwa data saham amazon tidak stasioner dalam ragam.

2.5 Partisi Data

2.5.1 Bagian 1

dt_amzn.ts1 <- amzn.ts[1:150] |> ts()
mean(dt_amzn.ts1 )
## [1] 136.7386
var(dt_amzn.ts1 )
## [1] 442.5268

2.5.1.1 Plot Time Series

install_load('tsibble')
dt_amzn.ts1 |> as_tsibble() |> 
  ggplot(aes(x = index, y = value)) +
  geom_line() + theme_bw() +
  xlab("Obs") + ylab("Nilai")

mean(amzn.ts)
## [1] 119.8065

Plot deret waktu di atas menunjukkan bahwa data tidak stasioner dalam rataan, ditandai dengan data yang tidak menyebar di sekitar nilai tengahnya (119.8065) dan stasioner dalam ragam, ditandai dengan lebar pita yang cenderung sama.

2.5.1.2 Plot ACF

acf(dt_amzn.ts1)

Berdasarkan plot ACF, terlihat bahwa plot ACF pada data tersebut cenderung tails off dan tidak membentuk gelombang sinus.

2.5.1.3 Uji ADF

tseries::adf.test(dt_amzn.ts1)
## 
##  Augmented Dickey-Fuller Test
## 
## data:  dt_amzn.ts1
## Dickey-Fuller = -1.5598, Lag order = 5, p-value = 0.7598
## alternative hypothesis: stationary

\(H_0\) : Data tidak stasioner dalam rataan

\(H_1\) : Data stasioner dalam rataan

Berdasarkan uji ADF tersebut, didapat p-value sebesar 0.7598 yang lebih besar dari taraf nyata 5% sehingga tak tolak \(H_0\) dan menandakan bahwa data tidak stasioner dalam rataan. Hal ini sesuai dengan hasil eksplorasi menggunakan plot time series dan plot ACF.

2.5.1.4 Plot Boxcox

index <- seq(1:150)
bc = boxcox(dt_amzn.ts1~index, lambda = seq(-2,6,by=1))

#Nilai Rounded Lambda
lambda <- bc$x[which.max(bc$y)]
lambda
## [1] 1.555556
#SK
bc$x[bc$y > max(bc$y) - 1/2 * qchisq(.95,1)]
##  [1] 0.5050505 0.5858586 0.6666667 0.7474747 0.8282828 0.9090909 0.9898990
##  [8] 1.0707071 1.1515152 1.2323232 1.3131313 1.3939394 1.4747475 1.5555556
## [15] 1.6363636 1.7171717 1.7979798 1.8787879 1.9595960 2.0404040 2.1212121
## [22] 2.2020202 2.2828283 2.3636364 2.4444444 2.5252525 2.6060606

Gambar di atas menunjukkan nilai rounded value (\(\lambda\)) optimum sebesar 1.555556 dan pada selang kepercayaan 95% nilai memiliki batas bawah 0.50505 dan batas atas 2.6060606. Selang tersebut memuat nilai satu sehingga dapat dikatakan bahwa data saham amazon stasioner dalam ragam.

2.5.2 Bagian 2

dt_stas2 <- amzn.ts[1:310] |> ts()
mean(dt_stas2)
## [1] 120.4367
var(dt_stas2)
## [1] 604.008

2.5.2.1 Plot Time Series

dt_stas2 |> as_tsibble() |> 
  ggplot(aes(x = index, y = value)) +
  geom_line() + theme_bw() +
  xlab("Obs") + ylab("Nilai")

mean(amzn.ts)
## [1] 119.8065

Plot deret waktu di atas menunjukkan bahwa data tidak stasioner dalam rataan, ditandai dengan data yang menyebar di sekitar nilai tengahnya (18) dan stasioner dalam ragam, ditandai dengan lebar pita yang cenderung sama.

2.5.2.2 Plot ACF

acf(dt_stas2)

Berdasarkan plot ACF, terlihat bahwa plot ACF pada data tersebut cenderung tails off dan membentuk gelombang sinus.

2.5.2.3 Uji ADF

adf.test(dt_stas2)
## 
##  Augmented Dickey-Fuller Test
## 
## data:  dt_stas2
## Dickey-Fuller = -2.7608, Lag order = 6, p-value = 0.2556
## alternative hypothesis: stationary

\(H_0\) : Data tidak stasioner dalam rataan

\(H_1\) : Data stasioner dalam rataan

Berdasarkan uji ADF tersebut, didapat p-value sebesar 0.01 yang lebih kecil dari taraf nyata 5% sehingga tolak \(H_0\) dan menandakan bahwa data stasioner dalam rataan. Hal ini sesuai dengan hasil eksplorasi menggunakan plot time series dan plot ACF.

2.5.2.4 Plot Boxcox

index <- seq(1:310)
bc = boxcox(dt_stas2~index, lambda = seq(0,6,by=1))

#Nilai Rounded Lambda
lambda <- bc$x[which.max(bc$y)]
lambda
## [1] 0
#SK
bc$x[bc$y > max(bc$y) - 1/2 * qchisq(.95,1)]
## [1] 0.00000000 0.06060606 0.12121212 0.18181818 0.24242424 0.30303030

Gambar di atas menunjukkan nilai rounded value (\(\lambda\)) optimum sebesar 0 dan pada selang kepercayaan 95% nilai memiliki batas bawah 0 dan batas atas 0.30303. Selang tersebut tidak memuat nilai satu sehingga dapat dikatakan bahwa data saham amazon tidak stasioner dalam ragam.

3 Kesimpulan

Data saham amazon partisi 1 tidak stasioner dalam rataan, namun stasioner dalam ragam. Sedangkan pada partisi 2 keduanya tidak stasioner baik rataan maupun ragam.