Rで連鎖律

連鎖律とは

2つの関数が合成されている場合、合成関数の微分は、それぞれの関数の微分の積として求められる。これを連鎖律(chain rule)と呼ぶ。

いま、次の2つの関数を考える。

\[
\begin{eqnarray}
y &=& f(u) \\
u &=& g(x)
\end{eqnarray}
\]

\(x\) : 独立変数
\(u\) : 中間変数(\(x\) の関数)
\(y\) : 従属変数(\(u\) の関数、すなわち \(x\) の合成関数)
\(f\) : \(u\) を \(y\) に対応させる関数
\(g\) : \(x\) を \(u\) に対応させる関数

このとき、\(y\) は \(x\) の合成関数 \(y = f(g(x))\) とみなせる。連鎖律は、この合成関数の \(x\) に関する微分を次のように与える。

\[
\frac{dy}{dx} = \frac{dy}{du} \cdot \frac{du}{dx}
\]

ここで、

\(\dfrac{dy}{du}\) : \(u\) が微小変化したときの \(y\) の変化率(\(f\) を \(u\) で微分したもの、\(f'(u)\))
\(\dfrac{du}{dx}\) : \(x\) が微小変化したときの \(u\) の変化率(\(g\) を \(x\) で微分したもの、\(g'(x)\))
\(\dfrac{dy}{dx}\) : \(x\) が微小変化したときの \(y\) の変化率

すなわち、

\[
\frac{dy}{dx} = f'(g(x)) \cdot g'(x)
\]

が成り立つ。直感的には、\(x\) が変化すると、まず \(u\) がその変化率 \(g'(x)\) に応じて変化し、続いて \(u\) の変化が \(y\) の変化率 \(f'(u)\) に応じて \(y\) に伝わる。この2段階の変化率の掛け算が全体の変化率になる、というのが連鎖律の意味である。

以下では、具体例として

\[
\begin{eqnarray}
u &=& g(x) = x^2 \\
y &=& f(u) = \sin(u)
\end{eqnarray}
\]

を取り上げ、連鎖律による解析的な微分と、数値微分による近似値が一致することを確認する。

Rコード

library(ggplot2)

# 内側の関数 g(x) = x^2 とその導関数 g'(x) = 2x
g <- function(x) x^2
g_prime <- function(x) 2 * x

# 外側の関数 f(u) = sin(u) とその導関数 f'(u) = cos(u)
f <- function(u) sin(u)
f_prime <- function(u) cos(u)

# 合成関数 y = f(g(x))
y_composed <- function(x) f(g(x))

# 連鎖律による解析的な導関数 dy/dx = f'(g(x)) * g'(x)
dy_dx_analytic <- function(x) f_prime(g(x)) * g_prime(x)

# 中心差分法による数値微分
# 中心差分法: (y(x + h) - y(x - h)) / (2h) で dy/dx を近似する
numerical_derivative <- function(fn, x, h = 1e-5) {
  (fn(x + h) - fn(x - h)) / (2 * h)
}

# 検証用のxの範囲
x_values <- seq(-2, 2, by = 0.1)

# 解析解(連鎖律)と数値微分の結果をデータフレームにまとめる
kekka <- data.frame(
  x = x_values,
  kaiseki = dy_dx_analytic(x_values),
  suuchi = numerical_derivative(y_composed, x_values)
) |>
  transform(gosa = abs(kaiseki - suuchi))

# 誤差の最大値を確認する
max(kekka$gosa)

[1] 9.865126e-10

誤差の最大値は 9.865126e-10、すなわち約 \(9.87 \times 10^{-10}\) でした。

これは、連鎖律から解析的に求めた導関数 \(\dfrac{dy}{dx} = f'(g(x)) \cdot g'(x)\) の値と、中心差分法による数値微分の値との差の絶対値を、\(x \in [-2, 2]\) の範囲(刻み幅 0.1)で計算し、その中で最大のものを取り出した結果です。

この値がほぼゼロに等しいことから、\(x\) の範囲全体にわたって解析解と数値微分の結果が事実上一致していることがわかります。

理論上、両者は完全な等式として成り立つはずですが、実際の計算では次の2つの要因によりごくわずかな誤差が生じます。

数値微分自体が近似計算であること
- 中心差分法 \(\dfrac{f(x+h) - f(x-h)}{2h}\) は、\(h\) を有限の値(ここでは \(h = 10^{-5}\))としているため、真の導関数とはテイラー展開の高次項に由来する誤差(打ち切り誤差)を持つ
浮動小数点演算に伴う丸め誤差
- コンピュータ上の実数計算は有限桁の精度で行われるため、加減算や除算のたびにごく小さな誤差が蓄積する

打ち切り誤差は \(h\) のオーダーの2乗(\(O(h^2)\))で減少する一方、丸め誤差は \(h\) を小さくしすぎるとかえって増大するという性質があり、両者のバランスが取れる \(h\) の付近で誤差が最小になります。

今回得られた \(10^{-9}\) 台という誤差の大きさは、\(h = 10^{-5}\) を用いた中心差分法としては適当な精度であり、この結果は連鎖律の関係式 \(\dfrac{dy}{dx} = f'(g(x)) \cdot g'(x)\) が数値的にも正しく成立していることを裏付けるものといえます。

# 解析解と数値微分の結果を重ねてプロットし、一致していることを可視化する
ggplot(kekka, aes(x = x)) +
  geom_line(aes(y = kaiseki, color = "連鎖律による解析解"), linewidth = 1) +
  geom_point(aes(y = suuchi, color = "数値微分による近似"), size = 1.5, shape = 1) +
  scale_color_manual(
    values = c(
      "連鎖律による解析解" = "steelblue",
      "数値微分による近似" = "firebrick"
    )
  ) +
  labs(
    title = "連鎖律の検証: y = sin(x^2) の導関数",
    x = "x",
    y = "dy/dx",
    color = "凡例"
  ) +
  theme_minimal()