NMRL

# Load necessary libraries
library(tidyverse)

## ── Attaching core tidyverse packages ──────────────────────── tidyverse 2.0.0 ──
## ✔ dplyr     1.1.4     ✔ readr     2.1.5
## ✔ forcats   1.0.0     ✔ stringr   1.5.1
## ✔ ggplot2   3.5.1     ✔ tibble    3.2.1
## ✔ lubridate 1.9.3     ✔ tidyr     1.3.1
## ✔ purrr     1.0.2     
## ── Conflicts ────────────────────────────────────────── tidyverse_conflicts() ──
## ✖ dplyr::filter() masks stats::filter()
## ✖ dplyr::lag()    masks stats::lag()
## ℹ Use the conflicted package (<http://conflicted.r-lib.org/>) to force all conflicts to become errors

library(readxl)
library(janitor)

## 
## Attaching package: 'janitor'
## 
## The following objects are masked from 'package:stats':
## 
##     chisq.test, fisher.test

# Load the dataset
df <- read_excel("~/project/R PROGRAM/NMRL/for data analysis.xlsx") %>% clean_names()
head(df)

## # A tibble: 6 × 13
##   sample_number villlage collection_date collection_time     source_of_water
##           <dbl> <chr>    <chr>           <dttm>              <chr>          
## 1             1 Makina   22/07/2024      1899-12-31 17:05:00 Water track    
## 2             2 Makina   22/07/2024      1899-12-31 17:06:00 Public supply  
## 3             3 Makina   22/07/2024      1899-12-31 17:09:00 Public supply  
## 4             4 Makina   22/07/2024      1899-12-31 17:10:00 Public supply  
## 5             5 Makina   22/07/2024      1899-12-31 17:12:00 Public supply  
## 6             6 Makina   22/07/2024      1899-12-31 17:15:00 Public supply  
## # ℹ 8 more variables: exact_site_of_collection <chr>,
## #   any_source_of_pollution <chr>, source_of_pollution <chr>,
## #   is_water_treated <chr>, treatment_method <chr>, organism <chr>,
## #   e_coli_cfu_in_100ml <dbl>, pathotypes <chr>

#change e_coli_cfu_in_100ml to water_contamination
df <- df %>% rename(water_contamination = e_coli_cfu_in_100ml)

#change columns to character
df <- df %>% mutate(across(where(is.character),as.factor))
#change to date
df$collection_date <- as.Date(df$collection_date, format = "%d/%m/%Y")

#water contamination by source of water
kruskal.test(water_contamination ~ source_of_water, data = df)

## 
##  Kruskal-Wallis rank sum test
## 
## data:  water_contamination by source_of_water
## Kruskal-Wallis chi-squared = 11.533, df = 3, p-value = 0.009166

#water contamination by source of pollution
kruskal.test(water_contamination ~ source_of_pollution, data = df)

## 
##  Kruskal-Wallis rank sum test
## 
## data:  water_contamination by source_of_pollution
## Kruskal-Wallis chi-squared = 8.1942, df = 3, p-value = 0.04216

#water contamination by is water treatment
kruskal.test(water_contamination ~ is_water_treated, data = df)

## 
##  Kruskal-Wallis rank sum test
## 
## data:  water_contamination by is_water_treated
## Kruskal-Wallis chi-squared = 21.928, df = 1, p-value = 2.831e-06

 # Chi-square test
chisq_res <- chisq.test(table(df$source_of_water, df$pathotypes))

## Warning in stats::chisq.test(x, y, ...): Chi-squared approximation may be
## incorrect

p_val <- signif(chisq_res$p.value, 3)  # Round to 3 significant figures

# Create a bar plot (position = "fill" gives proportions; use "stack" for counts)
ggplot(df, aes(x = source_of_water, fill = pathotypes)) +
  geom_bar(position = "fill") +
  scale_y_continuous(labels = scales::percent_format()) +
  labs(title = "Distribution of Pathotypes by Source of Water",
       x = "Source of Water",
       y = "Proportion",
       fill = "Pathotype") +
  theme_minimal() +
  annotate("text", x = Inf, y = Inf, label = paste("Chi-square p =", p_val),
           hjust = 1.1, vjust = 1.5, size = 5, color = "red")

#variation of phanotypes by source of water
df %>%
  filter(!is.na(pathotypes) & pathotypes != "Negative") %>%
  group_by(source_of_water, pathotypes) %>%
  summarise(count = n()) %>%
  ggplot(aes(x = source_of_water, y = count, fill = pathotypes)) +
  geom_bar(stat = "identity", position = "dodge") +
  labs(title = "Distribution of Pathotypes by Water Source",
       x = "Water Source", y = "Count") +
  theme_minimal()

## `summarise()` has grouped output by 'source_of_water'. You can override using
## the `.groups` argument.

#variation of phanotypes by source of pollution
df %>%
  filter(!is.na(pathotypes) & pathotypes != "Negative") %>%
  group_by(source_of_pollution, pathotypes) %>%
  summarise(count = n()) %>%
  ggplot(aes(x = source_of_pollution, y = count, fill = pathotypes)) +
  geom_bar(stat = "identity", position = "dodge") +
  labs(title = "Distribution of Pathotypes by Source of Pollution",
       x = "Source of Pollution", y = "Count") +
  theme_minimal()

## `summarise()` has grouped output by 'source_of_pollution'. You can override
## using the `.groups` argument.

# Treatment
ggplot(df, aes(x = is_water_treated, y = water_contamination)) +
  geom_boxplot(fill = "lightgreen") +
  labs(title = "E. coli by Water Treatment", x = "Treated?", y = "E. coli CFU in 100ml") +
  theme_minimal()