##==========================
## Importing required libraries
##==========================

library(readxl)
library(ggplot2)
library(corrplot)

##==========================
## Data Importing and Summary of the data
##==========================

# -----------
# DATA IMPORTING
# -----------

# Setting file path
file_path <- "D:/AirQualityUCI.xlsx"

# Importing dataset into a data frame
AirQuality <- read_excel(
  path = file_path
)

# ------------------
# CHECKING COLUMN NAMES
# ------------------

# Displaying column names
colnames(AirQuality)


# ----------------
# CHECKING DATA TYPES
# ----------------

# Inspectting structure and data types
str(AirQuality)


# --------------
# FIRST FEW ROWS
# --------------

# Viewing first few rows
head(AirQuality)


# ------------------
# DATASET DIMENSIONS
# ------------------

# Checking number of rows and columns
dim(AirQuality)


# ----------------------------
# DESCRIPTIVE STATISTICS
# ----------------------------

# Summary statistics for all variables
summary(AirQuality)

# Summary statistics for numeric variables only
summary(AirQuality[sapply(AirQuality, is.numeric)])

##=============================
## EXPLORATORY VISUALISATION
##=============================

# ---------------------------------------------------
# HANDLING MISSING VALUES (-200 indicates NA)
# ---------------------------------------------------

AirQuality[AirQuality == -200] <- NA


# -----------------------------------------
# CREATING PROPER DATETIME VARIABLE
# -----------------------------------------

time_only <- format(AirQuality$Time, "%H:%M:%S")

AirQuality$DateTime <- as.POSIXct(
  paste(format(AirQuality$Date, "%Y-%m-%d"), time_only),
  format = "%Y-%m-%d %H:%M:%S",
  tz = "UTC"
)

# Checking that DateTime is valid
summary(AirQuality$DateTime)


# ------------------------------
# TIME TRENDS OF TWO POLLUTANTS
# ------------------------------

# ---- Time Series: CO(GT) ----
ggplot(AirQuality, aes(x = DateTime, y = `CO(GT)`)) +
  geom_line(color = "blue", na.rm = TRUE) +
  labs(title = "Time Series of CO Concentration",
       x = "Date",
       y = "CO (mg/m^3)") +
  theme_minimal()


# ---- Time Series: NO2(GT) ----
ggplot(AirQuality, aes(x = DateTime, y = `NO2(GT)`)) +
  geom_line(color = "red", na.rm = TRUE) +
  labs(title = "Time Series of NO2 Concentration",
       x = "Date",
       y = "NO2 (µg/m^3)") +
  theme_minimal()



# --------------------------------------------
# RELATIONSHIPS BETWEEN POLLUTANTS AND WEATHER
# --------------------------------------------

# ---- Temperature vs CO(GT) ----
ggplot(AirQuality, aes(x = T, y = `CO(GT)`)) +
  geom_point(alpha = 0.4) +
  geom_smooth(method = "lm", color = "black", se = FALSE) +
  labs(title = "Relationship Between Temperature and CO",
       x = "Temperature (°C)",
       y = "CO (mg/m^3)") +
  theme_minimal()


# ---- Relative Humidity vs NO2(GT) ----
ggplot(AirQuality, aes(x = RH, y = `NO2(GT)`)) +
  geom_point(alpha = 0.4, color = "darkgreen") +
  geom_smooth(method = "lm", color = "black", se = FALSE) +
  labs(title = "Relationship Between Relative Humidity and NO2",
       x = "Relative Humidity (%)",
       y = "NO2 (µg/m^3)") +
  theme_minimal()

##========================
## CORRELATION ANALYSIS
##========================

# ------------------------------
# SELECTING POLLUTANTS FOR ANALYSIS
# ------------------------------
# Selecting four major pollutants for correlation analysis

pollutants <- AirQuality[, c("CO(GT)", "NO2(GT)", "NOx(GT)", "C6H6(GT)")]

# Removing rows with missing values
pollutants_clean <- na.omit(pollutants)


# ------------------------------
# CALCULATING CORRELATION MATRIX
# ------------------------------
# Pearson correlation is appropriate because variables are continuous

cor_matrix <- cor(pollutants_clean, method = "pearson")

# Displaying correlation values
cor_matrix


# -----------------------------------
# VISUALISING CORRELATION USING HEATMAP
# -----------------------------------

corrplot(cor_matrix,
         method = "color",
         type = "upper",
         addCoef.col = "black",
         tl.col = "black",
         tl.srt = 45)


# ------------------------------
# SCATTERPLOT MATRIX 
# ------------------------------

pairs(pollutants_clean,
      main = "Scatterplot Matrix of Selected Pollutants")


##===============================
## HYPOTHESIS TESTING (T-TEST)
##===============================

AirQuality[AirQuality == -200] <- NA

# Selecting relevant variables
data_test <- AirQuality[, c("CO(GT)", "T")]

# Removing missing values
data_test <- na.omit(data_test)


# ------------------------------------------
# CREATING TWO CONDITIONS (HIGH vs LOW TEMP)
# ------------------------------------------
# Using median split to classify temperature levels

median_temp <- median(data_test$T)

data_test$TempGroup <- ifelse(data_test$T > median_temp,
                              "High", "Low")

data_test$TempGroup <- as.factor(data_test$TempGroup)


# ----------
# HYPOTHESES
# ----------

# H0 (Null Hypothesis):
# There is no difference in mean CO levels between
# high-temperature and low-temperature conditions.

# H1 (Alternative Hypothesis):
# There is a significant difference in mean CO levels
# between high and low temperature conditions.


# ------------------------------------------
# CHECKING NORMALITY (ASSUMPTION FOR T-TEST)
# ------------------------------------------

shapiro_high <- shapiro.test(data_test$`CO(GT)`[data_test$TempGroup == "High"])
shapiro_low  <- shapiro.test(data_test$`CO(GT)`[data_test$TempGroup == "Low"])

shapiro_high
shapiro_low


# ---------------------------
# PERFORMING APPROPRIATE TEST
# ---------------------------

# If both groups approximately normal → need to use independent t-test
# Otherwise → need to use Wilcoxon rank-sum test

t_test_result <- t.test(`CO(GT)` ~ TempGroup,
                        data = data_test,
                        var.equal = FALSE)

t_test_result


# Non-parametric alternative
wilcox_test_result <- wilcox.test(`CO(GT)` ~ TempGroup,
                                  data = data_test)

wilcox_test_result


# ----------------------------
# VISUALISING GROUP DIFFERENCE
# ----------------------------


ggplot(data_test, aes(x = TempGroup, y = `CO(GT)`)) +
  geom_boxplot(fill = "lightblue") +
  labs(title = "Comparison of CO Levels Under High and Low Temperature",
       x = "Temperature Group",
       y = "CO (mg/m^3)") +
  theme_minimal()


##=============================
## SEASONAL GROUP COMPARISON
##=============================

# Selecting relevant variables
season_data <- AirQuality[, c("Date", "CO(GT)")]

# Removing missing values
season_data <- na.omit(season_data)
# -------------------------
# CREATING SEASON VARIABLE
# -------------------------

# Extracting month from Date
season_data$Month <- as.numeric(format(season_data$Date, "%m"))

# Assigning seasons
season_data$Season <- ifelse(season_data$Month %in% c(12,1,2), "Winter",
                             ifelse(season_data$Month %in% c(3,4,5), "Spring",
                                    ifelse(season_data$Month %in% c(6,7,8), "Summer",
                                           "Autumn")))

season_data$Season <- as.factor(season_data$Season)


# ----------------
# STATE HYPOTHESES
# ----------------

# H0 (Null Hypothesis):
# Mean CO concentrations are equal across all seasons.

# H1 (Alternative Hypothesis):
# At least one season has a significantly different
# mean CO concentration.


# ------------------
# CHECKING NORMALITY 
# ------------------

by(season_data$`CO(GT)`, season_data$Season, shapiro.test)


# ------------------------
# PERFORMING ONE-WAY ANOVA
# ------------------------

anova_model <- aov(`CO(GT)` ~ Season, data = season_data)

summary(anova_model)


# --------------
# POST-HOC TEST 
# --------------

TukeyHSD(anova_model)


# --------------------------
# NON-PARAMETRIC ALTERNATIVE 
# --------------------------

kruskal_result <- kruskal.test(`CO(GT)` ~ Season, data = season_data)

kruskal_result


# -------------------------
# VISUALISING USING BOXPLOT
# -------------------------

ggplot(season_data, aes(x = Season, y = `CO(GT)`)) +
  geom_boxplot(fill = "lightblue") +
  labs(title = "Seasonal Distribution of CO Concentration",
       x = "Season",
       y = "CO (mg/m^3)") +
  theme_minimal()
