# Fungsi interpolasi linier 
lin2 <- function(x0, y0, x1, y1, x) {
  if (x0 == x1) stop("x0 dan x1 tidak boleh sama (pembagi nol).")
  a1 <- (y1 - y0) / (x1 - x0)
  a0 <- (x1*y0 - x0*y1) / (x1 - x0)
  y  <- a0 + a1 * x 
  list(y = y, a0 = a0, a1 = a1)
}

# Berikut dua titik yang diketahui (x0, y0) dan (x1, y1)
x0 <- 1960 
y0 <- 179.3
x1 <- 1970 
y1 <- 203.2

# Input nilai xi untuk mencari hampiran yi
x  <- c(1967, 1968, 1969)             

# Kalkulasi
hasil <- lin2(x0, y0, x1, y1, x)

# Tampilkan fungsi interpolasi linier
cat("Persamaan interpolasi: y =", round(hasil$a0, 2), "+", round(hasil$a1, 2), "* x\n\n")

# Tampilkan hasil hampiran
data.frame(x = x, y_hampiran = hasil$y)

# Visualisasi
plot(c(x0, x1), c(y0, y1), type="p", pch=19, col="black",
     xlab="Tahun", ylab="Nilai", main="Interpolasi Linier")

# Tambahkan garis interpolasi
abline(hasil$a0, hasil$a1, col="black", lwd=2)

# Tambahkan titik interpolasi
points(x, hasil$y, pch=19, col="red")

# Tambahkan legenda
legend("topleft", legend=c("Titik asli", "Garis interpolasi", "Titik hasil interpolasi"),
       col=c("black", "black", "red"), pch=c(19, NA, 17), lty=c(NA,1,NA))

# Fungsi umum untuk banyak titik
poly_interp <- function(x, y) {
  stopifnot(length(x) == length(y), length(x) >= 2)
  if (any(duplicated(x))) stop("Semua x harus berbeda.")
  n <- length(x)
  V <- outer(x, 0:(n-1), `^`)       
  b <- solve(V, y)                 
  pred <- function(xnew) {
    X <- outer(xnew, 0:(n-1), `^`)
    as.vector(X %*% b)
  }
  list(coef = setNames(b, paste0("a", 0:(n-1))), degree = n-1, predict = pred)
}

# Berikut tiga titik yang diketahui (x0, x1, x2) dan (y0, y1, y2) 
x3 <- c(8,9,9.5)
y3 <- c(2.0794,2.1972,2.2513)   

# Membuat fungsi interpolasi dari tiga titik
poli3 <- poly_interp(x3, y3)
poli3$coef

# Tentukan nilai x yang ingin dihampiri y nya
x<-c(9.2,9.3,9.35)

# Tampilkan hasil hampiran
hampiran <- data.frame(x = x, y_hampiran = poli3$predict(x))
print(hampiran)

# Buat string persamaan polinom
make_poly_equation <- function(coef) {
  terms <- c()
  for (k in seq_along(coef)) {
    cval <- round(coef[k], 6)  # dibulatkan biar rapi
    if (k == 1) {
      terms <- c(terms, paste0(cval))
    } else if (k == 2) {
      terms <- c(terms, paste0(cval, " * x"))
    } else {
      terms <- c(terms, paste0(cval, " * x^", k-1))
    }
  }
  paste("y(x) =", paste(terms, collapse = " + "))
}

make_poly_equation(poli3$coef)

# Visualisasi
plot(x3, y3, pch=19, col="black",
     xlab="x", ylab="y", main="Interpolasi Polinomial")

# Kurva interpolasi
xx <- seq(min(x3), max(x3), length.out=200)
lines(xx, poli3$predict(xx), col="black", lwd=2)

# Titik hampiran
points(hampiran$x, hampiran$y_hampiran, pch=19, col="red")

# Legenda
legend("topleft",
       legend=c("Titik asli", "Kurva interpolasi", "Titik hampiran"),
       col=c("black", "black", "red"),
       pch=c(19, NA, 19), lty=c(NA,1,NA))

lagrange <- function(x, x_data, y_data, return_basis = FALSE) {
  stopifnot(length(x_data) == length(y_data))
  n <- length(x_data)
  
  # fungsi bantu untuk 1 titik
  lagrange_single <- function(x0) {
    L <- 0
    basis_list <- numeric(n)
    
    for (i in 1:n) {
      pi_val <- 1
      for (j in 1:n) {
        if (j != i) {
          if (x_data[i] == x_data[j]) {
            stop(paste("Error: Titik x_data duplikat pada indeks", i, "dan", j))
          } else {
            pi_val <- pi_val * ((x0 - x_data[j]) / (x_data[i] - x_data[j]))
          }
        }
      }
      basis_list[i] <- pi_val
      L <- L + y_data[i] * pi_val
    }
    return(list(interpolasi = L, basis = basis_list))
  }
  
  # kalau x vektor
  if (length(x) > 1) {
    hasil_list <- lapply(x, lagrange_single)
    if (return_basis) {
      interpolasi <- sapply(hasil_list, function(z) z$interpolasi)
      basis <- t(sapply(hasil_list, function(z) z$basis))
      colnames(basis) <- paste0("L", 0:(n-1))
      return(data.frame(x = x, interpolasi = interpolasi, basis))
    } else {
      return(sapply(hasil_list, function(z) z$interpolasi))
    }
  } else {
    # kasus x tunggal
    hasil <- lagrange_single(x)
    if (return_basis) {
      return(hasil)
    } else {
      return(hasil$interpolasi)
    }
  }
}

# Berikut empat titik yang diketahui (x0, x1, x2, x3) dan (y0, y1, y2, y3) 
x4 <- c(0, 0.4, 0.8, 1.2)
y4 <- c(1, 0.921061, 0.696707,0.362358)
 
# Tentukan nilai x yang ingin dihampiri y nya
x <- c(0.5, 1)

# Tampilkan hasil hampiran
hampiran <- data.frame(x = x, y_hampiran = lagrange(x, x4, y4))
print(hampiran)

lagrange(x, x4, y4,return_basis = TRUE)

# Visualisasi
plot(x4, y4, pch=19, col="black",
     xlab="x", ylab="y", main="Interpolasi Lagrange")

# Kurva interpolasi halus
xx <- seq(min(x4), max(x4), length.out=200)
lines(xx, lagrange(xx, x4, y4, return_basis = F), col="black", lwd=2)

# Titik hampiran
points(hampiran$x, hampiran$y_hampiran, pch=19, col="red")

# Legenda
legend("topright",
       legend=c("Titik asli", "Kurva interpolasi", "Titik hampiran"),
       col=c("black", "black", "red"),
       pch=c(19, NA, 19), lty=c(NA,1,NA))

newton<- function(x, y, x_eval) {
  n <- length(x)       # ini n+1 sebenarnya
  coef <- y
  for (j in 2:n) {
    for (i in n:j) {
      coef[i] <- (coef[i] - coef[i-1]) / (x[i] - x[i-j+1])
    }
  }
  
  sapply(x_eval, function(xx) {
    result <- coef[n]
    if (n > 1) {
      for (k in (n-1):1) {
        result <- result * (xx - x[k]) + coef[k]
      }
    }
    result
  })
}

x <- c(8, 9, 9.5, 11)
y <- c(2.079442, 2.197225, 2.251292, 2.397895)

# nilai hampiran di x=9.2
hampiran_x <- 9.2
hampiran_y <- newton(x, y, hampiran_x)
hampiran <- data.frame(x = hampiran_x, y_hampiran = hampiran_y)

# Plot titik asli
plot(x, y, pch=19, col="black",
     xlab="x", ylab="y", main="Interpolasi Newton")

# Kurva interpolasi halus
xx <- seq(min(x), max(x), length.out=200)
lines(xx, newton(x, y, xx), col="black", lwd=2)

# Titik hampiran
points(hampiran$x, hampiran$y_hampiran, pch=19, col="red")

# Legenda
legend("topleft",
       legend=c("Titik asli", "Kurva interpolasi", "Titik hampiran"),
       col=c("black", "black", "red"),
       pch=c(19, NA, 19), lty=c(NA,1,NA))