# ============================================================
# REGRESI SPLINE TRUNCATED MULTIVARIABEL 
# Adaptasi dari Skripsi Rewani Ega Tanaty
# Data: Tingkat Pengangguran Terbuka di Kabupaten/Kota Kalimantan
# ============================================================
# Perubahan dari versi awal:
#  1. Knot optimal diambil dengan which.min(GCV)
#  2. Parameter yang dicetak adalah parameter di knot OPTIMAL
#  3. Uji signifikansi memakai knot optimal (bukan knot ke-1)
#  4. MSE di uji signifikansi = SSE/(n - n1)
#  5. Nama file output diganti (akhiran _v2) supaya tidak menimpa file lama
# ============================================================

# ============================================================
# BAGIAN 1: STATISTIKA DESKRIPTIF & DETEKSI MULTIKOLINIERITAS
# ============================================================
library(lmtest)
## Loading required package: zoo
## 
## Attaching package: 'zoo'
## The following objects are masked from 'package:base':
## 
##     as.Date, as.Date.numeric
library(MASS)
library(car)
## Loading required package: carData
library(pastecs)
library(pracma)
## 
## Attaching package: 'pracma'
## The following object is masked from 'package:car':
## 
##     logit
library(ggplot2)
library(gridExtra)

data = read.table(file.choose(), header = TRUE)   # pilih file Data_Kabupaten_Kota_Kalimantan.txt
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
nrow(data)   # harus 56
## [1] 56
# Statistika deskriptif
stat.desc(data)   # ada SD, varians, dll
##                        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.2696148
library(ggplot2)
library(gridExtra)
sp1 = ggplot(data, aes(x=x1, y=y)) + geom_point() + stat_smooth()
sp2 = ggplot(data, aes(x=x2, y=y)) + geom_point() + stat_smooth()
sp3 = ggplot(data, aes(x=x3, y=y)) + geom_point() + stat_smooth()
sp4 = ggplot(data, aes(x=x4, y=y)) + geom_point() + stat_smooth()
grid.arrange(sp1, sp2, sp3, sp4, ncol=2)
## `geom_smooth()` using method = 'loess' and formula = 'y ~ x'
## `geom_smooth()` using method = 'loess' and formula = 'y ~ x'
## `geom_smooth()` using method = 'loess' and formula = 'y ~ x'
## `geom_smooth()` using method = 'loess' and formula = 'y ~ x'

# Deteksi multikolinieritas
model_lm = lm(formula = y ~ x1 + x2 + x3 + x4, data = data)
summary(model_lm)
## 
## 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-08
vif(model_lm)
##       x1       x2       x3       x4 
## 1.304895 1.316581 1.083064 1.099406
# ============================================================
# REGRESI SPLINE TRUNCATED MULTIVARIABEL - 1, 2, DAN 3 KNOT
# Data: Tingkat Pengangguran Terbuka di Kabupaten/Kota Kalimantan
# (y, x1, x2, x3, x4; n = 56)
# ============================================================
# Isi:
#   BAGIAN 1 : statistika deskriptif & deteksi multikolinieritas
#   BAGIAN A : persiapan (fungsi bantu)
#   BAGIAN B : fungsi pemilihan knot optimal (GCV) untuk k knot
#   BAGIAN C : fungsi uji signifikansi (simultan & parsial) untuk k knot
#   BAGIAN D : jalankan 1, 2, 3 knot sekaligus + tabel perbandingan
#
# Skema knot (sama untuk 1, 2, 3 knot):
#   - tiap variabel punya grid nk = 50 titik dari min sampai max,
#     titik min & max dibuang (jadi 48 kandidat)
#   - 1 knot : 48 kandidat                 -> C(48,1) = 48 kombinasi
#   - 2 knot : semua pasangan (i<j)        -> C(48,2) = 1128 kombinasi
#   - 3 knot : semua tripel (i<j<l)        -> C(48,3) = 17296 kombinasi
# ============================================================

# ============================================================
# BAGIAN A: FUNGSI BANTU
# ============================================================
# Membentuk matriks desain
#   knot : matriks k x m (baris = knot ke-s, kolom = variabel ke-j)
#   urut = "gcv" -> 1, x1..xm, lalu semua (x-k1)+, semua (x-k2)+, ...
#   urut = "uji" -> 1, x1, (x1-k1)+, (x1-k2)+, x2, (x2-k1)+, ...
buat_mx = function(X, knot, urut = "gcv")
{
  N = nrow(X); m = ncol(X); k = nrow(knot)
  if (urut == "gcv") {
    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 {
    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
}

# ============================================================
# BAGIAN B: PEMILIHAN KNOT OPTIMAL DENGAN GCV (k = 1, 2, atau 3)
# ============================================================
GCVk = function(data, k = 2, nk = 50)
{
  data = as.matrix(data)
  N = nrow(data)
  y = data[, 1]
  X = data[, 2:ncol(data), drop = FALSE]
  m = ncol(X)
  
  # grid kandidat knot per variabel (buang titik min & max)
  grid = sapply(1:m, function(j) seq(min(X[, j]), max(X[, j]), length.out = nk))
  grid = grid[2:(nk - 1), , drop = FALSE]
  ng = nrow(grid)
  
  # tiap baris = indeks kandidat (naik); matrix() menjaga bentuk saat k = 1
  komb = matrix(t(combn(ng, k)), ncol = k)
  nkomb = nrow(komb)
  cat("\n>>> Jumlah kombinasi knot (k =", k, "):", nkomb, "\n")
  
  GCV = rep(NA, nkomb)
  Rsq = rep(NA, nkomb)
  knotmat = matrix(NA, nkomb, k * m)
  nama.knot = as.vector(t(outer(1:k, 1:m, function(s, j) paste0("k", s, "_x", j))))
  
  for (i in 1:nkomb)
  {
    knot = matrix(0, k, m)
    for (s in 1:k) knot[s, ] = grid[komb[i, s], ]
    knotmat[i, ] = as.vector(t(knot))
    
    mx = buat_mx(X, knot, "gcv")
    C = pinv(t(mx) %*% mx)
    B = C %*% (t(mx) %*% y)
    yhat = mx %*% B
    SSE = sum((y - yhat)^2)
    SSR = sum((yhat - mean(y))^2)
    trA = sum((mx %*% C) * mx)               # trace matriks hat
    GCV[i] = (SSE / N) / (((N - trA) / N)^2)
    Rsq[i] = SSR / (SSR + SSE) * 100
    if (i %% 2000 == 0) cat("  selesai", i, "dari", nkomb, "\n")
  }
  
  dataAll = cbind(GCV = GCV, Rsq = Rsq, komb_ke = 1:nkomb, knotmat)
  colnames(dataAll)[4:ncol(dataAll)] = nama.knot
  file.out = paste0("dataAll_knot", k, "_v2.csv")
  write.csv(dataAll, file = file.out)
  
  dataG = dataAll[order(GCV), -2, drop = FALSE]   # GCV, komb_ke, knot...
  cat("\n==============================================\n")
  cat("HASIL GCV terkecil dengan", k, "knot\n")
  cat("==============================================\n")
  print(dataG[1, ])
  cat("\nNilai GCV 10 terkecil pertama\n")
  print(dataG[1:min(10, nkomb), ])
  
  # ---- estimasi parameter pada knot OPTIMAL ----
  best = which.min(GCV)
  knotopt = matrix(knotmat[best, ], nrow = k, byrow = TRUE,
                   dimnames = list(paste0("knot_", 1:k), paste0("x", 1:m)))
  mxopt = buat_mx(X, knotopt, "gcv")
  Bopt = pinv(t(mxopt) %*% mxopt) %*% t(mxopt) %*% y
  rownames(Bopt) = colnames(mxopt)
  
  cat("\nKombinasi knot optimal ke-", best, "  GCV =", min(GCV), "\n")
  print(knotopt)
  cat("\n==============================================\n")
  cat("HASIL ESTIMASI PARAMETER", k, "TITIK KNOT (OPTIMAL)\n")
  cat("==============================================\n")
  print(Bopt)
  cat("\nFile knot tersimpan di:", file.out, "\n")
  
  invisible(list(k = k, knotopt = knotopt, mingcv = min(GCV), B = Bopt, dataAll = dataAll))
}

# ============================================================
# BAGIAN C: UJI SIGNIFIKANSI SIMULTAN & PARSIAL (k = 1, 2, atau 3)
# ============================================================
uji_k = function(data, k = 2, alpha = 0.1)
{
  data = as.matrix(data)
  y = data[, 1]
  X = data[, 2:ncol(data), drop = FALSE]
  n = nrow(data); m = ncol(X)
  
  # ambil knot dengan GCV minimum dari file hasil GCVk
  file.knot = paste0("dataAll_knot", k, "_v2.csv")
  knot_all = read.csv(file.knot, header = TRUE, row.names = 1)
  best = which.min(knot_all$GCV)
  knot = matrix(as.numeric(knot_all[best, 4:ncol(knot_all)]), nrow = k, byrow = TRUE,
                dimnames = list(paste0("knot_", 1:k), paste0("x", 1:m)))
  cat("\n>>> UJI SIGNIFIKANSI,", k, "KNOT\n")
  cat("Kombinasi knot optimal ke-", best, "  GCV =", knot_all$GCV[best], "\n")
  print(knot)
  cat("\n")
  
  mx = buat_mx(X, knot, "uji")
  B = pinv(t(mx) %*% mx) %*% t(mx) %*% y
  p = nrow(B)
  yhat = mx %*% B
  ybar = mean(y)
  res = y - yhat
  SSE = sum(res^2)
  SSR = sum((yhat - ybar)^2)
  SST = sum((y - ybar)^2)
  MSE = SSE / (n - p)
  MSR = SSR / (p - 1)
  Rsq = SSR / (SSR + SSE) * 100
  
  # --- Uji simultan ---
  Fhit = MSR / MSE
  pvalue = pf(Fhit, p - 1, n - p, lower.tail = FALSE)
  cat('---------------------------------------\n')
  cat('Kesimpulan hasil uji simultan\n')
  cat('---------------------------------------\n')
  if (pvalue <= alpha) {
    cat('Tolak Ho yakni minimal terdapat 1 variabel bebas yang signifikan\n\n')
  } else {
    cat('Gagal Tolak Ho yakni semua variabel bebas tidak berpengaruh signifikan\n\n')
  }
  
  # --- Uji parsial ---
  SE = sqrt(diag(MSE * pinv(t(mx) %*% mx)))
  thit = B[, 1] / SE
  pval = 2 * pt(abs(thit), n - p, lower.tail = FALSE)
  cat('---------------------------------------------\n')
  cat('Kesimpulan hasil uji parsial\n')
  cat('---------------------------------------------\n')
  for (i in 1:p)
  {
    if (pval[i] <= alpha)
      cat('Parameter', i - 1, colnames(mx)[i], ': Tolak Ho, signifikan, pvalue =', pval[i], '\n') else
        cat('Parameter', i - 1, colnames(mx)[i], ': Gagal tolak Ho, tidak signifikan, pvalue =', pval[i], '\n')
  }
  
  cat('=============================================\n')
  cat('Estimasi parameter, t hitung, p-value\n')
  cat('=============================================\n')
  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('=============================================\n')
  cat('Sumber   df   SS   MS   Fhit\n')
  cat('Regresi ', p - 1, ' ', SSR, ' ', MSR, ' ', Fhit, '\n')
  cat('Error   ', n - p, ' ', SSE, ' ', MSE, '\n')
  cat('Total   ', n - 1, ' ', SST, '\n')
  cat('=============================================\n')
  cat('s =', sqrt(MSE), '  Rsq =', Rsq, '\n')
  cat('pvalue(F) =', pvalue, '\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(k = k, B = B, tab = tab, Rsq = Rsq, MSE = MSE, pvalueF = pvalue, p = p))
}

# ============================================================
# BAGIAN D: JALANKAN 1, 2, 3 KNOT SEKALIGUS + PERBANDINGAN
# ============================================================
K_set = 1:3
alpha = 0.1
m = ncol(data) - 1

hasil = list()
uji   = list()
for (k in K_set)
{
  cat("\n\n############################################################\n")
  cat("#                    ", k, "TITIK KNOT\n")
  cat("############################################################\n")
  hasil[[k]] = GCVk(data, k = k)                 # knot optimal + estimasi parameter
  uji[[k]]   = uji_k(data, k = k, alpha = alpha) # uji simultan & parsial
}
## 
## 
## ############################################################
## #                     1 TITIK KNOT
## ############################################################
## 
## >>> Jumlah kombinasi knot (k = 1 ): 48 
## 
## ==============================================
## HASIL GCV terkecil dengan 1 knot
## ==============================================
##       GCV   komb_ke     k1_x1     k1_x2     k1_x3     k1_x4 
##  1.508187 45.000000 76.286735 14.630612 69.183673 75.922653 
## 
## Nilai GCV 10 terkecil pertama
##            GCV komb_ke    k1_x1    k1_x2     k1_x3    k1_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
## 
## Kombinasi knot optimal ke- 45   GCV = 1.508187 
##              x1       x2       x3       x4
## knot_1 76.28673 14.63061 69.18367 75.92265
## 
## ==============================================
## HASIL ESTIMASI PARAMETER 1 TITIK KNOT (OPTIMAL)
## ==============================================
##                   [,1]
## 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.053106627
## 
## File knot tersimpan di: dataAll_knot1_v2.csv 
## 
## >>> UJI SIGNIFIKANSI, 1 KNOT
## Kombinasi knot optimal ke- 45   GCV = 1.508187 
##              x1       x2       x3       x4
## knot_1 76.28673 14.63061 69.18367 75.92265
## 
## ---------------------------------------
## Kesimpulan hasil uji simultan
## ---------------------------------------
## Tolak Ho yakni minimal terdapat 1 variabel bebas yang signifikan
## 
## ---------------------------------------------
## Kesimpulan hasil uji parsial
## ---------------------------------------------
## Parameter 0 b0 : Tolak Ho, signifikan, pvalue = 0.06500353 
## Parameter 1 x1 : Tolak Ho, signifikan, pvalue = 2.835411e-05 
## Parameter 2 (x1-k1)+ : Gagal tolak Ho, tidak signifikan, pvalue = 0.9156743 
## Parameter 3 x2 : Tolak Ho, signifikan, pvalue = 0.01250812 
## Parameter 4 (x2-k1)+ : Tolak Ho, signifikan, pvalue = 1.020449e-05 
## Parameter 5 x3 : Gagal tolak Ho, tidak signifikan, pvalue = 0.541134 
## Parameter 6 (x3-k1)+ : Tolak Ho, signifikan, pvalue = 0.09644074 
## Parameter 7 x4 : Gagal tolak Ho, tidak signifikan, pvalue = 0.1203336 
## Parameter 8 (x4-k1)+ : Tolak Ho, signifikan, pvalue = 2.335269e-06 
## =============================================
## Estimasi parameter, t hitung, p-value
## =============================================
##                   Beta         SE   t_hitung      p_value
## b0        10.122393144 5.35720770  1.8894905 6.500353e-02
## x1        -0.200697922 0.04328026 -4.6371698 2.835411e-05
## (x1-k1)+   0.072850412 0.68432960  0.1064552 9.156743e-01
## x2         0.762874047 0.29372705  2.5972210 1.250812e-02
## (x2-k1)+ -21.843749933 4.41930754 -4.9427992 1.020449e-05
## x3        -0.005244856 0.00852003 -0.6155913 5.411340e-01
## (x3-k1)+   0.415138942 0.24472964  1.6963165 9.644074e-02
## x4        -0.017306261 0.01093879 -1.5820999 1.203336e-01
## (x4-k1)+   1.053106627 0.19589216  5.3759508 2.335269e-06
## 
## Analysis of Variance
## =============================================
## Sumber   df   SS   MS   Fhit
## Regresi  8   167.6864   20.9608   16.55933 
## Error    47   59.49258   1.2658 
## Total    55   227.1789 
## =============================================
## s = 1.125078   Rsq = 73.81246 
## pvalue(F) = 2.451288e-11 
## 
## 
## ############################################################
## #                     2 TITIK KNOT
## ############################################################
## 
## >>> Jumlah kombinasi knot (k = 2 ): 1128 
## 
## ==============================================
## HASIL GCV terkecil dengan 2 knot
## ==============================================
##        GCV    komb_ke      k1_x1      k1_x2      k1_x3      k1_x4      k2_x1 
##   1.316068 135.000000  61.132449  11.382041   7.692245  25.025510  76.286735 
##      k2_x2      k2_x3      k2_x4 
##  14.630612  69.183673  75.922653 
## 
## Nilai GCV 10 terkecil pertama
##            GCV komb_ke    k1_x1    k1_x2     k1_x3    k1_x4    k2_x1    k2_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
##          k2_x3    k2_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
## 
## Kombinasi 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
## 
## ==============================================
## HASIL ESTIMASI PARAMETER 2 TITIK KNOT (OPTIMAL)
## ==============================================
##                 [,1]
## b0        95.3587177
## 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.9034632
## 
## File knot tersimpan di: dataAll_knot2_v2.csv 
## 
## >>> UJI SIGNIFIKANSI, 2 KNOT
## Kombinasi 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
## 
## ---------------------------------------
## Kesimpulan hasil uji simultan
## ---------------------------------------
## Tolak Ho yakni minimal terdapat 1 variabel bebas yang signifikan
## 
## ---------------------------------------------
## Kesimpulan hasil uji parsial
## ---------------------------------------------
## Parameter 0 b0 : Gagal tolak Ho, tidak signifikan, pvalue = 0.2671925 
## Parameter 1 x1 : Tolak Ho, signifikan, pvalue = 0.007544198 
## Parameter 2 (x1-k1)+ : Tolak Ho, signifikan, pvalue = 0.01250391 
## Parameter 3 (x1-k2)+ : Gagal tolak Ho, tidak signifikan, pvalue = 0.3045796 
## Parameter 4 x2 : Gagal tolak Ho, tidak signifikan, pvalue = 0.2127214 
## Parameter 5 (x2-k1)+ : Gagal tolak Ho, tidak signifikan, pvalue = 0.2845971 
## Parameter 6 (x2-k2)+ : Tolak Ho, signifikan, pvalue = 0.0001058266 
## Parameter 7 x3 : Gagal tolak Ho, tidak signifikan, pvalue = 0.2122837 
## Parameter 8 (x3-k1)+ : Gagal tolak Ho, tidak signifikan, pvalue = 0.2057346 
## Parameter 9 (x3-k2)+ : Gagal tolak Ho, tidak signifikan, pvalue = 0.2224808 
## Parameter 10 x4 : Tolak Ho, signifikan, pvalue = 0.05283963 
## Parameter 11 (x4-k1)+ : Tolak Ho, signifikan, pvalue = 0.0465475 
## Parameter 12 (x4-k2)+ : Tolak Ho, signifikan, pvalue = 2.97441e-05 
## =============================================
## Estimasi parameter, t hitung, p-value
## =============================================
##                 Beta         SE  t_hitung      p_value
## b0        95.3587180 84.8285642  1.124135 0.2671924572
## x1        -2.9213119  1.0418249 -2.804034 0.0075441984
## (x1-k1)+   2.7449916  1.0529268  2.607011 0.0125039144
## (x1-k2)+   0.7954129  0.7655043  1.039070 0.3045796349
## x2         6.1469387  4.8596443  1.264895 0.2127214060
## (x2-k1)+  -5.3197766  4.9095217 -1.083563 0.2845971101
## (x2-k2)+ -18.0705567  4.2319010 -4.270080 0.0001058266
## x3         0.2059982  0.1626993  1.266128 0.2122837024
## (x3-k1)+  -0.2144538  0.1669140 -1.284816 0.2057346472
## (x3-k2)+   0.2853691  0.2305324  1.237870 0.2224808174
## x4         0.7283636  0.3657962  1.991173 0.0528396322
## (x4-k1)+  -0.7619877  0.3717972 -2.049471 0.0465474966
## (x4-k2)+   0.9034632  0.1935247  4.668464 0.0000297441
## 
## Analysis of Variance
## =============================================
## Sumber   df   SS   MS   Fhit
## Regresi  12   183.7252   15.31043   15.15057 
## Error    43   43.45373   1.010552 
## Total    55   227.1789 
## =============================================
## s = 1.005262   Rsq = 80.87246 
## pvalue(F) = 9.568368e-12 
## 
## 
## ############################################################
## #                     3 TITIK KNOT
## ############################################################
## 
## >>> Jumlah kombinasi knot (k = 3 ): 17296 
##   selesai 2000 dari 17296 
##   selesai 4000 dari 17296 
##   selesai 6000 dari 17296 
##   selesai 8000 dari 17296 
##   selesai 10000 dari 17296 
##   selesai 12000 dari 17296 
##   selesai 14000 dari 17296 
##   selesai 16000 dari 17296 
## 
## ==============================================
## HASIL GCV terkecil dengan 3 knot
## ==============================================
##         GCV     komb_ke       k1_x1       k1_x2       k1_x3       k1_x4 
##    1.444246 2200.000000   61.132449   11.382041    7.692245   25.025510 
##       k2_x1       k2_x2       k2_x3       k2_x4       k3_x1       k3_x2 
##   61.854082   11.536735   10.620408   27.449184   76.286735   14.630612 
##       k3_x3       k3_x4 
##   69.183673   75.922653 
## 
## Nilai GCV 10 terkecil pertama
##            GCV komb_ke    k1_x1    k1_x2    k1_x3    k1_x4    k2_x1    k2_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
##           k2_x3    k2_x4    k3_x1    k3_x2    k3_x3    k3_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
## 
## Kombinasi 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
## 
## ==============================================
## HASIL ESTIMASI PARAMETER 3 TITIK KNOT (OPTIMAL)
## ==============================================
##                  [,1]
## b0        -4.87954771
## x1        -0.57070096
## x2         1.29861829
## x3         0.20162854
## x4         1.23897093
## (x1-k1)+  -2.60712938
## (x2-k1)+  13.80910089
## (x3-k1)+  -0.22960609
## (x4-k1)+  -1.94377922
## (x1-k2)+   3.03208075
## (x2-k2)+ -14.12188190
## (x3-k2)+   0.02361872
## (x4-k2)+   0.67438515
## (x1-k3)+   0.76491110
## (x2-k3)+ -19.59971920
## (x3-k3)+   0.25980650
## (x4-k3)+   0.92445503
## 
## File knot tersimpan di: dataAll_knot3_v2.csv 
## 
## >>> UJI SIGNIFIKANSI, 3 KNOT
## Kombinasi 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
## 
## ---------------------------------------
## Kesimpulan hasil uji simultan
## ---------------------------------------
## Tolak Ho yakni minimal terdapat 1 variabel bebas yang signifikan
## 
## ---------------------------------------------
## Kesimpulan hasil uji parsial
## ---------------------------------------------
## Parameter 0 b0 : Gagal tolak Ho, tidak signifikan, pvalue = 0.6691134 
## Parameter 1 x1 : Gagal tolak Ho, tidak signifikan, pvalue = 0.8101606 
## Parameter 2 (x1-k1)+ : Gagal tolak Ho, tidak signifikan, pvalue = 0.5935914 
## Parameter 3 (x1-k2)+ : Gagal tolak Ho, tidak signifikan, pvalue = 0.2569091 
## Parameter 4 (x1-k3)+ : Gagal tolak Ho, tidak signifikan, pvalue = 0.3875066 
## Parameter 5 x2 : Gagal tolak Ho, tidak signifikan, pvalue = 0.9145372 
## Parameter 6 (x2-k1)+ : Gagal tolak Ho, tidak signifikan, pvalue = 0.6780688 
## Parameter 7 (x2-k2)+ : Gagal tolak Ho, tidak signifikan, pvalue = 0.5260547 
## Parameter 8 (x2-k3)+ : Tolak Ho, signifikan, pvalue = 0.0001404269 
## Parameter 9 x3 : Gagal tolak Ho, tidak signifikan, pvalue = 0.5340848 
## Parameter 10 (x3-k1)+ : Gagal tolak Ho, tidak signifikan, pvalue = 0.7141141 
## Parameter 11 (x3-k2)+ : Gagal tolak Ho, tidak signifikan, pvalue = 0.942884 
## Parameter 12 (x3-k3)+ : Gagal tolak Ho, tidak signifikan, pvalue = 0.2870221 
## Parameter 13 x4 : Tolak Ho, signifikan, pvalue = 0.09910266 
## Parameter 14 (x4-k1)+ : Gagal tolak Ho, tidak signifikan, pvalue = 0.1848397 
## Parameter 15 (x4-k2)+ : Gagal tolak Ho, tidak signifikan, pvalue = 0.3857466 
## Parameter 16 (x4-k3)+ : Tolak Ho, signifikan, pvalue = 4.564589e-05 
## =============================================
## Estimasi parameter, t hitung, p-value
## =============================================
##                  Beta         SE    t_hitung      p_value
## b0        -4.87954767 11.3313796 -0.43062256 6.691134e-01
## x1        -0.57070097  2.3596928 -0.24185393 8.101606e-01
## (x1-k1)+  -2.60712937  4.8454015 -0.53806261 5.935914e-01
## (x1-k2)+   3.03208075  2.6352515  1.15058497 2.569091e-01
## (x1-k3)+   0.76491110  0.8752551  0.87392932 3.875066e-01
## x2         1.29861830 12.0225195  0.10801549 9.145372e-01
## (x2-k1)+  13.80910084 33.0174968  0.41823585 6.780688e-01
## (x2-k2)+ -14.12188187 22.0729359 -0.63978267 5.260547e-01
## (x2-k3)+ -19.59971920  4.6428330 -4.22149996 1.404269e-04
## x3         0.20162854  0.3213964  0.62735156 5.340848e-01
## (x3-k1)+  -0.22960609  0.6222171 -0.36901280 7.141141e-01
## (x3-k2)+   0.02361872  0.3275433  0.07210868 9.428840e-01
## (x3-k3)+   0.25980650  0.2406862  1.07944073 2.870221e-01
## x4         1.23897093  0.7333330  1.68950652 9.910266e-02
## (x4-k1)+  -1.94377922  1.4399711 -1.34987380 1.848397e-01
## (x4-k2)+   0.67438515  0.7687876  0.87720604 3.857466e-01
## (x4-k3)+   0.92445502  0.2015338  4.58709688 4.564589e-05
## 
## Analysis of Variance
## =============================================
## Sumber   df   SS   MS   Fhit
## Regresi  16   185.9148   11.61967   10.9821 
## Error    39   41.26416   1.058056 
## Total    55   227.1789 
## =============================================
## s = 1.028618   Rsq = 81.83627 
## pvalue(F) = 7.291165e-10
# ---- Tabel perbandingan 1, 2, 3 knot ----
perbandingan = data.frame(
  Knot       = K_set,
  Kombinasi  = sapply(K_set, function(k) choose(48, k)),
  Parameter  = sapply(K_set, function(k) uji[[k]]$p),
  GCV_min    = sapply(K_set, function(k) hasil[[k]]$mingcv),
  MSE        = sapply(K_set, function(k) uji[[k]]$MSE),
  Rsq        = sapply(K_set, function(k) uji[[k]]$Rsq),
  pvalue_F   = sapply(K_set, function(k) uji[[k]]$pvalueF)
)
cat("\n\n==============================================\n")
## 
## 
## ==============================================
cat("PERBANDINGAN MODEL 1, 2, 3 TITIK KNOT\n")
## PERBANDINGAN MODEL 1, 2, 3 TITIK KNOT
cat("==============================================\n")
## ==============================================
print(perbandingan)
##   Knot Kombinasi Parameter  GCV_min      MSE      Rsq     pvalue_F
## 1    1        48         9 1.508187 1.265800 73.81246 2.451288e-11
## 2    2      1128        13 1.316068 1.010552 80.87246 9.568368e-12
## 3    3     17296        17 1.444246 1.058056 81.83627 7.291165e-10
cat("\nModel terbaik (GCV minimum): knot =",
    perbandingan$Knot[which.min(perbandingan$GCV_min)], "\n")
## 
## Model terbaik (GCV minimum): knot = 2
write.csv(perbandingan, file = "perbandingan_knot_1_2_3.csv", row.names = FALSE)

R Markdown

This is an R Markdown document. Markdown is a simple formatting syntax for authoring HTML, PDF, and MS Word documents. For more details on using R Markdown see http://rmarkdown.rstudio.com.

When you click the Knit button a document will be generated that includes both content as well as the output of any embedded R code chunks within the document. You can embed an R code chunk like this:

summary(cars)
##      speed           dist       
##  Min.   : 4.0   Min.   :  2.00  
##  1st Qu.:12.0   1st Qu.: 26.00  
##  Median :15.0   Median : 36.00  
##  Mean   :15.4   Mean   : 42.98  
##  3rd Qu.:19.0   3rd Qu.: 56.00  
##  Max.   :25.0   Max.   :120.00

Including Plots

You can also embed plots, for example:

Note that the echo = FALSE parameter was added to the code chunk to prevent printing of the R code that generated the plot.