# TUGAS BIOSTATISTIKA INTERMEDIATE
# ANALISIS PENGUKURAN BERULANG HEMOGLOBIN (Hb)
# Disesuaikan untuk data long format: id, kelompok, usia, waktu, Hb
# Nama : Prisceilla Cindy Samosir
# NIM : 2611018010
# ========================= PARAMETER ANALISIS ==========================
FILE_DATA <- "data_hemoglobin_anemia_long.csv"
FOLDER_HASIL <- "hasil_analisis_hb"
KELOMPOK_FOKUS <- "Pemberian TTD + ANC"
INSTALL_PAKET_OTOMATIS <- TRUE
BUKA_HTML_OTOMATIS <- TRUE
TAMPILKAN_GRAFIK_RSTUDIO <- TRUE
# ============================== FUNGSI =================================
get_script_dir <- function() {
# Aman untuk source(), run-by-line di RStudio, dan eksekusi dari Rscript.
args <- commandArgs(trailingOnly = FALSE)
file_arg <- args[grepl("^--file=", args)]
if (length(file_arg) > 0L) {
path <- sub("^--file=", "", file_arg[1])
if (nzchar(path)) {
return(dirname(normalizePath(path, winslash = "/", mustWork = FALSE)))
}
}
if (requireNamespace("rstudioapi", quietly = TRUE) && rstudioapi::isAvailable()) {
p <- tryCatch(rstudioapi::getSourceEditorContext()$path, error = function(e) "")
if (nzchar(p)) return(dirname(normalizePath(p, winslash = "/", mustWork = FALSE)))
}
# Fallback terakhir: folder kerja saat ini.
getwd()
}
cek_paket <- function(paket) {
belum <- paket[!vapply(paket, requireNamespace, logical(1), quietly = TRUE)]
if (length(belum)) {
if (!INSTALL_PAKET_OTOMATIS) {
stop("Paket belum tersedia: ", paste(belum, collapse = ", "))
}
message("Menginstal paket: ", paste(belum, collapse = ", "))
install.packages(belum, repos = "https://cloud.r-project.org")
}
gagal <- paket[!vapply(paket, requireNamespace, logical(1), quietly = TRUE)]
if (length(gagal)) stop("Instalasi gagal: ", paste(gagal, collapse = ", "))
}
fp <- function(p) {
if (is.na(p)) return("tidak tersedia")
if (p < .001) "< 0,001" else paste0("= ", formatC(p, digits = 3, format = "f", decimal.mark = ","))
}
ci_mean <- function(x) {
x <- x[!is.na(x)]
n <- length(x)
if (n == 0L) {
return(data.frame(n = 0L, mean = NA_real_, sd = NA_real_, median = NA_real_,
min = NA_real_, max = NA_real_, se = NA_real_,
lower_ci = NA_real_, upper_ci = NA_real_))
}
mu <- mean(x)
sdv <- if (n > 1L) sd(x) else NA_real_
se <- if (n > 1L) sdv / sqrt(n) else NA_real_
ci <- if (n > 1L) qt(.975, df = n - 1L) * se else NA_real_
data.frame(
n = n,
mean = mu,
sd = sdv,
median = median(x),
min = min(x),
max = max(x),
se = se,
lower_ci = if (is.na(ci)) NA_real_ else mu - ci,
upper_ci = if (is.na(ci)) NA_real_ else mu + ci
)
}
safe_print <- function(x) {
teks <- capture.output(print(x))
cat(paste(teks, collapse = "\n"), "\n")
invisible(teks)
}
# ======================= JALANKAN ANALISIS ============================
jalankan_analisis <- function() {
waktu_mulai <- Sys.time()
paket <- c(
"dplyr", "tidyr", "ggplot2", "afex", "emmeans", "rstatix", "car",
"effectsize", "lme4", "lmerTest", "pbkrtest", "performance",
"ggpubr", "knitr", "htmltools", "base64enc"
)
cek_paket(paket)
suppressPackageStartupMessages({
library(dplyr)
library(tidyr)
library(ggplot2)
library(afex)
library(emmeans)
library(rstatix)
library(car)
library(effectsize)
library(lme4)
library(lmerTest)
library(pbkrtest)
library(performance)
library(ggpubr)
library(knitr)
library(htmltools)
library(base64enc)
})
opsi_lama <- options(contrasts = c("contr.sum", "contr.poly"))
on.exit(options(opsi_lama), add = TRUE)
script_dir <- tryCatch(get_script_dir(), error = function(e) getwd())
kandidat <- unique(c(
file.path(script_dir, FILE_DATA),
file.path(getwd(), FILE_DATA),
normalizePath(FILE_DATA, winslash = "/", mustWork = FALSE)
))
kandidat <- kandidat[file.exists(kandidat)]
if (!length(kandidat)) {
stop(
"Berkas data tidak ditemukan.\n",
"Nama file yang dicari: ", FILE_DATA, "\n",
"Letakkan file CSV di folder script atau folder kerja RStudio."
)
}
path_data <- normalizePath(kandidat[1], winslash = "/", mustWork = TRUE)
folder <- file.path(dirname(path_data), FOLDER_HASIL)
dir.create(folder, recursive = TRUE, showWarnings = FALSE)
folder <- normalizePath(folder, winslash = "/", mustWork = TRUE)
message("Data: ", path_data)
message("Hasil: ", folder)
# ---------- pembaca laporan ----------
isi <- character(0)
status <- data.frame()
tambah <- function(html) isi <<- c(isi, html)
escape <- function(x) htmltools::htmlEscape(paste(x, collapse = "\n"))
narasi <- function(x) tambah(paste0("<p>", escape(x), "</p>"))
output <- function(x) {
teks <- capture.output(print(x))
cat(paste(teks, collapse = "\n"), "\n")
tambah(paste0("<pre>", escape(teks), "</pre>"))
invisible(x)
}
tabel <- function(x, nama) {
x <- as.data.frame(x)
write.csv(x, file.path(folder, paste0(nama, ".csv")),
row.names = FALSE, na = "", fileEncoding = "UTF-8")
if (nrow(x) > 0L) {
tambah(paste0(
"<div class='table'>",
as.character(knitr::kable(x, format = "html", digits = 4)),
"</div>"
))
} else {
narasi("Tidak ada baris untuk tabel ini.")
}
print(x)
invisible(x)
}
gambar <- function(p, nama, lebar = 10, tinggi = 6) {
path <- file.path(folder, paste0(nama, ".png"))
ggplot2::ggsave(path, plot = p, width = lebar, height = tinggi, dpi = 180)
if (TAMPILKAN_GRAFIK_RSTUDIO) {
tryCatch(print(p), error = function(e) message("Plot gagal ditampilkan: ", conditionMessage(e)))
}
tambah(paste0(
"<figure><img src='", base64enc::dataURI(file = path, mime = "image/png"),
"' alt='", escape(nama), "'><figcaption>", escape(nama), "</figcaption></figure>"
))
message("Grafik tersimpan: ", path)
}
bagian <- function(judul, f) {
message("\n>>> ", judul)
tambah(paste0("<h2>", escape(judul), "</h2>"))
warnings <- character(0)
ok <- tryCatch(
withCallingHandlers({
f()
TRUE
}, warning = function(w) {
warnings <<- c(warnings, conditionMessage(w))
invokeRestart("muffleWarning")
}),
error = function(e) {
narasi(paste("BAGIAN GAGAL:", conditionMessage(e)))
message("Bagian gagal: ", conditionMessage(e))
FALSE
}
)
if (length(warnings)) {
narasi(paste("Peringatan estimasi:", paste(unique(warnings), collapse = "; ")))
}
status <<- rbind(
status,
data.frame(
bagian = judul,
status = if (ok) "Selesai" else "Gagal",
peringatan = paste(unique(warnings), collapse = "; "),
stringsAsFactors = FALSE
)
)
invisible(ok)
}
# ---------- 1. membaca data ----------
raw <- tryCatch(
read.csv(path_data, sep = ";", dec = ".", stringsAsFactors = FALSE,
check.names = FALSE, na.strings = c("", "NA", "NaN")),
error = function(e) stop("Gagal membaca CSV: ", conditionMessage(e))
)
names(raw) <- trimws(names(raw))
wajib <- c("id", "kelompok", "usia", "waktu", "Hb")
if (!all(wajib %in% names(raw))) {
stop("Kolom tidak ditemukan: ", paste(setdiff(wajib, names(raw)), collapse = ", "))
}
level_kelompok <- c("Kontrol", "Pemberian TTD", "Pemberian TTD + ANC")
level_waktu <- c("TM0", "TM1", "TM2", "TM3")
raw <- raw |>
mutate(
id = trimws(as.character(id)),
kelompok = trimws(as.character(kelompok)),
waktu = trimws(as.character(waktu)),
usia = suppressWarnings(as.numeric(usia)),
Hb = suppressWarnings(as.numeric(Hb))
)
if (anyNA(raw$id) || any(raw$id == "")) stop("Terdapat id kosong.")
if (anyNA(raw$kelompok) || any(raw$kelompok == "")) stop("Terdapat kelompok kosong.")
if (anyNA(raw$waktu) || any(raw$waktu == "")) stop("Terdapat waktu kosong.")
if (!all(raw$kelompok %in% level_kelompok)) {
stop("Label kelompok tidak dikenali: ",
paste(setdiff(unique(raw$kelompok), level_kelompok), collapse = ", "))
}
if (!all(raw$waktu %in% level_waktu)) {
stop("Label waktu tidak dikenali: ",
paste(setdiff(unique(raw$waktu), level_waktu), collapse = ", "))
}
if (any(!is.finite(raw$Hb) & !is.na(raw$Hb))) stop("Hb memiliki nilai tak hingga.")
if (any(!is.finite(raw$usia) & !is.na(raw$usia))) stop("Usia memiliki nilai tak hingga.")
if (anyDuplicated(raw[c("id", "waktu")])) {
stop("Terdapat duplikasi kombinasi id-waktu.")
}
if (any(tapply(raw$kelompok, raw$id, function(x) length(unique(x))) != 1L)) {
stop("Satu id memiliki lebih dari satu kelompok.")
}
dat_long <- raw |>
mutate(
kelompok = factor(kelompok, levels = level_kelompok),
waktu = factor(waktu, levels = level_waktu),
waktu_num = as.integer(waktu) - 1L
) |>
arrange(id, waktu)
dat_wide <- dat_long |>
select(id, kelompok, usia, waktu, Hb) |>
pivot_wider(names_from = waktu, values_from = Hb)
hb_cols <- level_waktu
complete_id <- complete.cases(dat_wide[, hb_cols])
wide_cc <- dat_wide[complete_id, , drop = FALSE]
long_cc <- dat_long |>
semi_join(wide_cc |> select(id), by = "id")
long_obs <- dat_long |>
filter(!is.na(Hb))
if (n_distinct(wide_cc$id) < 6L) stop("Peserta lengkap terlalu sedikit untuk analisis berulang.")
if (n_distinct(long_cc$id) < 6L) stop("Peserta lengkap terlalu sedikit untuk analisis berulang.")
if (any(table(dat_long$kelompok) < 3L)) stop("Setiap kelompok sebaiknya memiliki cukup observasi.")
write.csv(dat_wide, file.path(folder, "data_wide.csv"), row.names = FALSE, na = "")
write.csv(dat_long, file.path(folder, "data_long.csv"), row.names = FALSE, na = "")
write.csv(wide_cc, file.path(folder, "data_complete_case.csv"), row.names = FALSE, na = "")
# ---------- 2. deskriptif ----------
bagian("1. Sumber dan kelengkapan data", function() {
narasi(paste(
"Data hemoglobin anemia dianalisis dalam format long dengan variabel id, kelompok, usia, waktu, dan Hb.",
"Model utama menggunakan 3 kelompok dan 4 waktu pengukuran (TM0–TM3)."
))
narasi(paste(
n_distinct(dat_long$id), "subjek;",
nrow(dat_long), "baris observasi;",
nrow(wide_cc), "subjek lengkap untuk ANOVA;"
))
tabel(dat_long |>
group_by(kelompok, waktu) |>
summarise(
n_total = n(),
n_tersedia = sum(!is.na(Hb)),
n_hilang = sum(is.na(Hb)),
.groups = "drop"
),
"01_kelengkapan")
tabel(dat_wide |>
group_by(kelompok) |>
summarise(
n_total = n(),
n_lengkap = sum(complete_id),
n_tidak_lengkap = sum(!complete_id),
.groups = "drop"
),
"01_peserta")
narasi("ANOVA memakai peserta yang memiliki Hb lengkap pada seluruh waktu. LMM memakai seluruh pengukuran Hb yang tersedia.")
})
bagian("2. Statistik deskriptif dan visualisasi", function() {
desk <- function(d) {
d |>
group_by(kelompok, waktu) |>
summarise(ci_mean(Hb), .groups = "drop")
}
tabel(desk(dat_long), "02_deskriptif_semua_tersedia")
tabel(desk(long_cc), "02_deskriptif_complete_case")
gambar(
ggplot(dat_long, aes(waktu, Hb, colour = kelompok, group = kelompok)) +
stat_summary(fun = mean, geom = "line", linewidth = 1) +
stat_summary(fun = mean, geom = "point", size = 2) +
stat_summary(fun.data = mean_cl_normal, geom = "errorbar", width = .2) +
theme_bw(base_size = 12) +
theme(legend.position = "bottom") +
labs(title = "Profil Hb: semua pengukuran tersedia", x = "Waktu", y = "Hb (g/dL)", colour = "Kelompok"),
"02_profil_hb",
10, 6
)
gambar(
ggplot(long_cc, aes(waktu, Hb, group = id)) +
geom_line(alpha = .15) +
stat_summary(aes(group = kelompok), fun = mean, geom = "line", linewidth = 1.1) +
facet_wrap(~kelompok) +
theme_bw(base_size = 12) +
labs(title = "Lintasan individu: peserta lengkap", x = "Waktu", y = "Hb (g/dL)"),
"02_lintasan",
12, 7
)
for (g in level_kelompok) {
z <- wide_cc[wide_cc$kelompok == g, hb_cols, drop = FALSE]
narasi(paste("Kovarians dan korelasi dalam kelompok:", g))
output(round(cov(z), 3))
output(round(cor(z), 3))
}
})
# ---------- 3. asumsi ----------
bagian("3. Asumsi ANOVA pada peserta lengkap", function() {
tabel(long_cc |> group_by(kelompok, waktu) |> identify_outliers(Hb), "03_outlier")
tabel(long_cc |> group_by(kelompok, waktu) |> shapiro_test(Hb), "03_shapiro")
gambar(
ggpubr::ggqqplot(long_cc, "Hb", facet.by = c("kelompok", "waktu")),
"03_qq_per_sel",
12, 9
)
tabel(long_cc |> group_by(waktu) |> levene_test(Hb ~ kelompok), "03_levene")
tabel(rstatix::box_m(wide_cc[, hb_cols], wide_cc$kelompok), "03_box_m")
narasi("p > 0,05 tidak membuktikan normalitas atau kesamaan varians. Outlier statistik dipertahankan dalam analisis utama kecuali ada alasan substantif untuk mengecualikannya.")
})
# ---------- 4. RM ANOVA satu kelompok ----------
aov1 <- NULL
aov2 <- NULL
lmm <- NULL
lmm_ri <- NULL
lmm_rs <- NULL
bagian("4. RM ANOVA satu kelompok", function() {
if (!(KELOMPOK_FOKUS %in% level_kelompok)) {
stop("KELOMPOK_FOKUS tidak dikenali.")
}
d1 <- long_cc |>
filter(kelompok == KELOMPOK_FOKUS)
if (n_distinct(d1$id) < 4L) stop("Jumlah peserta lengkap pada kelompok fokus terlalu sedikit.")
narasi(paste("Kelompok fokus:", KELOMPOK_FOKUS, "—", n_distinct(d1$id), "peserta lengkap."))
gambar(
ggpubr::ggqqplot(d1, "Hb", facet.by = "waktu") +
labs(title = paste("Q-Q Hb per waktu:", KELOMPOK_FOKUS)),
"04_qq_kelompok_fokus",
10, 6
)
output(rstatix::anova_test(data = d1, dv = Hb, wid = id, within = waktu, effect.size = "pes"))
aov1 <<- afex::aov_ez(
id = "id",
dv = "Hb",
data = d1,
within = "waktu",
anova_table = list(es = c("ges", "pes"), correction = "GG")
)
output(aov1)
output(summary(aov1))
output(aov1$Anova)
tabel(as.data.frame(aov1$anova_table) |> tibble::rownames_to_column("efek"), "04_anova_gg")
output(effectsize::eta_squared(aov1, partial = TRUE))
p <- as.numeric(aov1$anova_table[1, "Pr(>F)"])
narasi(paste0("Efek waktu dengan koreksi GG: p ", fp(p), "."))
})
bagian("5. Post hoc satu kelompok dan tren waktu", function() {
if (is.null(aov1)) stop("Model RM ANOVA belum tersedia.")
em <- emmeans::emmeans(aov1, ~ waktu, model = "multivariate")
tabel(summary(em), "05_emmeans")
tabel(summary(pairs(em, adjust = "bonferroni"), infer = c(TRUE, TRUE)), "05_pasangan_bonferroni")
tabel(summary(contrast(em, "trt.vs.ctrl", ref = 1, adjust = "holm"), infer = c(TRUE, TRUE)), "05_vs_baseline_holm")
pol <- poly(0:3, degree = 3)
koef <- list(
linear = pol[, 1],
kuadratik = pol[, 2],
kubik = pol[, 3]
)
tabel(summary(contrast(em, method = koef, adjust = "holm"), infer = c(TRUE, TRUE)), "05_tren")
narasi("Tren waktu menggunakan urutan TM0–TM3. Koefisien ortogonal dinormalisasi sehingga estimate bukan kenaikan Hb per satuan waktu.")
})
# ---------- 6. nonparametrik ----------
bagian("6. Pembanding nonparametrik satu kelompok", function() {
if (is.null(aov1)) stop("RM ANOVA belum tersedia.")
d1 <- long_cc |> filter(kelompok == KELOMPOK_FOKUS)
tabel(rstatix::friedman_test(d1, Hb ~ waktu | id), "06_friedman")
tabel(rstatix::friedman_effsize(d1, Hb ~ waktu | id), "06_kendall_w")
w <- wide_cc[wide_cc$kelompok == KELOMPOK_FOKUS, ]
ij <- combn(seq_along(hb_cols), 2)
hasil <- lapply(seq_len(ncol(ij)), function(j) {
a <- ij[1, j]
b <- ij[2, j]
wt <- wilcox.test(w[[hb_cols[a]]], w[[hb_cols[b]]], paired = TRUE, exact = FALSE)
data.frame(
waktu_1 = level_waktu[a],
waktu_2 = level_waktu[b],
n_pasangan = nrow(w),
V = unname(wt$statistic),
p = wt$p.value,
stringsAsFactors = FALSE
)
}) |> bind_rows() |> mutate(p_bonferroni = p.adjust(p, "bonferroni"))
tabel(hasil, "06_wilcoxon")
})
# ---------- 7. mixed design ANOVA ----------
bagian("7. Mixed Design ANOVA kelompok × waktu", function() {
aov2 <<- afex::aov_ez(
id = "id",
dv = "Hb",
data = long_cc,
between = "kelompok",
within = "waktu",
anova_table = list(es = c("ges", "pes"), correction = "GG")
)
output(aov2)
output(summary(aov2))
output(aov2$Anova)
tabel(as.data.frame(aov2$anova_table) |> tibble::rownames_to_column("efek"), "07_anova_gg")
output(effectsize::eta_squared(aov2, partial = TRUE))
})
bagian("7b. Visualisasi interaksi kelompok dan waktu", function() {
if (is.null(aov2)) stop("Model Mixed Design ANOVA belum tersedia.")
gambar(
afex::afex_plot(aov2, x = "waktu", trace = "kelompok", error = "within",
mapping = c("colour", "shape", "linetype")) +
labs(title = "Mixed Design ANOVA: kelompok × waktu", x = "Waktu", y = "Hb (g/dL)"),
"07_interaksi_anova",
11, 7
)
})
# ---------- 8. efek sederhana ----------
bagian("8. Efek sederhana dan post hoc antarkelompok", function() {
if (is.null(aov2)) stop("Model Mixed Design ANOVA belum tersedia.")
jt1 <- as.data.frame(joint_tests(aov2, by = "kelompok", model = "multivariate"))
jt1$p_holm_antar_kelompok <- p.adjust(jt1$p.value, "holm")
tabel(jt1, "08_efek_waktu_per_kelompok")
jt2 <- as.data.frame(joint_tests(aov2, by = "waktu", model = "multivariate"))
jt2$p_holm_antar_waktu <- p.adjust(jt2$p.value, "holm")
tabel(jt2, "08_efek_kelompok_per_waktu")
ew <- emmeans(aov2, ~ waktu | kelompok, model = "multivariate")
tabel(summary(pairs(ew, adjust = "holm"), infer = c(TRUE, TRUE)), "08_pairwise_waktu_dalam_kelompok")
ek <- emmeans(aov2, ~ kelompok | waktu, model = "multivariate")
tabel(summary(pairs(ek, adjust = "tukey"), infer = c(TRUE, TRUE)), "08_kelompok_per_waktu")
})
# ---------- 9. lmm ----------
bagian("9. LMM dengan semua pengukuran tersedia", function() {
lmm_ri <<- lmerTest::lmer(
Hb ~ kelompok * waktu + (1 | id),
data = long_obs,
REML = TRUE,
control = lmerControl(optimizer = "bobyqa")
)
lmm_rs <<- tryCatch(
lmerTest::lmer(
Hb ~ kelompok * waktu + (1 + waktu_num | id),
data = long_obs,
REML = TRUE,
control = lmerControl(optimizer = "bobyqa", optCtrl = list(maxfun = 200000))
),
error = function(e) {
narasi(paste("Random slope gagal diestimasi; memakai random intercept:", conditionMessage(e)))
NULL
}
)
info_model <- function(m, nama) {
data.frame(
model = nama,
n_observasi = nobs(m),
AIC = AIC(m),
BIC = BIC(m),
singular = lme4::isSingular(m),
konvergensi = paste(m@optinfo$conv$lme4$messages, collapse = "; "),
stringsAsFactors = FALSE
)
}
info_rs <- if (!is.null(lmm_rs)) info_model(lmm_rs, "random_intercept_slope") else NULL
tabel(bind_rows(info_model(lmm_ri, "random_intercept"), info_rs), "09_perbandingan_model")
layak <- !is.null(lmm_rs) && !isSingular(lmm_rs) && is.null(lmm_rs@optinfo$conv$lme4$messages)
pakai_rs <- layak && AIC(lmm_rs) < AIC(lmm_ri)
lmm <<- if (pakai_rs) lmm_rs else lmm_ri
narasi(paste("Model yang dipakai:", if (pakai_rs) "random intercept + slope waktu." else "random intercept."))
output(summary(lmm))
output(VarCorr(lmm))
hasil <- anova(lmm, type = 3, ddf = "Kenward-Roger")
tabel(as.data.frame(hasil) |> tibble::rownames_to_column("efek"), "09_lmm_kr")
output(performance::icc(lmm_ri))
})
bagian("10. Diagnostik LMM", function() {
if (is.null(lmm)) stop("LMM belum tersedia.")
diag <- data.frame(
fitted = fitted(lmm),
residual = resid(lmm),
kelompok = long_obs$kelompok,
waktu = long_obs$waktu
)
gambar(
ggplot(diag, aes(sample = residual)) + stat_qq() + stat_qq_line() + theme_bw() +
labs(title = "Q-Q residual LMM"),
"10_qq_residual",
8, 6
)
gambar(
ggplot(diag, aes(fitted, residual, colour = kelompok)) +
geom_point(alpha = .35) +
geom_hline(yintercept = 0, linetype = 2) +
theme_bw() +
labs(x = "Prediksi", y = "Residual", title = "Residual terhadap prediksi"),
"10_residual_prediksi",
8, 6
)
gambar(
ggplot(diag, aes(waktu, residual, fill = kelompok)) +
geom_boxplot() + theme_bw() +
labs(title = "Residual menurut kelompok dan waktu"),
"10_residual_sel",
10, 6
)
output(performance::check_singularity(lmm))
output(performance::check_convergence(lmm))
})
bagian("11. Informasi sesi dan batas interpretasi", function() {
narasi("Taraf signifikansi ditetapkan pada α = 0,05. Interpretasi mempertimbangkan nilai p, ukuran efek, interval kepercayaan, dan arah perubahan.")
output(sessionInfo())
})
# ---------- simpan laporan ----------
tambah("<h2>Status setiap bagian</h2>")
tabel(status, "status_analisis")
css <- "body{font:16px/1.6 Arial,sans-serif;max-width:1100px;margin:40px auto;padding:0 24px;color:#18332e}h1,h2{color:#0b6b4f}h2{border-top:1px solid #d9e4e0;padding-top:24px}pre{background:#f0f5f3;padding:16px;overflow:auto;font-size:13px;color:#182822}.table{overflow:auto}table{border-collapse:collapse;font-size:13px;width:100%}td,th{border:1px solid #d9e4e0;padding:8px;text-align:left}th{background:#e7f2ee}img{max-width:100%}summary{cursor:pointer;font-weight:bold}"
selesai <- sum(status$status == "Selesai", na.rm = TRUE)
header <- paste0(
"<!doctype html><html lang='id'><head><meta charset='UTF-8'>",
"<meta name='viewport' content='width=device-width,initial-scale=1'>",
"<title>Analisis Hb — pengukuran berulang</title><style>", css, "</style></head><body>",
"<h1>Analisis Pengukuran Berulang Hemoglobin</h1>",
"<table class='identitas'><tbody>",
"<tr><th>File data</th><td>", escape(FILE_DATA), "</td></tr>",
"<tr><th>Kelompok fokus</th><td>", escape(KELOMPOK_FOKUS), "</td></tr>",
"<tr><th>Waktu</th><td>TM0, TM1, TM2, TM3</td></tr>",
"<tr><th>Subjek</th><td>", n_distinct(dat_long$id), "</td></tr>",
"</tbody></table>",
"<p>Dibuat: ", escape(format(Sys.time())), " · Bagian selesai: ", selesai, "/", nrow(status), "</p>"
)
path_html <- file.path(folder, "laporan_analisis_hb.html")
writeLines(c(header, isi, "</body></html>"), path_html, useBytes = TRUE)
saveRDS(
list(
data_wide = dat_wide,
data_long = dat_long,
aov_satu = aov1,
aov_campuran = aov2,
lmm = lmm,
status = status
),
file.path(folder, "objek_analisis.rds")
)
writeLines(capture.output(sessionInfo()), file.path(folder, "sessionInfo.txt"))
message("\nLaporan tersedia: ", path_html)
message("Durasi: ", round(as.numeric(difftime(Sys.time(), waktu_mulai, units = "secs")), 1), " detik.")
if (any(status$status == "Gagal")) {
message("Terdapat bagian analisis yang gagal; rinciannya tercatat dalam status_analisis.csv dan laporan HTML.")
}
if (BUKA_HTML_OTOMATIS && interactive()) {
try(utils::browseURL(path_html), silent = TRUE)
}
invisible(list(folder = folder, status = status))
}
# ============================ EKSEKUSI ================================
hasil_analisis <- jalankan_analisis()
## Registered S3 method overwritten by 'lme4':
## method from
## na.action.merMod car
## Data: D:/2026/S2 UNMUL/TUGAS/BIOSTAT/data_hemoglobin_anemia_long.csv
## Hasil: D:/2026/S2 UNMUL/TUGAS/BIOSTAT/hasil_analisis_hb
##
## >>> 1. Sumber dan kelengkapan data
## kelompok waktu n_total n_tersedia n_hilang
## 1 Kontrol TM0 30 30 0
## 2 Kontrol TM1 30 30 0
## 3 Kontrol TM2 30 30 0
## 4 Kontrol TM3 30 30 0
## 5 Pemberian TTD TM0 30 30 0
## 6 Pemberian TTD TM1 30 30 0
## 7 Pemberian TTD TM2 30 30 0
## 8 Pemberian TTD TM3 30 30 0
## 9 Pemberian TTD + ANC TM0 30 30 0
## 10 Pemberian TTD + ANC TM1 30 30 0
## 11 Pemberian TTD + ANC TM2 30 30 0
## 12 Pemberian TTD + ANC TM3 30 30 0
## kelompok n_total n_lengkap n_tidak_lengkap
## 1 Kontrol 30 90 0
## 2 Pemberian TTD 30 90 0
## 3 Pemberian TTD + ANC 30 90 0
##
## >>> 2. Statistik deskriptif dan visualisasi
## kelompok waktu n mean sd median min max se
## 1 Kontrol TM0 30 8.736667 0.5301160 8.70 8.0 9.8 0.09678550
## 2 Kontrol TM1 30 8.773333 0.5362278 8.65 7.9 9.9 0.09790135
## 3 Kontrol TM2 30 8.840000 0.5468405 8.85 8.0 10.0 0.09983895
## 4 Kontrol TM3 30 8.870000 0.5831839 8.80 7.9 10.1 0.10647432
## 5 Pemberian TTD TM0 30 8.763333 0.5792286 8.80 7.5 9.8 0.10575219
## 6 Pemberian TTD TM1 30 9.136667 0.6189442 9.15 7.9 10.2 0.11300324
## 7 Pemberian TTD TM2 30 9.503333 0.6122223 9.55 8.0 10.4 0.11177598
## 8 Pemberian TTD TM3 30 9.963333 0.6702976 9.95 8.6 10.9 0.12237904
## 9 Pemberian TTD + ANC TM0 30 8.780000 0.6509410 8.85 7.5 9.8 0.11884502
## 10 Pemberian TTD + ANC TM1 30 9.310000 0.6849868 9.35 8.0 10.3 0.12506090
## 11 Pemberian TTD + ANC TM2 30 9.930000 0.7066921 9.85 8.7 10.9 0.12902374
## 12 Pemberian TTD + ANC TM3 30 10.250000 0.5764038 10.20 9.3 10.9 0.10523646
## lower_ci upper_ci
## 1 8.538718 8.934615
## 2 8.573103 8.973564
## 3 8.635806 9.044194
## 4 8.652236 9.087764
## 5 8.547046 8.979621
## 6 8.905549 9.367784
## 7 9.274726 9.731941
## 8 9.713040 10.213627
## 9 8.536935 9.023065
## 10 9.054222 9.565778
## 11 9.666117 10.193883
## 12 10.034767 10.465233
## kelompok waktu n mean sd median min max se
## 1 Kontrol TM0 30 8.736667 0.5301160 8.70 8.0 9.8 0.09678550
## 2 Kontrol TM1 30 8.773333 0.5362278 8.65 7.9 9.9 0.09790135
## 3 Kontrol TM2 30 8.840000 0.5468405 8.85 8.0 10.0 0.09983895
## 4 Kontrol TM3 30 8.870000 0.5831839 8.80 7.9 10.1 0.10647432
## 5 Pemberian TTD TM0 30 8.763333 0.5792286 8.80 7.5 9.8 0.10575219
## 6 Pemberian TTD TM1 30 9.136667 0.6189442 9.15 7.9 10.2 0.11300324
## 7 Pemberian TTD TM2 30 9.503333 0.6122223 9.55 8.0 10.4 0.11177598
## 8 Pemberian TTD TM3 30 9.963333 0.6702976 9.95 8.6 10.9 0.12237904
## 9 Pemberian TTD + ANC TM0 30 8.780000 0.6509410 8.85 7.5 9.8 0.11884502
## 10 Pemberian TTD + ANC TM1 30 9.310000 0.6849868 9.35 8.0 10.3 0.12506090
## 11 Pemberian TTD + ANC TM2 30 9.930000 0.7066921 9.85 8.7 10.9 0.12902374
## 12 Pemberian TTD + ANC TM3 30 10.250000 0.5764038 10.20 9.3 10.9 0.10523646
## lower_ci upper_ci
## 1 8.538718 8.934615
## 2 8.573103 8.973564
## 3 8.635806 9.044194
## 4 8.652236 9.087764
## 5 8.547046 8.979621
## 6 8.905549 9.367784
## 7 9.274726 9.731941
## 8 9.713040 10.213627
## 9 8.536935 9.023065
## 10 9.054222 9.565778
## 11 9.666117 10.193883
## 12 10.034767 10.465233
## Grafik tersimpan: D:/2026/S2 UNMUL/TUGAS/BIOSTAT/hasil_analisis_hb/02_profil_hb.png
## Grafik tersimpan: D:/2026/S2 UNMUL/TUGAS/BIOSTAT/hasil_analisis_hb/02_lintasan.png
## TM0 TM1 TM2 TM3
## TM0 0.281 0.276 0.277 0.299
## TM1 0.276 0.288 0.273 0.290
## TM2 0.277 0.273 0.299 0.308
## TM3 0.299 0.290 0.308 0.340
## TM0 TM1 TM2 TM3
## TM0 1.000 0.970 0.956 0.967
## TM1 0.970 1.000 0.930 0.926
## TM2 0.956 0.930 1.000 0.966
## TM3 0.967 0.926 0.966 1.000
## TM0 TM1 TM2 TM3
## TM0 0.336 0.340 0.339 0.352
## TM1 0.340 0.383 0.345 0.359
## TM2 0.339 0.345 0.375 0.348
## TM3 0.352 0.359 0.348 0.449
## TM0 TM1 TM2 TM3
## TM0 1.000 0.948 0.957 0.908
## TM1 0.948 1.000 0.911 0.864
## TM2 0.957 0.911 1.000 0.848
## TM3 0.908 0.864 0.848 1.000
## TM0 TM1 TM2 TM3
## TM0 0.424 0.430 0.437 0.350
## TM1 0.430 0.469 0.452 0.365
## TM2 0.437 0.452 0.499 0.378
## TM3 0.350 0.365 0.378 0.332
## TM0 TM1 TM2 TM3
## TM0 1.000 0.965 0.950 0.934
## TM1 0.965 1.000 0.934 0.924
## TM2 0.950 0.934 1.000 0.927
## TM3 0.934 0.924 0.927 1.000
##
## >>> 3. Asumsi ANOVA pada peserta lengkap
## kelompok waktu id usia Hb waktu_num is.outlier is.extreme
## 1 Kontrol TM0 P007 24 9.8 0 TRUE FALSE
## 2 Kontrol TM0 P022 29 9.8 0 TRUE FALSE
## 3 Kontrol TM3 P007 24 10.1 3 TRUE FALSE
## 4 Pemberian TTD TM2 P045 39 8.0 2 TRUE FALSE
## 5 Pemberian TTD TM2 P048 30 8.1 2 TRUE FALSE
## kelompok waktu variable statistic p
## 1 Kontrol TM0 Hb 0.9256333 0.037641323
## 2 Kontrol TM1 Hb 0.9425960 0.106856847
## 3 Kontrol TM2 Hb 0.9572106 0.262429244
## 4 Kontrol TM3 Hb 0.9503149 0.172375376
## 5 Pemberian TTD TM0 Hb 0.9720809 0.597542969
## 6 Pemberian TTD TM1 Hb 0.9609214 0.326967390
## 7 Pemberian TTD TM2 Hb 0.9457828 0.130226887
## 8 Pemberian TTD TM3 Hb 0.9534881 0.209441459
## 9 Pemberian TTD + ANC TM0 Hb 0.9586031 0.285198854
## 10 Pemberian TTD + ANC TM1 Hb 0.9516830 0.187515319
## 11 Pemberian TTD + ANC TM2 Hb 0.9281162 0.043769627
## 12 Pemberian TTD + ANC TM3 Hb 0.8668000 0.001423043
## Grafik tersimpan: D:/2026/S2 UNMUL/TUGAS/BIOSTAT/hasil_analisis_hb/03_qq_per_sel.png
## waktu df1 df2 statistic p
## 1 TM0 2 87 1.1781739 0.3127031
## 2 TM1 2 87 1.7669237 0.1769335
## 3 TM2 2 87 2.0384721 0.1364017
## 4 TM3 2 87 0.8676822 0.4235254
## statistic p.value parameter
## 1 43.56665 0.00171942 20
## method
## 1 Box's M-test for Homogeneity of Covariance Matrices
##
## >>> 4. RM ANOVA satu kelompok
## Grafik tersimpan: D:/2026/S2 UNMUL/TUGAS/BIOSTAT/hasil_analisis_hb/04_qq_kelompok_fokus.png
## ANOVA Table (type III tests)
##
## $ANOVA
## Effect DFn DFd F p p<.05 pes
## 1 waktu 3 87 441.902 1.6e-52 * 0.938
##
## $`Mauchly's Test for Sphericity`
## Effect W p p<.05
## 1 waktu 0.776 0.219
##
## $`Sphericity Corrections`
## Effect GGe DF[GG] p[GG] p[GG]<.05 HFe DF[HF] p[HF]
## 1 waktu 0.884 2.65, 76.86 1.08e-46 * 0.981 2.94, 85.31 1.49e-51
## p[HF]<.05
## 1 *
##
## Anova Table (Type 3 tests)
##
## Response: Hb
## Effect df MSE F ges pes p.value
## 1 waktu 2.65, 76.86 0.03 441.90 *** .435 .938 <.001
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '+' 0.1 ' ' 1
##
## Sphericity correction method: GG
##
## Univariate Type III Repeated-Measures ANOVA Assuming Sphericity
##
## Sum Sq num Df Error SS den Df F value Pr(>F)
## (Intercept) 10984.4 1 47.486 29 6708.3 < 2.2e-16 ***
## waktu 38.5 3 2.527 87 441.9 < 2.2e-16 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
##
## Mauchly Tests for Sphericity
##
## Test statistic p-value
## waktu 0.77597 0.21857
##
##
## Greenhouse-Geisser and Huynh-Feldt Corrections
## for Departure from Sphericity
##
## GG eps Pr(>F[GG])
## waktu 0.8835 < 2.2e-16 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## HF eps Pr(>F[HF])
## waktu 0.9806045 1.49376e-51
##
## Type III Repeated Measures MANOVA Tests: Pillai test statistic
## Df test stat approx F num Df den Df Pr(>F)
## (Intercept) 1 0.99570 6708.3 1 29 < 2.2e-16 ***
## waktu 1 0.98236 501.3 3 27 < 2.2e-16 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## efek num Df den Df MSE F pes ges Pr(>F)
## 1 waktu 2.65051 76.86478 0.03287917 441.9022 0.9384161 0.4350298 1.081739e-46
## # Effect Size for ANOVA (Type III)
##
## Parameter | Eta2 (partial) | 95% CI
## -----------------------------------------
## waktu | 0.94 | [0.92, 1.00]
##
## - One-sided CIs: upper bound fixed at [1.00].
##
## >>> 5. Post hoc satu kelompok dan tren waktu
## waktu emmean SE df lower.CL upper.CL
## TM0 8.78 0.1188450 29 8.536935 9.023065
## TM1 9.31 0.1250609 29 9.054222 9.565778
## TM2 9.93 0.1290237 29 9.666117 10.193883
## TM3 10.25 0.1052365 29 10.034767 10.465233
##
## Confidence level used: 0.95
## contrast estimate SE df lower.CL upper.CL t.ratio p.value
## TM0 - TM1 -0.53 0.03292276 29 -0.6232225 -0.4367775 -16.098 <0.0001
## TM0 - TM2 -1.15 0.04032911 29 -1.2641940 -1.0358060 -28.515 <0.0001
## TM0 - TM3 -1.47 0.04292469 29 -1.5915435 -1.3484565 -34.246 <0.0001
## TM1 - TM2 -0.62 0.04633710 29 -0.7512059 -0.4887941 -13.380 <0.0001
## TM1 - TM3 -0.94 0.04880173 29 -1.0781847 -0.8018153 -19.262 <0.0001
## TM2 - TM3 -0.32 0.05037788 29 -0.4626476 -0.1773524 -6.352 <0.0001
##
## Confidence level used: 0.95
## Conf-level adjustment: bonferroni method for 6 estimates
## P value adjustment: bonferroni method for 6 tests
## contrast estimate SE df lower.CL upper.CL t.ratio p.value
## TM1 - TM0 0.53 0.03292276 29 0.4463463 0.6136537 16.098 <0.0001
## TM2 - TM0 1.15 0.04032911 29 1.0475274 1.2524726 28.515 <0.0001
## TM3 - TM0 1.47 0.04292469 29 1.3609323 1.5790677 34.246 <0.0001
##
## Confidence level used: 0.95
## Conf-level adjustment: bonferroni method for 3 estimates
## P value adjustment: holm method for 3 tests
## contrast estimate SE df lower.CL upper.CL t.ratio p.value
## linear 1.1247422 0.03153451 29 1.0446159 1.2048685 35.667 <0.0001
## kuadratik -0.1050000 0.03016716 29 -0.1816520 -0.0283480 -3.481 0.0032
## kubik -0.0872067 0.03162914 29 -0.1675734 -0.0068399 -2.757 0.0100
##
## Confidence level used: 0.95
## Conf-level adjustment: bonferroni method for 3 estimates
## P value adjustment: holm method for 3 tests
##
## >>> 6. Pembanding nonparametrik satu kelompok
## .y. n statistic df p method
## 1 Hb 30 86.41611 3 1.288722e-18 Friedman test
## .y. n effsize method magnitude
## 1 Hb 30 0.960179 Kendall W large
## waktu_1 waktu_2 n_pasangan V p p_bonferroni
## 1 TM0 TM1 30 0.0 1.735764e-06 1.041458e-05
## 2 TM0 TM2 30 0.0 1.781217e-06 1.068730e-05
## 3 TM0 TM3 30 0.0 1.758922e-06 1.055353e-05
## 4 TM1 TM2 30 0.0 1.734667e-06 1.040800e-05
## 5 TM1 TM3 30 0.0 1.780097e-06 1.068058e-05
## 6 TM2 TM3 30 8.5 9.781069e-06 5.868641e-05
##
## >>> 7. Mixed Design ANOVA kelompok × waktu
## Anova Table (Type 3 tests)
##
## Response: Hb
## Effect df MSE F ges pes p.value
## 1 kelompok 2, 87 1.41 13.06 *** .221 .231 <.001
## 2 waktu 2.44, 212.67 0.03 545.71 *** .257 .862 <.001
## 3 kelompok:waktu 4.89, 212.67 0.03 107.20 *** .120 .711 <.001
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '+' 0.1 ' ' 1
##
## Sphericity correction method: GG
##
## Univariate Type III Repeated-Measures ANOVA Assuming Sphericity
##
## Sum Sq num Df Error SS den Df F value Pr(>F)
## (Intercept) 30723.0 1 122.619 87 21798.351 < 2.2e-16 ***
## kelompok 36.8 2 122.619 87 13.061 1.096e-05 ***
## waktu 44.9 3 7.156 261 545.705 < 2.2e-16 ***
## kelompok:waktu 17.6 6 7.156 261 107.198 < 2.2e-16 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
##
## Mauchly Tests for Sphericity
##
## Test statistic p-value
## waktu 0.64119 3.6178e-07
## kelompok:waktu 0.64119 3.6178e-07
##
##
## Greenhouse-Geisser and Huynh-Feldt Corrections
## for Departure from Sphericity
##
## GG eps Pr(>F[GG])
## waktu 0.81483 < 2.2e-16 ***
## kelompok:waktu 0.81483 < 2.2e-16 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## HF eps Pr(>F[HF])
## waktu 0.8401417 1.216338e-94
## kelompok:waktu 0.8401417 3.923773e-57
##
## Type III Repeated Measures MANOVA Tests: Pillai test statistic
## Df test stat approx F num Df den Df Pr(>F)
## (Intercept) 1 0.99602 21798.4 1 87 < 2.2e-16 ***
## kelompok 2 0.23092 13.1 2 87 1.096e-05 ***
## waktu 1 0.96567 797.1 3 85 < 2.2e-16 ***
## kelompok:waktu 2 1.04668 31.5 6 172 < 2.2e-16 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## efek num Df den Df MSE F pes ges
## 1 kelompok 2.000000 87.0000 1.40941858 13.06096 0.2309183 0.2209997
## 2 waktu 2.444498 212.6713 0.03364621 545.70504 0.8624952 0.2569773
## 3 kelompok:waktu 4.888996 212.6713 0.03364621 107.19838 0.7113439 0.1196247
## Pr(>F)
## 1 1.095980e-05
## 2 7.064013e-92
## 3 1.720301e-55
## # Effect Size for ANOVA (Type III)
##
## Parameter | Eta2 (partial) | 95% CI
## ----------------------------------------------
## kelompok | 0.23 | [0.11, 1.00]
## waktu | 0.86 | [0.84, 1.00]
## kelompok:waktu | 0.71 | [0.67, 1.00]
##
## - One-sided CIs: upper bound fixed at [1.00].
##
## >>> 7b. Visualisasi interaksi kelompok dan waktu
## Grafik tersimpan: D:/2026/S2 UNMUL/TUGAS/BIOSTAT/hasil_analisis_hb/07_interaksi_anova.png
##
## >>> 8. Efek sederhana dan post hoc antarkelompok
## model term kelompok df1 df2 F.ratio p.value
## 1 waktu Kontrol 3 87 5.776 1.194158e-03
## 4 waktu Pemberian TTD 3 87 402.108 7.421211e-51
## 7 waktu Pemberian TTD + ANC 3 87 735.742 1.118461e-61
## p_holm_antar_kelompok
## 1 1.194158e-03
## 4 1.484242e-50
## 7 3.355382e-61
## model term waktu df1 df2 F.ratio p.value p_holm_antar_waktu
## 1 kelompok TM0 2 87 0.041 9.595254e-01 9.595254e-01
## 3 kelompok TM1 2 87 5.923 3.876646e-03 7.753292e-03
## 5 kelompok TM2 2 87 23.143 8.731459e-09 2.619438e-08
## 7 kelompok TM3 2 87 42.553 1.294062e-13 5.176249e-13
## kelompok = Kontrol:
## contrast estimate SE df lower.CL upper.CL t.ratio p.value
## TM0 - TM1 -0.0366667 0.03126535 87 -0.1210807 0.0477474 -1.173 0.4882
## TM0 - TM2 -0.1033333 0.03433365 87 -0.1960315 -0.0106351 -3.010 0.0171
## TM0 - TM3 -0.1333333 0.04193439 87 -0.2465530 -0.0201137 -3.180 0.0123
## TM1 - TM2 -0.0666667 0.04377465 87 -0.1848548 0.0515215 -1.523 0.3942
## TM1 - TM3 -0.0966667 0.05111902 87 -0.2346841 0.0413508 -1.891 0.2478
## TM2 - TM3 -0.0300000 0.05022173 87 -0.1655948 0.1055948 -0.597 0.5518
##
## kelompok = Pemberian TTD:
## contrast estimate SE df lower.CL upper.CL t.ratio p.value
## TM0 - TM1 -0.3733333 0.03126535 87 -0.4577474 -0.2889193 -11.941 <0.0001
## TM0 - TM2 -0.7400000 0.03433365 87 -0.8326982 -0.6473018 -21.553 <0.0001
## TM0 - TM3 -1.2000000 0.04193439 87 -1.3132196 -1.0867804 -28.616 <0.0001
## TM1 - TM2 -0.3666667 0.04377465 87 -0.4848548 -0.2484785 -8.376 <0.0001
## TM1 - TM3 -0.8266667 0.05111902 87 -0.9646841 -0.6886492 -16.171 <0.0001
## TM2 - TM3 -0.4600000 0.05022173 87 -0.5955948 -0.3244052 -9.159 <0.0001
##
## kelompok = Pemberian TTD + ANC:
## contrast estimate SE df lower.CL upper.CL t.ratio p.value
## TM0 - TM1 -0.5300000 0.03126535 87 -0.6144140 -0.4455860 -16.952 <0.0001
## TM0 - TM2 -1.1500000 0.03433365 87 -1.2426982 -1.0573018 -33.495 <0.0001
## TM0 - TM3 -1.4700000 0.04193439 87 -1.5832196 -1.3567804 -35.055 <0.0001
## TM1 - TM2 -0.6200000 0.04377465 87 -0.7381882 -0.5018118 -14.163 <0.0001
## TM1 - TM3 -0.9400000 0.05111902 87 -1.0780174 -0.8019826 -18.388 <0.0001
## TM2 - TM3 -0.3200000 0.05022173 87 -0.4555948 -0.1844052 -6.372 <0.0001
##
## Confidence level used: 0.95
## Conf-level adjustment: bonferroni method for 6 estimates
## P value adjustment: holm method for 6 tests
## waktu = TM0:
## contrast estimate SE df lower.CL
## Kontrol - Pemberian TTD -0.0266667 0.1520419 87 -0.3892074
## Kontrol - (Pemberian TTD + ANC) -0.0433333 0.1520419 87 -0.4058741
## Pemberian TTD - (Pemberian TTD + ANC) -0.0166667 0.1520419 87 -0.3792074
## upper.CL t.ratio p.value
## 0.3358741 -0.175 0.9832
## 0.3192074 -0.285 0.9562
## 0.3458741 -0.110 0.9934
##
## waktu = TM1:
## contrast estimate SE df lower.CL
## Kontrol - Pemberian TTD -0.3633333 0.1591532 87 -0.7428310
## Kontrol - (Pemberian TTD + ANC) -0.5366667 0.1591532 87 -0.9161643
## Pemberian TTD - (Pemberian TTD + ANC) -0.1733333 0.1591532 87 -0.5528310
## upper.CL t.ratio p.value
## 0.0161643 -2.283 0.0636
## -0.1571690 -3.372 0.0032
## 0.2061643 -1.089 0.5234
##
## waktu = TM2:
## contrast estimate SE df lower.CL
## Kontrol - Pemberian TTD -0.6633333 0.1614699 87 -1.0483551
## Kontrol - (Pemberian TTD + ANC) -1.0900000 0.1614699 87 -1.4750218
## Pemberian TTD - (Pemberian TTD + ANC) -0.4266667 0.1614699 87 -0.8116884
## upper.CL t.ratio p.value
## -0.2783116 -4.108 0.0003
## -0.7049782 -6.750 <0.0001
## -0.0416449 -2.642 0.0262
##
## waktu = TM3:
## contrast estimate SE df lower.CL
## Kontrol - Pemberian TTD -1.0933333 0.1578778 87 -1.4697898
## Kontrol - (Pemberian TTD + ANC) -1.3800000 0.1578778 87 -1.7564565
## Pemberian TTD - (Pemberian TTD + ANC) -0.2866667 0.1578778 87 -0.6631232
## upper.CL t.ratio p.value
## -0.7168768 -6.925 <0.0001
## -1.0035435 -8.741 <0.0001
## 0.0897898 -1.816 0.1705
##
## Confidence level used: 0.95
## Conf-level adjustment: tukey method for comparing a family of 3 estimates
## P value adjustment: tukey method for comparing a family of 3 estimates
##
## >>> 9. LMM dengan semua pengukuran tersedia
## model n_observasi AIC BIC singular konvergensi
## 1 random_intercept 360 164.6389 219.0444 FALSE
## 2 random_intercept_slope 360 168.5489 230.7266 FALSE
## Linear mixed model fit by REML. t-tests use Satterthwaite's method [
## lmerModLmerTest]
## Formula: Hb ~ kelompok * waktu + (1 | id)
## Data: long_obs
## Control: lmerControl(optimizer = "bobyqa")
##
## REML criterion at convergence: 136.6
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -2.17750 -0.56714 0.00529 0.54743 2.94272
##
## Random effects:
## Groups Name Variance Std.Dev.
## id (Intercept) 0.34550 0.5878
## Residual 0.02742 0.1656
## Number of obs: 360, groups: id, 90
##
## Fixed effects:
## Estimate Std. Error df t value Pr(>|t|)
## (Intercept) 9.23806 0.06257 87.00002 147.643 < 2e-16 ***
## kelompok1 -0.43306 0.08849 87.00002 -4.894 4.50e-06 ***
## kelompok2 0.10361 0.08849 87.00002 1.171 0.2448
## waktu1 -0.47806 0.01512 260.99999 -31.628 < 2e-16 ***
## waktu2 -0.16472 0.01512 260.99999 -10.898 < 2e-16 ***
## waktu3 0.18639 0.01512 260.99999 12.331 < 2e-16 ***
## kelompok1:waktu1 0.40972 0.02138 260.99999 19.167 < 2e-16 ***
## kelompok2:waktu1 -0.10028 0.02138 260.99999 -4.691 4.39e-06 ***
## kelompok1:waktu2 0.13306 0.02138 260.99999 6.225 1.92e-09 ***
## kelompok2:waktu2 -0.04028 0.02138 260.99999 -1.884 0.0606 .
## kelompok1:waktu3 -0.15139 0.02138 260.99999 -7.082 1.32e-11 ***
## kelompok2:waktu3 -0.02472 0.02138 260.99999 -1.157 0.2485
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Correlation of Fixed Effects:
## (Intr) klmpk1 klmpk2 waktu1 waktu2 waktu3 klm1:1 klm2:1 klm1:2
## kelompok1 0.000
## kelompok2 0.000 -0.500
## waktu1 0.000 0.000 0.000
## waktu2 0.000 0.000 0.000 -0.333
## waktu3 0.000 0.000 0.000 -0.333 -0.333
## klmpk1:wkt1 0.000 0.000 0.000 0.000 0.000 0.000
## klmpk2:wkt1 0.000 0.000 0.000 0.000 0.000 0.000 -0.500
## klmpk1:wkt2 0.000 0.000 0.000 0.000 0.000 0.000 -0.333 0.167
## klmpk2:wkt2 0.000 0.000 0.000 0.000 0.000 0.000 0.167 -0.333 -0.500
## klmpk1:wkt3 0.000 0.000 0.000 0.000 0.000 0.000 -0.333 0.167 -0.333
## klmpk2:wkt3 0.000 0.000 0.000 0.000 0.000 0.000 0.167 -0.333 0.167
## klm2:2 klm1:3
## kelompok1
## kelompok2
## waktu1
## waktu2
## waktu3
## klmpk1:wkt1
## klmpk2:wkt1
## klmpk1:wkt2
## klmpk2:wkt2
## klmpk1:wkt3 0.167
## klmpk2:wkt3 -0.333 -0.500
## Groups Name Std.Dev.
## id (Intercept) 0.58779
## Residual 0.16558
## efek Sum Sq Mean Sq NumDF DenDF F value Pr(>F)
## 1 kelompok 0.7161595 0.3580797 2 87 13.06096 1.095979e-05
## 2 waktu 44.8831868 14.9610623 3 261 545.70493 4.267270e-112
## 3 kelompok:waktu 17.6337222 2.9389537 6 261 107.19837 1.685124e-67
## # Intraclass Correlation Coefficient
##
## Adjusted ICC: 0.926
## Unadjusted ICC: 0.532
##
## >>> 10. Diagnostik LMM
## Grafik tersimpan: D:/2026/S2 UNMUL/TUGAS/BIOSTAT/hasil_analisis_hb/10_qq_residual.png
## Grafik tersimpan: D:/2026/S2 UNMUL/TUGAS/BIOSTAT/hasil_analisis_hb/10_residual_prediksi.png
## Grafik tersimpan: D:/2026/S2 UNMUL/TUGAS/BIOSTAT/hasil_analisis_hb/10_residual_sel.png
## [1] FALSE
## [1] TRUE
## attr(,"gradient")
## [1] 3.025428e-07
##
## >>> 11. Informasi sesi dan batas interpretasi
## R version 4.6.1 (2026-06-24 ucrt)
## Platform: x86_64-w64-mingw32/x64
## Running under: Windows 11 x64 (build 26200)
##
## Matrix products: default
## LAPACK version 3.12.1
##
## locale:
## [1] LC_COLLATE=English_Indonesia.utf8 LC_CTYPE=English_Indonesia.utf8
## [3] LC_MONETARY=English_Indonesia.utf8 LC_NUMERIC=C
## [5] LC_TIME=English_Indonesia.utf8
##
## time zone: Asia/Makassar
## tzcode source: internal
##
## attached base packages:
## [1] stats graphics grDevices utils datasets methods base
##
## other attached packages:
## [1] base64enc_0.1-6 htmltools_0.5.9 knitr_1.51 ggpubr_1.0.0
## [5] performance_0.18.2 pbkrtest_0.5.5 lmerTest_3.2-1 effectsize_1.0.3
## [9] car_3.1-5 carData_3.0-6 rstatix_1.1.0 emmeans_2.0.4
## [13] afex_1.5-1 lme4_2.0-6 Matrix_1.7-5 ggplot2_4.0.3
## [17] tidyr_1.3.2 dplyr_1.2.1
##
## loaded via a namespace (and not attached):
## [1] gtable_0.3.6 xfun_0.60 bslib_0.12.0
## [4] bayestestR_0.19.0 insight_1.5.4 lattice_0.22-9
## [7] numDeriv_2016.8-1.1 vctrs_0.7.3 tools_4.6.1
## [10] Rdpack_2.6.6 generics_0.1.4 stats4_4.6.1
## [13] datawizard_1.4.0 parallel_4.6.1 tibble_3.3.1
## [16] pkgconfig_2.0.3 RColorBrewer_1.1-3 S7_0.2.2
## [19] lifecycle_1.0.5 compiler_4.6.1 farver_2.1.2
## [22] stringr_1.6.0 textshaping_1.0.5 sass_0.4.10
## [25] yaml_2.3.12 Formula_1.2-6 pillar_1.11.1
## [28] nloptr_2.2.1 jquerylib_0.1.4 MASS_7.3-65
## [31] cachem_1.1.0 reformulas_0.4.4 boot_1.3-32
## [34] abind_1.4-8 nlme_3.1-169 tidyselect_1.2.1
## [37] digest_0.6.39 mvtnorm_1.4-2 stringi_1.8.9
## [40] reshape2_1.4.5 purrr_1.2.2 labeling_0.4.3
## [43] splines_4.6.1 fastmap_1.2.0 grid_4.6.1
## [46] cli_3.6.6 magrittr_2.0.5 broom_1.0.13
## [49] withr_3.0.3 scales_1.4.0 backports_1.5.1
## [52] estimability_2.0.0 rmarkdown_2.31 otel_0.2.0
## [55] ggsignif_0.6.4 ragg_1.5.2 evaluate_1.0.5
## [58] parameters_0.29.3 rbibutils_2.4.1 rlang_1.3.0
## [61] Rcpp_1.1.2 glue_1.8.1 rstudioapi_0.19.0
## [64] minqa_1.2.8 jsonlite_2.0.0 R6_2.6.1
## [67] plyr_1.8.9 systemfonts_1.3.2
## bagian status
## 1 1. Sumber dan kelengkapan data Selesai
## 2 2. Statistik deskriptif dan visualisasi Selesai
## 3 3. Asumsi ANOVA pada peserta lengkap Selesai
## 4 4. RM ANOVA satu kelompok Selesai
## 5 5. Post hoc satu kelompok dan tren waktu Selesai
## 6 6. Pembanding nonparametrik satu kelompok Selesai
## 7 7. Mixed Design ANOVA kelompok × waktu Selesai
## 8 7b. Visualisasi interaksi kelompok dan waktu Selesai
## 9 8. Efek sederhana dan post hoc antarkelompok Selesai
## 10 9. LMM dengan semua pengukuran tersedia Selesai
## 11 10. Diagnostik LMM Selesai
## 12 11. Informasi sesi dan batas interpretasi Selesai
## peringatan
## 1
## 2 Computation failed in `stat_summary()`.\nCaused by error in `fun.data()`:\n! The package "Hmisc" is required.
## 3
## 4
## 5
## 6
## 7
## 8 Panel(s) show a mixed within-between-design.\nError bars do not allow comparisons across all means.\nSuppress error bars with: error = "none"
## 9 NaNs produced
## 10
## 11
## 12
##
## Laporan tersedia: D:/2026/S2 UNMUL/TUGAS/BIOSTAT/hasil_analisis_hb/laporan_analisis_hb.html
## Durasi: 15.5 detik.