Code

library(ggplot2)
library(rpart)
library(rpart.plot)
library(caret)

agri_data <- read.csv("agricultural_data.csv",
                      stringsAsFactors = FALSE)

required_cols <- c(
  "State", "Crop", "N", "P", "K",
  "Temperature_C", "Humidity_pct", "Rainfall_mm",
  "Soil_pH", "Fertilizer_kg_ha", "Yield_tonnes_ha"
)

missing_cols <- setdiff(required_cols, names(agri_data))

if (length(missing_cols) > 0) {
  stop(paste("Missing columns:",
             paste(missing_cols, collapse = ", ")))
}

agri_data$Crop <- factor(agri_data$Crop)

numeric_cols <- setdiff(required_cols, c("State", "Crop"))
agri_data[numeric_cols] <- lapply(
  agri_data[numeric_cols], as.numeric
)

agri_data <- na.omit(agri_data)

# TASK 1: ANALYSE CROPS WITH HIGH AVERAGE YIELDS

avg_yields <- aggregate(
  Yield_tonnes_ha ~ Crop,
  data = agri_data,
  FUN = mean
)

avg_yields <- avg_yields[
  order(-avg_yields$Yield_tonnes_ha), ]

print(avg_yields)

ggplot(avg_yields,
       aes(x = reorder(Crop, Yield_tonnes_ha),
           y = Yield_tonnes_ha)) +
  geom_col() +
  coord_flip() +
  labs(
    title = "Average Yield by Crop",
    x = "Crop",
    y = "Average Yield (tonnes/ha)"
  ) +
  theme_minimal()

# TASK 2: RAINFALL AND CROP YIELD RELATIONSHIP

correlation <- cor(
  agri_data$Rainfall_mm,
  agri_data$Yield_tonnes_ha,
  use = "complete.obs",
  method = "pearson"
)

cat("Correlation between rainfall and yield:",
    round(correlation, 4), "\n")

ggplot(agri_data,
       aes(x = Rainfall_mm,
           y = Yield_tonnes_ha,
           color = Crop)) +
  geom_point(size = 3) +
  geom_smooth(
    aes(group = 1),
    method = "lm",
    se = FALSE,
    color = "black"
  ) +
  labs(
    title = "Rainfall vs Crop Yield",
    x = "Rainfall (mm)",
    y = "Yield (tonnes/ha)"
  ) +
  theme_minimal()

# TASK 3: REGRESSION MODEL TO PREDICT CROP YIELD

set.seed(123)

train_rows <- createDataPartition(
  agri_data$Crop,
  p = 0.75,
  list = FALSE
)

train_data <- agri_data[train_rows, ]
test_data <- agri_data[-train_rows, ]

regression_model <- lm(
  Yield_tonnes_ha ~ N + P + K +
    Temperature_C + Humidity_pct +
    Rainfall_mm + Soil_pH + Fertilizer_kg_ha,
  data = train_data
)

print(summary(regression_model))

# TASK 4: DECISION TREE FOR CROP RECOMMENDATION

classification_model <- rpart(
  Crop ~ N + P + K +
    Temperature_C + Humidity_pct +
    Rainfall_mm + Soil_pH + Fertilizer_kg_ha,
  data = train_data,
  method = "class",
  control = rpart.control(
    cp = 0.001,
    minsplit = 2,
    minbucket = 1,
    maxdepth = 5,
    xval = 5
  )
)

rpart.plot(
  classification_model,
  main = "Crop Recommendation Decision Tree",
  type = 3,
  extra = 101,
  fallen.leaves = TRUE
)

# TASK 5: EVALUATE REGRESSION AND CLASSIFICATION MODELS

pred_yields <- predict(
  regression_model,
  newdata = test_data
)

actual_yields <- test_data$Yield_tonnes_ha

mae_value <- mean(
  abs(actual_yields - pred_yields)
)

rmse_value <- sqrt(
  mean((actual_yields - pred_yields)^2)
)

sse <- sum((actual_yields - pred_yields)^2)
sst <- sum((actual_yields - mean(actual_yields))^2)

r2_value <- if (sst == 0) {
  NA_real_
} else {
  1 - sse / sst
}

cat("Test MAE:", round(mae_value, 4), "\n")
cat("Test RMSE:", round(rmse_value, 4), "\n")
cat("Test R-squared:", round(r2_value, 4), "\n")

pred_crops <- predict(
  classification_model,
  newdata = test_data,
  type = "class"
)

pred_crops <- factor(
  pred_crops,
  levels = levels(agri_data$Crop)
)

actual_crops <- factor(
  test_data$Crop,
  levels = levels(agri_data$Crop)
)

conf_matrix <- confusionMatrix(
  data = pred_crops,
  reference = actual_crops
)

print(conf_matrix$table)

cat(
  "Test Classification Accuracy:",
  round(conf_matrix$overall["Accuracy"] * 100, 2),
  "%\n"
)