# =====================================================================
# Ch02 常態分配、標準誤與信賴區間 (Nielsen 4.3.1 後半)
# 學習目標：68-95-99.7 法則 / Z 值與 t 值 / 信賴區間 CI 的計算
# =====================================================================

# ---------- 1. 常態分配：68-95-99.7 法則 ----------
# 大多數隨機量測誤差服從常態分配（鐘形曲線）
#   ±1 SD 內約 68% 的數據
#   ±2 SD 內約 95%
#   ±3 SD 內約 99.7%
# R 可以用 pnorm() 直接算「小於某值的機率」

pnorm(1) - pnorm(-1)    # P(-1SD < X < +1SD)  -> 0.6827 (68%)
pnorm(2) - pnorm(-2)    # -> 0.9545 (95%)
pnorm(3) - pnorm(-3)    # -> 0.9973 (99.7%)

# ---------- 2. Z 值：大樣本(n>30)時的信賴區間 ----------
# CI = x_bar ± Z * SD / sqrt(n)
# Nielsen Table 4.2：信心水準對應 Z 值
#   80%->1.29, 90%->1.64, 95%->1.96, 99%->2.58, 99.9%->3.29
# 注意：課本 Table 4.2 的 80% 那一列是 1.29，但 R 算出來的精確值是
# qnorm(0.90) = 1.2816（四捨五入其實比較接近 1.28）；課本的 1.29 是查表捨入的結果，
# 兩者差異純粹來自捨入方式，不是計算錯誤。
# R 用 qnorm(p) 可反查（p = 中央面積）
qnorm(0.80 + 0.10)      # 1.2816 ~ 1.29（課本捨入值，R 精確值為 1.28）
qnorm(0.90 + 0.05)      # 1.6449 ~ 1.64
qnorm(0.95 + 0.025)     # 1.9600 = 1.96
qnorm(0.99 + 0.005)     # 2.5758 ~ 2.58

# 課本範例：假設水分測了 25 次，x_bar=64.72, SD=0.2927
n_big <- 25; xbar <- 64.72; s <- 0.2927
ci_z <- xbar + c(-1, 1) * qnorm(0.975) * s / sqrt(n_big)
round(ci_z, 2)          # -> 64.61  64.83，即 64.72 ± 0.115%

# 標準誤 (Standard Error of the Mean, SEM)
sem <- s / sqrt(n_big)
sem                     # 平均數的「不確定度」，n 越大 SEM 越小

# ---------- 3. t 值：小樣本(n<30)才是食品分析的日常 ----------
# 分析通常只做 3~5 次重複，必須用 t 分布代替常態
# CI = x_bar ± t * SD / sqrt(n)，df = n - 1
# R 用 qt(p, df) 反查 t 值表（Nielsen Table 4.3）
qt(0.975, 1)   # df=1, 95% -> 12.706（課本 12.7）
qt(0.975, 2)   # -> 4.303（課本 4.30）
qt(0.975, 3)   # -> 3.182（課本 3.18）
qt(0.975, 9)   # -> 2.262（課本 2.26）

# 回到真正的 4 次水分測定
moisture <- c(64.53, 64.45, 65.10, 64.78)
xbar <- mean(moisture); s <- sd(moisture); n <- length(moisture)
t_crit <- qt(0.975, df = n - 1)         # 3.182
half_width <- t_crit * s / sqrt(n)
half_width                              # -> 0.4656（課本用查表值 t=3.18 手算得 0.465）
c(xbar - half_width, xbar + half_width) # -> 64.249 ~ 65.181

cat(sprintf("95%% CI = %.2f ± %.3f %%", xbar, half_width))

# ---------- 4. 模擬實驗：親眼看看「重複抽樣」 ----------
# 我們模擬「每次抽 n=4 個樣本算平均」1000 次，
# 看看 95% CI 是否真的蓋住真實平均值約 95 次/100 次
set.seed(2024)                    # 設定亂數種子讓結果可重現
true_mu <- 65.05                  # 假設的真值 (%)
hits <- 0                         # 計數：CI 有蓋住真值的次數
for (i in 1:1000) {
  samp <- rnorm(4, mean = true_mu, sd = 0.293)  # 抽 4 個樣本
  ci <- mean(samp) + c(-1, 1) * qt(0.975, 3) * sd(samp)/sqrt(4)
  if (true_mu >= ci[1] && true_mu <= ci[2]) hits <- hits + 1
}
hits / 1000        # 大約 0.94~0.96，證明「95% 信心」不是口號！

# n 越大，CI 越窄（更精確）— 試著把 n 改成 16 再跑一次：
# （用 n16 這個新變數名，避免覆寫掉前面 n=4 水分範例的 n，
#   不然檔尾的小結會誤印成 n=16 但數值其實還是 n=4 算出來的）
n16 <- 16
samp <- rnorm(n16, mean = true_mu, sd = 0.293)
mean(samp) + c(-1,1) * qt(0.975, n16-1) * sd(samp)/sqrt(n16)

# ---------- 5. 小結 ----------
# 報告結果的專業寫法：平均值 ± 半寬 (信心水準)
cat(sprintf("水分含量 = %.2f ± %.2f %% (95%%, n=%d)",
            xbar, half_width, n))