正規分布に従う確率変数をつくる

確率変数\(original\)\(N(5, 9)\)にしたがう。確率変数\(standarized\)\(original\)を標準化した変数である。

  1. これら2つのベクトル
  2. 2つのベクトルを列として含むワイド型データフレーム
  3. 2をtidyにした縦に長いデータフレーム

をつくる。

original <- rnorm(n = 1000000, mean = 5, sd = 3)
standarized <- (original - mean(original)) / sd(original)

dfWide <- data.frame(original, standarized)
dfLong <- dfWide %>% pivot_longer(cols = 1:2, names_to = "type")

分布を確認する

dfLong %>% 
  group_by(type) %>% 
  summarise(mean(value), sd(value), min(value), max(value))
## # A tibble: 2 x 5
##   type        `mean(value)` `sd(value)` `min(value)` `max(value)`
## * <chr>               <dbl>       <dbl>        <dbl>        <dbl>
## 1 original         5.00e+ 0        3.00        -9.99        19.9 
## 2 standarized      1.95e-17        1           -5.01         4.96

\(standarized\)がほぼほぼ\(N(0, 1)\)に従っていることがわかる。

描画する

dfLong %>% 
  ggplot(mapping = aes(x = value, y = ..density.., fill = type)) +
  geom_histogram(position = "dodge", bins = 50) +
  scale_x_continuous(name = "x")

変換前と後でかなり形が変わっている。

標準化のプロセスを確認する

  1. まず\(N(0, 1)\)に従う\(standarized\)における-3の下側確率\[F(x) = \int_{-\infty}^{-3} f(x) dx = p\]を求める。ただし\(f(x) = P(standarized = x), F(x) = P(standarized \leq x)\)
  2. 1で求めた下側確率\(p\)を用いて、\[G(q) = \int_{-\infty}^q g(x)dx = p\]を満たす\(q\)を求める。ただし\(g(x) = P(original = x), G(x) = P(standarized \leq x)\)
  3. \(q\)を標準化した値\(q\prime\)を求める。すなわち\[q\prime = \frac{q-\mu}{\sigma}\]
  4. 理論的にいえば\(q\prime = -3\)となっているはずなので、それを確認する。
p <- pnorm(q = -3, mean = 0, sd = 1)
p
## [1] 0.001349898
q <- qnorm(p = p, mean = 5, sd = 3)
q
## [1] -4
qPrime <- (q - 5)/3
qPrime
## [1] -3

ここで、上の図を見ながら考えてみると、\[|F(q\prime) - E(standarized)| = |-3 -0| = 3 \\ |G(q) - E(original)| = |-4 -3| = 7\]である。同じく標準正規分布に従う確率変数といえども、標準化したほうが期待値と下側確率が\(p\)となる値との差は狭まっている。つまりキュッと締まった形になっている。