ANALISIS DATA SKRIPSI
Ramzi
2026-10-03
#MEMANGGIL LIBRARY YANG DIGUNAKAN#
library(lmtest) ## Loading required package: zoo##
## Attaching package: 'zoo'## The following objects are masked from 'package:base':
##
## as.Date, as.Date.numericlibrary(MASS)
library(car)## Loading required package: carDatalibrary(pastecs)
library(pracma)##
## Attaching package: 'pracma'## The following object is masked from 'package:car':
##
## logitlibrary(Matrix)##
## Attaching package: 'Matrix'## The following objects are masked from 'package:pracma':
##
## expm, lu, tril, triulibrary(ggplot2)
library(gridExtra)#MEMANGGIL DATA#
data=read.table(file.choose(), header=TRUE)
data## Y X1 X2 X3 X4
## 1 4.52 67.88 13.00 14.86 72.04
## 2 4.97 71.02 12.89 14.38 51.19
## 3 5.70 61.98 13.58 18.56 73.59
## 4 5.45 68.96 12.78 29.61 43.00
## 5 5.08 67.40 13.31 14.89 74.71
## 6 6.22 69.04 12.55 55.70 71.41
## 7 3.49 76.22 12.50 10.46 67.09
## 8 9.00 62.90 14.13 15.59 80.01
## 9 8.26 65.16 14.70 75.04 30.11
## 10 9.46 69.24 12.90 31.21 80.02
## 11 3.71 74.28 12.60 20.67 67.03
## 12 3.91 75.81 12.08 8.67 47.87
## 13 3.38 71.78 12.39 10.74 65.98
## 14 7.55 64.14 12.33 8.54 25.74
## 15 3.52 70.38 11.56 19.95 65.77
## 16 7.30 60.75 11.79 28.13 67.17
## 17 4.50 75.57 12.02 14.71 36.88
## 18 4.02 74.09 12.04 10.28 65.69
## 19 3.39 77.53 11.57 6.56 64.76
## 20 2.70 73.93 11.15 25.23 25.55
## 21 3.71 65.53 11.81 4.21 22.68
## 22 7.14 67.71 13.64 28.93 67.95
## 23 12.36 60.05 13.64 37.68 79.44
## 24 8.78 63.84 12.89 10.14 41.94
## 25 3.57 72.03 11.96 10.16 69.38
## 26 4.96 64.68 11.92 17.31 68.86
## 27 3.87 72.55 12.28 11.73 59.18
## 28 2.93 74.61 12.38 5.76 66.22
## 29 3.73 70.17 11.86 6.35 70.11
## 30 2.24 73.15 12.10 4.65 48.85
## 31 3.90 71.15 12.19 4.83 68.84
## 32 4.49 70.08 12.88 3.30 65.59
## 33 3.07 69.27 12.59 14.48 72.19
## 34 6.95 70.16 12.36 15.39 50.71
## 35 2.46 76.50 12.37 9.17 68.82
## 36 8.32 62.07 13.92 21.93 77.10
## 37 5.54 66.82 14.80 6.11 79.10
## 38 5.08 66.44 13.29 11.21 71.94
## 39 4.45 67.38 12.99 18.72 41.10
## 40 4.83 67.81 12.20 5.75 66.97
## 41 4.14 66.91 12.63 26.30 45.79
## 42 5.86 65.65 13.73 38.11 75.83
## 43 4.76 73.01 12.71 33.79 72.87
## 44 5.25 67.41 12.69 56.63 71.31
## 45 4.98 70.04 12.90 45.81 69.48
## 46 4.21 64.68 12.54 45.53 70.22
## 47 5.29 71.54 12.48 51.50 30.59
## 48 4.70 65.50 12.11 66.27 68.03
## 49 2.83 70.50 12.47 67.99 70.51
## 50 4.30 65.04 11.98 40.96 37.58
## 51 5.69 64.55 12.51 48.22 68.68
## 52 2.63 72.77 12.40 43.84 68.45
## 53 2.49 71.22 11.77 51.25 70.81
## 54 2.91 77.73 12.82 54.56 21.39
## 55 3.10 63.93 11.74 62.92 67.98
## 56 5.95 62.71 14.94 61.01 80.77#ANALISIS STATISTIKA DESKRIPTIF#
stat.desc(data)## Y X1 X2 X3 X4
## nbr.val 56.0000000 5.600000e+01 56.00000000 56.0000000 56.0000000
## nbr.null 0.0000000 0.000000e+00 0.00000000 0.0000000 0.0000000
## nbr.na 0.0000000 0.000000e+00 0.00000000 0.0000000 0.0000000
## min 2.2400000 6.005000e+01 11.15000000 3.3000000 21.3900000
## max 12.3600000 7.773000e+01 14.94000000 75.0400000 80.7700000
## range 10.1200000 1.768000e+01 3.79000000 71.7400000 59.3800000
## sum 277.6000000 3.863250e+03 708.36000000 1476.2800000 3422.8700000
## median 4.5100000 6.914000e+01 12.50500000 18.6400000 67.9650000
## mean 4.9571429 6.898661e+01 12.64928571 26.3621429 61.1226786
## SE.mean 0.2715868 6.008087e-01 0.10790865 2.6750203 2.2021767
## CI.mean.0.95 0.5442721 1.204048e+00 0.21625376 5.3608606 4.4132607
## var 4.1305262 2.021438e+01 0.65207948 400.7210935 271.5766018
## std.dev 2.0323696 4.496041e+00 0.80751438 20.0180192 16.4795814
## coef.var 0.4099881 6.517266e-02 0.06383873 0.7593472 0.2696148model=(lm(formula=Y~X1+X2+X3+X4,data=data))
model##
## Call:
## lm(formula = Y ~ X1 + X2 + X3 + X4, data = data)
##
## Coefficients:
## (Intercept) X1 X2 X3 X4
## 10.498612 -0.237422 0.919839 -0.008344 -0.009455summary(model)##
## Call:
## lm(formula = Y ~ X1 + X2 + X3 + X4, data = data)
##
## Residuals:
## Min 1Q Median 3Q Max
## -2.1296 -0.9922 -0.1119 0.7639 4.6374
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 10.498612 5.796395 1.811 0.07600 .
## X1 -0.237422 0.049782 -4.769 1.59e-05 ***
## X2 0.919839 0.278410 3.304 0.00175 **
## X3 -0.008344 0.010186 -0.819 0.41651
## X4 -0.009455 0.012466 -0.758 0.45170
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 1.453 on 51 degrees of freedom
## Multiple R-squared: 0.526, Adjusted R-squared: 0.4888
## F-statistic: 14.15 on 4 and 51 DF, p-value: 7.786e-08vif(model)## X1 X2 X3 X4
## 1.304895 1.316581 1.083064 1.099406#SCATTERPLOT#
par(mfrow = c(2, 2))
for (j in 2:5) {
plot(
data[, j], data[, 1],
pch = 19,
xlab = names(data)[j],
ylab = names(data)[1],
main = paste(names(data)[1], "vs", names(data)[j])
)
}
par(mfrow = c(1, 1))###Matriks desain spline linear truncated
buat_mx <- function(X, knot, urut = "gcv")
{
X <- as.matrix(X)
knot <- as.matrix(knot)
N <- nrow(X)
m <- ncol(X)
k <- nrow(knot)
if (urut == "gcv")
{
# urutan kolom: b0, x1..xm, lalu per knot: (x1-ks)+ ... (xm-ks)+
trunc <- do.call(
cbind,
lapply(
1:k,
function(s)
sapply(
1:m,
function(j)
pmax(X[, j] - knot[s, j], 0)
)
)
)
mx <- cbind(1, X, trunc)
nama <- c(
"b0",
paste0("x", 1:m),
unlist(
lapply(
1:k,
function(s) paste0("(x", 1:m, "-k", s, ")+")
)
)
)
}
else
{
# urutan kolom: b0, lalu per variabel: xj, (xj-k1)+, (xj-k2)+, ...
kol <- list(rep(1, N))
nama <- "b0"
for (j in 1:m)
{
kol[[length(kol) + 1]] <- X[, j]
nama <- c(nama, paste0("x", j))
for (s in 1:k)
{
kol[[length(kol) + 1]] <- pmax(X[, j] - knot[s, j], 0)
nama <- c(nama, paste0("(x", j, "-k", s, ")+"))
}
}
mx <- do.call(cbind, kol)
}
colnames(mx) <- nama
mx
}GCV untuk k titik knot (k = 1, 2, 3, …)
Fungsi ini menggantikan GCV1 dan GCVk yang sebelumnya tidak ada. Kandidat knot tiap variabel adalah nk titik sama jarak antara min dan max (titik ujung dibuang). Kombinasi knot dipilih berdasarkan indeks yang sama untuk semua variabel, dan knot pada tiap variabel diurutkan (knot1 < knot2 < …).
GCVk <- function(data, k = 1, nk = 50, simpan = TRUE)
{
data <- as.matrix(data)
N <- nrow(data)
M <- ncol(data)
y <- data[, 1]
X <- data[, 2:M, drop = FALSE]
m <- ncol(X)
# kandidat knot tiap variabel
grid <- sapply(
1:m,
function(j) seq(min(X[, j]), max(X[, j]), length.out = nk)
)
grid <- grid[2:(nk - 1), , drop = FALSE]
nk1 <- nrow(grid)
# semua kombinasi indeks knot (k baris x nkomb kolom)
kombinasi <- combn(nk1, k)
nkomb <- ncol(kombinasi)
GCV <- rep(NA_real_, nkomb)
MSE <- rep(NA_real_, nkomb)
Rsq <- rep(NA_real_, nkomb)
for (i in 1:nkomb)
{
knot <- grid[kombinasi[, i], , drop = FALSE]
mx <- buat_mx(X, knot, "gcv")
C <- pinv(crossprod(mx))
B <- C %*% crossprod(mx, y)
yhat <- mx %*% B
SSE <- sum((y - yhat)^2)
SSR <- sum((yhat - mean(y))^2)
MSE[i] <- SSE / N
Rsq[i] <- SSR / (SSR + SSE) * 100
# trace matriks hat = trace(mx C mx')
trA <- sum((mx %*% C) * mx)
GCV[i] <- MSE[i] / ((N - trA) / N)^2
}
# tabel semua kombinasi knot
knotmat <- t(
apply(
kombinasi, 2,
function(idx) as.vector(t(grid[idx, , drop = FALSE]))
)
)
if (nrow(knotmat) != nkomb) knotmat <- t(knotmat)
colnames(knotmat) <- paste0(
"knot", rep(1:k, each = m),
"_x", rep(1:m, times = k)
)
dataAll <- cbind(
GCV = GCV,
Rsq = Rsq,
knot_ke = 1:nkomb,
knotmat
)
if (simpan)
{
write.csv(
dataAll,
file = paste0("dataAll_knot", k, "_v2.csv")
)
}
cat("\n=== ", k, " titik knot: 10 GCV terkecil ===\n", sep = "")
print(head(dataAll[order(GCV), -2], 10))
# knot optimal
best <- which.min(GCV)
knotopt <- grid[kombinasi[, best], , drop = FALSE]
dimnames(knotopt) <- list(paste0("knot_", 1:k), paste0("x", 1:m))
mxopt <- buat_mx(X, knotopt, "gcv")
Bopt <- pinv(crossprod(mxopt)) %*% crossprod(mxopt, y)
rownames(Bopt) <- colnames(mxopt)
colnames(Bopt) <- "Estimasi"
cat("\nKnot optimal ke-", best, " GCV = ", min(GCV), "\n", sep = "")
print(knotopt)
cat("\nEstimasi parameter\n")
print(Bopt)
invisible(
list(
knotopt = knotopt,
mingcv = min(GCV),
Rsq = Rsq[best],
B = Bopt,
dataAll = dataAll
)
)
}Uji signifikansi parameter
uji_k <- function(data, hasil, alpha = 0.1)
{
data <- as.matrix(data)
y <- data[, 1]
X <- data[, 2:ncol(data), drop = FALSE]
n <- nrow(data)
knot <- hasil$knotopt
k <- nrow(knot)
print(knot)
mx <- buat_mx(X, knot, "gcv")
p <- ncol(mx)
C <- pinv(crossprod(mx))
B <- C %*% crossprod(mx, y)
yhat <- mx %*% B
res <- y - yhat
SSE <- sum(res^2)
SST <- sum((y - mean(y))^2)
SSR <- SST - SSE
df_reg <- p - 1
df_error <- n - p
MSR <- SSR / df_reg
MSE <- SSE / df_error
Fhit <- MSR / MSE
pvalue <- pf(Fhit, df_reg, df_error, lower.tail = FALSE)
Rsq <- SSR / SST * 100
SE <- sqrt(diag(MSE * C))
thit <- B[, 1] / SE
pval <- 2 * pt(abs(thit), df_error, lower.tail = FALSE)
tab <- cbind(
Beta = B[, 1],
SE = SE,
t_hitung = thit,
p_value = pval
)
rownames(tab) <- colnames(mx)
print(tab)
cat("\nAnalysis of Variance\n")
cat("Regresi:", df_reg, SSR, MSR, Fhit, "\n")
cat("Error :", df_error, SSE, MSE, "\n")
cat("Total :", n - 1, SST, "\n")
cat("s =", sqrt(MSE), "\n")
cat("R-squared =", Rsq, "\n")
cat("p-value F =", pvalue, "\n")
cat("\nUji simultan:\n")
if (pvalue <= alpha) {
cat("Tolak H0: minimal terdapat satu parameter yang signifikan.\n")
} else {
cat("Gagal menolak H0.\n")
}
cat("\nUji parsial:\n")
for (i in 1:p)
{
if (pval[i] <= alpha) {
cat(colnames(mx)[i], ": signifikan, p-value =", pval[i], "\n")
} else {
cat(colnames(mx)[i], ": tidak signifikan, p-value =", pval[i], "\n")
}
}
write.csv(res, file = paste0("output_uji_residual_knot", k, "_v2.csv"))
write.csv(mx, file = paste0("output_uji_mx_knot", k, "_v2.csv"))
write.csv(yhat, file = paste0("output_uji_yhat_knot", k, "_v2.csv"))
invisible(
list(
B = B,
tab = tab,
Rsq = Rsq,
pvalue = pvalue,
residual = res,
yhat = yhat
)
)
}Pemilihan knot optimal dengan GCV
Catatan: 3 knot memerlukan waktu lebih lama karena jumlah kombinasinya jauh lebih banyak. Kurangi nk (mis. nk = 30) bila terlalu lambat.
hasil1 <- GCVk(data, k = 1)##
## === 1 titik knot: 10 GCV terkecil ===
## GCV knot_ke knot1_x1 knot1_x2 knot1_x3 knot1_x4
## [1,] 1.508187 45 76.28673 14.63061 69.183673 75.92265
## [2,] 1.515931 44 75.92592 14.55327 67.719592 74.71082
## [3,] 1.603414 46 76.64755 14.70796 70.647755 77.13449
## [4,] 1.615275 43 75.56510 14.47592 66.255510 73.49898
## [5,] 1.734461 42 75.20429 14.39857 64.791429 72.28714
## [6,] 1.763515 47 77.00837 14.78531 72.111837 78.34633
## [7,] 1.830095 3 61.13245 11.38204 7.692245 25.02551
## [8,] 1.861406 4 61.49327 11.45939 9.156327 26.23735
## [9,] 1.869770 41 74.84347 14.32122 63.327347 71.07531
## [10,] 1.912402 5 61.85408 11.53673 10.620408 27.44918
##
## Knot optimal ke-45 GCV = 1.508187
## x1 x2 x3 x4
## knot_1 76.28673 14.63061 69.18367 75.92265
##
## Estimasi parameter
## Estimasi
## b0 10.122393144
## x1 -0.200697922
## x2 0.762874047
## x3 -0.005244856
## x4 -0.017306261
## (x1-k1)+ 0.072850412
## (x2-k1)+ -21.843749933
## (x3-k1)+ 0.415138942
## (x4-k1)+ 1.053106627hasil2 <- GCVk(data, k = 2)##
## === 2 titik knot: 10 GCV terkecil ===
## GCV knot_ke knot1_x1 knot1_x2 knot1_x3 knot1_x4 knot2_x1 knot2_x2
## [1,] 1.316068 135 61.13245 11.38204 7.692245 25.02551 76.28673 14.63061
## [2,] 1.317894 179 61.49327 11.45939 9.156327 26.23735 76.28673 14.63061
## [3,] 1.321114 178 61.49327 11.45939 9.156327 26.23735 75.92592 14.55327
## [4,] 1.321402 134 61.13245 11.38204 7.692245 25.02551 75.92592 14.55327
## [5,] 1.358738 221 61.85408 11.53673 10.620408 27.44918 75.92592 14.55327
## [6,] 1.359630 222 61.85408 11.53673 10.620408 27.44918 76.28673 14.63061
## [7,] 1.359989 138 61.13245 11.38204 7.692245 25.02551 77.36918 14.86265
## [8,] 1.388799 136 61.13245 11.38204 7.692245 25.02551 76.64755 14.70796
## [9,] 1.393180 180 61.49327 11.45939 9.156327 26.23735 76.64755 14.70796
## [10,] 1.396087 177 61.49327 11.45939 9.156327 26.23735 75.56510 14.47592
## knot2_x3 knot2_x4
## [1,] 69.18367 75.92265
## [2,] 69.18367 75.92265
## [3,] 67.71959 74.71082
## [4,] 67.71959 74.71082
## [5,] 67.71959 74.71082
## [6,] 69.18367 75.92265
## [7,] 73.57592 79.55816
## [8,] 70.64776 77.13449
## [9,] 70.64776 77.13449
## [10,] 66.25551 73.49898
##
## Knot optimal ke-135 GCV = 1.316068
## x1 x2 x3 x4
## knot_1 61.13245 11.38204 7.692245 25.02551
## knot_2 76.28673 14.63061 69.183673 75.92265
##
## Estimasi parameter
## Estimasi
## b0 95.3587174
## x1 -2.9213119
## x2 6.1469387
## x3 0.2059982
## x4 0.7283636
## (x1-k1)+ 2.7449916
## (x2-k1)+ -5.3197766
## (x3-k1)+ -0.2144538
## (x4-k1)+ -0.7619877
## (x1-k2)+ 0.7954129
## (x2-k2)+ -18.0705567
## (x3-k2)+ 0.2853691
## (x4-k2)+ 0.9034632hasil3 <- GCVk(data, k = 3)##
## === 3 titik knot: 10 GCV terkecil ===
## GCV knot_ke knot1_x1 knot1_x2 knot1_x3 knot1_x4 knot2_x1 knot2_x2
## [1,] 1.444246 2200 61.13245 11.38204 7.692245 25.02551 61.85408 11.53673
## [2,] 1.445793 3146 61.49327 11.45939 9.156327 26.23735 61.85408 11.53673
## [3,] 1.447104 2157 61.13245 11.38204 7.692245 25.02551 61.49327 11.45939
## [4,] 1.452451 3106 61.13245 11.38204 7.692245 25.02551 77.00837 14.78531
## [5,] 1.453903 2199 61.13245 11.38204 7.692245 25.02551 61.85408 11.53673
## [6,] 1.454858 3145 61.49327 11.45939 9.156327 26.23735 61.85408 11.53673
## [7,] 1.458180 2156 61.13245 11.38204 7.692245 25.02551 61.49327 11.45939
## [8,] 1.458931 2242 61.13245 11.38204 7.692245 25.02551 62.21490 11.61408
## [9,] 1.464584 3037 61.13245 11.38204 7.692245 25.02551 73.03939 13.93449
## [10,] 1.465580 3101 61.13245 11.38204 7.692245 25.02551 76.28673 14.63061
## knot2_x3 knot2_x4 knot3_x1 knot3_x2 knot3_x3 knot3_x4
## [1,] 10.620408 27.44918 76.28673 14.63061 69.18367 75.92265
## [2,] 10.620408 27.44918 76.28673 14.63061 69.18367 75.92265
## [3,] 9.156327 26.23735 76.28673 14.63061 69.18367 75.92265
## [4,] 72.111837 78.34633 77.36918 14.86265 73.57592 79.55816
## [5,] 10.620408 27.44918 75.92592 14.55327 67.71959 74.71082
## [6,] 10.620408 27.44918 75.92592 14.55327 67.71959 74.71082
## [7,] 9.156327 26.23735 75.92592 14.55327 67.71959 74.71082
## [8,] 12.084490 28.66102 76.28673 14.63061 69.18367 75.92265
## [9,] 56.006939 65.01612 76.28673 14.63061 69.18367 75.92265
## [10,] 69.183673 75.92265 76.64755 14.70796 70.64776 77.13449
##
## Knot optimal ke-2200 GCV = 1.444246
## x1 x2 x3 x4
## knot_1 61.13245 11.38204 7.692245 25.02551
## knot_2 61.85408 11.53673 10.620408 27.44918
## knot_3 76.28673 14.63061 69.183673 75.92265
##
## Estimasi parameter
## Estimasi
## b0 -4.87954771
## x1 -0.57070096
## x2 1.29861828
## x3 0.20162854
## x4 1.23897093
## (x1-k1)+ -2.60712938
## (x2-k1)+ 13.80910091
## (x3-k1)+ -0.22960609
## (x4-k1)+ -1.94377923
## (x1-k2)+ 3.03208076
## (x2-k2)+ -14.12188191
## (x3-k2)+ 0.02361872
## (x4-k2)+ 0.67438515
## (x1-k3)+ 0.76491110
## (x2-k3)+ -19.59971920
## (x3-k3)+ 0.25980650
## (x4-k3)+ 0.92445503Perbandingan GCV
gcv_perbandingan <-
data.frame(
Model = c("1 Knot", "2 Knot", "3 Knot"),
GCV = c(hasil1$mingcv, hasil2$mingcv, hasil3$mingcv)
)
gcv_perbandingan## Model GCV
## 1 1 Knot 1.508187
## 2 2 Knot 1.316068
## 3 3 Knot 1.444246ggplot(
gcv_perbandingan,
aes(x = Model, y = GCV, group = 1)
) +
geom_line(linewidth = 1) +
geom_point(size = 3) +
geom_text(aes(label = round(GCV, 5)), vjust = -0.8) +
labs(
title = "Perbandingan Nilai GCV",
x = "Jumlah Titik Knot",
y = "GCV"
) +
theme_minimal()
gcv_minimum <- min(hasil1$mingcv, hasil2$mingcv, hasil3$mingcv)
model_terbaik <-
c("1 Knot", "2 Knot", "3 Knot")[
which.min(c(hasil1$mingcv, hasil2$mingcv, hasil3$mingcv))
]
cat("Model dengan GCV minimum:", model_terbaik, "\n")## Model dengan GCV minimum: 2 Knotcat("Nilai GCV:", gcv_minimum, "\n")## Nilai GCV: 1.316068knot_optimal <- hasil2$knotopt
knot_optimal## x1 x2 x3 x4
## knot_1 61.13245 11.38204 7.692245 25.02551
## knot_2 76.28673 14.63061 69.183673 75.92265B_optimal <- hasil2$B
B_optimal## Estimasi
## b0 95.3587174
## x1 -2.9213119
## x2 6.1469387
## x3 0.2059982
## x4 0.7283636
## (x1-k1)+ 2.7449916
## (x2-k1)+ -5.3197766
## (x3-k1)+ -0.2144538
## (x4-k1)+ -0.7619877
## (x1-k2)+ 0.7954129
## (x2-k2)+ -18.0705567
## (x3-k2)+ 0.2853691
## (x4-k2)+ 0.9034632Uji signifikansi parameter
uji2 <- uji_k(data, hasil2, alpha = 0.1)## x1 x2 x3 x4
## knot_1 61.13245 11.38204 7.692245 25.02551
## knot_2 76.28673 14.63061 69.183673 75.92265
## Beta SE t_hitung p_value
## b0 95.3587174 84.8285642 1.124135 0.2671924603
## x1 -2.9213119 1.0418249 -2.804034 0.0075441985
## x2 6.1469387 4.8596443 1.264895 0.2127214042
## x3 0.2059982 0.1626993 1.266128 0.2122837024
## x4 0.7283636 0.3657962 1.991173 0.0528396322
## (x1-k1)+ 2.7449916 1.0529268 2.607011 0.0125039146
## (x2-k1)+ -5.3197766 4.9095217 -1.083563 0.2845971079
## (x3-k1)+ -0.2144538 0.1669140 -1.284816 0.2057346472
## (x4-k1)+ -0.7619877 0.3717972 -2.049471 0.0465474966
## (x1-k2)+ 0.7954129 0.7655043 1.039070 0.3045796347
## (x2-k2)+ -18.0705567 4.2319010 -4.270080 0.0001058266
## (x3-k2)+ 0.2853691 0.2305324 1.237870 0.2224808174
## (x4-k2)+ 0.9034632 0.1935247 4.668464 0.0000297441
##
## Analysis of Variance
## Regresi: 12 183.7252 15.31043 15.15057
## Error : 43 43.45373 1.010552
## Total : 55 227.1789
## s = 1.005262
## R-squared = 80.87246
## p-value F = 9.568368e-12
##
## Uji simultan:
## Tolak H0: minimal terdapat satu parameter yang signifikan.
##
## Uji parsial:
## b0 : tidak signifikan, p-value = 0.2671925
## x1 : signifikan, p-value = 0.007544199
## x2 : tidak signifikan, p-value = 0.2127214
## x3 : tidak signifikan, p-value = 0.2122837
## x4 : signifikan, p-value = 0.05283963
## (x1-k1)+ : signifikan, p-value = 0.01250391
## (x2-k1)+ : tidak signifikan, p-value = 0.2845971
## (x3-k1)+ : tidak signifikan, p-value = 0.2057346
## (x4-k1)+ : signifikan, p-value = 0.0465475
## (x1-k2)+ : tidak signifikan, p-value = 0.3045796
## (x2-k2)+ : signifikan, p-value = 0.0001058266
## (x3-k2)+ : tidak signifikan, p-value = 0.2224808
## (x4-k2)+ : signifikan, p-value = 2.97441e-05uji3 <- uji_k(data, hasil3, alpha = 0.1)## x1 x2 x3 x4
## knot_1 61.13245 11.38204 7.692245 25.02551
## knot_2 61.85408 11.53673 10.620408 27.44918
## knot_3 76.28673 14.63061 69.183673 75.92265
## Beta SE t_hitung p_value
## b0 -4.87954771 11.3313795 -0.43062257 6.691134e-01
## x1 -0.57070096 2.3596928 -0.24185393 8.101606e-01
## x2 1.29861828 12.0225196 0.10801548 9.145372e-01
## x3 0.20162854 0.3213964 0.62735155 5.340848e-01
## x4 1.23897093 0.7333330 1.68950652 9.910266e-02
## (x1-k1)+ -2.60712938 4.8454015 -0.53806261 5.935914e-01
## (x2-k1)+ 13.80910091 33.0174969 0.41823585 6.780688e-01
## (x3-k1)+ -0.22960609 0.6222171 -0.36901280 7.141141e-01
## (x4-k1)+ -1.94377923 1.4399711 -1.34987380 1.848397e-01
## (x1-k2)+ 3.03208076 2.6352515 1.15058497 2.569091e-01
## (x2-k2)+ -14.12188191 22.0729360 -0.63978267 5.260547e-01
## (x3-k2)+ 0.02361872 0.3275433 0.07210868 9.428840e-01
## (x4-k2)+ 0.67438515 0.7687876 0.87720605 3.857466e-01
## (x1-k3)+ 0.76491110 0.8752551 0.87392932 3.875066e-01
## (x2-k3)+ -19.59971920 4.6428330 -4.22149996 1.404269e-04
## (x3-k3)+ 0.25980650 0.2406862 1.07944073 2.870221e-01
## (x4-k3)+ 0.92445503 0.2015338 4.58709688 4.564589e-05
##
## Analysis of Variance
## Regresi: 16 185.9148 11.61967 10.9821
## Error : 39 41.26416 1.058056
## Total : 55 227.1789
## s = 1.028618
## R-squared = 81.83627
## p-value F = 7.291165e-10
##
## Uji simultan:
## Tolak H0: minimal terdapat satu parameter yang signifikan.
##
## Uji parsial:
## b0 : tidak signifikan, p-value = 0.6691134
## x1 : tidak signifikan, p-value = 0.8101606
## x2 : tidak signifikan, p-value = 0.9145372
## x3 : tidak signifikan, p-value = 0.5340848
## x4 : signifikan, p-value = 0.09910266
## (x1-k1)+ : tidak signifikan, p-value = 0.5935914
## (x2-k1)+ : tidak signifikan, p-value = 0.6780688
## (x3-k1)+ : tidak signifikan, p-value = 0.7141141
## (x4-k1)+ : tidak signifikan, p-value = 0.1848397
## (x1-k2)+ : tidak signifikan, p-value = 0.2569091
## (x2-k2)+ : tidak signifikan, p-value = 0.5260547
## (x3-k2)+ : tidak signifikan, p-value = 0.942884
## (x4-k2)+ : tidak signifikan, p-value = 0.3857466
## (x1-k3)+ : tidak signifikan, p-value = 0.3875066
## (x2-k3)+ : signifikan, p-value = 0.0001404269
## (x3-k3)+ : tidak signifikan, p-value = 0.2870221
## (x4-k3)+ : signifikan, p-value = 4.564589e-05Plot nilai aktual vs prediksi (model 2 knot)
X <- as.matrix(data[, 2:5])
y <- data[, 1]
knot <- hasil2$knotopt
mx <- buat_mx(X, knot, "gcv")
yhat <- as.vector(mx %*% hasil2$B)
data_plot <-
data.frame(
Observasi = 1:nrow(data),
Aktual = y,
Prediksi = yhat
)
ggplot(data_plot, aes(x = Observasi)) +
geom_point(aes(y = Aktual), size = 2) +
geom_line(aes(y = Aktual, colour = "Aktual"), linewidth = 0.7) +
geom_line(aes(y = Prediksi, colour = "Prediksi"), linewidth = 1) +
scale_colour_manual(
name = NULL,
values = c(Aktual = "black", Prediksi = "red")
) +
labs(
title = "Nilai Aktual dan Prediksi Model Spline 2 Knot",
x = "Pengamatan",
y = "Nilai"
) +
theme_minimal()
Scatterplot dengan kurva loess
buat_sp <- function(j)
{
ggplot(
data,
aes(
x = .data[[names(data)[j]]],
y = .data[[names(data)[1]]]
)
) +
geom_point(size = 2) +
geom_smooth(method = "loess", se = FALSE) +
labs(x = names(data)[j], y = names(data)[1]) +
theme_minimal()
}
grid.arrange(
buat_sp(2),
buat_sp(3),
buat_sp(4),
buat_sp(5),
ncol = 2
)## `geom_smooth()` using formula = 'y ~ x'
## `geom_smooth()` using formula = 'y ~ x'
## `geom_smooth()` using formula = 'y ~ x'
## `geom_smooth()` using formula = 'y ~ x'
Ringkasan GCV
hasil_gcv <-
data.frame(
Jumlah_Knot = c(1, 2, 3),
GCV = c(hasil1$mingcv, hasil2$mingcv, hasil3$mingcv)
)
hasil_gcv## Jumlah_Knot GCV
## 1 1 1.508187
## 2 2 1.316068
## 3 3 1.444246hasil_gcv[which.min(hasil_gcv$GCV), ]## Jumlah_Knot GCV
## 2 2 1.316068