# Load libraries
library(sf)
## Linking to GEOS 3.8.0, GDAL 3.0.4, PROJ 6.3.1; sf_use_s2() is TRUE
library(ggplot2)
library(dplyr)
## 
## 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
library(readr)

#Load kml file
map_data <- st_read("kerala.kml")
## Reading layer `Kerala Districts' from data source `/cloud/project/kerala.kml' using driver `LIBKML'
## Simple feature collection with 14 features and 18 fields
## Geometry type: POLYGON
## Dimension:     XYZ
## Bounding box:  xmin: 74.86434 ymin: 8.293018 xmax: 77.4124 ymax: 12.79553
## z_range:       zmin: 0 zmax: 0
## Geodetic CRS:  WGS 84
print(map_data)
## Simple feature collection with 14 features and 18 fields
## Geometry type: POLYGON
## Dimension:     XYZ
## Bounding box:  xmin: 74.86434 ymin: 8.293018 xmax: 77.4124 ymax: 12.79553
## z_range:       zmin: 0 zmax: 0
## Geodetic CRS:  WGS 84
## First 10 features:
##          Name description timestamp begin  end altitudeMode tessellate extrude
## 1   ALAPPUZHA        <NA>      <NA>  <NA> <NA>         <NA>         -1       0
## 2   ERNAKULAM        <NA>      <NA>  <NA> <NA>         <NA>         -1       0
## 3      IDUKKI        <NA>      <NA>  <NA> <NA>         <NA>         -1       0
## 4      KANNUR        <NA>      <NA>  <NA> <NA>         <NA>         -1       0
## 5   KASARAGOD        <NA>      <NA>  <NA> <NA>         <NA>         -1       0
## 6      KOLLAM        <NA>      <NA>  <NA> <NA>         <NA>         -1       0
## 7    KOTTAYAM        <NA>      <NA>  <NA> <NA>         <NA>         -1       0
## 8   KOZHIKODE        <NA>      <NA>  <NA> <NA>         <NA>         -1       0
## 9  MALAPPURAM        <NA>      <NA>  <NA> <NA>         <NA>         -1       0
## 10   PALAKKAD        <NA>      <NA>  <NA> <NA>         <NA>         -1       0
##    visibility drawOrder icon   District  STATE REMARKS State_LGD DISTRICT_L
## 1          -1        NA <NA>  ALAPPUZHA KERALA                32        554
## 2          -1        NA <NA>  ERNAKULAM KERALA                32        555
## 3          -1        NA <NA>     IDUKKI KERALA                32        556
## 4          -1        NA <NA>     KANNUR KERALA                32        557
## 5          -1        NA <NA>  KASARAGOD KERALA                32        558
## 6          -1        NA <NA>     KOLLAM KERALA                32        559
## 7          -1        NA <NA>   KOTTAYAM KERALA                32        560
## 8          -1        NA <NA>  KOZHIKODE KERALA                32        561
## 9          -1        NA <NA> MALAPPURAM KERALA                32        562
## 10         -1        NA <NA>   PALAKKAD KERALA                32        563
##    Shape_Leng Shape_Area                       geometry
## 1      295138 1453750000 POLYGON Z ((76.33289 9.8769...
## 2      414199 3115070000 POLYGON Z ((76.58487 10.298...
## 3      578311 4444600000 POLYGON Z ((77.17662 10.358...
## 4      339137 2976760000 POLYGON Z ((75.4706 12.3004...
## 5      344681 1993000000 POLYGON Z ((75.41595 12.500...
## 6      351984 2552850000 POLYGON Z ((76.96151 9.1529...
## 7      322326 2253780000 POLYGON Z ((76.62856 9.8545...
## 8      296316 2364870000 POLYGON Z ((75.78204 11.797...
## 9      405137 3589220000 POLYGON Z ((76.23719 11.527...
## 10     552163 4537720000 POLYGON Z ((76.70614 11.238...
district_column <- "Name"
#Load data
datafile <- read.csv("datafile.csv")
head(datafile)
##    District Geographical.Area X2011.Assessment.Very.Dense.Forest
## 1 ALAPPUZHA              1414                                  0
## 2 ERNAKULAM              2407                                 12
## 3    IDUKKI              5019                                350
## 4    KANNUR              2966                                 21
## 5 KASARAGOD              1992                                  0
## 6    KOLLAM              2491                                 75
##   X2011.Assessment...Mod..Dense.Forest X2011.Assessment...Open.Forest
## 1                                   12                             26
## 2                                  298                            385
## 3                                 2159                           1421
## 4                                  351                            269
## 5                                  307                            285
## 6                                  632                            623
##   X2011.Assessment...Total Percent.of.GA Change Scrub
## 1                       38          2.69      0     0
## 2                      695         28.87     -1     1
## 3                     3930         78.30     -2     5
## 4                      641         21.61      0     0
## 5                      592         29.72      0     1
## 6                     1330         53.39     -7     0
# Prepare the forest density data frame for merging
forest_density_data <- data.frame(
  District = c(datafile$District),
  Forest_Density_Percent = c(datafile$Percent.of.GA)) 
#  Merge the data
merged_data <- map_data %>%
  left_join(forest_density_data, by = setNames("District", district_column)) 


# Create the Choropleth Map using ggplot2
kerala_choropleth <- ggplot(data = merged_data) +
  geom_sf(aes(fill = Forest_Density_Percent), 
          color = "white",  
          size = 0.5) +   
  
  # Set a color scale for forest density (e.g., a green gradient)
  scale_fill_viridis_c(
    option = "G", 
    direction = -1, 
    name = "Forest Density (%)", 
    guide = guide_colorbar(barheight = 10) 
  ) +
  
  # Add labels and title
  labs(
    title = "Kerala District-wise Forest Density",
    subtitle = "Visualizing the percentage of forest cover by district",
    caption = "Data Source: data.gov.in"
  ) +
  
  # Customize the theme 
  theme_minimal() +
  theme(
    plot.title = element_text(hjust = 0.5, face = "bold"),
    plot.subtitle = element_text(hjust = 0.5),
    axis.title = element_blank(),
    axis.text = element_blank(),
    panel.grid = element_blank()
  )

#  Display the map
print(kerala_choropleth)