Test 01: Z-test for a Population Mean (Variance Known)

# remove all variable
rm(list = ls())

# load required functions
source(file = "../r/func.R")

# load required packages
wants = c("pacman", "dplyr", "psych", "ggpubr", "ggplot2", "DescTools")
has = wants %in% rownames(x = installed.packages())
if (any(!has)) install.packages(wants[!has])
pacman::p_load(char = wants)

# load required data
data = read.csv(file = "../data/Data_Test_01.csv", header = TRUE)

# show data
psych::headTail(x = data, top = 2, bottom = 2)

# visualizing sample data:

# normality plot
nPlot = ggpubr::ggqqplot(data = data$Yield, ylab = "Yield (kg/ha)")

# histogram plot
hPlot = ggpubr::gghistogram(data = data$Yield, binwidth = 1000, add = "mean", 
                    fill = "gray", xlab = "Yield (kg/ha)",
                    add.params = list(size = 1.2, linetype = 2))

## Warning in (function (mapping = NULL, data = NULL, ..., xintercept, na.rm = FALSE, : Using both `xintercept` and `mapping` may not have the desired result as mapping is overwritten if `xintercept` is specified

# box plot
bPlot = ggpubr::ggboxplot(data = data$Yield,
                          add = c("mean"), fill = "gray", width = 0.5, xlab = "Year",
                          ylab = "Yield (kg/ha)", orientation = "horizontal") +
  ggplot2::geom_hline(yintercept = 2843.95, linetype = "dashed", color = "red", size = 1.5) +
  ggplot2::annotate(geom = "text", x = 1.3, y = 2700, label = "P. Mean", color = "red", angle = 90)

# arrange on one page
ggpubr::ggarrange(ggpubr::ggarrange(hPlot, nPlot, nrow = 1, ncol = 2),
                  bPlot, nrow = 2, ncol = 1, heights = c(4, 3))

# method 1
# setting initial parameter values: population
mu.p = 2843.95
sd.p = 587.01

# setting initial parameter values: sample
n.x = nrow(x = data)
mean.x = mean(x = data$Yield)

# test statistic
Z = (mean.x - mu.p) / (sd.p / sqrt(n.x))

# compute the critical values and p-value
alpha = 0.05
alternative = "two.sided"

if (alternative != "two.sided") {
  lower.critical = qnorm(p = (1 - alpha), lower.tail = FALSE)
  upper.critical = qnorm(p = (1 - alpha), lower.tail = TRUE)
} else {
  lower.critical = qnorm(p = (1 - (alpha / 2)), lower.tail = FALSE)
  upper.critical = qnorm(p = (1 - (alpha / 2)), lower.tail = TRUE)
}

p_value = pnorm(q = Z, lower.tail = TRUE)
p_value = switch(EXPR = alternative,
                 two.sided = 2 * min(p_value, 1 - p_value),
                 less = p_value,
                 greater = 1 - p_value)

cat("Test Statistic:", Z, "\n")

## Test Statistic: 7.445936

cat("Critical Values:", c(lower.critical, upper.critical), "\n")

## Critical Values: -1.959964 1.959964

cat("P-value:", p_value, "\n")

## P-value: 9.636736e-14

plotDistStat(distribution = "norm", statistic = Z, alpha = 0.05,
             alternative = "two.sided")

# method 2
# setting initial parameter values: population
mu.p = 2843.95
sd.p = 587.01

# setting initial parameter values:  sample
x = data$Yield

# use "ZTest" function from "DescTools" package
DescTools::ZTest(x = x, mu = mu.p, sd_pop = sd.p, 
                 alternative = "two.sided", conf.level = 0.95)

## 
##  One Sample z-test
## 
## data:  x
## z = 7.4459, Std. Dev. Population = 587.01, p-value = 9.626e-14
## alternative hypothesis: true mean is not equal to 2843.95
## 95 percent confidence interval:
##  3299.372 3624.788
## sample estimates:
## mean of x 
##   3462.08