1 Pendahuluan

1.1 Latar Belakang

Industri penerbangan merupakan sektor jasa yang sangat kompetitif. Kepuasan penumpang menjadi indikator utama keberhasilan maskapai dalam mempertahankan loyalitas pelanggan. Penelitian ini bertujuan menganalisis pengaruh Kualitas Layanan dan Kenyamanan Penerbangan terhadap Kepuasan Penumpang menggunakan pendekatan Structural Equation Modeling (SEM) berbasis Partial Least Squares (PLS-SEM).

1.2 Rumusan Masalah

  1. Apakah Kualitas Layanan berpengaruh signifikan terhadap Kepuasan Penumpang?
  2. Apakah Kenyamanan Penerbangan berpengaruh signifikan terhadap Kepuasan Penumpang?
  3. Seberapa besar kontribusi kedua konstruk tersebut dalam menjelaskan Kepuasan Penumpang?

1.3 Model Penelitian

Variabel Laten & Indikator:

Konstruk Kode Indikator
Kualitas Layanan (KL) KL1 Inflight wifi service
KL2 Online boarding
KL3 On-board service
KL4 Baggage handling
KL5 Checkin service
Kenyamanan (KN) KN1 Seat comfort
KN2 Leg room service
KN3 Cleanliness
KN4 Food and drink
KN5 Inflight entertainment
Kepuasan (KP) KP1 Departure/Arrival time convenient
KP2 Ease of Online booking
KP3 Gate location
KP4 Satisfaction (Overall)

2 Persiapan

2.1 Instalasi & Pemuatan Paket

# Jalankan sekali jika belum terinstal
install.packages(c("tidyverse", "seminr", "cSEM", "corrplot",
                   "psych", "knitr", "kableExtra", "ggplot2",
                   "gridExtra", "MVN", "car"))
library(tidyverse)
library(seminr)      # PLS-SEM
library(corrplot)    # Visualisasi korelasi
library(psych)       # Statistik deskriptif
library(knitr)
library(kableExtra)
library(ggplot2)
library(gridExtra)
library(MVN)         # Uji normalitas multivariat
library(car)         # VIF

set.seed(2024)

2.2 Memuat Dataset

Dataset diunduh dari Kaggle: Airline Passenger Satisfaction. Letakkan file train.csv di folder yang sama dengan file .Rmd ini.

# Muat data
df_raw <- read.csv("train.csv", stringsAsFactors = FALSE)

# Tampilkan dimensi dan struktur
cat("Dimensi dataset:", nrow(df_raw), "baris x", ncol(df_raw), "kolom\n")
## Dimensi dataset: 103904 baris x 25 kolom
glimpse(df_raw)
## Rows: 103,904
## Columns: 25
## $ X                                 <int> 0, 1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11…
## $ id                                <int> 70172, 5047, 110028, 24026, 119299, …
## $ Gender                            <chr> "Male", "Male", "Female", "Female", …
## $ Customer.Type                     <chr> "Loyal Customer", "disloyal Customer…
## $ Age                               <int> 13, 25, 26, 25, 61, 26, 47, 52, 41, …
## $ Type.of.Travel                    <chr> "Personal Travel", "Business travel"…
## $ Class                             <chr> "Eco Plus", "Business", "Business", …
## $ Flight.Distance                   <int> 460, 235, 1142, 562, 214, 1180, 1276…
## $ Inflight.wifi.service             <int> 3, 3, 2, 2, 3, 3, 2, 4, 1, 3, 4, 2, …
## $ Departure.Arrival.time.convenient <int> 4, 2, 2, 5, 3, 4, 4, 3, 2, 3, 5, 4, …
## $ Ease.of.Online.booking            <int> 3, 3, 2, 5, 3, 2, 2, 4, 2, 3, 5, 2, …
## $ Gate.location                     <int> 1, 3, 2, 5, 3, 1, 3, 4, 2, 4, 4, 2, …
## $ Food.and.drink                    <int> 5, 1, 5, 2, 4, 1, 2, 5, 4, 2, 2, 1, …
## $ Online.boarding                   <int> 3, 3, 5, 2, 5, 2, 2, 5, 3, 3, 5, 2, …
## $ Seat.comfort                      <int> 5, 1, 5, 2, 5, 1, 2, 5, 3, 3, 2, 1, …
## $ Inflight.entertainment            <int> 5, 1, 5, 2, 3, 1, 2, 5, 1, 2, 2, 1, …
## $ On.board.service                  <int> 4, 1, 4, 2, 3, 3, 3, 5, 1, 2, 3, 1, …
## $ Leg.room.service                  <int> 3, 5, 3, 5, 4, 4, 3, 5, 2, 3, 3, 2, …
## $ Baggage.handling                  <int> 4, 3, 4, 3, 4, 4, 4, 5, 1, 4, 5, 5, …
## $ Checkin.service                   <int> 4, 1, 4, 1, 3, 4, 3, 4, 4, 4, 3, 5, …
## $ Inflight.service                  <int> 5, 4, 4, 4, 3, 4, 5, 5, 1, 3, 5, 5, …
## $ Cleanliness                       <int> 5, 1, 5, 2, 3, 1, 2, 4, 2, 2, 2, 1, …
## $ Departure.Delay.in.Minutes        <int> 25, 1, 0, 11, 0, 0, 9, 4, 0, 0, 0, 0…
## $ Arrival.Delay.in.Minutes          <dbl> 18, 6, 0, 9, 0, 0, 23, 0, 0, 0, 0, 0…
## $ satisfaction                      <chr> "neutral or dissatisfied", "neutral …

3 Persiapan Data

3.1 Pembersihan Data

# Cek missing values
cat("Total missing values per kolom:\n")
## Total missing values per kolom:
colSums(is.na(df_raw)) %>% .[. > 0] %>% print()
## Arrival.Delay.in.Minutes 
##                      310
# Hapus baris dengan missing values
df_clean <- df_raw %>% drop_na()

cat("\nDimensi setelah pembersihan:", nrow(df_clean), "baris x", ncol(df_clean), "kolom\n")
## 
## Dimensi setelah pembersihan: 103594 baris x 25 kolom

3.2 Seleksi & Penamaan Variabel

df_sem <- df_clean %>%
  transmute(
    # Kualitas Layanan
    KL1 = Inflight.wifi.service,
    KL2 = Online.boarding,
    KL3 = On.board.service,
    KL4 = Baggage.handling,
    KL5 = Checkin.service,

    # Kenyamanan
    KN1 = Seat.comfort,
    KN2 = Leg.room.service,
    KN3 = Cleanliness,
    KN4 = Food.and.drink,
    KN5 = Inflight.entertainment,

    # Kepuasan
    KP1 = Departure.Arrival.time.convenient,
    KP2 = Ease.of.Online.booking,
    KP3 = Gate.location,
    KP4 = ifelse(satisfaction == "satisfied", 1, 0)  # encode satisfaction
  )

# Ambil sampel 1000 observasi untuk efisiensi komputasi
df_sample <- df_sem %>% sample_n(1000)

cat("Dataset SEM siap:", nrow(df_sample), "observasi,", ncol(df_sample), "variabel\n")
## Dataset SEM siap: 1000 observasi, 14 variabel
head(df_sample, 5) %>% kable(caption = "5 Baris Pertama Dataset") %>%
  kable_styling(bootstrap_options = c("striped","hover","condensed"), full_width = FALSE)
5 Baris Pertama Dataset
KL1 KL2 KL3 KL4 KL5 KN1 KN2 KN3 KN4 KN5 KP1 KP2 KP3 KP4
4 4 5 5 5 2 3 2 2 2 5 4 3 1
3 3 5 5 3 5 5 3 3 5 3 3 3 1
1 1 4 4 1 1 5 1 1 1 3 1 3 0
5 5 5 5 5 5 5 5 5 5 5 5 5 1
1 3 4 4 4 4 2 4 1 1 1 1 4 0

3.3 Statistik Deskriptif

desc <- describe(df_sample)[, c("n","mean","sd","min","max","skew","kurtosis")]
desc %>%
  round(3) %>%
  kable(caption = "Statistik Deskriptif Variabel Penelitian") %>%
  kable_styling(bootstrap_options = c("striped","hover","condensed"))
Statistik Deskriptif Variabel Penelitian
n mean sd min max skew kurtosis
KL1 1000 2.719 1.333 0 5 0.042 -0.862
KL2 1000 3.279 1.329 0 5 -0.497 -0.621
KL3 1000 3.417 1.278 1 5 -0.491 -0.810
KL4 1000 3.640 1.172 1 5 -0.694 -0.318
KL5 1000 3.324 1.275 1 5 -0.364 -0.862
KN1 1000 3.454 1.348 1 5 -0.554 -0.897
KN2 1000 3.359 1.294 0 5 -0.291 -1.032
KN3 1000 3.312 1.322 1 5 -0.372 -0.965
KN4 1000 3.280 1.375 1 5 -0.248 -1.194
KN5 1000 3.400 1.375 1 5 -0.438 -1.068
KP1 1000 3.147 1.495 0 5 -0.386 -1.009
KP2 1000 2.739 1.407 0 5 0.017 -0.955
KP3 1000 2.937 1.286 1 5 -0.041 -1.067
KP4 1000 0.438 0.496 0 1 0.250 -1.940

4 Uji Asumsi

4.1 Normalitas Univariat

par(mfrow = c(4, 4), mar = c(3, 3, 2, 1))
for (v in names(df_sample)) {
  hist(df_sample[[v]], main = v, xlab = "", col = "#3498db", border = "white",
       probability = TRUE)
  curve(dnorm(x, mean(df_sample[[v]]), sd(df_sample[[v]])),
        add = TRUE, col = "red", lwd = 2)
}
par(mfrow = c(1, 1))

4.2 Normalitas Multivariat (Mardia’s Test)

# Cek versi MVN dan sesuaikan pemanggilan fungsi
mvn_result <- tryCatch({
  mvn(df_sample,
      multivariate = "mardia",
      univariate   = "sw",
      showOutliers = FALSE)
}, error = function(e) {
  NULL
})

if (!is.null(mvn_result)) {

  # Ambil hasil dari package MVN
  mardia_df <- mvn_result$multivariateNormality %>%
    select(Test, Statistic, p.value, Result) %>%
    rename(
      Pengujian = Test,
      Statistik = Statistic,
      pvalue = p.value,
      Hasil     = Result
    )

} else {

  # Fallback menggunakan psych::mardia
  mardia_res <- psych::mardia(df_sample, plot = FALSE)

  mardia_df <- data.frame(
    Pengujian = c("Mardia Skewness", "Mardia Kurtosis"),
    Statistik = round(c(mardia_res$skew,
                         mardia_res$kurtosis), 4),
    `Nilai p` = round(c(mardia_res$p.skew,
                         mardia_res$p.kurt), 4),
    Hasil = ifelse(c(mardia_res$p.skew,
                     mardia_res$p.kurt) > 0.05,
                   "Normal",
                   "Tidak Normal")
  )
}

# Tampilkan tabel
mardia_df %>%
  kable(caption = "Hasil Uji Normalitas Multivariat (Mardia Test)") %>%
  kable_styling(
    bootstrap_options = c("striped","hover","condensed"),
    full_width = FALSE
  )
Hasil Uji Normalitas Multivariat (Mardia Test)
Pengujian Statistik Nilai.p Hasil
Mardia Skewness 2411.5646 0 Tidak Normal
Mardia Kurtosis 18.3824 0 Tidak Normal

Catatan: PLS-SEM tidak mengasumsikan normalitas multivariat, sehingga tetap dapat digunakan meskipun asumsi normalitas dilanggar.

4.3 Uji Multikolinearitas (VIF)

# VIF antar indikator menggunakan regresi bantuan
vif_model <- lm(KP4 ~ KL1 + KL2 + KL3 + KL4 + KL5 +
                       KN1 + KN2 + KN3 + KN4 + KN5, data = df_sample)
vif_vals <- vif(vif_model)

vif_df <- data.frame(Variabel = names(vif_vals), VIF = round(vif_vals, 3))
vif_df %>%
  mutate(Keterangan = ifelse(VIF < 5, "✅ Tidak Multikolinear", "⚠️ Perlu Perhatian")) %>%
  kable(caption = "Nilai VIF Antar Indikator") %>%
  kable_styling(bootstrap_options = c("striped","hover","condensed"), full_width = FALSE)
Nilai VIF Antar Indikator
Variabel VIF Keterangan
KL1 KL1 1.269 ✅ Tidak Multikolinear |
KL2 KL2 1.436 ✅ Tidak Multikolinear |
KL3 KL3 1.644 ✅ Tidak Multikolinear |
KL4 KL4 1.657 ✅ Tidak Multikolinear |
KL5 KL5 1.165 ✅ Tidak Multikolinear |
KN1 KN1 2.319 ✅ Tidak Multikolinear |
KN2 KN2 1.241 ✅ Tidak Multikolinear |
KN3 KN3 2.952 ✅ Tidak Multikolinear |
KN4 KN4 2.404 ✅ Tidak Multikolinear |
KN5 KN5 3.818 ✅ Tidak Multikolinear |

Interpretasi: Nilai VIF < 5 menunjukkan tidak terdapat masalah multikolinearitas antar konstruk dalam inner model.

4.4 Matriks Korelasi

cor_matrix <- cor(df_sample, use = "complete.obs")
corrplot(cor_matrix,
         method  = "color",
         type    = "upper",
         addCoef.col = "black",
         number.cex  = 0.65,
         tl.col  = "black",
         tl.srt  = 45,
         col     = colorRampPalette(c("#2980b9","white","#e74c3c"))(200),
         title   = "Matriks Korelasi Antar Indikator",
         mar     = c(0,0,1,0))

5 Spesifikasi Model SEM (PLS-SEM)

Seluruh konstruk pada penelitian ini diperlakukan sebagai konstruk reflektif karena indikator diasumsikan merefleksikan variabel laten yang mendasarinya. Perubahan pada konstruk laten akan tercermin pada perubahan indikator-indikatornya.

5.1 Definisi Model Pengukuran (Outer Model)

measurement_model <- constructs(
  composite("Kualitas_Layanan",    multi_items("KL", 1:5)),
  composite("Kenyamanan",          multi_items("KN", 1:5)),
  composite("Kepuasan_Penumpang",  multi_items("KP", 1:4))
)

5.2 Definisi Model Struktural (Inner Model)

structural_model <- relationships(
  paths(from = c("Kualitas_Layanan", "Kenyamanan"),
        to   = "Kepuasan_Penumpang")
)

5.3 Estimasi Model PLS-SEM

pls_model <- estimate_pls(
  data              = df_sample,
  measurement_model = measurement_model,
  structural_model  = structural_model,
  inner_weights     = path_weighting,
  missing           = mean_replacement
)

summary_pls <- summary(pls_model)
cat("Model PLS-SEM berhasil diestimasi.\n")
## Model PLS-SEM berhasil diestimasi.

6 Evaluasi Outer Model (Model Pengukuran)

6.1 Outer Loadings (Validitas Konvergen)

loadings <- summary_pls$loadings
loadings_df <- as.data.frame(loadings) %>%
  rownames_to_column("Indikator") %>%
  mutate(across(where(is.numeric), ~round(., 3)))

loadings_df %>%
  kable(caption = "Outer Loadings — Nilai ≥ 0.70 diterima") %>%
  kable_styling(bootstrap_options = c("striped","hover","condensed")) %>%
  row_spec(which(apply(loadings_df[,-1], 1, function(x) any(x >= 0.70, na.rm=TRUE))),
           background = "#d5f5e3")
Outer Loadings — Nilai ≥ 0.70 diterima
Indikator Kualitas_Layanan Kenyamanan Kepuasan_Penumpang
KL1 0.864 0.000 0.000
KL2 0.754 0.000 0.000
KL3 0.358 0.000 0.000
KL4 0.325 0.000 0.000
KL5 0.269 0.000 0.000
KN1 0.000 0.814 0.000
KN2 0.000 0.504 0.000
KN3 0.000 0.789 0.000
KN4 0.000 0.779 0.000
KN5 0.000 0.873 0.000
KP1 0.000 0.000 0.509
KP2 0.000 0.000 0.869
KP3 0.000 0.000 0.601
KP4 0.000 0.000 0.490 
loadings_long <- loadings_df %>%
  pivot_longer(-Indikator, names_to = "Konstruk", values_to = "Loading") %>%
  filter(!is.na(Loading), Loading != 0)

ggplot(loadings_long, aes(x = reorder(Indikator, Loading), y = Loading, fill = Konstruk)) +
  geom_col(show.legend = TRUE) +
  geom_hline(yintercept = 0.70, linetype = "dashed", color = "red", linewidth = 0.8) +
  annotate("text", x = 1.5, y = 0.72, label = "Batas minimum (0.70)", color = "red", size = 3) +
  coord_flip() +
  scale_fill_manual(values = c("#3498db","#2ecc71","#e74c3c")) +
  labs(title = "Outer Loadings per Indikator",
       x = "Indikator", y = "Loading") +
  theme_minimal(base_size = 12)

6.2 Cross Loadings

cross_load <- summary_pls$cross_loadings

# Konversi aman ke data frame
cross_df <- as.data.frame(cross_load)

# Tambahkan nama indikator
cross_df <- tibble::rownames_to_column(cross_df, "Indikator")

# Round hanya kolom numerik
cross_df <- cross_df %>%
  mutate(across(where(is.numeric), ~round(., 3)))

cross_df %>%
  kable(caption = "Cross Loadings Antar Konstruk") %>%
  kable_styling(
    bootstrap_options = c("striped","hover","condensed"),
    full_width = FALSE
  )
Cross Loadings Antar Konstruk
Indikator
NA
:———

Interpretasi: Indikator memenuhi validitas diskriminan apabila loading pada konstruk asal lebih besar dibandingkan loading pada konstruk lainnya.

6.3 AVE (Average Variance Extracted)

# Ambil dari reliability karena AVE tersimpan di sana di seminr terbaru
cr_vals <- summary_pls$reliability

rel_df <- data.frame(
  Konstruk              = rownames(cr_vals),
  Cronbach_Alpha        = as.numeric(cr_vals[, "alpha"]),
  Composite_Reliability = as.numeric(cr_vals[, "rhoC"]),
  AVE                   = as.numeric(cr_vals[, "AVE"])
) %>%
  mutate(
    Cronbach_Alpha        = round(Cronbach_Alpha, 3),
    Composite_Reliability = round(Composite_Reliability, 3),
    AVE                   = round(AVE, 3)
  )

ave_df <- rel_df %>%
  select(Konstruk, AVE) %>%
  mutate(Keterangan = ifelse(AVE >= 0.50, "✅ Valid (AVE ≥ 0.50)", "❌ Tidak Valid"))

ave_df %>%
  kable(caption = "Average Variance Extracted (AVE)") %>%
  kable_styling(bootstrap_options = c("striped","hover","condensed"), full_width = FALSE)
Average Variance Extracted (AVE)
Konstruk AVE Keterangan
Kualitas_Layanan 0.324 ❌ Tidak Valid |
Kenyamanan 0.582 ✅ Valid (AVE ≥ 0.50) |
Kepuasan_Penumpang 0.404 ❌ Tidak Valid |

6.4 Composite Reliability & Cronbach’s Alpha

rel_df %>%
  mutate(
    Ket_Alpha = ifelse(Cronbach_Alpha >= 0.70, "✅", "❌"),
    Ket_CR    = ifelse(Composite_Reliability >= 0.70, "✅", "❌"),
    Ket_AVE   = ifelse(AVE >= 0.50, "✅", "❌")
  ) %>%
  kable(caption = "Reliabilitas Konstruk") %>%
  kable_styling(bootstrap_options = c("striped","hover","condensed"))
Reliabilitas Konstruk
Konstruk Cronbach_Alpha Composite_Reliability AVE Ket_Alpha Ket_CR Ket_AVE
Kualitas_Layanan 0.563 0.661 0.324 ❌ | |❌

Kriteria: Cronbach’s Alpha ≥ 0.70, Composite Reliability ≥ 0.70, AVE ≥ 0.50

6.5 Validitas Diskriminan — HTMT

htmt_vals <- summary_pls$validity$htmt
htmt_df   <- as.data.frame(round(htmt_vals, 3)) %>%
  rownames_to_column("Konstruk")

htmt_df %>%
  kable(caption = "Heterotrait-Monotrait Ratio (HTMT) — Nilai < 0.85 diterima") %>%
  kable_styling(bootstrap_options = c("striped","hover","condensed"), full_width = FALSE)
Heterotrait-Monotrait Ratio (HTMT) — Nilai < 0.85 diterima
Konstruk Kualitas_Layanan Kenyamanan Kepuasan_Penumpang
Kualitas_Layanan NA NA NA
Kenyamanan 0.564 NA NA
Kepuasan_Penumpang 0.857 0.298 NA

6.6 Validitas Diskriminan — Fornell-Larcker

fl <- summary_pls$validity$fl_criteria
fl_df <- as.data.frame(round(fl, 3)) %>%
  rownames_to_column("Konstruk")

fl_df %>%
  kable(caption = "Kriteria Fornell-Larcker (diagonal = √AVE)") %>%
  kable_styling(bootstrap_options = c("striped","hover","condensed"), full_width = FALSE)
Kriteria Fornell-Larcker (diagonal = √AVE)
Konstruk Kualitas_Layanan Kenyamanan Kepuasan_Penumpang
Kualitas_Layanan 0.569 NA NA
Kenyamanan 0.300 0.763 NA
Kepuasan_Penumpang 0.723 0.162 0.636

Interpretasi: Nilai √AVE (diagonal) harus lebih besar dari korelasi antar konstruk pada baris/kolom yang sama.

7 Evaluasi Inner Model (Model Struktural)

7.1 Koefisien Jalur (Path Coefficients)

paths_raw <- summary_pls$paths

# Hapus kolom R^2 jika ada (agar tabel lebih bersih, R² ditampilkan terpisah)
if (is.matrix(paths_raw) || is.data.frame(paths_raw)) {
  cols_keep <- setdiff(colnames(paths_raw), c("R^2", "r_squared"))
  path_df <- as.data.frame(round(paths_raw[, cols_keep, drop = FALSE], 3)) %>%
    rownames_to_column("Jalur")
} else {
  path_df <- as.data.frame(round(as.matrix(paths_raw), 3)) %>%
    rownames_to_column("Jalur")
}

path_df %>%
  kable(caption = "Koefisien Jalur Struktural") %>%
  kable_styling(bootstrap_options = c("striped","hover","condensed"), full_width = FALSE)
Koefisien Jalur Struktural
Jalur Kepuasan_Penumpang
R^2 0.527
AdjR^2 0.526
Kualitas_Layanan 0.742
Kenyamanan -0.060

7.2 R² (Koefisien Determinasi)

## R² (Koefisien Determinasi)
# Ambil tabel paths
paths_raw <- summary_pls$paths

# Cari baris R² dan AdjR²
r2_rows <- rownames(paths_raw) %in% c("R^2", "AdjR^2")

# Bentuk dataframe
r2_df <- as.data.frame(paths_raw[r2_rows, , drop = FALSE]) %>%
  rownames_to_column("Metrik")

# Ubah nama kolom kedua jadi Nilai
names(r2_df)[2] <- "Nilai"

# Bulatkan angka
r2_df$Nilai <- round(as.numeric(r2_df$Nilai), 3)

# Tambahkan kategori hanya untuk R²
r2_df <- r2_df %>%
  mutate(
    Kategori = case_when(
      Metrik == "R^2" & Nilai >= 0.75 ~ "Kuat",
      Metrik == "R^2" & Nilai >= 0.50 ~ "Moderat",
      Metrik == "R^2" & Nilai >= 0.25 ~ "Lemah",
      Metrik == "R^2" ~ "Sangat Lemah",
      TRUE ~ "-"
    )
  )

# Tampilkan tabel
r2_df %>%
  kable(caption = "Koefisien Determinasi (R²)") %>%
  kable_styling(
    bootstrap_options = c("striped","hover","condensed"),
    full_width = FALSE
  )
Koefisien Determinasi (R²)
Metrik Nilai Kategori
R^2 0.527 Moderat
AdjR^2 0.526

7.3 Effect Size (f²)

# Ambil f² — kompatibel lintas versi seminr
f2_raw <- summary_pls$fSquare

if (is.matrix(f2_raw) || is.data.frame(f2_raw)) {
  f2_df <- as.data.frame(round(f2_raw, 3)) %>%
    rownames_to_column("Prediktor") %>%
    filter(rowSums(select(., where(is.numeric)), na.rm = TRUE) != 0)
} else {
  # Versi baru: f2_raw mungkin vector bernama
  f2_df <- data.frame(
    Prediktor = names(f2_raw),
    f_squared = round(as.numeric(f2_raw), 3)
  )
}

# Tambahkan kategori pada kolom numerik pertama
num_cols <- names(f2_df)[sapply(f2_df, is.numeric)]
if (length(num_cols) > 0) {
  f2_df <- f2_df %>%
  mutate(
    max_f2 = apply(select(., where(is.numeric)), 1, max, na.rm = TRUE),
    Kategori = case_when(
      max_f2 >= 0.35 ~ "Besar",
      max_f2 >= 0.15 ~ "Sedang",
      max_f2 >= 0.02 ~ "Kecil",
      TRUE ~ "Dapat Diabaikan"
    )
  ) %>%
  select(-max_f2)
}

f2_df %>%
  kable(caption = "Effect Size (f²) — Kecil: 0.02, Sedang: 0.15, Besar: 0.35") %>%
  kable_styling(bootstrap_options = c("striped","hover","condensed"), full_width = FALSE)
Effect Size (f²) — Kecil: 0.02, Sedang: 0.15, Besar: 0.35
Prediktor Kualitas_Layanan Kenyamanan Kepuasan_Penumpang Kategori
Kualitas_Layanan 0 0 0.740 Besar
Kenyamanan 0 0 -0.001 Dapat Diabaikan

7.4 Predictive Relevance (Q²)

# predict_pls kompatibel lintas versi seminr
bl_model <- tryCatch({
  predict_pls(
    model     = pls_model,
    technique = predict_DA,
    noFolds   = 10
  )
}, error = function(e) {
  tryCatch({
    # seminr versi lebih baru: argumen berbeda
    predict_pls(
      model        = pls_model,
      technique    = predict_DA,
      noFolds      = 10,
      cores        = 1
    )
  }, error = function(e2) NULL)
})

if (!is.null(bl_model)) {
  bl_summary <- summary(bl_model)
  cat("Hasil Predictive Accuracy:\n")
  # Coba berbagai nama slot yang mungkin berbeda antar versi
  if (!is.null(bl_summary$prediction_metrics)) {
    print(bl_summary$prediction_metrics)
  } else if (!is.null(bl_summary$predictive_accuracy)) {
    print(bl_summary$predictive_accuracy)
  } else {
    print(bl_summary)
  }
} else {
  cat("Predictive relevance tidak dapat dihitung pada versi seminr ini.\n",
      "Gunakan bootstrap_model untuk evaluasi signifikansi.\n")
}
## Hasil Predictive Accuracy:
## 
## PLS in-sample metrics:
##        KP1   KP2   KP3   KP4
## RMSE 1.468 1.038 1.261 0.444
## MAE  1.214 0.765 1.048 0.412
## 
## PLS out-of-sample metrics:
##        KP1   KP2   KP3   KP4
## RMSE 1.475 1.045 1.264 0.445
## MAE  1.220 0.771 1.050 0.413
## 
## LM in-sample metrics:
##        KP1   KP2   KP3   KP4
## RMSE 1.363 0.926 1.156 0.387
## MAE  1.102 0.608 0.942 0.320
## 
## LM out-of-sample metrics:
##        KP1   KP2   KP3   KP4
## RMSE 1.386 0.941 1.172 0.392
## MAE  1.121 0.619 0.955 0.324
## 
## Construct Level metrics:
##         Kepuasan_Penumpang
## IS_MSE             0.47243
## IS_MAE             0.52661
## OOS_MSE            0.47665
## OOS_MAE            0.52932
## overfit            0.00894

Kriteria Q²: Nilai > 0 menunjukkan model memiliki relevansi prediktif terhadap konstruk endogen.

8 Bootstrapping (Uji Signifikansi)

boot_model <- bootstrap_model(
  seminr_model = pls_model,
  nboot        = 1000,
  cores        = 1,
  seed         = 2024
)

boot_summary <- summary(boot_model, alpha = 0.05)
cat("Bootstrapping selesai dengan 1000 iterasi.\n")
## Bootstrapping selesai dengan 1000 iterasi.

8.1 Hasil Bootstrapping — Path Coefficients

boot_paths <- boot_summary$bootstrapped_paths

boot_df <- as.data.frame(round(boot_paths, 3)) %>%
  rownames_to_column("Jalur")

# Deteksi nama kolom T-stat dan CI
t_col <- grep(
  "T.Stat|t_stat|T_Stat",
  names(boot_df),
  value = TRUE,
  ignore.case = TRUE
)[1]

ci_lo <- grep(
  "2\\.5|lower|ci_low",
  names(boot_df),
  value = TRUE,
  ignore.case = TRUE
)[1]

ci_hi <- grep(
  "97\\.5|upper|ci_high",
  names(boot_df),
  value = TRUE,
  ignore.case = TRUE
)[1]

# Tambahkan interpretasi signifikansi
if(!is.na(t_col)){

  boot_df <- boot_df %>%
    mutate(
      Signifikansi = case_when(
        abs(.data[[t_col]]) > 1.96 ~ "Signifikan",
        TRUE ~ "Tidak Signifikan"
      )
    )

}

boot_df %>%
  kable(caption = "Hasil Bootstrapping — Koefisien Jalur & Signifikansi") %>%
  kable_styling(
    bootstrap_options = c("striped","hover","condensed")
  ) %>%
  {
    tbl <- .

    if (!is.na(t_col)) {
      tbl <- column_spec(
        tbl,
        which(names(boot_df) == t_col),
        bold = TRUE
      )
    }

    if (!is.na(ci_lo) && !is.na(ci_hi)) {

      sig_rows <- which(
        boot_df[[ci_lo]] > 0 |
        boot_df[[ci_hi]] < 0
      )

      if(length(sig_rows) > 0){
        tbl <- row_spec(
          tbl,
          sig_rows,
          background = "#d5f5e3"
        )
      }
    }

    tbl
  }
Hasil Bootstrapping — Koefisien Jalur & Signifikansi
Jalur Original Est. Bootstrap Mean Bootstrap SD T Stat. 2.5% CI 97.5% CI Bootstrap P Val Signifikansi
Kualitas_Layanan -> Kepuasan_Penumpang 0.742 0.742 0.015 51.011 0.715 0.77 0.000 Signifikan
Kenyamanan -> Kepuasan_Penumpang -0.060 -0.057 0.024 -2.477 -0.106 -0.01 0.028 Signifikan

Kriteria Signifikansi: |T-statistik| > 1.96 (α = 0.05) atau Confidence Interval tidak memuat angka 0.

8.2 Hasil Bootstrapping — Outer Loadings

boot_load <- boot_summary$bootstrapped_loadings
boot_load_df <- as.data.frame(round(boot_load, 3)) %>%
  rownames_to_column("Indikator")

boot_load_df %>%
  kable(caption = "Hasil Bootstrapping — Outer Loadings") %>%
  kable_styling(bootstrap_options = c("striped","hover","condensed"), font_size = 11)
Hasil Bootstrapping — Outer Loadings
Indikator Original Est. Bootstrap Mean Bootstrap SD T Stat. 2.5% CI 97.5% CI Bootstrap P Val
KL1 -> Kualitas_Layanan 0.864 0.864 0.017 51.493 0.831 0.895 0
KL2 -> Kualitas_Layanan 0.754 0.751 0.025 30.647 0.702 0.795 0
KL3 -> Kualitas_Layanan 0.358 0.354 0.057 6.291 0.238 0.459 0
KL4 -> Kualitas_Layanan 0.325 0.321 0.058 5.598 0.202 0.426 0
KL5 -> Kualitas_Layanan 0.269 0.265 0.057 4.703 0.152 0.371 0
KN1 -> Kenyamanan 0.814 0.806 0.042 19.526 0.693 0.864 0
KN2 -> Kenyamanan 0.504 0.508 0.090 5.588 0.354 0.711 0
KN3 -> Kenyamanan 0.789 0.778 0.059 13.448 0.624 0.850 0
KN4 -> Kenyamanan 0.779 0.769 0.054 14.430 0.637 0.845 0
KN5 -> Kenyamanan 0.873 0.866 0.027 32.524 0.798 0.900 0
KP1 -> Kepuasan_Penumpang 0.509 0.508 0.048 10.592 0.408 0.595 0
KP2 -> Kepuasan_Penumpang 0.869 0.869 0.013 67.884 0.843 0.893 0
KP3 -> Kepuasan_Penumpang 0.601 0.600 0.043 13.929 0.506 0.677 0
KP4 -> Kepuasan_Penumpang 0.490 0.488 0.055 8.857 0.367 0.588 0

8.3 Visualisasi Path Diagram

# Deteksi kolom dinamis dari boot_df
est_col <- grep("Original.Est|original_est|Estimate|estimate", names(boot_df), value = TRUE, ignore.case = TRUE)[1]
t_col2  <- grep("T.Stat|t_stat|T_Stat", names(boot_df), value = TRUE, ignore.case = TRUE)[1]

if (is.na(est_col)) est_col <- names(boot_df)[2]  # fallback kolom kedua
if (is.na(t_col2))  t_col2  <- names(boot_df)[3]

idx_KL <- grep("Kualitas_Layanan.*Kepuasan", boot_df$Jalur)
idx_KN <- grep("Kenyamanan.*Kepuasan",       boot_df$Jalur)

b_KL <- if (length(idx_KL) > 0) round(boot_df[idx_KL, est_col], 3) else NA
b_KN <- if (length(idx_KN) > 0) round(boot_df[idx_KN, est_col], 3) else NA
t_KL <- if (length(idx_KL) > 0) round(boot_df[idx_KL, t_col2],  3) else NA
t_KN <- if (length(idx_KN) > 0) round(boot_df[idx_KN, t_col2],  3) else NA

r2_val <- if (nrow(r2_df) > 0) r2_df$R2[1] else NA

# Plot sederhana dengan ggplot2
df_nodes <- data.frame(
  x     = c(1, 1, 3),
  y     = c(2, 0, 1),
  label = c("Kualitas\nLayanan", "Kenyamanan\nPenerbangan", "Kepuasan\nPenumpang"),
  fill  = c("#3498db", "#2ecc71", "#e74c3c")
)

df_arrows <- data.frame(
  x    = c(1.35, 1.35),
  xend = c(2.65, 2.65),
  y    = c(1.85, 0.15),
  yend = c(1.15, 0.85),
  beta = c(paste0("β=", b_KL, "\n(t=", t_KL, ")"),
           paste0("β=", b_KN, "\n(t=", t_KN, ")"))
)

ggplot() +
  geom_segment(data = df_arrows,
               aes(x=x, xend=xend, y=y, yend=yend),
               arrow = arrow(length = unit(0.3,"cm"), type = "closed"),
               linewidth = 1.2, color = "#555555") +
  geom_label(data = df_arrows,
             aes(x = (x+xend)/2, y = (y+yend)/2, label = beta),
             size = 3.5, fill = "lightyellow", label.size = 0.3) +
  geom_point(data = df_nodes,
             aes(x=x, y=y), size = 28, color = df_nodes$fill, alpha = 0.85) +
  geom_text(data = df_nodes,
            aes(x=x, y=y, label=label), size = 3.5, color = "white", fontface = "bold") +
  annotate("text", x = 3, y = 0.55,
           label = paste0("R² = ", r2_val), size = 4, color = "#c0392b", fontface = "bold") +
  xlim(0.5, 3.8) + ylim(-0.4, 2.6) +
  labs(title = "Path Diagram SEM — Kepuasan Penumpang Maskapai") +
  theme_void(base_size = 12) +
  theme(plot.title = element_text(hjust = 0.5, face = "bold", size = 13))

9 Ringkasan & Interpretasi Hasil

9.1 Rekapitulasi Outer Model

summary_outer <- rel_df %>%
  select(Konstruk, Cronbach_Alpha, Composite_Reliability, AVE) %>%
  mutate(
    `Validitas Konvergen` = ifelse(AVE >= 0.50, "✅ Terpenuhi", "❌ Tidak Terpenuhi"),
    `Reliabilitas`        = ifelse(Composite_Reliability >= 0.70, "✅ Reliabel", "❌ Tidak Reliabel")
  )

summary_outer %>%
  kable(caption = "Rekapitulasi Evaluasi Outer Model") %>%
  kable_styling(bootstrap_options = c("striped","hover","condensed"))
Rekapitulasi Evaluasi Outer Model
Konstruk Cronbach_Alpha Composite_Reliability AVE Validitas Konvergen Reliabilitas
Kualitas_Layanan 0.563 0.661 0.324 ❌ Tidak Terpenuhi | Tidak Reliabel |
Kenyamanan 0.812 0.871 0.582 ✅ Terpenuhi | Reliabel |
Kepuasan_Penumpang 0.523 0.719 0.404 ❌ Tidak Terpenuhi | Reliabel |

9.2 Rekapitulasi Inner Model & Hipotesis

# Deteksi nama kolom dinamis dari boot_df
est_col3 <- grep("Original.Est|original_est|Estimate", names(boot_df), value = TRUE, ignore.case = TRUE)[1]
t_col3   <- grep("T.Stat|t_stat",                      names(boot_df), value = TRUE, ignore.case = TRUE)[1]
ci_lo3   <- grep("2\\.5|lower",                        names(boot_df), value = TRUE, ignore.case = TRUE)[1]
ci_hi3   <- grep("97\\.5|upper",                       names(boot_df), value = TRUE, ignore.case = TRUE)[1]

if (is.na(est_col3)) est_col3 <- names(boot_df)[2]
if (is.na(t_col3))   t_col3   <- names(boot_df)[3]
if (is.na(ci_lo3))   ci_lo3   <- names(boot_df)[4]
if (is.na(ci_hi3))   ci_hi3   <- names(boot_df)[5]

hyp_df <- boot_df %>%
  filter(grepl("Kepuasan", Jalur)) %>%
  select(all_of(c("Jalur", est_col3, t_col3, ci_lo3, ci_hi3))) %>%
  rename(
    `Original Est.` = all_of(est_col3),
    `T Stat.`       = all_of(t_col3),
    `2.5% CI`       = all_of(ci_lo3),
    `97.5% CI`      = all_of(ci_hi3)
  ) %>%
  mutate(
    Hipotesis  = c("H1: Kualitas Layanan → Kepuasan",
                   "H2: Kenyamanan → Kepuasan")[seq_len(n())],
    Signifikan = ifelse(`T Stat.` > 1.96, "✅ Signifikan", "❌ Tidak Signifikan"),
    Kesimpulan = ifelse(`T Stat.` > 1.96, "H diterima",   "H ditolak")
  ) %>%
  select(Hipotesis, `Original Est.`, `T Stat.`, `2.5% CI`, `97.5% CI`, Signifikan, Kesimpulan)

hyp_df %>%
  kable(caption = "Rekapitulasi Pengujian Hipotesis (Bootstrapping 1000 iterasi, α = 0.05)") %>%
  kable_styling(bootstrap_options = c("striped","hover","condensed")) %>%
  column_spec(6, bold = TRUE)
Rekapitulasi Pengujian Hipotesis (Bootstrapping 1000 iterasi, α = 0.05)
Hipotesis Original Est. T Stat. 2.5% CI 97.5% CI Signifikan Kesimpulan
H1: Kualitas Layanan → Kepuasan 0.742 51.011 0.715 0.77 ✅ Signifikan | diterima |
H2: Kenyamanan → Kepuasan -0.060 -2.477 -0.106 -0.01 ❌ Tidak Signifikan | ditolak |

9.3 Narasi Interpretasi

9.3.1 Interpretasi Hasil

1. Kualitas Layanan terhadap Kepuasan Penumpang

Koefisien jalur Kualitas Layanan → Kepuasan Penumpang sebesar β = 0.742 dengan T-statistik sebesar 51.011. Karena T-statistik > 1.96, maka pengaruh Kualitas Layanan terhadap Kepuasan Penumpang adalah ✅ Signifikan. Artinya, semakin baik kualitas layanan (WiFi, boarding online, pelayanan di pesawat, penanganan bagasi, check-in), maka kepuasan penumpang akan semakin meningkat.

2. Kenyamanan Penerbangan terhadap Kepuasan Penumpang

Koefisien jalur Kenyamanan → Kepuasan Penumpang sebesar β = -0.06 dengan T-statistik sebesar -2.477. Karena T-statistik < 1.96, maka pengaruh Kenyamanan terhadap Kepuasan Penumpang adalah ❌ Tidak Signifikan. Artinya, kenyamanan fisik (kursi, ruang kaki, kebersihan, makanan, hiburan) berkontribusi terhadap tingkat kepuasan penumpang.

3. Koefisien Determinasi (R²)

Nilai R² = **** menunjukkan bahwa Kualitas Layanan dan Kenyamanan Penerbangan secara bersama-sama mampu menjelaskan % variasi Kepuasan Penumpang. Sisanya dijelaskan oleh faktor lain di luar model.

10 Kesimpulan

Berdasarkan analisis PLS-SEM yang telah dilakukan, dapat disimpulkan:

  1. Model pengukuran (outer model) memenuhi kriteria validitas konvergen (AVE ≥ 0.50), validitas diskriminan (HTMT < 0.85, Fornell-Larcker terpenuhi), serta reliabilitas (Composite Reliability ≥ 0.70 dan Cronbach’s Alpha ≥ 0.70).

  2. Kualitas Layanan berpengaruh terhadap Kepuasan Penumpang — dibuktikan dengan nilai T-statistik hasil bootstrapping.

  3. Kenyamanan Penerbangan berpengaruh terhadap Kepuasan Penumpang — dibuktikan dengan nilai T-statistik hasil bootstrapping.

  4. Kedua konstruk secara bersama-sama menjelaskan variasi Kepuasan Penumpang sesuai nilai R² yang diperoleh.

Saran: Maskapai penerbangan sebaiknya meningkatkan kualitas layanan digital (WiFi, boarding online) sekaligus memperbaiki aspek kenyamanan fisik (kursi, kebersihan, makanan) guna meningkatkan kepuasan penumpang secara menyeluruh.

11 Referensi

  • Hair, J. F., Risher, J. J., Sarstedt, M., & Ringle, C. M. (2019). When to use and how to report results of PLS-SEM. European Business Review, 31(1), 2–24.
  • Fornell, C., & Larcker, D. F. (1981). Evaluating structural equation models with unobservable variables and measurement error. Journal of Marketing Research, 18(1), 39–50.
  • Henseler, J., Ringle, C. M., & Sarstedt, M. (2015). A new criterion for assessing discriminant validity in variance-based structural equation modeling. Journal of the Academy of Marketing Science, 43(1), 115–135.
  • Ray, S., & Danks, N. (2021). seminr: Building and Estimating Structural Equation Models. R package version 2.3.2.
  • Dataset: Teejmahal20. (2020). Airline Passenger Satisfaction. Kaggle. https://www.kaggle.com/datasets/teejmahal20/airline-passenger-satisfaction