==============================================================================

Dashboard Indeks Kerentanan Sosial (SoVI) - Versi Premium

Author: Meldiro_222313204

Date: 2025-07-22

==============================================================================

# ===== 1. INSTALASI DAN LOAD PACKAGES =====
required_packages <- c(
  "shiny", "shinydashboard", "shinydashboardPlus", "shinyWidgets",
  "dplyr", "ggplot2", "plotly", "DT", "readr",
  "car", "lmtest", "psych", "broom", "nortest",
  "shinycssloaders", "shinyjs", "leaflet", 
  "rmarkdown", "knitr", "htmltools", "tinytex"
)

# Install paket yang belum ada
missing_packages <- required_packages[!required_packages %in% installed.packages()[,"Package"]]
if(length(missing_packages) > 0) {
  install.packages(missing_packages, dependencies = TRUE)
}

# Install tinytex jika belum ada (untuk PDF generation)
if (!requireNamespace("tinytex", quietly = TRUE)) {
  install.packages("tinytex")
  tinytex::install_tinytex()
}

# Load semua packages
invisible(lapply(required_packages, library, character.only = TRUE))
## 
## Attaching package: 'shinydashboard'
## The following object is masked from 'package:graphics':
## 
##     box
## 
## Attaching package: 'shinydashboardPlus'
## The following objects are masked from 'package:shinydashboard':
## 
##     box, dashboardHeader, dashboardPage, dashboardSidebar, messageItem,
##     notificationItem, taskItem
## The following object is masked from 'package:graphics':
## 
##     box
## 
## Attaching package: 'shinyWidgets'
## The following object is masked from 'package:shinydashboardPlus':
## 
##     progressBar
## 
## Attaching package: 'dplyr'
## The following objects are masked from 'package:stats':
## 
##     filter, lag
## The following objects are masked from 'package:base':
## 
##     intersect, setdiff, setequal, union
## 
## Attaching package: 'plotly'
## The following object is masked from 'package:ggplot2':
## 
##     last_plot
## The following object is masked from 'package:stats':
## 
##     filter
## The following object is masked from 'package:graphics':
## 
##     layout
## 
## Attaching package: 'DT'
## The following objects are masked from 'package:shiny':
## 
##     dataTableOutput, renderDataTable
## Loading required package: carData
## 
## Attaching package: 'car'
## The following object is masked from 'package:dplyr':
## 
##     recode
## Loading required package: zoo
## 
## Attaching package: 'zoo'
## The following objects are masked from 'package:base':
## 
##     as.Date, as.Date.numeric
## 
## Attaching package: 'psych'
## The following object is masked from 'package:car':
## 
##     logit
## The following objects are masked from 'package:ggplot2':
## 
##     %+%, alpha
## The following object is masked from 'package:shinyWidgets':
## 
##     progressBar
## The following object is masked from 'package:shinydashboardPlus':
## 
##     progressBar
## 
## Attaching package: 'shinyjs'
## The following object is masked from 'package:lmtest':
## 
##     reset
## The following object is masked from 'package:shinyWidgets':
## 
##     alert
## The following object is masked from 'package:shiny':
## 
##     runExample
## The following objects are masked from 'package:methods':
## 
##     removeClass, show
# ===== 2. LOAD DATA =====
# Fungsi untuk load data dengan error handling
load_sovi_data <- function() {
  tryCatch({
    # Coba load dari URL
    sovi_data <- read_csv("https://raw.githubusercontent.com/bmlmcmc/naspaclust/main/data/sovi_data.csv", 
                          show_col_types = FALSE)
    distance <- read_csv("https://raw.githubusercontent.com/bmlmcmc/naspaclust/main/data/distance.csv", 
                         show_col_types = FALSE)
    
    # Convert DISTRICTCODE ke factor
    sovi_data$DISTRICTCODE <- as.factor(sovi_data$DISTRICTCODE)
    
    return(list(sovi_data = sovi_data, distance = distance))
  }, error = function(e) {
    # Jika file tidak ditemukan, buat sample data
    warning("File data tidak ditemukan. Membuat sample data untuk demo.")
    
    set.seed(123)
    n <- 100
    
    sovi_data <- data.frame(
      DISTRICTCODE = factor(sample(1:10, n, replace = TRUE)),
      CHILDREN = runif(n, 0, 0.4),
      FEMALE = runif(n, 0.45, 0.55),
      ELDERLY = runif(n, 0, 0.3),
      FHEAD = runif(n, 0, 0.4),
      FAMILYSIZE = rnorm(n, 4, 1.5),
      NOELECTRIC = runif(n, 0, 0.2),
      LOWEDU = runif(n, 0, 0.6),
      GROWTH = rnorm(n, 0.02, 0.05),
      POVERTY = runif(n, 0, 0.4),
      ILLITERATE = runif(n, 0, 0.3),
      NOTRAINING = runif(n, 0, 0.5),
      DPRONE = runif(n, 0, 0.3),
      RENTED = runif(n, 0, 0.4),
      NOSEWER = runif(n, 0, 0.5),
      TAPWATER = runif(n, 0.5, 1),
      POPULATION = sample(1000:50000, n, replace = TRUE),
      # Tambahkan koordinat untuk peta
      LATITUDE = runif(n, -8.5, -8.0),
      LONGITUDE = runif(n, 115.0, 115.5)
    )
    
    distance <- data.frame(
      DISTRICTCODE = factor(1:10),
      DISTANCE_TO_CITY = runif(10, 5, 100)
    )
    
    return(list(sovi_data = sovi_data, distance = distance))
  })
}

# Load data
data_list <- load_sovi_data()
## New names:
## • `` -> `...1`
sovi_data <- data_list$sovi_data
distance <- data_list$distance

# Tambahkan koordinat jika tidak ada
if(!"LATITUDE" %in% names(sovi_data)) {
  set.seed(123)
  sovi_data$LATITUDE <- runif(nrow(sovi_data), -8.5, -8.0)
  sovi_data$LONGITUDE <- runif(nrow(sovi_data), 115.0, 115.5)
}

# ===== 3. PREMIUM BLUE STYLING =====
custom_css <- tags$head(
  tags$style(HTML("
    @import url('https://fonts.googleapis.com/css2?family=Poppins:wght@300;400;500;600;700&display=swap');
    
    /* Global Styles */
    * {
      font-family: 'Poppins', sans-serif !important;
    }
    
    body {
      background: linear-gradient(135deg, #E3F2FD 0%, #BBDEFB 100%) !important;
      min-height: 100vh;
    }
    
    /* Header Styles */
    .main-header .navbar {
      background: linear-gradient(135deg, #1976D2 0%, #2196F3 100%) !important;
      border: none !important;
      box-shadow: 0 4px 20px rgba(0,0,0,0.1) !important;
    }
    
    .main-header .logo {
      background: linear-gradient(135deg, #0D47A1 0%, #1565C0 100%) !important;
      color: white !important;
      font-weight: 700 !important;
      font-size: 18px !important;
      text-shadow: 2px 2px 4px rgba(0,0,0,0.3) !important;
      border: none !important;
    }
    
    .main-header .logo:hover {
      background: linear-gradient(135deg, #1565C0 0%, #0D47A1 100%) !important;
      transition: all 0.3s ease !important;
    }
    
    /* Sidebar Styles */
    .main-sidebar {
      background: linear-gradient(180deg, #1A237E 0%, #283593 100%) !important;
      box-shadow: 4px 0 20px rgba(0,0,0,0.1) !important;
    }
    
    .sidebar-menu > li > a {
      color: #E3F2FD !important;
      font-weight: 500 !important;
      transition: all 0.3s ease !important;
      border-left: 3px solid transparent !important;
    }
    
    .sidebar-menu > li > a:hover {
      background: linear-gradient(90deg, #1976D2, #2196F3) !important;
      color: white !important;
      border-left: 3px solid #64B5F6 !important;
      transform: translateX(5px) !important;
    }
    
    .sidebar-menu > li.active > a {
      background: linear-gradient(90deg, #1976D2, #2196F3) !important;
      color: white !important;
      border-left: 3px solid #64B5F6 !important;
      box-shadow: 0 4px 15px rgba(33, 150, 243, 0.3) !important;
    }
    
    .sidebar-menu .treeview-menu > li > a {
      color: #BBDEFB !important;
      padding-left: 35px !important;
    }
    
    .sidebar-menu .treeview-menu > li > a:hover {
      background: rgba(33, 150, 243, 0.2) !important;
      color: #E3F2FD !important;
    }
    
    /* Content Area */
    .content-wrapper {
      background: linear-gradient(135deg, #E3F2FD 0%, #BBDEFB 100%) !important;
      min-height: 100vh !important;
    }
    
    /* Box Styles */
    .box {
      border-radius: 15px !important;
      box-shadow: 0 10px 30px rgba(0,0,0,0.1) !important;
      border: none !important;
      overflow: hidden !important;
      transition: all 0.3s ease !important;
      background: white !important;
    }
    
    .box:hover {
      transform: translateY(-5px) !important;
      box-shadow: 0 20px 40px rgba(0,0,0,0.15) !important;
    }
    
    .box-header {
      border-radius: 15px 15px 0 0 !important;
      padding: 20px !important;
    }
    
    .box.box-solid.box-primary > .box-header {
      background: linear-gradient(135deg, #1976D2, #2196F3) !important;
      color: white !important;
    }
    
    .box.box-solid.box-info > .box-header {
      background: linear-gradient(135deg, #0288D1, #03A9F4) !important;
      color: white !important;
    }
    
    .box.box-solid.box-success > .box-header {
      background: linear-gradient(135deg, #388E3C, #4CAF50) !important;
      color: white !important;
    }
    
    .box.box-solid.box-warning > .box-header {
      background: linear-gradient(135deg, #FFA000, #FFC107) !important;
      color: white !important;
    }
    
    .box.box-solid.box-danger > .box-header {
      background: linear-gradient(135deg, #D32F2F, #F44336) !important;
      color: white !important;
    }
    
    .box-title {
      font-weight: 600 !important;
      font-size: 18px !important;
      text-shadow: 1px 1px 2px rgba(0,0,0,0.1) !important;
    }
    
    /* Value Boxes */
    .small-box {
      border-radius: 20px !important;
      box-shadow: 0 15px 35px rgba(0,0,0,0.1) !important;
      transition: all 0.3s ease !important;
      overflow: hidden !important;
      position: relative !important;
    }
    
    .small-box::before {
      content: '';
      position: absolute;
      top: 0;
      left: 0;
      right: 0;
      bottom: 0;
      background: linear-gradient(45deg, rgba(255,255,255,0.1), transparent);
      pointer-events: none;
    }
    
    .small-box:hover {
      transform: translateY(-10px) scale(1.02) !important;
      box-shadow: 0 25px 50px rgba(0,0,0,0.2) !important;
    }
    
    .small-box h3 {
      font-weight: 700 !important;
      font-size: 2.5em !important;
      text-shadow: 2px 2px 4px rgba(0,0,0,0.3) !important;
    }
    
    .small-box p {
      font-weight: 500 !important;
      text-shadow: 1px 1px 2px rgba(0,0,0,0.2) !important;
    }
    
    /* Buttons */
    .btn {
      border-radius: 25px !important;
      font-weight: 600 !important;
      padding: 10px 25px !important;
      transition: all 0.3s ease !important;
      border: none !important;
      text-transform: uppercase !important;
      letter-spacing: 1px !important;
      position: relative !important;
      overflow: hidden !important;
    }
    
    .btn::before {
      content: '';
      position: absolute;
      top: 0;
      left: -100%;
      width: 100%;
      height: 100%;
      background: linear-gradient(90deg, transparent, rgba(255,255,255,0.2), transparent);
      transition: left 0.5s;
    }
    
    .btn:hover::before {
      left: 100%;
    }
    
    .btn-primary {
      background: linear-gradient(135deg, #1976D2, #2196F3) !important;
      box-shadow: 0 5px 15px rgba(33, 150, 243, 0.3) !important;
    }
    
    .btn-primary:hover {
      background: linear-gradient(135deg, #2196F3, #1976D2) !important;
      transform: translateY(-2px) !important;
      box-shadow: 0 8px 25px rgba(33, 150, 243, 0.4) !important;
    }
    
    .btn-success {
      background: linear-gradient(135deg, #388E3C, #4CAF50) !important;
      box-shadow: 0 5px 15px rgba(76, 175, 80, 0.3) !important;
    }
    
    .btn-success:hover {
      background: linear-gradient(135deg, #4CAF50, #388E3C) !important;
      transform: translateY(-2px) !important;
      box-shadow: 0 8px 25px rgba(76, 175, 80, 0.4) !important;
    }
    
    /* Form Controls */
    .form-control {
      border-radius: 10px !important;
      border: 2px solid #BBDEFB !important;
      padding: 12px 15px !important;
      transition: all 0.3s ease !important;
      font-weight: 500 !important;
    }
    
    .form-control:focus {
      border-color: #2196F3 !important;
      box-shadow: 0 0 0 0.2rem rgba(33, 150, 243, 0.25) !important;
      transform: translateY(-1px) !important;
    }
    
    /* DataTables */
    .dataTables_wrapper {
      background: white !important;
      border-radius: 15px !important;
      padding: 20px !important;
      box-shadow: 0 5px 15px rgba(0,0,0,0.08) !important;
    }
    
    .dataTables_wrapper .dataTables_length,
    .dataTables_wrapper .dataTables_filter {
      margin-bottom: 20px !important;
    }
    
    .dataTables_wrapper .dataTables_filter input {
      border-radius: 20px !important;
      border: 2px solid #BBDEFB !important;
      padding: 8px 15px !important;
    }
    
    /* Leaflet Map */
    .leaflet-container {
      height: 500px !important;
      border-radius: 15px !important;
      box-shadow: 0 10px 30px rgba(0,0,0,0.1) !important;
    }
    
    /* Tab Panels */
    .nav-tabs-custom > .nav-tabs {
      border-bottom: 3px solid #2196F3 !important;
      background: linear-gradient(90deg, #E3F2FD, #BBDEFB) !important;
      border-radius: 10px 10px 0 0 !important;
    }
    
    .nav-tabs-custom > .nav-tabs > li.active {
      border-top-color: #1976D2 !important;
    }
    
    .nav-tabs-custom > .nav-tabs > li.active > a {
      background: linear-gradient(135deg, #1976D2, #2196F3) !important;
      color: white !important;
      font-weight: 600 !important;
      border-radius: 10px 10px 0 0 !important;
    }
    
    /* Animations */
    @keyframes fadeInUp {
      from {
        opacity: 0;
        transform: translateY(30px);
      }
      to {
        opacity: 1;
        transform: translateY(0);
      }
    }
    
    .box, .small-box {
      animation: fadeInUp 0.6s ease-out !important;
    }
    
    /* Loading Spinners */
    .shiny-spinner-output-container {
      background: rgba(255,255,255,0.9) !important;
      border-radius: 15px !important;
    }
    
    /* Responsive Design */
    @media (max-width: 768px) {
      .box {
        margin: 10px 5px !important;
      }
      
      .small-box h3 {
        font-size: 2em !important;
      }
    }
    
    /* Custom Scrollbar */
    ::-webkit-scrollbar {
      width: 8px;
    }
    
    ::-webkit-scrollbar-track {
      background: #E3F2FD;
      border-radius: 10px;
    }
    
    ::-webkit-scrollbar-thumb {
      background: linear-gradient(135deg, #1976D2, #2196F3);
      border-radius: 10px;
    }
    
    ::-webkit-scrollbar-thumb:hover {
      background: linear-gradient(135deg, #2196F3, #1976D2);
    }
    
    /* Header Title Animation */
    .main-header .logo {
      position: relative;
      overflow: hidden;
    }
    
    .main-header .logo::after {
      content: '';
      position: absolute;
      top: 0;
      left: -100%;
      width: 100%;
      height: 100%;
      background: linear-gradient(90deg, transparent, rgba(255,255,255,0.2), transparent);
      transition: left 2s;
    }
    
    .main-header .logo:hover::after {
      left: 100%;
    }
    
    /* Custom interpretation box */
    .interpretation-box {
      background-color: #E3F2FD;
      border: 1px solid #BBDEFB;
      border-radius: 8px;
      padding: 15px;
      margin-top: 15px;
      box-shadow: 0 2px 5px rgba(0,0,0,0.1);
    }

    .interpretation-box h4 {
      color: #1565C0;
      margin-top: 0;
      border-bottom: 1px solid #BBDEFB;
      padding-bottom: 8px;
    }

    .interpretation-box p {
      color: #37474F;
      margin-bottom: 8px;
    }
    
    /* Status messages */
    .status-success {
      background-color: #E8F5E9;
      color: #2E7D32;
      padding: 10px;
      border-radius: 5px;
      border: 1px solid #4CAF50;
      margin-top: 10px;
      box-shadow: 0 2px 4px rgba(0,0,0,0.1);
    }

    .status-error {
      background-color: #FFEBEE;
      color: #C62828;
      padding: 10px;
      border-radius: 5px;
      border: 1px solid #F44336;
      margin-top: 10px;
      box-shadow: 0 2px 4px rgba(0,0,0,0.1);
    }
    
    /* Welcome header */
    .welcome-header {
      color: #1565C0;
      text-align: center;
      font-family: 'Poppins', sans-serif;
      font-weight: bold;
      margin-bottom: 20px;
      text-shadow: 1px 1px 2px rgba(0,0,0,0.1);
    }
  "))
)

# ===== 4. USER INTERFACE (UI) =====
ui <- dashboardPage(
  header = dashboardHeader(
    title = "🎯 SoVI Dashboard Premium - Meldiro"
  ),
  
  sidebar = dashboardSidebar(
    width = 300,
    sidebarMenu(
      id = "sidebar_menu",
      menuItem("🏠 Beranda", tabName = "beranda", icon = icon("home")),
      menuItem("⚙️ Manajemen Data", tabName = "manajemen_data", icon = icon("cogs")),
      menuItem("📊 Eksplorasi Data", tabName = "eksplorasi_data", icon = icon("chart-bar")),
      menuItem("🗺️ Peta Interaktif", tabName = "peta", icon = icon("map")),
      menuItem("✅ Uji Asumsi Data", tabName = "uji_asumsi", icon = icon("check-circle")),
      menuItem("🔬 Statistik Inferensia", icon = icon("flask"),
               menuSubItem("📈 Uji Beda Rata-rata", tabName = "uji_beda_rata", icon = icon("balance-scale")),
               menuSubItem("📊 Uji Proporsi & Varians", tabName = "uji_proporsi_varians", icon = icon("percent")),
               menuSubItem("🔍 Uji ANOVA", tabName = "uji_anova", icon = icon("search"))
      ),
      menuItem("📈 Regresi Linear", tabName = "regresi_linear", icon = icon("chart-line")),
      menuItem("📥 Unduh Laporan", tabName = "unduh_laporan", icon = icon("download"))
    )
  ),
  
  body = dashboardBody(
    custom_css,
    useShinyjs(),
    
    tabItems(
      # ===== TAB BERANDA =====
      tabItem(tabName = "beranda",
        fluidRow(
          box(
            title = "🎉 Selamat Datang di Dashboard Indeks Kerentanan Sosial (SoVI)",
            width = 12, solidHeader = TRUE, status = "primary",
            div(class = "text-center",
              h2("🚀 Dashboard Analisis Komprehensif SoVI", style = "color: #1565C0; margin-bottom: 20px; font-weight: 700;"),
              h3("👨‍💻 Author: Meldiro_222313204", style = "color: #1976D2; margin-bottom: 20px; font-weight: 500;"),
              p("Dashboard interaktif premium ini dirancang untuk melakukan analisis mendalam terhadap data Indeks Kerentanan Sosial (SoVI) dengan teknologi terdepan dan visualisasi yang memukau.",
                 style = "font-size: 18px; line-height: 1.8; color: #37474F; font-weight: 400;")
            )
          )
        ),
        
        fluidRow(
          valueBoxOutput("total_districts", width = 3),
          valueBoxOutput("avg_population", width = 3),
          valueBoxOutput("poverty_rate", width = 3),
          valueBoxOutput("data_completeness", width = 3)
        ),
        
        fluidRow(
          box(
            title = "✨ Fitur Utama Dashboard Premium", width = 6, solidHeader = TRUE, status = "info",
            tags$div(
              tags$ul(
                tags$li("🎨 Eksplorasi data interaktif dengan visualisasi modern"),
                tags$li("🗺️ Peta interaktif untuk visualisasi spasial yang memukau"),
                tags$li("🔄 Transformasi data (binning) untuk analisis kategorikal"),
                tags$li("✅ Pengujian asumsi statistik (normalitas, homogenitas)"),
                tags$li("🔬 Uji statistik inferensial lengkap dan komprehensif"),
                tags$li("📈 Analisis regresi linear berganda dengan diagnostik"),
                tags$li("📄 Ekspor laporan PDF formal dengan desain profesional")
              ),
              style = "font-size: 16px; line-height: 2; color: #37474F; font-weight: 500;"
            )
          ),
          
          box(
            title = "📋 Informasi Dataset", width = 6, solidHeader = TRUE, status = "success",
            withSpinner(DTOutput("dataset_info"), type = 6, color = "#2196F3")
          )
        )
      ),
      
      # ===== TAB MANAJEMEN DATA =====
      tabItem(tabName = "manajemen_data",
        h2("⚙️ Manajemen Data: Transformasi Variabel", style = "color: #1565C0; margin-bottom: 30px; font-weight: 700;"),
        
        fluidRow(
          box(
            title = "🎛️ Pengaturan Binning", width = 4, solidHeader = TRUE, status = "info",
            selectInput("bin_variable", "📊 Pilih Variabel untuk Binning:",
                       choices = names(select_if(sovi_data, is.numeric)),
                       selected = "POPULATION"),
            
            sliderInput("num_bins", "🔢 Jumlah Bin (Kategori):",
                        value = 4, min = 2, max = 10, step = 1),
            
            radioButtons("bin_method", "⚡ Metode Binning:",
                        choices = list("Equal Width" = "width", "Equal Frequency" = "frequency"),
                        selected = "width"),
            
            actionButton("apply_binning", "🚀 Terapkan Binning", 
                         class = "btn-primary", icon = icon("play"))
          ),
          
          box(
            title = "📈 Hasil Binning", width = 8, solidHeader = TRUE, status = "success",
            tabsetPanel(
              tabPanel("📋 Ringkasan",
                       withSpinner(verbatimTextOutput("binning_summary"), type = 6)),
              tabPanel("📊 Visualisasi",
                       withSpinner(plotlyOutput("binning_plot"), type = 6)),
              tabPanel("📄 Data",
                       withSpinner(DTOutput("binned_data_table"), type = 6))
            ),
            br(),
            downloadButton("download_binned_data", "💾 Unduh Data Hasil Binning",
                           class = "btn-success")
          )
        )
      ),
      
      # ===== TAB EKSPLORASI DATA =====
      tabItem(tabName = "eksplorasi_data",
        h2("📊 Eksplorasi Data Komprehensif", style = "color: #1565C0; margin-bottom: 30px; font-weight: 700;"),
        
        fluidRow(
          box(
            title = "📈 Statistik Deskriptif", width = 12, solidHeader = TRUE, status = "primary",
            fluidRow(
              column(6,
                selectInput("desc_var", "📊 Pilih Variabel:",
                            choices = names(sovi_data),
                            selected = names(select_if(sovi_data, is.numeric))[1:min(5, length(names(select_if(sovi_data, is.numeric))))],
                           multiple = TRUE)
              ),
              column(6,
                actionButton("update_desc", "🔄 Update Statistik",
                            class = "btn-primary")
              )
            ),
            br(),
            withSpinner(DTOutput("descriptive_stats"), type = 6)
          )
        ),
        
        fluidRow(
          box(
            title = "📊 Visualisasi Univariat", width = 6, solidHeader = TRUE, status = "warning",
            selectInput("plot_type", "🎨 Jenis Visualisasi:",
                        choices = c("Histogram" = "hist", "Density Plot" = "density",
                                  "Boxplot" = "box", "Violin Plot" = "violin")),
            selectInput("hist_var", "📊 Variabel:",
                        choices = names(select_if(sovi_data, is.numeric))),
            withSpinner(plotlyOutput("univariate_plot"), type = 6),
            br(),
            downloadButton("download_univariate_png", "📥 Unduh Plot (PNG)", class = "btn-primary"),
            downloadButton("download_univariate_pdf", "📥 Unduh Plot (PDF)", class = "btn-success")
          ),
          
          box(
            title = "📈 Visualisasi Bivariat", width = 6, solidHeader = TRUE, status = "warning",
            selectInput("scatter_x", "📊 Variabel X:",
                        choices = names(select_if(sovi_data, is.numeric)),
                        selected = names(select_if(sovi_data, is.numeric))[1]),
            selectInput("scatter_y", "📊 Variabel Y:",
                        choices = names(select_if(sovi_data, is.numeric)),
                        selected = names(select_if(sovi_data, is.numeric))[min(2, length(names(select_if(sovi_data, is.numeric))))]),
            checkboxInput("add_regression", "📈 Tambahkan Garis Regresi", value = TRUE),
            withSpinner(plotlyOutput("bivariate_plot"), type = 6),
            br(),
            downloadButton("download_bivariate_png", "📥 Unduh Plot (PNG)", class = "btn-primary"),
            downloadButton("download_bivariate_pdf", "📥 Unduh Plot (PDF)", class = "btn-success")
          )
        )
      ),
      
      # ===== TAB PETA INTERAKTIF =====
      tabItem(tabName = "peta",
        h2("🗺️ Peta Interaktif SoVI Premium", style = "color: #1565C0; margin-bottom: 30px; font-weight: 700;"),
        
        fluidRow(
          box(
            title = "🎛️ Pengaturan Peta", width = 4, solidHeader = TRUE, status = "info",
            selectInput("map_variable", "📊 Pilih Variabel untuk Visualisasi:",
                       choices = names(select_if(sovi_data, is.numeric)),
                       selected = "POVERTY"),
            
            selectInput("map_color", "🎨 Skema Warna:",
                       choices = c("Blues", "Reds", "Greens", "Oranges", "Purples", "Viridis"),
                       selected = "Blues"),
            
            checkboxInput("show_districts", "🏷️ Tampilkan Label District", value = TRUE),
            
            actionButton("update_map", "🔄 Update Peta", 
                         class = "btn-primary", icon = icon("refresh"))
          ),
          
          box(
            title = "🗺️ Peta SoVI Interaktif", width = 8, solidHeader = TRUE, status = "primary",
            withSpinner(leafletOutput("sovi_map", height = 500), type = 6),
            br(),
            downloadButton("download_map_png", "📥 Unduh Peta (PNG)", class = "btn-primary")
          )
        ),
        
        fluidRow(
          box(
            title = "📊 Statistik Spasial", width = 12, solidHeader = TRUE, status = "success",
            withSpinner(DTOutput("spatial_stats"), type = 6)
          )
        )
      ),
      
      # ===== TAB UJI ASUMSI =====
      tabItem(tabName = "uji_asumsi",
        h2("✅ Pengujian Asumsi Statistik", style = "color: #1565C0; margin-bottom: 30px; font-weight: 700;"),
        
        fluidRow(
          box(
            title = "📊 Uji Normalitas", width = 6, solidHeader = TRUE, status = "info",
            selectInput("normality_var", "📊 Pilih Variabel:",
                        choices = names(select_if(sovi_data, is.numeric))),
            actionButton("run_normality", "🚀 Jalankan Uji", class = "btn-primary"),
            br(), br(),
            withSpinner(verbatimTextOutput("normality_results"), type = 6),
            withSpinner(plotOutput("normality_plot"), type = 6),
            br(),
            downloadButton("download_normality_png", "📥 Unduh Plot (PNG)", class = "btn-primary"),
            downloadButton("download_normality_pdf", "📥 Unduh Plot (PDF)", class = "btn-success")
          ),
          
          box(
            title = "⚖️ Uji Homogenitas Varians", width = 6, solidHeader = TRUE, status = "success",
            selectInput("homogeneity_var", "📊 Variabel Numerik:",
                        choices = names(select_if(sovi_data, is.numeric))),
            selectInput("homogeneity_group", "👥 Variabel Kelompok:",
                        choices = names(select_if(sovi_data, is.factor))),
            actionButton("run_homogeneity", "🚀 Jalankan Uji", class = "btn-success"),
            br(), br(),
            withSpinner(verbatimTextOutput("homogeneity_results"), type = 6),
            withSpinner(plotOutput("homogeneity_plot"), type = 6),
            br(),
            downloadButton("download_homogeneity_png", "📥 Unduh Plot (PNG)", class = "btn-primary"),
            downloadButton("download_homogeneity_pdf", "📥 Unduh Plot (PDF)", class = "btn-success")
          )
        )
      ),
      
      # ===== TAB UJI BEDA RATA-RATA =====
      tabItem(tabName = "uji_beda_rata",
        h2("📈 Uji Beda Rata-rata", style = "color: #1565C0; margin-bottom: 30px; font-weight: 700;"),
        
        fluidRow(
          box(
            title = "🎛️ Pengaturan Uji", width = 4, solidHeader = TRUE, status = "primary",
            selectInput("ttest_var", "📊 Variabel Numerik:",
                        choices = names(select_if(sovi_data, is.numeric))),
            radioButtons("ttest_type", "🔬 Jenis Uji:",
                        choices = list("One Sample t-test" = "one",
                                     "Two Sample t-test" = "two")),
            conditionalPanel(
              condition = "input.ttest_type == 'one'",
              numericInput("mu_value", "🎯 Nilai μ₀:", value = 0)
            ),
            conditionalPanel(
              condition = "input.ttest_type == 'two'",
              selectInput("ttest_group", "👥 Variabel Kelompok:",
                          choices = names(select_if(sovi_data, is.factor)))
            ),
            numericInput("alpha_level", "⚡ Tingkat Signifikansi (α):",
                         value = 0.05, min = 0.01, max = 0.1, step = 0.01),
            actionButton("run_ttest", "🚀 Jalankan Uji t", class = "btn-primary")
          ),
          
          box(
            title = "📊 Hasil Uji t", width = 8, solidHeader = TRUE, status = "success",
            withSpinner(verbatimTextOutput("ttest_results"), type = 6),
            withSpinner(plotOutput("ttest_plot"), type = 6),
            br(),
            downloadButton("download_ttest_png", "📥 Unduh Plot (PNG)", class = "btn-primary"),
            downloadButton("download_ttest_pdf", "📥 Unduh Plot (PDF)", class = "btn-success")
          )
        )
      ),
      
      # ===== TAB UJI PROPORSI & VARIANS =====
      tabItem(tabName = "uji_proporsi_varians",
        h2("📊 Uji Proporsi & Varians", style = "color: #1565C0; margin-bottom: 30px; font-weight: 700;"),
        
        fluidRow(
          box(
            title = "📊 Uji Proporsi", width = 6, solidHeader = TRUE, status = "info",
            selectInput("prop_var", "👥 Variabel Kategorikal:",
                        choices = names(select_if(sovi_data, is.factor))),
            numericInput("prop_value", "🎯 Proporsi yang diuji (p₀):",
                         value = 0.5, min = 0, max = 1, step = 0.01),
            actionButton("run_prop_test", "🚀 Jalankan Uji Proporsi",
                         class = "btn-info"),
            br(), br(),
            withSpinner(verbatimTextOutput("prop_test_results"), type = 6)
          ),
          
          box(
            title = "⚖️ Uji Varians", width = 6, solidHeader = TRUE, status = "warning",
            selectInput("var_test_var1", "📊 Variabel 1:",
                        choices = names(select_if(sovi_data, is.numeric))),
            selectInput("var_test_var2", "📊 Variabel 2:",
                        choices = names(select_if(sovi_data, is.numeric))),
            actionButton("run_var_test", "🚀 Jalankan Uji Varians",
                         class = "btn-warning"),
            br(), br(),
            withSpinner(verbatimTextOutput("var_test_results"), type = 6),
            withSpinner(plotOutput("var_test_plot"), type = 6),
            br(),
            downloadButton("download_var_test_png", "📥 Unduh Plot (PNG)", class = "btn-primary"),
            downloadButton("download_var_test_pdf", "📥 Unduh Plot (PDF)", class = "btn-success")
          )
        )
      ),
      
      # ===== TAB UJI ANOVA =====
      tabItem(tabName = "uji_anova",
        h2("🔍 Uji ANOVA", style = "color: #1565C0; margin-bottom: 30px; font-weight: 700;"),
        
        fluidRow(
          box(
            title = "🎛️ Pengaturan ANOVA", width = 4, solidHeader = TRUE, status = "primary",
            selectInput("anova_dependent", "📊 Variabel Dependen:",
                        choices = names(select_if(sovi_data, is.numeric))),
            selectInput("anova_independent", "👥 Variabel Independen (Faktor):",
                        choices = names(select_if(sovi_data, is.factor))),
            actionButton("run_anova", "🚀 Jalankan ANOVA",
                         class = "btn-primary")
          ),
          
          box(
            title = "📊 Hasil ANOVA", width = 8, solidHeader = TRUE, status = "success",
            withSpinner(verbatimTextOutput("anova_results"), type = 6),
            withSpinner(plotOutput("anova_plot"), type = 6),
            br(),
            downloadButton("download_anova_png", "📥 Unduh Plot (PNG)", class = "btn-primary"),
            downloadButton("download_anova_pdf", "📥 Unduh Plot (PDF)", class = "btn-success")
          )
        )
      ),
      
      # ===== TAB REGRESI LINEAR =====
      tabItem(tabName = "regresi_linear",
        h2("📈 Analisis Regresi Linear Berganda", style = "color: #1565C0; margin-bottom: 30px; font-weight: 700;"),
        
        fluidRow(
          box(
            title = "🎛️ Spesifikasi Model", width = 4, solidHeader = TRUE, status = "primary",
            selectInput("dependent_var", "🎯 Variabel Dependen (Y):",
                        choices = names(select_if(sovi_data, is.numeric))),
            selectInput("independent_vars", "📊 Variabel Independen (X):",
                        choices = names(select_if(sovi_data, is.numeric)),
                       multiple = TRUE),
            actionButton("run_regression", "🚀 Jalankan Regresi",
                         class = "btn-primary", icon = icon("play"))
          ),
          
          box(
            title = "📊 Hasil Regresi", width = 8, solidHeader = TRUE, status = "success",
            tabsetPanel(
              tabPanel("📋 Ringkasan Model", withSpinner(verbatimTextOutput("regression_summary"), type = 6)),
              tabPanel("🔍 Diagnostik", withSpinner(plotOutput("regression_diagnostics"), type = 6)),
              tabPanel("🎯 Prediksi", withSpinner(plotlyOutput("regression_prediction"), type = 6))
            ),
            br(),
            downloadButton("download_regression_png", "📥 Unduh Plot Diagnostik (PNG)", class = "btn-primary"),
            downloadButton("download_regression_pdf", "📥 Unduh Plot Diagnostik (PDF)", class = "btn-success")
          )
        )
      ),
      
      # ===== TAB UNDUH LAPORAN =====
      tabItem(tabName = "unduh_laporan",
        h2("📥 Unduh Laporan Analisis Premium", style = "color: #1565C0; margin-bottom: 30px; font-weight: 700;"),
        
        fluidRow(
          box(
            title = "📄 Generator Laporan PDF Premium", width = 12, solidHeader = TRUE, status = "primary",
            fluidRow(
              column(6,
                h4("✅ Pilih Komponen Laporan:", style = "color: #1565C0; font-weight: 600;"),
                checkboxGroupInput("report_sections", "",
                                 choices = list(
                                   "📊 Statistik Deskriptif" = "descriptive",
                                   "✅ Uji Normalitas" = "normality",
                                   "⚖️ Uji Homogenitas" = "homogeneity", 
                                   "📈 Uji t-test" = "ttest",
                                   "🔍 Uji ANOVA" = "anova",
                                   "📈 Analisis Regresi" = "regression"
                                 ),
                                 selected = c("descriptive", "normality", "ttest"))
              ),
              column(6,
                h4("📝 Informasi Laporan:", style = "color: #1565C0; font-weight: 600;"),
                textInput("report_title", "📋 Judul Laporan:",
                          value = "Analisis SoVI - Social Vulnerability Index"),
                textInput("report_author", "👨‍💻 Nama Penulis:",
                          value = "Meldiro_222313204"),
                textInput("report_institution", "🏛️ Institusi:",
                          value = "Universitas"),
                br(),
                downloadButton("download_report", "📄 Generate & Download PDF Report",
                              class = "btn-success btn-lg", icon = icon("file-pdf"))
              )
            )
          )
        )
      )
    )
  )
)

# ===== 5. SERVER LOGIC =====
server <- function(input, output, session) {
  
  # Reactive values untuk menyimpan data dan hasil analisis
  values <- reactiveValues(
    binned_data = NULL,
    current_model = NULL,
    normality_result = NULL,
    homogeneity_result = NULL,
    ttest_result = NULL,
    anova_result = NULL
  )
  
  # ===== BERANDA OUTPUTS =====
  output$total_districts <- renderValueBox({
    valueBox(
      value = length(unique(sovi_data$DISTRICTCODE)),
      subtitle = "Total Districts",
      icon = icon("map"),
      color = "blue"
    )
  })
  
  output$avg_population <- renderValueBox({
    valueBox(
      value = format(round(mean(sovi_data$POPULATION, na.rm = TRUE)), big.mark = ","),
      subtitle = "Avg Population",
      icon = icon("users"),
      color = "green"
    )
  })
  
  output$poverty_rate <- renderValueBox({
    valueBox(
      value = paste0(round(mean(sovi_data$POVERTY, na.rm = TRUE) * 100, 1), "%"),
      subtitle = "Avg Poverty Rate",
      icon = icon("chart-line"),
      color = "yellow"
    )
  })
  
  output$data_completeness <- renderValueBox({
    completeness <- round((1 - sum(is.na(sovi_data)) / (nrow(sovi_data) * ncol(sovi_data))) * 100, 1)
    valueBox(
      value = paste0(completeness, "%"),
      subtitle = "Data Completeness",
      icon = icon("check-circle"),
      color = "purple"
    )
  })
  
  output$dataset_info <- renderDT({
    info_df <- data.frame(
      Variable = names(sovi_data),
      Type = sapply(sovi_data, class),
      Missing = sapply(sovi_data, function(x) sum(is.na(x))),
      Missing_Percent = round(sapply(sovi_data, function(x) sum(is.na(x))/length(x) * 100), 2),
      stringsAsFactors = FALSE
    )
    
    datatable(info_df,
               options = list(pageLength = 10, scrollX = TRUE),
              rownames = FALSE) %>%
      formatStyle("Missing_Percent",
                   backgroundColor = styleInterval(c(5, 15), c("lightgreen", "yellow", "lightcoral")))
  })
  
  # ===== MANAJEMEN DATA =====
  observeEvent(input$apply_binning, {
    req(input$bin_variable, input$num_bins)
    
    var_to_bin <- input$bin_variable
    num_bins <- input$num_bins
    
    if (!is.numeric(sovi_data[[var_to_bin]])) {
      showNotification("Variabel yang dipilih bukan numerik!", type = "error")
      return()
    }
    
    # Lakukan binning
    tryCatch({
      if (input$bin_method == "width") {
        binned_col <- cut(sovi_data[[var_to_bin]],
                          breaks = num_bins,
                          include.lowest = TRUE,
                          dig.lab = 4)
      } else {
        binned_col <- cut(sovi_data[[var_to_bin]],
                          breaks = quantile(sovi_data[[var_to_bin]],
                                          probs = seq(0, 1, length.out = num_bins + 1),
                                          na.rm = TRUE),
                         include.lowest = TRUE,
                          dig.lab = 4)
      }
      
      values$binned_data <- sovi_data %>%
        mutate(!!paste0(var_to_bin, "_Binned") := binned_col)
      
      showNotification("Binning berhasil diterapkan!", type = "success")
    }, error = function(e) {
      showNotification(paste("Error dalam binning:", e$message), type = "error")
    })
  })
  
  output$binning_summary <- renderPrint({
    req(values$binned_data)
    var_to_bin <- input$bin_variable
    binned_col_name <- paste0(var_to_bin, "_Binned")
    
    if (binned_col_name %in% names(values$binned_data)) {
      cat("=== RINGKASAN BINNING ===\n")
      cat("Variabel:", var_to_bin, "\n")
      cat("Metode:", ifelse(input$bin_method == "width", "Equal Width", "Equal Frequency"), "\n")
      cat("Jumlah Bin:", input$num_bins, "\n\n")
      
      summary_table <- table(values$binned_data[[binned_col_name]], useNA = "ifany")
      print(summary_table)
      
      cat("\n=== STATISTIK DESKRIPTIF PER KATEGORI ===\n")
      desc_by_bin <- values$binned_data %>%
        group_by(!!sym(binned_col_name)) %>%
        summarise(
          Count = n(),
          Mean = round(mean(!!sym(var_to_bin), na.rm = TRUE), 3),
          Median = round(median(!!sym(var_to_bin), na.rm = TRUE), 3),
          SD = round(sd(!!sym(var_to_bin), na.rm = TRUE), 3),
          .groups = 'drop'
        )
      print(desc_by_bin)
    }
  })
  
  output$binning_plot <- renderPlotly({
    req(values$binned_data)
    var_to_bin <- input$bin_variable
    binned_col_name <- paste0(var_to_bin, "_Binned")
    
    if (binned_col_name %in% names(values$binned_data)) {
      p <- ggplot(values$binned_data, aes_string(x = binned_col_name)) +
        geom_bar(fill = "#2196F3", alpha = 0.7, color = "white") +
        labs(title = paste("Distribusi Kategori:", var_to_bin),
             x = "Kategori", y = "Frekuensi") +
        theme_minimal() +
        theme(axis.text.x = element_text(angle = 45, hjust = 1))
      
      ggplotly(p)
    }
  })
  
  output$binned_data_table <- renderDT({
    req(values$binned_data)
    datatable(values$binned_data,
               options = list(pageLength = 10, scrollX = TRUE),
              filter = 'top')
  })
  
  # ===== EKSPLORASI DATA =====
  output$descriptive_stats <- renderDT({
    req(input$desc_var)
    
    data_for_desc <- sovi_data %>% select(all_of(input$desc_var))
    numeric_cols <- select_if(data_for_desc, is.numeric)
    
    if (ncol(numeric_cols) > 0) {
      desc_stats <- describe(numeric_cols) %>%
        as.data.frame() %>%
        select(n, mean, sd, median, min, max, skew, kurtosis) %>%
        mutate_if(is.numeric, round, 3)
      
      datatable(desc_stats, 
                options = list(pageLength = 15, scrollX = TRUE),
                rownames = TRUE) %>%
        formatRound(columns = c('mean', 'sd', 'median', 'min', 'max', 'skew', 'kurtosis'), digits = 3)
    }
  })
  
  output$univariate_plot <- renderPlotly({
    req(input$hist_var, input$plot_type)
    
    p <- switch(input$plot_type,
      "hist" = ggplot(sovi_data, aes_string(x = input$hist_var)) +
        geom_histogram(bins = 30, fill = "#2196F3", color = "white", alpha = 0.7) +
        labs(title = paste("Histogram:", input$hist_var), x = input$hist_var, y = "Frekuensi"),
      
      "density" = ggplot(sovi_data, aes_string(x = input$hist_var)) +
        geom_density(fill = "#2196F3", alpha = 0.7) +
        labs(title = paste("Density Plot:", input$hist_var), x = input$hist_var, y = "Density"),
      
      "box" = ggplot(sovi_data, aes_string(y = input$hist_var)) +
        geom_boxplot(fill = "#2196F3", alpha = 0.7) +
        labs(title = paste("Boxplot:", input$hist_var), y = input$hist_var),
      
      "violin" = ggplot(sovi_data, aes_string(x = "''", y = input$hist_var)) +
        geom_violin(fill = "#2196F3", alpha = 0.7) +
        labs(title = paste("Violin Plot:", input$hist_var), y = input$hist_var, x = "")
    )
    
    p <- p + theme_minimal()
    ggplotly(p)
  })
  
  output$bivariate_plot <- renderPlotly({
    req(input$scatter_x, input$scatter_y)
    
    p <- ggplot(sovi_data, aes_string(x = input$scatter_x, y = input$scatter_y)) +
      geom_point(color = "#1976D2", alpha = 0.6, size = 2) +
      labs(title = paste("Scatter Plot:", input$scatter_x, "vs", input$scatter_y),
           x = input$scatter_x, y = input$scatter_y) +
      theme_minimal()
    
    if (input$add_regression) {
      p <- p + geom_smooth(method = "lm", color = "#F44336", se = TRUE, alpha = 0.3)
    }
    
    ggplotly(p)
  })
  
  # Download handlers for plots
  output$download_univariate_png <- downloadHandler(
    filename = function() {
      paste0("univariate_plot_", input$hist_var, "_", Sys.Date(), ".png")
    },
    content = function(file) {
      p <- switch(input$plot_type,
        "hist" = ggplot(sovi_data, aes_string(x = input$hist_var)) +
          geom_histogram(bins = 30, fill = "#2196F3", color = "white", alpha = 0.7) +
          labs(title = paste("Histogram:", input$hist_var), x = input$hist_var, y = "Frekuensi") +
          theme_minimal(),
        
        "density" = ggplot(sovi_data, aes_string(x = input$hist_var)) +
          geom_density(fill = "#2196F3", alpha = 0.7) +
          labs(title = paste("Density Plot:", input$hist_var), x = input$hist_var, y = "Density") +
          theme_minimal(),
        
        "box" = ggplot(sovi_data, aes_string(y = input$hist_var)) +
          geom_boxplot(fill = "#2196F3", alpha = 0.7) +
          labs(title = paste("Boxplot:", input$hist_var), y = input$hist_var) +
          theme_minimal(),
        
        "violin" = ggplot(sovi_data, aes_string(x = "''", y = input$hist_var)) +
          geom_violin(fill = "#2196F3", alpha = 0.7) +
          labs(title = paste("Violin Plot:", input$hist_var), y = input$hist_var, x = "") +
          theme_minimal()
      )
      
      ggsave(file, plot = p, width = 10, height = 6, dpi = 300)
    }
  )
  
  output$download_univariate_pdf <- downloadHandler(
    filename = function() {
      paste0("univariate_plot_", input$hist_var, "_", Sys.Date(), ".pdf")
    },
    content = function(file) {
      p <- switch(input$plot_type,
        "hist" = ggplot(sovi_data, aes_string(x = input$hist_var)) +
          geom_histogram(bins = 30, fill = "#2196F3", color = "white", alpha = 0.7) +
          labs(title = paste("Histogram:", input$hist_var), x = input$hist_var, y = "Frekuensi") +
          theme_minimal(),
        
        "density" = ggplot(sovi_data, aes_string(x = input$hist_var)) +
          geom_density(fill = "#2196F3", alpha = 0.7) +
          labs(title = paste("Density Plot:", input$hist_var), x = input$hist_var, y = "Density") +
          theme_minimal(),
        
        "box" = ggplot(sovi_data, aes_string(y = input$hist_var)) +
          geom_boxplot(fill = "#2196F3", alpha = 0.7) +
          labs(title = paste("Boxplot:", input$hist_var), y = input$hist_var) +
          theme_minimal(),
        
        "violin" = ggplot(sovi_data, aes_string(x = "''", y = input$hist_var)) +
          geom_violin(fill = "#2196F3", alpha = 0.7) +
          labs(title = paste("Violin Plot:", input$hist_var), y = input$hist_var, x = "") +
          theme_minimal()
      )
      
      ggsave(file, plot = p, width = 10, height = 6, device = "pdf")
    }
  )
  
  output$download_bivariate_png <- downloadHandler(
    filename = function() {
      paste0("bivariate_plot_", input$scatter_x, "_vs_", input$scatter_y, "_", Sys.Date(), ".png")
    },
    content = function(file) {
      p <- ggplot(sovi_data, aes_string(x = input$scatter_x, y = input$scatter_y)) +
        geom_point(color = "#1976D2", alpha = 0.6, size = 2) +
        labs(title = paste("Scatter Plot:", input$scatter_x, "vs", input$scatter_y),
             x = input$scatter_x, y = input$scatter_y) +
        theme_minimal()
      
      if (input$add_regression) {
        p <- p + geom_smooth(method = "lm", color = "#F44336", se = TRUE, alpha = 0.3)
      }
      
      ggsave(file, plot = p, width = 10, height = 6, dpi = 300)
    }
  )
  
  output$download_bivariate_pdf <- downloadHandler(
    filename = function() {
      paste0("bivariate_plot_", input$scatter_x, "_vs_", input$scatter_y, "_", Sys.Date(), ".pdf")
    },
    content = function(file) {
      p <- ggplot(sovi_data, aes_string(x = input$scatter_x, y = input$scatter_y)) +
        geom_point(color = "#1976D2", alpha = 0.6, size = 2) +
        labs(title = paste("Scatter Plot:", input$scatter_x, "vs", input$scatter_y),
             x = input$scatter_x, y = input$scatter_y) +
        theme_minimal()
      
      if (input$add_regression) {
        p <- p + geom_smooth(method = "lm", color = "#F44336", se = TRUE, alpha = 0.3)
      }
      
      ggsave(file, plot = p, width = 10, height = 6, device = "pdf")
    }
  )
  
  # ===== PETA INTERAKTIF =====
  output$sovi_map <- renderLeaflet({
    leaflet(sovi_data) %>%
      addTiles() %>%
      setView(lng = mean(sovi_data$LONGITUDE, na.rm = TRUE),
               lat = mean(sovi_data$LATITUDE, na.rm = TRUE),
               zoom = 10)
  })
  
  observeEvent(input$update_map, {
    req(input$map_variable)
    
    # Buat color palette
    pal <- colorNumeric(
      palette = input$map_color,
      domain = sovi_data[[input$map_variable]]
    )
    
    # Update peta
    leafletProxy("sovi_map", data = sovi_data) %>%
      clearMarkers() %>%
      clearControls() %>%
      addCircleMarkers(
        lng = ~LONGITUDE, 
        lat = ~LATITUDE,
        radius = ~sqrt(get(input$map_variable)) * 10,
        color = ~pal(get(input$map_variable)),
        fillOpacity = 0.7,
        popup = ~paste(
          "<b>District:</b>", DISTRICTCODE, "<br>",
          "<b>", input$map_variable, ":</b>", round(get(input$map_variable), 3), "<br>",
          "<b>Population:</b>", format(POPULATION, big.mark = ",")
        ),
        label = if(input$show_districts) ~as.character(DISTRICTCODE) else NULL
      ) %>%
      addLegend(
        pal = pal, 
        values = ~get(input$map_variable),
        title = input$map_variable,
        position = "bottomright"
      )
  })
  
  output$spatial_stats <- renderDT({
    req(input$map_variable)
    
    spatial_summary <- sovi_data %>%
      group_by(DISTRICTCODE) %>%
      summarise(
        Count = n(),
        Mean_Value = round(mean(get(input$map_variable), na.rm = TRUE), 3),
        Min_Value = round(min(get(input$map_variable), na.rm = TRUE), 3),
        Max_Value = round(max(get(input$map_variable), na.rm = TRUE), 3),
        Avg_Lat = round(mean(LATITUDE, na.rm = TRUE), 4),
        Avg_Lng = round(mean(LONGITUDE, na.rm = TRUE), 4),
        .groups = 'drop'
      )
    
    datatable(spatial_summary,
               options = list(pageLength = 10, scrollX = TRUE),
              rownames = FALSE)
  })
  
  # Download map as PNG
  output$download_map_png <- downloadHandler(
    filename = function() {
      paste0("sovi_map_", input$map_variable, "_", Sys.Date(), ".png")
    },
    content = function(file) {
      # Create a static map using ggplot for download
      map_data <- sovi_data
      
      p <- ggplot(map_data, aes(x = LONGITUDE, y = LATITUDE, color = get(input$map_variable))) +
        geom_point(size = 3, alpha = 0.7) +
        scale_color_gradientn(colors = colorRampPalette(brewer.pal(9, input$map_color))(100)) +
        labs(title = paste("SoVI Map:", input$map_variable),
             color = input$map_variable) +
        theme_minimal() +
        coord_fixed(ratio = 1.3)
      
      ggsave(file, plot = p, width = 10, height = 8, dpi = 300)
    }
  )
  
  # ===== UJI ASUMSI =====
  observeEvent(input$run_normality, {
    req(input$normality_var)
    
    var_data <- sovi_data[[input$normality_var]]
    var_data <- var_data[!is.na(var_data)]
    
    # Simpan hasil untuk laporan
    values$normality_result <- list(
      variable = input$normality_var,
      n = length(var_data),
      shapiro = if(length(var_data) <= 5000) shapiro.test(var_data) else NULL,
      ad = ad.test(var_data)
    )
    
    output$normality_results <- renderPrint({
      cat("=== UJI NORMALITAS ===\n")
      cat("Variabel:", input$normality_var, "\n")
      cat("Jumlah observasi:", length(var_data), "\n\n")
      
      # Shapiro-Wilk test
      if (length(var_data) <= 5000) {
        shapiro_result <- shapiro.test(var_data)
        cat("Shapiro-Wilk Test:\n")
        cat("W =", round(shapiro_result$statistic, 4), "\n")
        cat("p-value =", format.pval(shapiro_result$p.value), "\n")
        cat("Interpretasi:", ifelse(shapiro_result$p.value > 0.05, "Data berdistribusi normal", "Data tidak berdistribusi normal"), "\n\n")
      }
      
      # Anderson-Darling test
      ad_result <- ad.test(var_data)
      cat("Anderson-Darling Test:\n")
      cat("A =", round(ad_result$statistic, 4), "\n")
      cat("p-value =", format.pval(ad_result$p.value), "\n")
      cat("Interpretasi:", ifelse(ad_result$p.value > 0.05, "Data berdistribusi normal", "Data tidak berdistribusi normal"), "\n")
    })
    
    output$normality_plot <- renderPlot({
      par(mfrow = c(1, 2))
      
      # Q-Q plot
      qqnorm(var_data, main = paste("Q-Q Plot:", input$normality_var), 
             col = "#1976D2", pch = 16)
      qqline(var_data, col = "#F44336", lwd = 2)
      
      # Histogram with normal curve
      hist(var_data, freq = FALSE, main = paste("Histogram:", input$normality_var),
           xlab = input$normality_var, col = "#2196F3", border = "white")
      curve(dnorm(x, mean = mean(var_data), sd = sd(var_data)),
             add = TRUE, col = "#F44336", lwd = 2)
    })
  })
  
  # Download normality plots
  output$download_normality_png <- downloadHandler(
    filename = function() {
      paste0("normality_", input$normality_var, "_", Sys.Date(), ".png")
    },
    content = function(file) {
      var_data <- sovi_data[[input$normality_var]]
      var_data <- var_data[!is.na(var_data)]
      
      png(file, width = 1000, height = 500)
      par(mfrow = c(1, 2))
      
      # Q-Q plot
      qqnorm(var_data, main = paste("Q-Q Plot:", input$normality_var), 
             col = "#1976D2", pch = 16)
      qqline(var_data, col = "#F44336", lwd = 2)
      
      # Histogram with normal curve
      hist(var_data, freq = FALSE, main = paste("Histogram:", input$normality_var),
           xlab = input$normality_var, col = "#2196F3", border = "white")
      curve(dnorm(x, mean = mean(var_data), sd = sd(var_data)),
             add = TRUE, col = "#F44336", lwd = 2)
      
      dev.off()
    }
  )
  
  output$download_normality_pdf <- downloadHandler(
    filename = function() {
      paste0("normality_", input$normality_var, "_", Sys.Date(), ".pdf")
    },
    content = function(file) {
      var_data <- sovi_data[[input$normality_var]]
      var_data <- var_data[!is.na(var_data)]
      
      pdf(file, width = 10, height = 5)
      par(mfrow = c(1, 2))
      
      # Q-Q plot
      qqnorm(var_data, main = paste("Q-Q Plot:", input$normality_var), 
             col = "#1976D2", pch = 16)
      qqline(var_data, col = "#F44336", lwd = 2)
      
      # Histogram with normal curve
      hist(var_data, freq = FALSE, main = paste("Histogram:", input$normality_var),
           xlab = input$normality_var, col = "#2196F3", border = "white")
      curve(dnorm(x, mean = mean(var_data), sd = sd(var_data)),
             add = TRUE, col = "#F44336", lwd = 2)
      
      dev.off()
    }
  )
  
  # Uji Homogenitas
  observeEvent(input$run_homogeneity, {
    req(input$homogeneity_var, input$homogeneity_group)
    
    # Ensure the group variable is a factor
    sovi_data[[input$homogeneity_group]] <- as.factor(sovi_data[[input$homogeneity_group]])
    
    # Create formula for leveneTest
    formula_str <- paste(input$homogeneity_var, "~", input$homogeneity_group)
    
    # Run the test
    tryCatch({
      levene_result <- leveneTest(as.formula(formula_str), data = sovi_data)
      
      # Simpan hasil untuk laporan
      values$homogeneity_result <- list(
        variable = input$homogeneity_var,
        group = input$homogeneity_group,
        levene = levene_result
      )
      
      output$homogeneity_results <- renderPrint({
        cat("=== UJI HOMOGENITAS VARIANS ===\n")
        cat("Variabel:", input$homogeneity_var, "\n")
        cat("Kelompok:", input$homogeneity_group, "\n\n")
        
        # Levene's test
        cat("Levene's Test:\n")
        print(levene_result)
        
        p_value <- levene_result$`Pr(>F)`[1]
        cat("\nInterpretasi:", ifelse(p_value > 0.05, "Varians homogen", "Varians tidak homogen"), "\n")
      })
      
      # Create boxplot for visualization
      output$homogeneity_plot <- renderPlot({
        ggplot(sovi_data, aes_string(x = input$homogeneity_group, y = input$homogeneity_var)) +
          geom_boxplot(fill = "#2196F3", alpha = 0.7) +
          labs(title = paste("Boxplot by Group:", input$homogeneity_var, "by", input$homogeneity_group),
               x = input$homogeneity_group, y = input$homogeneity_var) +
          theme_minimal()
      })
      
    }, error = function(e) {
      output$homogeneity_results <- renderPrint({
        cat("Error dalam uji homogenitas:", e$message, "\n")
        cat("Pastikan variabel grup memiliki setidaknya dua level dan tidak ada nilai NA.")
      })
    })
  })
  
  # Download homogeneity plots
  output$download_homogeneity_png <- downloadHandler(
    filename = function() {
      paste0("homogeneity_", input$homogeneity_var, "_by_", input$homogeneity_group, "_", Sys.Date(), ".png")
    },
    content = function(file) {
      p <- ggplot(sovi_data, aes_string(x = input$homogeneity_group, y = input$homogeneity_var)) +
        geom_boxplot(fill = "#2196F3", alpha = 0.7) +
        labs(title = paste("Boxplot by Group:", input$homogeneity_var, "by", input$homogeneity_group),
             x = input$homogeneity_group, y = input$homogeneity_var) +
        theme_minimal()
      
      ggsave(file, plot = p, width = 10, height = 6, dpi = 300)
    }
  )
  
  output$download_homogeneity_pdf <- downloadHandler(
    filename = function() {
      paste0("homogeneity_", input$homogeneity_var, "_by_", input$homogeneity_group, "_", Sys.Date(), ".pdf")
    },
    content = function(file) {
      p <- ggplot(sovi_data, aes_string(x = input$homogeneity_group, y = input$homogeneity_var)) +
        geom_boxplot(fill = "#2196F3", alpha = 0.7) +
        labs(title = paste("Boxplot by Group:", input$homogeneity_var, "by", input$homogeneity_group),
             x = input$homogeneity_group, y = input$homogeneity_var) +
        theme_minimal()
      
      ggsave(file, plot = p, width = 10, height = 6, device = "pdf")
    }
  )
  
  # ===== UJI BEDA RATA-RATA =====
  observeEvent(input$run_ttest, {
    req(input$ttest_var, input$ttest_type)
    
    if (input$ttest_type == "one") {
      # One sample t-test
      req(input$mu_value)
      
      tryCatch({
        t_test_result <- t.test(sovi_data[[input$ttest_var]], mu = input$mu_value)
        
        # Save result for report
        values$ttest_result <- list(
          type = "one",
          variable = input$ttest_var,
          mu = input$mu_value,
          result = t_test_result,
          alpha = input$alpha_level
        )
        
        output$ttest_results <- renderPrint({
          cat("=== UJI T-TEST SATU SAMPEL ===\n")
          cat("Variabel:", input$ttest_var, "\n")
          cat("Nilai μ₀:", input$mu_value, "\n")
          cat("Alpha:", input$alpha_level, "\n\n")
          
          print(t_test_result)
          
          cat("\nInterpretasi:\n")
          if (t_test_result$p.value < input$alpha_level) {
            cat("Tolak H₀: Terdapat perbedaan signifikan antara rata-rata sampel dan nilai hipotesis μ₀.\n")
          } else {
            cat("Terima H₀: Tidak terdapat perbedaan signifikan antara rata-rata sampel dan nilai hipotesis μ₀.\n")
          }
        })
        
        output$ttest_plot <- renderPlot({
          var_data <- sovi_data[[input$ttest_var]]
          var_data <- var_data[!is.na(var_data)]
          
          # Create histogram with vertical lines for mean and mu
          ggplot(data.frame(x = var_data), aes(x = x)) +
            geom_histogram(bins = 30, fill = "#2196F3", color = "white", alpha = 0.7) +
            geom_vline(xintercept = mean(var_data), color = "#F44336", size = 1.5) +
            geom_vline(xintercept = input$mu_value, color = "#4CAF50", size = 1.5, linetype = "dashed") +
            labs(title = paste("Histogram with Mean and μ₀:", input$ttest_var),
                 x = input$ttest_var, y = "Frequency") +
            annotate("text", x = mean(var_data), y = 0, label = "Sample Mean", 
                     color = "#F44336", vjust = -1, hjust = 1.1) +
            annotate("text", x = input$mu_value, y = 0, label = "μ₀", 
                     color = "#4CAF50", vjust = -1, hjust = -0.1) +
            theme_minimal()
        })
        
      }, error = function(e) {
        output$ttest_results <- renderPrint({
          cat("Error dalam uji t satu sampel:", e$message, "\n")
        })
      })
      
    } else {
      # Two sample t-test
      req(input$ttest_group)
      
      # Ensure the group variable is a factor
      sovi_data[[input$ttest_group]] <- as.factor(sovi_data[[input$ttest_group]])
      
      tryCatch({
        # Check if there are exactly two levels
        if (length(levels(sovi_data[[input$ttest_group]])) != 2) {
          output$ttest_results <- renderPrint({
            cat("Error: Variabel grup harus memiliki tepat dua level.\n")
            cat("Level yang ditemukan:", paste(levels(sovi_data[[input$ttest_group]]), collapse = ", "))
          })
          return()
        }
        
        # Create formula for t.test
        formula_str <- paste(input$ttest_var, "~", input$ttest_group)
        
        # Run the test
        t_test_result <- t.test(as.formula(formula_str), 
                               data = sovi_data,
                               var.equal = TRUE)
        
        # Save result for report
        values$ttest_result <- list(
          type = "two",
          variable = input$ttest_var,
          group = input$ttest_group,
          result = t_test_result,
          alpha = input$alpha_level
        )
        
        output$ttest_results <- renderPrint({
          cat("=== UJI T-TEST DUA SAMPEL ===\n")
          cat("Variabel:", input$ttest_var, "\n")
          cat("Grup:", input$ttest_group, "\n")
          cat("Alpha:", input$alpha_level, "\n\n")
          
          print(t_test_result)
          
          cat("\nInterpretasi:\n")
          if (t_test_result$p.value < input$alpha_level) {
            cat("Tolak H₀: Terdapat perbedaan signifikan antara rata-rata kedua grup.\n")
          } else {
            cat("Terima H₀: Tidak terdapat perbedaan signifikan antara rata-rata kedua grup.\n")
          }
        })
        
        output$ttest_plot <- renderPlot({
          # Create boxplot for visualization
          ggplot(sovi_data, aes_string(x = input$ttest_group, y = input$ttest_var, fill = input$ttest_group)) +
            geom_boxplot(alpha = 0.7) +
            scale_fill_manual(values = c("#2196F3", "#F44336")) +
            labs(title = paste("Boxplot by Group:", input$ttest_var, "by", input$ttest_group),
                 x = input$ttest_group, y = input$ttest_var) +
            theme_minimal() +
            theme(legend.position = "none")
        })
        
      }, error = function(e) {
        output$ttest_results <- renderPrint({
          cat("Error dalam uji t dua sampel:", e$message, "\n")
        })
      })
    }
  })
  
  # Download t-test plots
  output$download_ttest_png <- downloadHandler(
    filename = function() {
      if (input$ttest_type == "one") {
        paste0("ttest_one_sample_", input$ttest_var, "_", Sys.Date(), ".png")
      } else {
        paste0("ttest_two_sample_", input$ttest_var, "_by_", input$ttest_group, "_", Sys.Date(), ".png")
      }
    },
    content = function(file) {
      if (input$ttest_type == "one") {
        var_data <- sovi_data[[input$ttest_var]]
        var_data <- var_data[!is.na(var_data)]
        
        p <- ggplot(data.frame(x = var_data), aes(x = x)) +
          geom_histogram(bins = 30, fill = "#2196F3", color = "white", alpha = 0.7) +
          geom_vline(xintercept = mean(var_data), color = "#F44336", size = 1.5) +
          geom_vline(xintercept = input$mu_value, color = "#4CAF50", size = 1.5, linetype = "dashed") +
          labs(title = paste("Histogram with Mean and μ₀:", input$ttest_var),
               x = input$ttest_var, y = "Frequency") +
          annotate("text", x = mean(var_data), y = 0, label = "Sample Mean", 
                   color = "#F44336", vjust = -1, hjust = 1.1) +
          annotate("text", x = input$mu_value, y = 0, label = "μ₀", 
                   color = "#4CAF50", vjust = -1, hjust = -0.1) +
          theme_minimal()
      } else {
        p <- ggplot(sovi_data, aes_string(x = input$ttest_group, y = input$ttest_var, fill = input$ttest_group)) +
          geom_boxplot(alpha = 0.7) +
          scale_fill_manual(values = c("#2196F3", "#F44336")) +
          labs(title = paste("Boxplot by Group:", input$ttest_var, "by", input$ttest_group),
               x = input$ttest_group, y = input$ttest_var) +
          theme_minimal() +
          theme(legend.position = "none")
      }
      
      ggsave(file, plot = p, width = 10, height = 6, dpi = 300)
    }
  )
  
  output$download_ttest_pdf <- downloadHandler(
    filename = function() {
      if (input$ttest_type == "one") {
        paste0("ttest_one_sample_", input$ttest_var, "_", Sys.Date(), ".pdf")
      } else {
        paste0("ttest_two_sample_", input$ttest_var, "_by_", input$ttest_group, "_", Sys.Date(), ".pdf")
      }
    },
    content = function(file) {
      if (input$ttest_type == "one") {
        var_data <- sovi_data[[input$ttest_var]]
        var_data <- var_data[!is.na(var_data)]
        
        p <- ggplot(data.frame(x = var_data), aes(x = x)) +
          geom_histogram(bins = 30, fill = "#2196F3", color = "white", alpha = 0.7) +
          geom_vline(xintercept = mean(var_data), color = "#F44336", size = 1.5) +
          geom_vline(xintercept = input$mu_value, color = "#4CAF50", size = 1.5, linetype = "dashed") +
          labs(title = paste("Histogram with Mean and μ₀:", input$ttest_var),
               x = input$ttest_var, y = "Frequency") +
          annotate("text", x = mean(var_data), y = 0, label = "Sample Mean", 
                   color = "#F44336", vjust = -1, hjust = 1.1) +
          annotate("text", x = input$mu_value, y = 0, label = "μ₀", 
                   color = "#4CAF50", vjust = -1, hjust = -0.1) +
          theme_minimal()
      } else {
        p <- ggplot(sovi_data, aes_string(x = input$ttest_group, y = input$ttest_var, fill = input$ttest_group)) +
          geom_boxplot(alpha = 0.7) +
          scale_fill_manual(values = c("#2196F3", "#F44336")) +
          labs(title = paste("Boxplot by Group:", input$ttest_var, "by", input$ttest_group),
               x = input$ttest_group, y = input$ttest_var) +
          theme_minimal() +
          theme(legend.position = "none")
      }
      
      ggsave(file, plot = p, width = 10, height = 6, device = "pdf")
    }
  )
  
  # ===== UJI PROPORSI & VARIANS =====
  observeEvent(input$run_prop_test, {
    req(input$prop_var, input$prop_value)
    
    # Ensure the variable is a factor
    sovi_data[[input$prop_var]] <- as.factor(sovi_data[[input$prop_var]])
    
    tryCatch({
      # Get counts for each level
      table_counts <- table(sovi_data[[input$prop_var]])
      
      # Use the first level as "success"
      success_count <- table_counts[1]
      total_count <- sum(table_counts)
      
      # Run the test
      prop_test_result <- prop.test(success_count, total_count, p = input$prop_value)
      
      output$prop_test_results <- renderPrint({
        cat("=== UJI PROPORSI SATU SAMPEL ===\n")
        cat("Variabel:", input$prop_var, "\n")
        cat("Level yang diuji:", names(table_counts)[1], "\n")
        cat("Jumlah sukses:", success_count, "\n")
        cat("Total observasi:", total_count, "\n")
        cat("Proporsi sampel:", round(success_count/total_count, 4), "\n")
        cat("Proporsi hipotesis (p₀):", input$prop_value, "\n\n")
        
        print(prop_test_result)
        
        cat("\nInterpretasi:\n")
        if (prop_test_result$p.value < 0.05) {
          cat("Tolak H₀: Proporsi sampel berbeda signifikan dari proporsi hipotesis.\n")
        } else {
          cat("Terima H₀: Proporsi sampel tidak berbeda signifikan dari proporsi hipotesis.\n")
        }
      })
    }, error = function(e) {
      output$prop_test_results <- renderPrint({
        cat("Error dalam uji proporsi:", e$message, "\n")
      })
    })
  })
  
  observeEvent(input$run_var_test, {
    req(input$var_test_var1, input$var_test_var2)
    
    tryCatch({
      # Get the data
      var1_data <- sovi_data[[input$var_test_var1]]
      var2_data <- sovi_data[[input$var_test_var2]]
      
      # Remove NA values
      var1_data <- var1_data[!is.na(var1_data)]
      var2_data <- var2_data[!is.na(var2_data)]
      
      # Run the test
      var_test_result <- var.test(var1_data, var2_data)
      
      output$var_test_results <- renderPrint({
        cat("=== UJI VARIANS DUA SAMPEL (F-TEST) ===\n")
        cat("Variabel 1:", input$var_test_var1, "\n")
        cat("Variabel 2:", input$var_test_var2, "\n")
        cat("Varians 1:", round(var(var1_data, na.rm = TRUE), 4), "\n")
        cat("Varians 2:", round(var(var2_data, na.rm = TRUE), 4), "\n\n")
        
        print(var_test_result)
        
        cat("\nInterpretasi:\n")
        if (var_test_result$p.value < 0.05) {
          cat("Tolak H₀: Varians kedua variabel berbeda signifikan.\n")
        } else {
          cat("Terima H₀: Varians kedua variabel tidak berbeda signifikan.\n")
        }
      })
      
      # Create boxplot for visualization
      output$var_test_plot <- renderPlot({
        # Combine data for plotting
        plot_data <- data.frame(
          value = c(var1_data, var2_data),
          variable = c(rep(input$var_test_var1, length(var1_data)),
                      rep(input$var_test_var2, length(var2_data)))
        )
        
        ggplot(plot_data, aes(x = variable, y = value, fill = variable)) +
          geom_boxplot(alpha = 0.7) +
          scale_fill_manual(values = c("#2196F3", "#F44336")) +
          labs(title = "Boxplot Comparison of Variances",
               x = "Variable", y = "Value") +
          theme_minimal() +
          theme(legend.position = "none")
      })
      
    }, error = function(e) {
      output$var_test_results <- renderPrint({
        cat("Error dalam uji varians:", e$message, "\n")
      })
    })
  })
  
  # Download variance test plots
  output$download_var_test_png <- downloadHandler(
    filename = function() {
      paste0("var_test_", input$var_test_var1, "_vs_", input$var_test_var2, "_", Sys.Date(), ".png")
    },
    content = function(file) {
      var1_data <- sovi_data[[input$var_test_var1]]
      var2_data <- sovi_data[[input$var_test_var2]]
      
      # Remove NA values
      var1_data <- var1_data[!is.na(var1_data)]
      var2_data <- var2_data[!is.na(var2_data)]
      
      # Combine data for plotting
      plot_data <- data.frame(
        value = c(var1_data, var2_data),
        variable = c(rep(input$var_test_var1, length(var1_data)),
                    rep(input$var_test_var2, length(var2_data)))
      )
      
      p <- ggplot(plot_data, aes(x = variable, y = value, fill = variable)) +
        geom_boxplot(alpha = 0.7) +
        scale_fill_manual(values = c("#2196F3", "#F44336")) +
        labs(title = "Boxplot Comparison of Variances",
             x = "Variable", y = "Value") +
        theme_minimal() +
        theme(legend.position = "none")
      
      ggsave(file, plot = p, width = 10, height = 6, dpi = 300)
    }
  )
  
  output$download_var_test_pdf <- downloadHandler(
    filename = function() {
      paste0("var_test_", input$var_test_var1, "_vs_", input$var_test_var2, "_", Sys.Date(), ".pdf")
    },
    content = function(file) {
      var1_data <- sovi_data[[input$var_test_var1]]
      var2_data <- sovi_data[[input$var_test_var2]]
      
      # Remove NA values
      var1_data <- var1_data[!is.na(var1_data)]
      var2_data <- var2_data[!is.na(var2_data)]
      
      # Combine data for plotting
      plot_data <- data.frame(
        value = c(var1_data, var2_data),
        variable = c(rep(input$var_test_var1, length(var1_data)),
                    rep(input$var_test_var2, length(var2_data)))
      )
      
      p <- ggplot(plot_data, aes(x = variable, y = value, fill = variable)) +
        geom_boxplot(alpha = 0.7) +
        scale_fill_manual(values = c("#2196F3", "#F44336")) +
        labs(title = "Boxplot Comparison of Variances",
             x = "Variable", y = "Value") +
        theme_minimal() +
        theme(legend.position = "none")
      
      ggsave(file, plot = p, width = 10, height = 6, device = "pdf")
    }
  )
  
  # ===== UJI ANOVA =====
  observeEvent(input$run_anova, {
    req(input$anova_dependent, input$anova_independent)
    
    # Ensure the independent variable is a factor
    sovi_data[[input$anova_independent]] <- as.factor(sovi_data[[input$anova_independent]])
    
    tryCatch({
      # Create formula for aov
      formula_str <- paste(input$anova_dependent, "~", input$anova_independent)
      
      # Run ANOVA
      anova_model <- aov(as.formula(formula_str), data = sovi_data)
      anova_summary <- summary(anova_model)
      
      # Save result for report
      values$anova_result <- list(
        dependent = input$anova_dependent,
        independent = input$anova_independent,
        model = anova_model,
        summary = anova_summary
      )
      
      output$anova_results <- renderPrint({
        cat("=== UJI ANOVA SATU ARAH ===\n")
        cat("Variabel Dependen:", input$anova_dependent, "\n")
        cat("Variabel Independen (Faktor):", input$anova_independent, "\n\n")
        
        print(anova_summary)
        
        # Get p-value
        p_value <- anova_summary[[1]][["Pr(>F)"]][1]
        
        cat("\nInterpretasi:\n")
        if (!is.na(p_value) && p_value < 0.05) {
          cat("Tolak H₀: Terdapat perbedaan signifikan antara rata-rata grup.\n")
          
          # Add Tukey HSD test if ANOVA is significant
          cat("\nUji Post-hoc (Tukey HSD):\n")
          tukey_result <- TukeyHSD(anova_model)
          print(tukey_result)
        } else {
          cat("Terima H₀: Tidak terdapat perbedaan signifikan antara rata-rata grup.\n")
        }
      })
      
      output$anova_plot <- renderPlot({
        # Create boxplot for visualization
        ggplot(sovi_data, aes_string(x = input$anova_independent, y = input$anova_dependent, fill = input$anova_independent)) +
          geom_boxplot(alpha = 0.7) +
          scale_fill_brewer(palette = "Blues") +
          labs(title = paste("Boxplot by Group:", input$anova_dependent, "by", input$anova_independent),
               x = input$anova_independent, y = input$anova_dependent) +
          theme_minimal() +
          theme(legend.position = "none")
      })
      
    }, error = function(e) {
      output$anova_results <- renderPrint({
        cat("Error dalam uji ANOVA:", e$message, "\n")
      })
    })
  })
  
  # Download ANOVA plots
  output$download_anova_png <- downloadHandler(
    filename = function() {
      paste0("anova_", input$anova_dependent, "_by_", input$anova_independent, "_", Sys.Date(), ".png")
    },
    content = function(file) {
      p <- ggplot(sovi_data, aes_string(x = input$anova_independent, y = input$anova_dependent, fill = input$anova_independent)) +
        geom_boxplot(alpha = 0.7) +
        scale_fill_brewer(palette = "Blues") +
        labs(title = paste("Boxplot by Group:", input$anova_dependent, "by", input$anova_independent),
             x = input$anova_independent, y = input$anova_dependent) +
        theme_minimal() +
        theme(legend.position = "none")
      
      ggsave(file, plot = p, width = 10, height = 6, dpi = 300)
    }
  )
  
  output$download_anova_pdf <- downloadHandler(
    filename = function() {
      paste0("anova_", input$anova_dependent, "_by_", input$anova_independent, "_", Sys.Date(), ".pdf")
    },
    content = function(file) {
      p <- ggplot(sovi_data, aes_string(x = input$anova_independent, y = input$anova_dependent, fill = input$anova_independent)) +
        geom_boxplot(alpha = 0.7) +
        scale_fill_brewer(palette = "Blues") +
        labs(title = paste("Boxplot by Group:", input$anova_dependent, "by", input$anova_independent),
             x = input$anova_independent, y = input$anova_dependent) +
        theme_minimal() +
        theme(legend.position = "none")
      
      ggsave(file, plot = p, width = 10, height = 6, device = "pdf")
    }
  )
  
  # ===== REGRESI LINEAR =====
  observeEvent(input$run_regression, {
    req(input$dependent_var, input$independent_vars)
    
    tryCatch({
      # Create formula
      formula_str <- paste(input$dependent_var, "~", paste(input$independent_vars, collapse = " + "))
      
      # Fit model
      model <- lm(as.formula(formula_str), data = sovi_data)
      
      # Save model for later use
      values$current_model <- model
      
      output$regression_summary <- renderPrint({
        cat("=== HASIL REGRESI LINEAR BERGANDA ===\n\n")
        print(summary(model))
        
        # Additional statistics
        cat("\n=== STATISTIK TAMBAHAN ===\n")
        cat("R-squared:", round(summary(model)$r.squared, 4), "\n")
        cat("Adjusted R-squared:", round(summary(model)$adj.r.squared, 4), "\n")
        cat("F-statistic:", round(summary(model)$fstatistic[1], 4), 
            "on", summary(model)$fstatistic[2], "and", summary(model)$fstatistic[3], "DF\n")
        
        # Calculate p-value for F-statistic
        f_p_value <- pf(summary(model)$fstatistic[1], 
                        summary(model)$fstatistic[2], 
                        summary(model)$fstatistic[3], 
                        lower.tail = FALSE)
        cat("p-value:", format.pval(f_p_value), "\n")
        
        # Multicollinearity check
        if (length(input$independent_vars) > 1) {
          cat("\n=== UJI MULTIKOLINEARITAS ===\n")
          vif_values <- vif(model)
          print(vif_values)
          
          cat("\nInterpretasi VIF:\n")
          cat("- VIF < 5: Tidak ada multikolinearitas yang signifikan\n")
          cat("- 5 <= VIF < 10: Multikolinearitas moderat\n")
          cat("- VIF >= 10: Multikolinearitas tinggi\n")
        }
        
        # Durbin-Watson test for autocorrelation
        cat("\n=== UJI AUTOKORELASI ===\n")
        dw_test <- durbinWatsonTest(model)
        cat("Durbin-Watson statistic:", round(dw_test$dw, 4), "\n")
        cat("p-value:", format.pval(dw_test$p), "\n")
        
        if (dw_test$p < 0.05) {
          cat("Interpretasi: Terdapat autokorelasi dalam residual.\n")
        } else {
          cat("Interpretasi: Tidak terdapat autokorelasi dalam residual.\n")
        }
        
        # Breusch-Pagan test for heteroscedasticity
        cat("\n=== UJI HETEROSKEDASTISITAS ===\n")
        bp_test <- bptest(model)
        cat("Breusch-Pagan test statistic:", round(bp_test$statistic, 4), "\n")
        cat("p-value:", format.pval(bp_test$p.value), "\n")
        
        if (bp_test$p.value < 0.05) {
          cat("Interpretasi: Terdapat heteroskedastisitas dalam residual.\n")
        } else {
          cat("Interpretasi: Tidak terdapat heteroskedastisitas dalam residual (homoskedastisitas).\n")
        }
      })
      
      output$regression_diagnostics <- renderPlot({
        par(mfrow = c(2, 2))
        plot(model, col = "#1976D2", pch = 16)
      })
      
      output$regression_prediction <- renderPlotly({
        # Create data frame with actual and predicted values
        pred_data <- data.frame(
          Actual = sovi_data[[input$dependent_var]],
          Predicted = fitted(model),
          Residuals = residuals(model)
        )
        
        # Create scatter plot
        p <- ggplot(pred_data, aes(x = Actual, y = Predicted)) +
          geom_point(color = "#1976D2", alpha = 0.6) +
          geom_abline(intercept = 0, slope = 1, color = "#F44336", linetype = "dashed") +
          labs(title = "Actual vs Predicted Values",
               x = "Actual Values", y = "Predicted Values") +
          theme_minimal()
        
        ggplotly(p)
      })
      
    }, error = function(e) {
      output$regression_summary <- renderPrint({
        cat("Error dalam analisis regresi:", e$message, "\n")
      })
    })
  })
  
  # Download regression plots
  output$download_regression_png <- downloadHandler(
    filename = function() {
      paste0("regression_diagnostics_", input$dependent_var, "_", Sys.Date(), ".png")
    },
    content = function(file) {
      req(values$current_model)
      
      png(file, width = 1000, height = 1000)
      par(mfrow = c(2, 2))
      plot(values$current_model, col = "#1976D2", pch = 16)
      dev.off()
    }
  )
  
  output$download_regression_pdf <- downloadHandler(
    filename = function() {
      paste0("regression_diagnostics_", input$dependent_var, "_", Sys.Date(), ".pdf")
    },
    content = function(file) {
      req(values$current_model)
      
      pdf(file, width = 10, height = 10)
      par(mfrow = c(2, 2))
      plot(values$current_model, col = "#1976D2", pch = 16)
      dev.off()
    }
  )
  
  # ===== DOWNLOAD BINNED DATA =====
  output$download_binned_data <- downloadHandler(
    filename = function() {
      paste0("binned_sovi_data_", input$bin_variable, "_", Sys.Date(), ".csv")
    },
    content = function(file) {
      req(values$binned_data)
      write_csv(values$binned_data, file)
    }
  )
  
  # ===== GENERATE REPORT =====
  output$download_report <- downloadHandler(
    filename = function() {
      paste0("SoVI_Analysis_Report_", Sys.Date(), ".pdf")
    },
    content = function(file) {
      # Create temporary Rmd file
      tempReport <- file.path(tempdir(), "report.Rmd")
      
      # Generate report content
      report_content <- generate_report_content(input, values, sovi_data)
      
      # Write Rmd file
      writeLines(report_content, tempReport)
      
      # Render PDF dengan output format yang eksplisit
      tryCatch({
        rmarkdown::render(
          tempReport, 
          output_format = rmarkdown::pdf_document(
            toc = TRUE,
            toc_depth = 3,
            number_sections = TRUE,
            fig_caption = TRUE,
            latex_engine = "pdflatex"
          ),
          output_file = file,
          params = list(
            title = input$report_title,
            author = input$report_author,
            institution = input$report_institution,
            sections = input$report_sections
          ),
          envir = new.env(parent = globalenv()),
          quiet = TRUE
        )
        showNotification("✅ Laporan PDF berhasil dibuat!", type = "success")
      }, error = function(e) {
        showNotification(paste("❌ Error membuat PDF:", e$message), type = "error")
        # Fallback ke HTML jika PDF gagal
        rmarkdown::render(
          tempReport,
          output_format = rmarkdown::html_document(),
          output_file = gsub("\\.pdf$", ".html", file),
          envir = new.env(parent = globalenv()),
          quiet = TRUE
        )
      })
    }
  )
}

# ===== FUNGSI GENERATE REPORT =====
generate_report_content <- function(input, values, sovi_data) {
  content <- c(
    "---",
    paste("title: '", input$report_title, "'", sep = ""),
    paste("author: '", input$report_author, "'", sep = ""),
    paste("date: '", Sys.Date(), "'", sep = ""),
    "output:",
    "  pdf_document:",
    "    toc: true",
    "    toc_depth: 3",
    "    number_sections: true",
    "    fig_caption: true",
    "    latex_engine: pdflatex",
    "header-includes:",
    "  - \\usepackage{float}",
    "  - \\usepackage{booktabs}",
    "  - \\usepackage{longtable}",
    "  - \\usepackage{array}",
    "  - \\usepackage{multirow}",
    "  - \\usepackage{wrapfig}",
    "  - \\usepackage{colortbl}",
    "  - \\usepackage{pdflscape}",
    "  - \\usepackage{tabu}",
    "  - \\usepackage{threeparttable}",
    "  - \\usepackage{threeparttablex}",
    "  - \\usepackage[normalem]{ulem}",
    "  - \\usepackage{makecell}",
    "  - \\usepackage{xcolor}",
    "---",
    "",
    "```{r setup, include=FALSE}",
    "knitr::opts_chunk$set(echo = FALSE, warning = FALSE, message = FALSE, fig.pos = 'H')",
    "library(knitr)",
    "library(ggplot2)",
    "library(dplyr)",
    "library(psych)",
    "```",
    "",
    "\\newpage",
    "",
    "# Pendahuluan",
    "",
    "Laporan ini menyajikan analisis komprehensif terhadap data **Social Vulnerability Index (SoVI)**. Analisis dilakukan menggunakan dashboard interaktif premium yang telah dikembangkan untuk memahami pola kerentanan sosial di berbagai wilayah dengan teknologi terdepan.",
    "",
    paste("**Institusi:** ", input$report_institution),
    paste("**Tanggal Analisis:** ", Sys.Date()),
    paste("**Total Observasi:** ", nrow(sovi_data)),
    paste("**Total Variabel:** ", ncol(sovi_data)),
    "",
    "\\newpage",
    "",
    "# Ringkasan Dataset",
    "",
    "Dataset SoVI yang dianalisis memiliki karakteristik sebagai berikut:",
    "",
    paste("- **Jumlah District:** ", length(unique(sovi_data$DISTRICTCODE))),
    paste("- **Rata-rata Populasi:** ", format(round(mean(sovi_data$POPULATION, na.rm = TRUE)), big.mark = ",")),
    paste("- **Rata-rata Tingkat Kemiskinan:** ", paste0(round(mean(sovi_data$POVERTY, na.rm = TRUE) * 100, 1), "%")),
    ""
  )
  
  # Add sections based on user selection
  if ("descriptive" %in% input$report_sections) {
    content <- c(content,
      "\\newpage",
      "",
      "# Statistik Deskriptif",
      "",
      "Berikut adalah statistik deskriptif untuk variabel-variabel utama dalam dataset:",
      "",
      "```{r descriptive-stats, results='asis'}",
      "library(psych)",
      "numeric_vars <- sovi_data[sapply(sovi_data, is.numeric)]",
      "desc_stats <- describe(numeric_vars)",
      "desc_table <- desc_stats[,c('n', 'mean', 'sd', 'median', 'min', 'max')]",
      "kable(desc_table, digits = 3, caption = 'Statistik Deskriptif Variabel Numerik',",
      "      booktabs = TRUE, longtable = TRUE) %>%",
      "  kable_styling(latex_options = c('striped', 'hold_position'))",
      "```",
      ""
    )
  }
  
  if ("normality" %in% input$report_sections && !is.null(values$normality_result)) {
    content <- c(content,
      "\\newpage",
      "",
      "# Uji Normalitas",
      "",
      paste("Uji normalitas dilakukan terhadap variabel **", values$normality_result$variable, "** dengan jumlah observasi ", values$normality_result$n, "."),
      "",
      "## Hasil Uji Shapiro-Wilk",
      ""
    )
    
    if (!is.null(values$normality_result$shapiro)) {
      shapiro <- values$normality_result$shapiro
      content <- c(content,
        paste("- **Statistik W:** ", round(shapiro$statistic, 4)),
        paste("- **p-value:** ", format.pval(shapiro$p.value)),
        paste("- **Interpretasi:** ", ifelse(shapiro$p.value > 0.05, "Data berdistribusi normal", "Data tidak berdistribusi normal")),
        ""
      )
    }
    
    content <- c(content,
      "## Hasil Uji Anderson-Darling",
      ""
    )
    
    if (!is.null(values$normality_result$ad)) {
      ad <- values$normality_result$ad
      content <- c(content,
        paste("- **Statistik A:** ", round(ad$statistic, 4)),
        paste("- **p-value:** ", format.pval(ad$p.value)),
        paste("- **Interpretasi:** ", ifelse(ad$p.value > 0.05, "Data berdistribusi normal", "Data tidak berdistribusi normal")),
        ""
      )
    }
  }
  
  if ("homogeneity" %in% input$report_sections && !is.null(values$homogeneity_result)) {
    content <- c(content,
      "\\newpage",
      "",
      "# Uji Homogenitas Varians",
      "",
      paste("Uji homogenitas varians dilakukan untuk menguji kesamaan varians variabel **",
             values$homogeneity_result$variable, "** antar kelompok **",
             values$homogeneity_result$group, "**."),
      "",
      "## Hasil Uji Levene",
      ""
    )
    
    levene <- values$homogeneity_result$levene
    p_value <- levene$`Pr(>F)`[1]
    content <- c(content,
      paste("- **F-statistik:** ", round(levene$`F value`[1], 4)),
      paste("- **p-value:** ", format.pval(p_value)),
      paste("- **Interpretasi:** ", ifelse(p_value > 0.05, "Varians homogen", "Varians tidak homogen")),
      ""
    )
  }
  
  if ("ttest" %in% input$report_sections && !is.null(values$ttest_result)) {
    content <- c(content,
      "\\newpage",
      "",
      "# Uji t-test",
      "",
      paste("Uji t dilakukan terhadap variabel **", values$ttest_result$variable, "**."),
      ""
    )
    
    if (values$ttest_result$type == "one") {
      content <- c(content,
        "## One Sample t-test",
        "",
        paste("- **Hipotesis:** H₀: μ =", values$ttest_result$mu),
        paste("- **t-statistik:** ", round(values$ttest_result$result$statistic, 4)),
        paste("- **p-value:** ", format.pval(values$ttest_result$result$p.value)),
        paste("- **Kesimpulan:** ", ifelse(values$ttest_result$result$p.value < values$ttest_result$alpha,
                                           "Tolak H₀ (ada perbedaan signifikan)",
                                           "Terima H₀ (tidak ada perbedaan signifikan)")),
        ""
      )
    } else {
      content <- c(content,
        "## Two Sample t-test",
        "",
        paste("- **Kelompok:** ", values$ttest_result$group),
        paste("- **t-statistik:** ", round(values$ttest_result$result$statistic, 4)),
        paste("- **p-value:** ", format.pval(values$ttest_result$result$p.value)),
        paste("- **Kesimpulan:** ", ifelse(values$ttest_result$result$p.value < values$ttest_result$alpha,
                                           "Tolak H₀ (ada perbedaan signifikan)",
                                           "Terima H₀ (tidak ada perbedaan signifikan)")),
        ""
      )
    }
  }
  
  if ("anova" %in% input$report_sections && !is.null(values$anova_result)) {
    content <- c(content,
      "\\newpage",
      "",
      "# Uji ANOVA",
      "",
      paste("Analisis varians (ANOVA) dilakukan untuk menguji perbedaan rata-rata variabel **",
             values$anova_result$dependent, "** antar kelompok **",
             values$anova_result$independent, "**."),
      "",
      "## Hasil ANOVA",
      "",
      paste("- **F-statistik:** ", round(values$anova_result$summary[[1]]$`F value`[1], 4)),
      paste("- **p-value:** ", format.pval(values$anova_result$summary[[1]]$`Pr(>F)`[1])),
      paste("- **Kesimpulan:** ", ifelse(values$anova_result$summary[[1]]$`Pr(>F)`[1] < 0.05,
                                         "Tolak H₀ (ada perbedaan rata-rata antar kelompok)",
                                         "Terima H₀ (tidak ada perbedaan rata-rata antar kelompok)")),
      ""
    )
    
    if (values$anova_result$summary[[1]]$`Pr(>F)`[1] < 0.05) {
      content <- c(content,
        "## Post-hoc Test (Tukey HSD)",
        "",
        "Karena hasil ANOVA signifikan, dilakukan uji lanjut Tukey HSD untuk mengetahui kelompok mana yang berbeda secara signifikan.",
        ""
      )
    }
  }
  
  if ("regression" %in% input$report_sections && !is.null(values$current_model)) {
    content <- c(content,
      "\\newpage",
      "",
      "# Analisis Regresi Linear",
      "",
      "Analisis regresi linear berganda dilakukan untuk memodelkan hubungan antar variabel.",
      "",
      "```{r regression-summary, results='asis'}",
      "if (!is.null(values$current_model)) {",
      "  summary_model <- summary(values$current_model)",
      "  cat('**R-squared:** ', round(summary_model$r.squared, 4), '\\n\\n')",
      "  cat('**Adjusted R-squared:** ', round(summary_model$adj.r.squared, 4), '\\n\\n')",
      "  cat('**F-statistic:** ', round(summary_model$fstatistic[1], 4), '\\n\\n')",
      "  f_pvalue <- pf(summary_model$fstatistic[1], summary_model$fstatistic[2], summary_model$fstatistic[3], lower.tail = FALSE)",
      "  cat('**p-value:** ', format.pval(f_pvalue), '\\n\\n')",
      "}",
      "```",
      ""
    )
  }
  
  content <- c(content,
    "\\newpage",
    "",
    "# Kesimpulan",
    "",
    "Berdasarkan analisis yang telah dilakukan terhadap data Social Vulnerability Index (SoVI), dapat disimpulkan bahwa:",
    "",
    "1. Dataset memiliki karakteristik yang beragam dengan berbagai tingkat kerentanan sosial",
    "2. Uji statistik yang dilakukan memberikan wawasan mendalam tentang distribusi dan hubungan antar variabel",
    "3. Hasil analisis dapat digunakan sebagai dasar untuk pengambilan kebijakan terkait pengurangan kerentanan sosial",
    "4. Dashboard premium yang dikembangkan memberikan kemudahan dalam eksplorasi dan analisis data",
    "",
    "## Rekomendasi",
    "",
    "Berdasarkan hasil analisis, disarankan untuk:",
    "",
    "- Melakukan monitoring berkelanjutan terhadap indikator kerentanan sosial",
    "- Mengembangkan program intervensi yang tepat sasaran",
    "- Meningkatkan kualitas data dan sistem informasi",
    "",
    "---",
    "",
    "\\begin{center}",
    "\\textit{Laporan ini dibuat secara otomatis menggunakan Dashboard SoVI Premium}",
    "\\end{center}",
    "",
    "\\begin{center}",
    paste("\\textit{Kontak: ", input$report_author, " - ", input$report_institution, "}"),
    "\\end{center}"
  )
  
  return(content)
}

# ===== 6. JALANKAN APLIKASI =====
shinyApp(ui = ui, server = server)
## 
## Listening on http://127.0.0.1:8616