## ----include = FALSE---------------------------------------------------------- knitr::opts_chunk$set( collapse = TRUE, comment = "#>", warning = FALSE, message = FALSE, fig.width = 10, fig.height = 8 ) ## ----enhanced-ndvi, eval=FALSE------------------------------------------------ # # Enhanced NDVI with quality filtering # ndvi_enhanced <- calculate_ndvi_enhanced( # red_data = red_band, # nir_data = nir_band, # quality_filter = TRUE, # temporal_smoothing = FALSE, # Single date # verbose = TRUE # ) # # # Visualize enhanced NDVI # create_ndvi_map( # ndvi_data = ndvi_enhanced, # title = "Enhanced NDVI with Quality Filtering", # ndvi_classes = "none" # ) # # print("Enhanced NDVI Analysis:") # ndvi_values <- terra::values(ndvi_enhanced, mat = FALSE) # valid_ndvi <- ndvi_values[!is.na(ndvi_values)] # print(paste("Mean NDVI:", round(mean(valid_ndvi), 3))) # print(paste("NDVI range:", paste(round(range(valid_ndvi), 3), collapse = " - "))) ## ----multi-index-analysis, eval=FALSE----------------------------------------- # # Calculate multiple indices relevant to agriculture # agricultural_indices <- calculate_multiple_indices( # red = red_band, # nir = nir_band, # indices = c("NDVI", "SAVI", "MSAVI", "DVI", "RVI"), # output_stack = TRUE, # verbose = TRUE # ) # # print("Agricultural Indices Calculated:") # print(names(agricultural_indices)) # # # Visualize key indices # if (terra::nlyr(agricultural_indices) >= 2) { # # NDVI visualization # create_spatial_map( # spatial_data = agricultural_indices[[1]], # title = "NDVI for Agricultural Assessment", # color_scheme = "ndvi" # ) # # # SAVI visualization (soil-adjusted) # if ("SAVI" %in% names(agricultural_indices)) { # plot_raster_fast( # agricultural_indices[["SAVI"]], # title = "SAVI (Soil Adjusted Vegetation Index)", # color_scheme = "viridis" # ) # } # } ## ----yield-assessment, eval=FALSE--------------------------------------------- # # Extract yield analysis if available # if (!is.null(crop_analysis$analysis_results$yield_analysis)) { # yield_results <- crop_analysis$analysis_results$yield_analysis # # print("Yield Potential Analysis:") # if (!"error" %in% names(yield_results)) { # print(paste("Composite yield index:", round(yield_results$composite_yield_index, 3))) # print(paste("Yield potential class:", yield_results$yield_potential_class)) # print(paste("Classification confidence:", round(yield_results$classification_confidence, 3))) # # print("Index contributions to yield assessment:") # for (idx in names(yield_results$index_contributions)) { # contrib <- yield_results$index_contributions[[idx]] # print(paste(" ", idx, "- Normalized:", round(contrib$mean_normalized, 3), # "Raw mean:", round(contrib$raw_mean, 3))) # } # } # } # # # Create yield potential map visualization # yield_potential_raster <- agricultural_indices[[1]] # Use NDVI as proxy # terra::values(yield_potential_raster) <- ifelse(terra::values(yield_potential_raster) > 0.6, 3, # ifelse(terra::values(yield_potential_raster) > 0.4, 2, 1)) # names(yield_potential_raster) <- "Yield_Potential" # # plot_raster_fast( # yield_potential_raster, # title = "Yield Potential Classification", # color_scheme = "terrain" # ) ## ----crop-performance, eval=FALSE--------------------------------------------- # # Calculate comprehensive vegetation statistics # if (inherits(agricultural_indices, "SpatRaster")) { # veg_stats <- list() # # for (i in 1:terra::nlyr(agricultural_indices)) { # index_name <- names(agricultural_indices)[i] # values <- terra::values(agricultural_indices[[i]], mat = FALSE) # valid_values <- values[!is.na(values)] # # if (length(valid_values) > 0) { # veg_stats[[index_name]] <- list( # mean = mean(valid_values), # median = median(valid_values), # sd = sd(valid_values), # cv = sd(valid_values) / mean(valid_values), # min = min(valid_values), # max = max(valid_values), # q25 = quantile(valid_values, 0.25), # q75 = quantile(valid_values, 0.75) # ) # } # } # # print("Crop Performance Metrics:") # for (index in names(veg_stats)) { # stats <- veg_stats[[index]] # print(paste(index, "- Mean:", round(stats$mean, 3), # "CV:", round(stats$cv, 3), # "Range:", paste(round(c(stats$min, stats$max), 3), collapse = "-"))) # } # } ## ----field-analysis, eval=FALSE----------------------------------------------- # # Create sample field boundaries # field_1 <- sf::st_polygon(list(matrix(c( # -83.5, 40.0, -83.3, 40.0, -83.3, 40.2, -83.5, 40.2, -83.5, 40.0 # ), ncol = 2, byrow = TRUE))) # # field_2 <- sf::st_polygon(list(matrix(c( # -83.3, 40.0, -83.1, 40.0, -83.1, 40.2, -83.3, 40.2, -83.3, 40.0 # ), ncol = 2, byrow = TRUE))) # # fields_sf <- sf::st_sf( # field_id = c("Field_1", "Field_2"), # crop_type = c("Corn", "Soybeans"), # geometry = sf::st_sfc(field_1, field_2, crs = 4326) # ) # # # Extract vegetation statistics by field using spatial join # field_analysis <- universal_spatial_join( # source_data = agricultural_indices, # target_data = fields_sf, # method = "zonal", # summary_function = "mean", # verbose = TRUE # ) # # print("Field-Level Analysis Results:") # if (inherits(field_analysis, "sf")) { # print(sf::st_drop_geometry(field_analysis)) # } ## ----variable-rate, eval=FALSE------------------------------------------------ # # Calculate coefficient of variation for management zones # if (inherits(agricultural_indices, "SpatRaster") && terra::nlyr(agricultural_indices) > 0) { # # Use NDVI for variability analysis # ndvi_layer <- agricultural_indices[[1]] # # # Calculate local variability using focal statistics # local_cv <- terra::focal(ndvi_layer, w = matrix(1, 3, 3), # fun = function(x) sd(x, na.rm = TRUE) / mean(x, na.rm = TRUE)) # names(local_cv) <- "Local_CV" # # # Classify variability zones # variability_zones <- terra::classify(local_cv, # matrix(c(0, 0.1, 1, # Low variability # 0.1, 0.2, 2, # Medium variability # 0.2, 1, 3), # High variability # ncol = 3, byrow = TRUE)) # names(variability_zones) <- "Variability_Zones" # # # Visualize variability # plot_raster_fast( # local_cv, # title = "Spatial Variability (Coefficient of Variation)", # color_scheme = "plasma" # ) # # plot_raster_fast( # variability_zones, # title = "Management Zones (1=Low, 2=Medium, 3=High Variability)", # color_scheme = "categorical" # ) # } ## ----stress-monitoring, eval=FALSE-------------------------------------------- # # Create stress detection workflow # stress_indices <- c("NDVI", "SAVI", "DVI") # # # Calculate stress thresholds based on crop type # corn_thresholds <- list( # healthy = list(NDVI = c(0.6, 1.0), SAVI = c(0.4, 0.8)), # stressed = list(NDVI = c(0.3, 0.6), SAVI = c(0.2, 0.4)), # severely_stressed = list(NDVI = c(0.0, 0.3), SAVI = c(0.0, 0.2)) # ) # # # Apply stress classification # if (inherits(agricultural_indices, "SpatRaster") && "NDVI" %in% names(agricultural_indices)) { # ndvi_values <- terra::values(agricultural_indices[["NDVI"]], mat = FALSE) # # # Classify stress levels # stress_classification <- ifelse(ndvi_values >= 0.6, "Healthy", # ifelse(ndvi_values >= 0.3, "Moderate Stress", "Severe Stress")) # # # Create stress map # stress_raster <- agricultural_indices[["NDVI"]] # stress_numeric <- ifelse(ndvi_values >= 0.6, 3, # ifelse(ndvi_values >= 0.3, 2, 1)) # terra::values(stress_raster) <- stress_numeric # names(stress_raster) <- "Stress_Level" # # plot_raster_fast( # stress_raster, # title = "Crop Stress Classification (1=Severe, 2=Moderate, 3=Healthy)", # color_scheme = "terrain" # ) # # # Calculate stress statistics # stress_table <- table(stress_classification) # stress_percent <- round(prop.table(stress_table) * 100, 1) # # print("Crop Stress Distribution:") # for (level in names(stress_table)) { # print(paste(level, ":", stress_table[[level]], "pixels (", stress_percent[[level]], "%)")) # } # } ## ----early-warning, eval=FALSE------------------------------------------------ # # Define warning thresholds # warning_thresholds <- list( # drought_stress = 0.4, # NDVI below this indicates drought stress # disease_risk = 0.3, # NDVI below this indicates disease risk # optimal_growth = 0.7 # NDVI above this indicates optimal conditions # ) # # # Generate alerts # if (exists("ndvi_values")) { # alerts <- list() # # # Calculate percentages in each category # drought_pixels <- sum(ndvi_values < warning_thresholds$drought_stress, na.rm = TRUE) # disease_pixels <- sum(ndvi_values < warning_thresholds$disease_risk, na.rm = TRUE) # optimal_pixels <- sum(ndvi_values > warning_thresholds$optimal_growth, na.rm = TRUE) # total_pixels <- sum(!is.na(ndvi_values)) # # alerts$drought_risk <- round((drought_pixels / total_pixels) * 100, 1) # alerts$disease_risk <- round((disease_pixels / total_pixels) * 100, 1) # alerts$optimal_conditions <- round((optimal_pixels / total_pixels) * 100, 1) # # print("Agricultural Alert System:") # print(paste("Drought stress risk:", alerts$drought_risk, "% of field")) # print(paste("Disease risk:", alerts$disease_risk, "% of field")) # print(paste("Optimal conditions:", alerts$optimal_conditions, "% of field")) # # # Generate recommendations # if (alerts$drought_risk > 20) { # print("RECOMMENDATION: Consider irrigation scheduling") # } # if (alerts$disease_risk > 10) { # print("RECOMMENDATION: Increase disease monitoring") # } # if (alerts$optimal_conditions > 70) { # print("STATUS: Crop conditions are generally favorable") # } # } ## ----seasonal-monitoring, eval=FALSE------------------------------------------ # # Simulate time series data for growing season # growth_stages <- c("Planting", "V6", "V12", "R1", "R3", "R6") # ndvi_progression <- c(0.2, 0.4, 0.7, 0.8, 0.75, 0.6) # # # Create time series visualization concept # seasonal_data <- data.frame( # Stage = growth_stages, # NDVI = ndvi_progression, # DOY = c(120, 150, 180, 200, 220, 260) # Day of year # ) # # print("Seasonal NDVI Progression:") # print(seasonal_data) # # # Create conceptual growth curve # if (requireNamespace("ggplot2", quietly = TRUE)) { # library(ggplot2) # # growth_plot <- ggplot(seasonal_data, aes(x = DOY, y = NDVI)) + # geom_line(color = "darkgreen", size = 1.2) + # geom_point(color = "red", size = 3) + # geom_text(aes(label = Stage), vjust = -0.5) + # labs(title = "Typical Corn NDVI Progression", # x = "Day of Year", # y = "NDVI Value") + # theme_minimal() + # ylim(0, 1) # # print(growth_plot) # } ## ----harvest-timing, eval=FALSE----------------------------------------------- # # Define maturity indicators based on NDVI decline # maturity_thresholds <- list( # corn = list( # early_maturity = 0.7, # NDVI starts declining # harvest_ready = 0.5, # Optimal harvest window # post_harvest = 0.3 # Past optimal timing # ), # soybeans = list( # early_maturity = 0.6, # harvest_ready = 0.4, # post_harvest = 0.25 # ) # ) # # # Calculate harvest readiness # if (exists("ndvi_values")) { # # Use corn thresholds for demonstration # thresholds <- maturity_thresholds$corn # # harvest_assessment <- list( # early_maturity = sum(ndvi_values <= thresholds$early_maturity & # ndvi_values > thresholds$harvest_ready, na.rm = TRUE), # harvest_ready = sum(ndvi_values <= thresholds$harvest_ready & # ndvi_values > thresholds$post_harvest, na.rm = TRUE), # post_optimal = sum(ndvi_values <= thresholds$post_harvest, na.rm = TRUE), # still_growing = sum(ndvi_values > thresholds$early_maturity, na.rm = TRUE) # ) # # total_pixels <- sum(!is.na(ndvi_values)) # # print("Harvest Timing Assessment:") # for (stage in names(harvest_assessment)) { # percentage <- round((harvest_assessment[[stage]] / total_pixels) * 100, 1) # print(paste(stage, ":", harvest_assessment[[stage]], "pixels (", percentage, "%)")) # } # } ## ----irrigation-assessment, eval=FALSE---------------------------------------- # # Calculate water-related indices # green_band <- red_band # terra::values(green_band) <- runif(2500, 0.1, 0.2) # names(green_band) <- "Green" # # # Calculate NDWI for water stress detection # ndwi <- calculate_water_index( # green = green_band, # nir = nir_band, # index_type = "NDWI", # verbose = TRUE # ) # # # Irrigation needs based on NDWI # ndwi_values <- terra::values(ndwi, mat = FALSE) # valid_ndwi <- ndwi_values[!is.na(ndwi_values)] # # if (length(valid_ndwi) > 0) { # irrigation_zones <- ifelse(valid_ndwi < -0.3, "High Need", # ifelse(valid_ndwi < -0.1, "Moderate Need", "Adequate")) # # irrigation_table <- table(irrigation_zones) # irrigation_percent <- round(prop.table(irrigation_table) * 100, 1) # # print("Irrigation Needs Assessment:") # for (zone in names(irrigation_table)) { # print(paste(zone, ":", irrigation_table[[zone]], "pixels (", irrigation_percent[[zone]], "%)")) # } # # # Visualize irrigation zones # irrigation_raster <- ndwi # terra::values(irrigation_raster) <- as.numeric(as.factor(irrigation_zones)) # names(irrigation_raster) <- "Irrigation_Needs" # # plot_raster_fast( # irrigation_raster, # title = "Irrigation Needs Assessment (1=Adequate, 2=High, 3=Moderate)", # color_scheme = "water" # ) # } ## ----crop-rotation, eval=FALSE------------------------------------------------ # # Simulate multi-year CDL data # create_rotation_analysis <- function() { # # Create sample rotation data # rotation_data <- data.frame( # year = rep(2021:2023, each = 4), # field = rep(c("Field_A", "Field_B", "Field_C", "Field_D"), 3), # crop = c( # "Corn", "Soybeans", "Corn", "Wheat", # 2021 # "Soybeans", "Corn", "Soybeans", "Corn", # 2022 # "Corn", "Soybeans", "Corn", "Soybeans" # 2023 # ), # yield = runif(12, 80, 200) # ) # # return(rotation_data) # } # # rotation_analysis <- create_rotation_analysis() # print("Crop Rotation Analysis:") # print(rotation_analysis) # # # Analyze rotation patterns # rotation_patterns <- list() # fields <- unique(rotation_analysis$field) # # for (field in fields) { # field_data <- rotation_analysis[rotation_analysis$field == field, ] # rotation_patterns[[field]] <- paste(field_data$crop, collapse = " → ") # } # # print("Rotation Patterns:") # for (field in names(rotation_patterns)) { # print(paste(field, ":", rotation_patterns[[field]])) # } ## ----sustainability-metrics, eval=FALSE--------------------------------------- # # Calculate diversity index for crop rotation # calculate_crop_diversity <- function(rotation_data) { # diversity_scores <- list() # # for (field in unique(rotation_data$field)) { # field_crops <- rotation_data$crop[rotation_data$field == field] # unique_crops <- length(unique(field_crops)) # total_years <- length(field_crops) # # # Simple diversity score (0-1, higher = more diverse) # diversity_scores[[field]] <- unique_crops / total_years # } # # return(diversity_scores) # } # # diversity_scores <- calculate_crop_diversity(rotation_analysis) # # print("Crop Diversity Scores (Higher = More Diverse):") # for (field in names(diversity_scores)) { # print(paste(field, ":", round(diversity_scores[[field]], 2))) # } # # # Sustainability recommendations # avg_diversity <- mean(unlist(diversity_scores)) # print(paste("Average field diversity:", round(avg_diversity, 2))) # # if (avg_diversity < 0.5) { # print("RECOMMENDATION: Consider more diverse crop rotations") # } else if (avg_diversity > 0.8) { # print("STATUS: Good crop diversity maintained") # } ## ----farm-integration, eval=FALSE--------------------------------------------- # # Create farm-ready data export # create_farm_export <- function(analysis_results, field_boundaries) { # # Compile key metrics for farm management # farm_data <- list( # summary_statistics = list(), # recommendations = list(), # alerts = list() # ) # # # Extract key metrics # if (exists("veg_stats") && length(veg_stats) > 0) { # farm_data$summary_statistics <- veg_stats # } # # # Generate management recommendations # if (exists("alerts")) { # if (alerts$drought_risk > 20) { # farm_data$recommendations <- c(farm_data$recommendations, # "Increase irrigation frequency") # } # if (alerts$disease_risk > 15) { # farm_data$recommendations <- c(farm_data$recommendations, # "Monitor for disease symptoms") # } # } # # return(farm_data) # } # # # Export field-level results # if (exists("field_analysis") && inherits(field_analysis, "sf")) { # field_summary <- sf::st_drop_geometry(field_analysis) # # # Add coordinates for GPS guidance # field_centroids <- sf::st_centroid(field_analysis) # field_coords <- sf::st_coordinates(field_centroids) # field_summary$centroid_lon <- field_coords[, 1] # field_summary$centroid_lat <- field_coords[, 2] # # print("Field Summary for Farm Management:") # print(field_summary) # } ## ----economic-analysis, eval=FALSE-------------------------------------------- # # Estimate economic value based on vegetation health # calculate_economic_metrics <- function(ndvi_data, crop_prices) { # if (!exists("ndvi_values")) return(NULL) # # # Simplified yield prediction based on NDVI # # (In practice, would use crop-specific models) # predicted_yield <- ifelse(ndvi_values > 0.7, "High", # ifelse(ndvi_values > 0.5, "Medium", "Low")) # # yield_table <- table(predicted_yield) # # # Economic projections (simplified) # economic_zones <- list( # high_yield = sum(predicted_yield == "High", na.rm = TRUE), # medium_yield = sum(predicted_yield == "Medium", na.rm = TRUE), # low_yield = sum(predicted_yield == "Low", na.rm = TRUE) # ) # # return(economic_zones) # } # # # Calculate economic zones # if (exists("ndvi_values")) { # economic_analysis <- calculate_economic_metrics(ndvi_values, # list(corn = 5.50, soybeans = 13.00)) # # if (!is.null(economic_analysis)) { # total_pixels <- sum(unlist(economic_analysis)) # # print("Economic Potential Analysis:") # for (zone in names(economic_analysis)) { # pixels <- economic_analysis[[zone]] # percentage <- round((pixels / total_pixels) * 100, 1) # print(paste(zone, ":", pixels, "pixels (", percentage, "%)")) # } # } # } ## ----quality-assurance, eval=FALSE-------------------------------------------- # # Vegetation index quality assessment # quality_check <- function(vegetation_indices) { # qc_results <- list() # # if (inherits(vegetation_indices, "SpatRaster")) { # for (i in 1:terra::nlyr(vegetation_indices)) { # index_name <- names(vegetation_indices)[i] # values <- terra::values(vegetation_indices[[i]], mat = FALSE) # # qc_results[[index_name]] <- list( # total_pixels = length(values), # valid_pixels = sum(!is.na(values)), # coverage_percent = round((sum(!is.na(values)) / length(values)) * 100, 1), # value_range = range(values, na.rm = TRUE), # outliers = sum(values < -1 | values > 1, na.rm = TRUE) # For normalized indices # ) # } # } # # return(qc_results) # } # # # Perform quality check # if (exists("agricultural_indices")) { # qc_results <- quality_check(agricultural_indices) # # print("Data Quality Assessment:") # for (index in names(qc_results)) { # qc <- qc_results[[index]] # print(paste(index, "- Coverage:", qc$coverage_percent, "%, Outliers:", qc$outliers)) # } # } ## ----field-validation, eval=FALSE--------------------------------------------- # # Create sampling points for field validation # create_validation_points <- function(field_boundary, n_points = 10) { # # Generate random points within field boundary # if (inherits(field_boundary, "sf")) { # sample_points <- sf::st_sample(field_boundary, n_points) # validation_sf <- sf::st_sf( # point_id = paste0("VP_", 1:length(sample_points)), # validation_type = "Ground_Truth", # geometry = sample_points # ) # return(validation_sf) # } # return(NULL) # } # # # Create validation points for first field # if (exists("fields_sf")) { # validation_points <- create_validation_points(fields_sf[1, ], n_points = 5) # # if (!is.null(validation_points)) { # print("Validation Points Created:") # print(sf::st_drop_geometry(validation_points)) # # # Extract vegetation index values at validation points # if (exists("agricultural_indices")) { # validation_extracted <- universal_spatial_join( # source_data = validation_points, # target_data = agricultural_indices, # method = "extract", # verbose = FALSE # ) # # print("Extracted Values at Validation Points:") # validation_summary <- sf::st_drop_geometry(validation_extracted) # print(head(validation_summary)) # } # } # } ## ----complete-workflow, eval=FALSE-------------------------------------------- # # Complete agricultural analysis workflow # run_agricultural_workflow <- function(spectral_data, cdl_data = NULL, # region_boundary = NULL) { # workflow_results <- list() # # # Step 1: Calculate vegetation indices # message("Step 1: Calculating vegetation indices...") # indices <- calculate_multiple_indices( # spectral_data = spectral_data, # indices = c("NDVI", "SAVI", "DVI"), # output_stack = TRUE # ) # workflow_results$vegetation_indices <- indices # # # Step 2: Crop classification (if CDL available) # if (!is.null(cdl_data)) { # message("Step 2: Analyzing crop distribution...") # crop_mask <- create_crop_mask(cdl_data, "corn") # workflow_results$crop_mask <- crop_mask # } # # # Step 3: Stress assessment # message("Step 3: Assessing crop stress...") # if (inherits(indices, "SpatRaster") && "NDVI" %in% names(indices)) { # ndvi_vals <- terra::values(indices[["NDVI"]], mat = FALSE) # stress_percent <- sum(ndvi_vals < 0.5, na.rm = TRUE) / sum(!is.na(ndvi_vals)) * 100 # workflow_results$stress_assessment <- round(stress_percent, 1) # } # # # Step 4: Generate recommendations # message("Step 4: Generating management recommendations...") # recommendations <- character() # # if (!is.null(workflow_results$stress_assessment)) { # if (workflow_results$stress_assessment > 25) { # recommendations <- c(recommendations, "High stress detected - investigate causes") # } # if (workflow_results$stress_assessment < 10) { # recommendations <- c(recommendations, "Crop conditions appear favorable") # } # } # # workflow_results$recommendations <- recommendations # # return(workflow_results) # } # # # Run the workflow # workflow_output <- run_agricultural_workflow( # spectral_data = spectral_stack # ) # # print("Complete Agricultural Workflow Results:") # print(paste("Stress assessment:", workflow_output$stress_assessment, "% of area")) # print("Recommendations:") # for (rec in workflow_output$recommendations) { # print(paste("-", rec)) # } ## ----best-practices, eval=FALSE----------------------------------------------- # # 1. Always validate your data # print("Data Validation Checklist:") # print("✓ Check coordinate reference systems") # print("✓ Verify date ranges match growing season") # print("✓ Validate vegetation index ranges") # print("✓ Confirm crop mask accuracy") # # # 2. Use appropriate indices for your crop type # crop_index_recommendations <- list( # corn = c("NDVI", "EVI2", "SAVI", "DVI"), # soybeans = c("NDVI", "SAVI", "GNDVI", "DVI"), # wheat = c("NDVI", "SAVI", "MSAVI"), # general = c("NDVI", "SAVI", "DVI", "RVI") # ) # # print("Recommended indices by crop type:") # for (crop in names(crop_index_recommendations)) { # indices <- crop_index_recommendations[[crop]] # print(paste(crop, ":", paste(indices, collapse = ", "))) # } # # # 3. Monitor throughout growing season # print("Seasonal Monitoring Schedule:") # print("- Planting: Establish baseline measurements") # print("- Early season: Monitor emergence and establishment") # print("- Mid-season: Peak growth assessment and stress detection") # print("- Late season: Maturity evaluation and harvest timing") # print("- Post-harvest: Residue analysis and planning") ## ----precision-ag-integration, eval=FALSE------------------------------------- # # Precision agriculture workflow components # precision_ag_components <- list( # data_collection = c("Satellite imagery", "UAV surveys", "Ground sensors"), # analysis_methods = c("Vegetation indices", "Stress detection", "Yield mapping"), # decision_support = c("Variable rate applications", "Irrigation scheduling", "Harvest timing"), # validation = c("Ground truth sampling", "Yield monitoring", "Economic analysis") # ) # # print("Precision Agriculture Integration:") # for (component in names(precision_ag_components)) { # print(paste(component, ":", paste(precision_ag_components[[component]], collapse = ", "))) # }