Pendahuluan

Sebagian besar sirup maple dunia diproduksi di Kanada Tenggara dan Amerika Serikat bagian timur laut dengan sirup maple diproduksi dalam jumlah kecil di beberapa negara lain - seperti Jepang dan Korea Selatan. Hingga tahun 1930-an, Amerika Serikat memproduksi sebagian besar sirup maple dunia dengan Vermont menjadi produsen dalam negeri terbesar. Namun, setelah pertumbuhan industri yang pesat pada 1990-an, Kanada sekarang memproduksi lebih dari 80 persen sirup maple dunia sementara produsen Amerika berjuang untuk mengimbanginya.

Kumpulan data ini menggunakan data dari laporan tahunan Layanan Statistik Pertanian Nasional Departemen Pertanian Amerika Serikat tentang produksi sirup maple dalam negeri yang diterbitkan setiap tahun pada bulan Juni sebagai bagian dari laporan produksi tanaman bulanan departemen. Kumpulan data dimulai pada tahun 2000 saat pelaporan sirup maple mengalami perubahan besar saat pergantian abad. (D) = ditahan untuk menghindari pengungkapan data untuk operasi individu.

Eksplorasi Data

Garis Besar Data

State : Data produksi negara dilaporkan
Num_of_Taps : Jumlah keran yang digunakan untuk mengekstrak sirup maple (dalam 1000 keran)
Yield_per_Tap : Jumlah sirup maple yang dikumpulkan per keran (dalam galon)
Production : Jumlah galon sirup maple yang diproduksi (dalam 1000 galon)
Avg_Price : Harga rata-rata per galon sirup maple (dalam dolar)
Value : Nilai produksi keseluruhan (dalam 1000 dolar)
Retail_Price : Harga eceran per galon (dalam dolar)
Wholesale_Price : Harga grosir per galon (dalam rupiah)
Bulk_P_Price : Harga massal untuk semua kadar sirup (dalam dolar per pon)
Bulk_G_Price : Harga massal untuk semua kadar sirup (dalam dolar per galon)
Date_Open : Kira-kira hari pertama getah dikumpulkan
Date_Closed : Kira-kira hari terakhir getah dikumpulkan
Year : Tahun Panen

Load Library

library(dplyr) # perselisihan data

## 
## Attaching package: 'dplyr'

## The following objects are masked from 'package:stats':
## 
##     filter, lag

## The following objects are masked from 'package:base':
## 
##     intersect, setdiff, setequal, union

library(ggplot2) # memvisualisasikan data
library(gridExtra) # menampilkan banyak grafik

## 
## Attaching package: 'gridExtra'

## The following object is masked from 'package:dplyr':
## 
##     combine

library(inspectdf) # untuk EDA

## Warning: package 'inspectdf' was built under R version 4.3.1

library(tidymodels) # membangun model yang rapi

## Warning: package 'tidymodels' was built under R version 4.3.1

## ── Attaching packages ────────────────────────────────────── tidymodels 1.1.0 ──

## ✔ broom        1.0.5     ✔ rsample      1.1.1
## ✔ dials        1.2.0     ✔ tibble       3.2.1
## ✔ infer        1.0.4     ✔ tidyr        1.3.0
## ✔ modeldata    1.1.0     ✔ tune         1.1.1
## ✔ parsnip      1.1.0     ✔ workflows    1.1.3
## ✔ purrr        1.0.1     ✔ workflowsets 1.0.1
## ✔ recipes      1.0.6     ✔ yardstick    1.2.0

## Warning: package 'broom' was built under R version 4.3.1

## Warning: package 'dials' was built under R version 4.3.1

## Warning: package 'infer' was built under R version 4.3.1

## Warning: package 'modeldata' was built under R version 4.3.1

## Warning: package 'parsnip' was built under R version 4.3.1

## Warning: package 'purrr' was built under R version 4.3.1

## Warning: package 'recipes' was built under R version 4.3.1

## Warning: package 'rsample' was built under R version 4.3.1

## Warning: package 'tidyr' was built under R version 4.3.1

## Warning: package 'tune' was built under R version 4.3.1

## Warning: package 'workflows' was built under R version 4.3.1

## Warning: package 'workflowsets' was built under R version 4.3.1

## Warning: package 'yardstick' was built under R version 4.3.1

## ── Conflicts ───────────────────────────────────────── tidymodels_conflicts() ──
## ✖ gridExtra::combine() masks dplyr::combine()
## ✖ purrr::discard()     masks scales::discard()
## ✖ dplyr::filter()      masks stats::filter()
## ✖ dplyr::lag()         masks stats::lag()
## ✖ recipes::step()      masks stats::step()
## • Learn how to get started at https://www.tidymodels.org/start/

library(caret) # melakukan pra-proses data

## Warning: package 'caret' was built under R version 4.3.1

## Loading required package: lattice

## 
## Attaching package: 'caret'

## The following objects are masked from 'package:yardstick':
## 
##     precision, recall, sensitivity, specificity

## The following object is masked from 'package:purrr':
## 
##     lift

library(rsample)

Load Data Set

Dataset sirup maple diperoleh dari laporan tahunan Layanan Statistik Pertanian Nasional Departemen Pertanian Amerika Serikat.

maple <- read.csv("data_input/Maple_Syrup.csv", stringsAsFactors = T)
head(maple)

Eksplorasi Data Set

Data Wrangling

Sama seperti artikel sebelumnya, Kami akan membuat model prediksi untuk mengklasifikasikan produsen Vermont, dimana skor kualitas 7-10 dianggap “Sangat Baik”.Oleh karena itu, saya membuat subset data Vermont dan menghapus tipe untuk analisis data yang lebih bersih.

# filter hanya data produsen Vermont dan beberapa variable yang dibuang
maple <- maple %>% 
  filter(State == "Vermont") %>%
  select(-Date_Open, 
         -State,
         -Date_Closed,
         -Year) %>% 
  mutate(Production = as.integer(Production),
         Value = as.numeric(Value),
         Wholesale_Price = as.numeric(Wholesale_Price),
         Bulk_P_Price = as.numeric(Bulk_P_Price),
         Bulk_G_Price = as.numeric(Bulk_G_Price))
maple

maple <- maple %>%
  mutate(Value = as.factor(ifelse(Value > 100, "up", "down")))

# Periksa missing value
round(prop.table(table(is.na(maple))),3)

## 
## FALSE 
##     1

summary(maple)

##   Num_of_Taps   Yield_per_Tap      Production       Avg_Price      Value   
##  Min.   :2030   Min.   :0.1320   Min.   :  1.00   Min.   :27.00   down: 9  
##  1st Qu.:2170   1st Qu.:0.2180   1st Qu.: 41.25   1st Qu.:27.85   up  :13  
##  Median :3225   Median :0.2760   Median : 74.00   Median :30.10            
##  Mean   :3690   Mean   :0.2748   Mean   : 66.68   Mean   :30.93            
##  3rd Qu.:4775   3rd Qu.:0.3412   3rd Qu.: 88.00   3rd Qu.:33.30            
##  Max.   :6500   Max.   :0.4100   Max.   :145.00   Max.   :39.50            
##   Retail_Price   Wholesale_Price   Bulk_P_Price    Bulk_G_Price   
##  Min.   :31.40   Min.   :  8.00   Min.   :12.00   Min.   :  6.00  
##  1st Qu.:32.90   1st Qu.: 33.00   1st Qu.:17.25   1st Qu.: 40.50  
##  Median :44.20   Median : 83.50   Median :30.00   Median : 60.00  
##  Mean   :40.55   Mean   : 70.05   Mean   :34.82   Mean   : 67.27  
##  3rd Qu.:45.45   3rd Qu.: 95.25   3rd Qu.:51.50   3rd Qu.: 98.25  
##  Max.   :47.40   Max.   :138.00   Max.   :66.00   Max.   :120.00

# memeriksa distribusi data variabel numerik
inspect_num(maple)

inspect_num(maple) %>% 
  show_plot()

Masih ada beberapa nilai yang hilang dalam data kami dan saya akan mencoba menghubungkannya dengan nilai median karena data kami tidak terdistribusi secara normal.

# NA median impute
prevalues <- preProcess(maple, method=c("medianImpute"))
maple <- predict(prevalues, maple)
colSums(is.na(maple))

##     Num_of_Taps   Yield_per_Tap      Production       Avg_Price           Value 
##               0               0               0               0               0 
##    Retail_Price Wholesale_Price    Bulk_P_Price    Bulk_G_Price 
##               0               0               0               0

Exploratory Data Analysis

Proporsi Variabel Target

inspectdf::inspect_cat(maple)

Untuk menghindari hilangnya varian, saya akan menggunakan upsampling (bukan downsampling) untuk menyeimbangkan proporsi.

# inspect corelation between predictors
GGally::ggcorr(maple[,-12], hjust = 1, layout.exp = 2, label = T, label_size = 2.9)

## Registered S3 method overwritten by 'GGally':
##   method from   
##   +.gg   ggplot2

## Warning in GGally::ggcorr(maple[, -12], hjust = 1, layout.exp = 2, label = T, :
## data in column(s) 'Value' are not numeric and were ignored

Berdasarkan plot di atas, terdapat beberapa variabel prediktor yang memiliki korelasi tinggi satu sama lain.

Cross Validation

set.seed(417)
split <- initial_split(data = maple, prop = 0.8, strata = "Value")

train <- training(split)
test <- testing(split)

Data Pre-processing

maple_recipe <- recipe(. ~ ., data = train) %>%
  step_nzv(all_predictors()) %>%
  prep()

maple_train <- juice(maple_recipe)
maple_test <- bake(maple_recipe, test)

Pembuatan Model

Naive Bayes

Berdasarkan karakteristiknya kami akan mencoba menggunakan Naive Bayes untuk membangun model nya.

library(e1071)

## Warning: package 'e1071' was built under R version 4.3.1

## 
## Attaching package: 'e1071'

## The following object is masked from 'package:tune':
## 
##     tune

## The following object is masked from 'package:rsample':
## 
##     permutations

## The following object is masked from 'package:parsnip':
## 
##     tune

# model building
naive <- naiveBayes(maple_train[-5], maple_train$Value, laplace = 1)

# Model fitting
naive_pred <- predict(naive, maple_test, type = "class") # for the class prediction
naive_prob <- predict(naive, maple_test, type = "raw") # for the probability

# Result
naive_table <- select(maple_test, Value) %>%
  bind_cols(Value_pred = naive_pred) %>% 
  bind_cols(Value_eprob = round(naive_prob[,1],4)) %>% 
  bind_cols(Value_pprob = round(naive_prob[,2],4))

# Performance evaluation - confusion matrix
naive_table %>% 
  conf_mat(Value, Value_pred) %>% 
  autoplot(type = "heatmap")

naive_table %>%
  summarise(
    accuracy = accuracy_vec(Value, Value_pred),
    sensitivity = sens_vec(Value, Value_pred),
    specificity = spec_vec(Value, Value_pred),
    precision = precision_vec(Value, Value_pred)
  )

# ROC
naive_roc <- data.frame(prediction=round(naive_prob[,1],4),
                      trueclass=as.numeric(naive_table$Value=="up"))
head(naive_roc)

library(ROCR)

## Warning: package 'ROCR' was built under R version 4.3.1

naive_roc <- ROCR::prediction(naive_roc$prediction, naive_roc$trueclass) 

# ROC curve
plot(performance(naive_roc, "fpr", "tpr"),
     main = "ROC")
abline(a = 0, b = 1)

# AUC
auc_ROCR_n <- performance(naive_roc, measure = "auc")
auc_ROCR_n <- auc_ROCR_n@y.values[[1]]
auc_ROCR_n

## [1] 0

Untuk informasi Anda, ini adalah metrik yang akan kami evaluasi dari model: - Accuracy:kemampuan untuk memprediksi dengan benar kedua kelas dari total pengamatan. - Precision:kemampuan untuk memprediksi kelas positif dengan benar dari total kelas prediksi-positif (positif palsu rendah). - Recall:kemampuan untuk memprediksi dengan benar kelas positif dari total kelas aktual-positif (negatif palsu rendah). - Specificity:kemampuan untuk memprediksi dengan benar kelas negatif dari total kelas negatif aktual.

Selain itu, ada juga kurva Receiver Operating Characteristics (ROC) dan Area Under Curve (AUC) untuk mengevaluasi model kami.ROC memplot proporsi rasio positif sebenarnya (TPR atau Sensitivitas) ke proporsi rasio negatif palsu (FNR atau 1-Specificity). ROC adalah kurva probabilitas dan AUC mewakili tingkat atau ukuran keterpisahan.Ini memberi tahu seberapa banyak model mampu membedakan antara kelas.Semakin dekat kurva mencapai kiri atas plot (Positif benar tinggi sedangkan negatif palsu rendah), semakin baik model kita. Semakin tinggi skor AUC, semakin baik model kita memisahkan kelas target kita.Untuk memudahkan pemahaman Anda, Anda dapat melihat ilustrasi di bawah ini.

Berdasarkan hasil tersebut, akurasi, sensitivitas, presisi dan spesifisitas model kami dapat diterima karena model tersebut sudah akurat (100%).Kurva ROC kami juga menunjukkan pemisahan yang cukup baik dengan skor AUC 1.

Dalam kasus kami, Data produksi untuk “Vermont” diklasifikasikan untuk rekomendasi kategorisasi produksi paling bagus.

Kesimpulan

Model klasifikasi yang di bangun telah menunjukkan performa yang luar biasa dengan akurasi, sensitivitas, presisi, dan spesifisitas mencapai 100%. Kurva ROC juga menunjukkan bahwa model dengan sempurna memisahkan kategori produksi paling bagus dari kategori lainnya. Oleh karena itu, model ini dapat diandalkan untuk memberikan rekomendasi kategorisasi produksi paling bagus untuk data produksi di wilayah “Vermont”.

Penting untuk selalu melakukan validasi tambahan pada model, seperti menggunakan data uji yang belum pernah dilihat sebelumnya, untuk memastikan bahwa kinerja yang luar biasa ini tidak disebabkan oleh overfitting atau faktor lain. Juga, pastikan bahwa dataset yang digunakan mewakili populasi atau konteks yang sesuai untuk penerapan model.

LBB: Classification in Machine Learning II

Carli

2023-07-27