library(readr)
## Warning: package 'readr' was built under R version 4.5.3
library(dplyr)
## Warning: package 'dplyr' was built under R version 4.5.3
##
## 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(car)
## Warning: package 'car' was built under R version 4.5.3
## Loading required package: carData
##
## Attaching package: 'car'
## The following object is masked from 'package:dplyr':
##
## recode
library(MVN)
## Warning: package 'MVN' was built under R version 4.5.3
library(biotools)
## Warning: package 'biotools' was built under R version 4.5.3
## Loading required package: MASS
##
## Attaching package: 'MASS'
## The following object is masked from 'package:dplyr':
##
## select
## ---
## biotools version 4.3
library(heplots)
## Warning: package 'heplots' was built under R version 4.5.3
## Loading required package: broom
##
## Attaching package: 'heplots'
## The following object is masked from 'package:biotools':
##
## boxM
data2 <- read_csv("manufacturing_iot_quality_10000.csv")
## Rows: 10000 Columns: 20
## ── Column specification ────────────────────────────────────────────────────────
## Delimiter: ","
## chr (6): plant, line_id, machine_id, shift, material_batch, defect_type
## dbl (14): record_id, operator_experience_years, ambient_temp_c, ambient_humi...
##
## ℹ Use `spec()` to retrieve the full column specification for this data.
## ℹ Specify the column types or set `show_col_types = FALSE` to quiet this message.
# Lihat struktur data
str(data2)
## spc_tbl_ [10,000 × 20] (S3: spec_tbl_df/tbl_df/tbl/data.frame)
## $ record_id : num [1:10000] 1 2 3 4 5 6 7 8 9 10 ...
## $ plant : chr [1:10000] "Plant_A" "Plant_A" "Plant_A" "Plant_C" ...
## $ line_id : chr [1:10000] "L3" "L1" "L4" "L3" ...
## $ machine_id : chr [1:10000] "M07" "M18" "M10" "M01" ...
## $ shift : chr [1:10000] "evening" "night" "evening" "evening" ...
## $ operator_experience_years: num [1:10000] 5.27 1.61 8.35 4.66 6.77 0 9.39 7.05 7.68 1.8 ...
## $ material_batch : chr [1:10000] "B090" "B122" "B071" "B069" ...
## $ ambient_temp_c : num [1:10000] 23.7 19 19.9 21.3 24.9 ...
## $ ambient_humidity_pct : num [1:10000] 67.9 34.8 50.4 34.8 32.7 90 48 55.1 41.9 44.7 ...
## $ vibration_rms : num [1:10000] 0.981 1.29 0.871 0.995 0.697 ...
## $ pressure_bar : num [1:10000] 5.08 4.8 5.64 5.28 5.05 ...
## $ feed_rate_mm_s : num [1:10000] 19.5 26.6 23.9 21.4 15 ...
## $ laser_power_w : num [1:10000] 304 377 329 323 341 ...
## $ coolant_flow_l_min : num [1:10000] 7.63 6.53 6.22 7.47 4.32 5.75 6.7 6.62 8.89 8.05 ...
## $ thickness_mm : num [1:10000] 1.22 1.28 1.22 1.2 1.22 ...
## $ tensile_strength_mpa : num [1:10000] 517 537 539 542 594 ...
## $ surface_roughness_um : num [1:10000] 2.68 2.81 2.14 2.2 1.34 ...
## $ defect_flag : num [1:10000] 0 0 0 0 0 0 0 1 0 0 ...
## $ defect_type : chr [1:10000] "none" "none" "none" "none" ...
## $ rework_needed : num [1:10000] 0 0 0 0 0 0 0 1 0 0 ...
## - attr(*, "spec")=
## .. cols(
## .. record_id = col_double(),
## .. plant = col_character(),
## .. line_id = col_character(),
## .. machine_id = col_character(),
## .. shift = col_character(),
## .. operator_experience_years = col_double(),
## .. material_batch = col_character(),
## .. ambient_temp_c = col_double(),
## .. ambient_humidity_pct = col_double(),
## .. vibration_rms = col_double(),
## .. pressure_bar = col_double(),
## .. feed_rate_mm_s = col_double(),
## .. laser_power_w = col_double(),
## .. coolant_flow_l_min = col_double(),
## .. thickness_mm = col_double(),
## .. tensile_strength_mpa = col_double(),
## .. surface_roughness_um = col_double(),
## .. defect_flag = col_double(),
## .. defect_type = col_character(),
## .. rework_needed = col_double()
## .. )
## - attr(*, "problems")=<externalptr>
head(data2)
## # A tibble: 6 × 20
## record_id plant line_id machine_id shift operator_experience_…¹ material_batch
## <dbl> <chr> <chr> <chr> <chr> <dbl> <chr>
## 1 1 Plan… L3 M07 even… 5.27 B090
## 2 2 Plan… L1 M18 night 1.61 B122
## 3 3 Plan… L4 M10 even… 8.35 B071
## 4 4 Plan… L3 M01 even… 4.66 B069
## 5 5 Plan… L1 M08 even… 6.77 B141
## 6 6 Plan… L3 M13 day 0 B049
## # ℹ abbreviated name: ¹​operator_experience_years
## # ℹ 13 more variables: ambient_temp_c <dbl>, ambient_humidity_pct <dbl>,
## # vibration_rms <dbl>, pressure_bar <dbl>, feed_rate_mm_s <dbl>,
## # laser_power_w <dbl>, coolant_flow_l_min <dbl>, thickness_mm <dbl>,
## # tensile_strength_mpa <dbl>, surface_roughness_um <dbl>, defect_flag <dbl>,
## # defect_type <chr>, rework_needed <dbl>
colSums(is.na(data2))
## record_id plant line_id
## 0 0 0
## machine_id shift operator_experience_years
## 0 0 0
## material_batch ambient_temp_c ambient_humidity_pct
## 0 0 0
## vibration_rms pressure_bar feed_rate_mm_s
## 0 0 0
## laser_power_w coolant_flow_l_min thickness_mm
## 0 0 0
## tensile_strength_mpa surface_roughness_um defect_flag
## 0 0 0
## defect_type rework_needed
## 0 0
# Hapus missing value jika ada
data2_clean <- na.omit(data2)
data2_clean$shift <- as.factor(data2_clean$shift)
data2_clean$tensile_strength_mpa <- as.numeric(data2_clean$tensile_strength_mpa)
data2_clean$surface_roughness_um <- as.numeric(data2_clean$surface_roughness_um)
data2_clean$ambient_humidity_pct <- as.numeric(data2_clean$ambient_humidity_pct)
Y2 <- data2_clean[, c("tensile_strength_mpa",
"surface_roughness_um")]
# Uji Normalitas Univariat
# install.packages("nortest")
library(nortest)
ad.test(data2_clean$tensile_strength_mpa)
##
## Anderson-Darling normality test
##
## data: data2_clean$tensile_strength_mpa
## A = 0.27246, p-value = 0.6695
ad.test(data2_clean$surface_roughness_um)
##
## Anderson-Darling normality test
##
## data: data2_clean$surface_roughness_um
## A = 0.719, p-value = 0.0607
library(MVN)
Y2_clean <- data2_clean[, c("tensile_strength_mpa",
"surface_roughness_um")]
mvn_result <- mvn(
data = Y2_clean,
mvn_test = "mardia"
)
mvn_result$multivariate_normality
## Test Statistic p.value Method MVN
## 1 Mardia Skewness 7.448 0.114 asymptotic ✓ Normal
## 2 Mardia Kurtosis 0.920 0.358 asymptotic ✓ Normal
# Homogenitas Matriks Kovarians
boxM(Y2, data2_clean$shift)
##
## Box's M-test for Homogeneity of Covariance Matrices
##
## data: Y2 by data2_clean$shift
## Chi-Sq (approx.) = 3.3034, df = 6, p-value = 0.7699
# Homogenitas Varians per DV
#--------------------------------------------------------
leveneTest(tensile_strength_mpa ~ shift, data = data2_clean)
## Levene's Test for Homogeneity of Variance (center = median)
## Df F value Pr(>F)
## group 2 0.1823 0.8333
## 9997
leveneTest(surface_roughness_um ~ shift, data = data2_clean)
## Levene's Test for Homogeneity of Variance (center = median)
## Df F value Pr(>F)
## group 2 1.8579 0.156
## 9997
plot(data2_clean$ambient_humidity_pct,
data2_clean$tensile_strength_mpa)

plot(data2_clean$ambient_humidity_pct,
data2_clean$surface_roughness_um)

model_mancova <- manova(
cbind(tensile_strength_mpa,
surface_roughness_um)
~ shift + ambient_humidity_pct,
data = data2_clean
)
summary(model_mancova, test = "Pillai")
## Df Pillai approx F num Df den Df Pr(>F)
## shift 2 0.007486 18.78 4 19992 2.019e-15 ***
## ambient_humidity_pct 1 0.086195 471.39 2 9995 < 2.2e-16 ***
## Residuals 9996
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
# Bisa juga:
summary(model_mancova, test = "Wilks")
## Df Wilks approx F num Df den Df Pr(>F)
## shift 2 0.99251 18.81 4 19990 1.894e-15 ***
## ambient_humidity_pct 1 0.91381 471.39 2 9995 < 2.2e-16 ***
## Residuals 9996
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
summary.aov(model_mancova)
## Response tensile_strength_mpa :
## Df Sum Sq Mean Sq F value Pr(>F)
## shift 2 375 187 0.1744 0.84
## ambient_humidity_pct 1 441559 441559 411.1148 <2e-16 ***
## Residuals 9996 10736240 1074
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Response surface_roughness_um :
## Df Sum Sq Mean Sq F value Pr(>F)
## shift 2 12.49 6.245 37.523 < 2.2e-16 ***
## ambient_humidity_pct 1 89.45 89.455 537.521 < 2.2e-16 ***
## Residuals 9996 1663.54 0.166
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1