Rの関数:density {stats}

# 必要なライブラリを読み込みます
library(ggplot2)
library(purrr)
library(dplyr)

# --- 1. データの準備 ---
# 再現性を確保するために乱数シードを設定します
seed <- 20250902
set.seed(seed)
# 標準正規分布に従う100個の乱数を生成します
sim_data <- rnorm(100, mean = 0, sd = 1)


# --- 2. 各カーネルで密度推定を実行 ---
# 比較するカーネルの種類のリストを作成します
kernel_types <- c(
  "gaussian", "epanechnikov", "rectangular", "triangular",
  "biweight", "cosine", "optcosine"
)

# 公平な比較のため、全てのカーネルで同じバンド幅を使用します
# ここではRのデフォルトのアルゴリズム "nrd0" でバンド幅を計算します
common_bw <- bw.nrd0(sim_data)

# purrr::map_dfr を使い、各カーネルでdensity()を実行し、結果を一つのデータフレームにまとめます
density_df <- map_dfr(kernel_types, ~ {
  # density関数で推定を実行
  d <- density(sim_data, kernel = .x, bw = common_bw, n = 1024)

  # ggplotで扱いやすいように、x座標、y座標、カーネルの種類を持つデータフレームに変換
  data.frame(
    x = d$x,
    y = d$y,
    kernel_type = .x
  )
})

# --- 3. ggplot2で結果をプロット ---
# カーネルの種類ごとに色分けして、密度曲線を重ねて描画します
g <- ggplot(density_df, aes(x = x, y = y, color = kernel_type)) +
  geom_line(linewidth = 0.8, alpha = 0.9) +
  # 元のデータの分布が分かりやすいように、ヒストグラムを背景に薄く描画します
  geom_histogram(
    data = data.frame(x = sim_data),
    aes(x = x, y = after_stat(density)),
    inherit.aes = FALSE,
    bins = 15,
    fill = "gray",
    alpha = 0.3
  ) +
  labs(
    title = "カーネル関数の種類による密度推定曲線の比較",
    subtitle = paste0("同じデータ・同じバンド幅(bw = ", round(common_bw, 4), ")で7種類のカーネルを比較"),
    x = "値",
    y = "密度",
    color = "カーネルの種類"
  ) +
  theme_minimal() +
  theme(legend.position = "bottom")
g
g + facet_wrap(. ~ kernel_type) + theme(legend.position = "none")

# --- 1. データの準備 ---
set.seed(seed)
sim_data <- rnorm(100, mean = 0, sd = 1)


# --- 2. 比較するバンド幅を設定 ---
# Rのデフォルトアルゴリズム("nrd0")で計算されるバンド幅を基準とします
bw_default <- bw.nrd0(sim_data)

# 比較のために3種類のバンド幅を用意します
bw_list <- list(
  "小 (過剰適合)" = bw_default * 0.2,
  "中 (デフォルト)" = bw_default,
  "大 (過小適合)" = bw_default * 2.5
)

cat("比較に使用するバンド幅の値:\n")
# 各バンド幅の具体的な値を表示します
walk2(names(bw_list), bw_list, ~ cat(paste0(.x, ": bw = ", round(.y, 4), "\n")))


# --- 3. 各バンド幅で密度推定を実行 ---
# purrr::map_dfrを使い、各バンド幅でdensity()を実行し、結果を一つのデータフレームにまとめます
density_df <- map_dfr(names(bw_list), ~ {
  # density関数を実行。kernelは"gaussian"に固定
  d <- density(sim_data, bw = bw_list[[.x]], kernel = "gaussian", n = 1024)

  data.frame(
    x = d$x,
    y = d$y,
    bw_label = .x # バンド幅のラベル
  )
})

# 凡例の順序を意図通りにするため、factor型に変換します
density_df$bw_label <- factor(density_df$bw_label, levels = names(bw_list))


# --- 4. ggplot2で結果をプロット ---
# バンド幅のラベルごとに色分けして、密度曲線を重ねて描画します
ggplot(density_df, aes(x = x, y = y, color = bw_label)) +
  geom_line(linewidth = 0.8, alpha = 0.9) +
  # 元のデータの分布が分かりやすいように、ヒストグラムを背景に薄く描画します
  geom_histogram(
    data = data.frame(x = sim_data),
    aes(x = x, y = after_stat(density)),
    inherit.aes = FALSE,
    bins = 15,
    fill = "gray",
    alpha = 0.3
  ) +
  labs(
    title = "バンド幅(bw)の変更による密度推定曲線の比較",
    subtitle = "カーネルを'gaussian'に固定し、バンド幅の大きさを3パターンで比較",
    x = "値",
    y = "密度",
    color = "バンド幅の大きさ"
  ) +
  theme_minimal() +
  theme(legend.position = "bottom")