estimated_standard_deviation
This is an old revision of the document!
Table of Contents
Why n-1
문제.
우리는 모집단의 (population) 평균값을 알고 있다면 샘플의 분산값을 다음과 같이 구할 수 있다. 그리고, 이것을 모집단의 분산값으로 추정할 수 있다.
\begin{eqnarray*} \hat{\sigma^{2}} = \dfrac {\displaystyle\sum_{i=1}^{n}{(X_{i}-\mu)}} {n} \end{eqnarray*}
그러나, 현재 우리가 가지고 있는 것은 샘플 밖에 없다. 즉, 모집단의 평균은 알지 못하는 상태이기에 모집단 분산을 추정하는 계산에 사용할 수 없다. 따라서 샘플의 평균을 사용한다. 그런데, 샘플의 평균을 사용할 때는 분모에 N 대신에 n-1을 사용해야 한다. 왜 n-1을 사용하는것이 모집단의 분산값 추정에 도움이 되는가가 문제이다.
\begin{eqnarray*} \hat{\sigma}^{2} \neq \frac {\displaystyle\sum_{i=1}^{n}{(X_{i}-\overline{X})}} {n} \end{eqnarray*}
\begin{eqnarray*} \hat{\sigma}^{2} = \frac {\displaystyle\sum_{i=1}^{n}{(X_{i}-\overline{X})}} {n-1} \end{eqnarray*}
이 모집단의 분산값을 ($\sigma^2$) 대표함을 알아보는 것이 문제이다.
직관적 이해
분산은 $ \text{SS}/ \text{df} $ 라고 배웠다. SS = Sum of Something Square라고 설명하고, 여기서 Something 은 error (개인의 점수를 평균으로 추측했을 때 틀린 만큼의 에러값) 하였다. 이렇게 평균을 가지고 개인점수를 예측하는 것이 가장 작은 오차를 갖는 방법이라고 하였다 (개인점수 중에서 평균이 제일 많이 나오므로, 평균으로 개인점수를 예측하면 제일 들 틀린다). 따라서 어느 한 집합에서 (샘플에서) 개인점수에서 평균을 빼고 이를 제곱하여 모두 더한 값은 (평균이 아닌) 다른 값으로 구한 값보다 항상 작게 된다 (최소값을 갖는다).
이를 그림으로도 설명할 수 있다. 아래에서 녹색의 세로선은 모집단의 평균값이고, 붉은색의 세로선은 3개로 이루어진 샘플의 평균값이다. 그리고 녹색 가로선은 3개의 샘플요소와 모집단평균과의 ($\mu$) 차이값들이고, 적색가로선은 3개의 샘플요소와 샘플평균과의 ($\overline{X}$) 차이값이다. 이 차이값들을 모아서 길이를 비교한 것이 그래프의 하단이다. 적색가로선 세개의 합이 녹색가로선 세개의 합보다 작다. 이는 샘플평균을 사용했을 때와 모집단의 평균을 사용했을 때를 비교하는 것이지만 모집단 평균외에 다른 값을 썼어도 마찬가지이다.
> m.p1 <- 12 > s1 <- c(6, 9, 15) > hist(s1) > s1 [1] 6 9 15 > abline(v=s1, lwd=3, lty=2) > abline(v=10, lwd=3, lty=1, col="red") > abline(v=12, lwd=3, lty=1, col="green") > > m.s1 <- mean(s1) > m.p1 <- 12 > abs(s1-m.s1) [1] 4 1 5 > abs(s1-m.p1) [1] 6 3 3 >
실험적, R에서 시뮬레이션으로 이해
아래 output의 코멘트를 읽을 것
rm(list=ls())
rnorm2 <- function(n,mean,sd){
mean+sd*scale(rnorm(n))
}
# set.seed(191)
nx <- 1000
mx <- 50
sdx <- mx * 0.1
sdx # 5
x <- rnorm2(nx, mx, sdx)
# x <- rnorm2(1000, 50, 5) 와 동일
mean(x)
sd(x)
length(x)
hist(x)
x.span <- seq(from = mean(x)-3*sd(x),
to = mean(x)+3*sd(x),
by = .1)
x.span
residuals <- function(x, v) {
return(x - v)
}
# sum of square residual 값을
# 구하는 펑션
ssr <- function(x, v) {
residuals <- (x - v)
return(sum(residuals^2))
}
# mean square residual 값을
# 구하는 펑션 (mean square
# residual = variance)
msr <- function(x, v) {
residuals <- (x - v)
return((mean(residuals^2)))
}
ssrs <- c() # sum of square residuals
msrs <- c() # mean square residuals = variance
vs <- c() # the value of v in (x - v)
# x.span의 값들을 v값으로 삼아 sum(x-x.span)^2 처럼 구하면
# SS값을 구한 것이 된다. 우리가 배운 SS값은 x.span의 값으로
# 샘플의 평균을 사용했을 때의 residual 값이다. x.span은
# 샘플의 평균을 중심으로 여러가지 값을 사용하는 것을 가정한다.
for (i in x.span) {
res.x <- residuals(x,i)
msr.x <- msr(x,i)
msrs <- append(msrs, msr.x)
vs <- append(vs, i)
}
# 아래 plot은 SS값들이나 (두번째는) MS값들을 v값이 변화함에
# 따라서 (x.span의 범위에 따라서 변화) 어떻게 변화하는지
# 구한 것
plot(msrs)
# v값이 x.span에 따라서 변화하여 대입되었을 때의
# MS값들을 (msr 펑션으로 구한 mean square값)
# 모아 놓은 값이 msrs
msrs
# 아래는 위에서 계산한 msr 값들을 저장한 msrs값들 중에서 최소값이
# 되는 것을 찾은 것. 우리는 이것이 샘플의 평균임을 안다.
min(msrs)
# 최소값일 때의 위치 (msrs에서 몇번째인지)
min.pos.msrs <- which(msrs == min(msrs))
min.pos.msrs
# msr 최소값이 구해졌을 때 사용된 v값
vs[min.pos.msrs]
output
> rm(list=ls())
>
> rnorm2 <- function(n,mean,sd){
+ mean+sd*scale(rnorm(n))
+ }
>
> # set.seed(191)
> nx <- 1000
> mx <- 50
> sdx <- mx * 0.1
> sdx # 5
[1] 5
> x <- rnorm2(nx, mx, sdx)
> # x <- rnorm2(1000, 50, 5) 와 동일
>
> mean(x)
[1] 50
> sd(x)
[1] 5
> length(x)
[1] 1000
> hist(x)
>
SS = sum(x-mean(x))^2 인데, mean(x)을 즉, x집합의 평균을, x 원소값을 예측하는데 (빼는데) 사용하면 SS값이 최소값이 된다고 하였다. 이것을 R에서 simulation으로 알아보기 위해서 mean(x) 대신에 다른 숫자들을 넣어보려고 한다. 이를 v라고 하면 sum(x-v)^2이라는 SS값들을 구해서 비교하려는 것이다. 대입할 숫자들은 (v) mean(x) +- 3 sd(x) 를 범위로 하고, 그 범위의 시작 숫자에서 (시작은 mean(x)-3sd(x)가 된다) 0.1씩 증가시키면서 대입하고, 각 숫자마다 (처음 숫자는 35, 다음 숫자는 35.1 . . . ) SS값을 구해서 저장하여 그것을 그래프로 그려보고 최소값이 어떤 것인지 보는 것이 진행하려는 작업이다.
단, 이 코드에서 SS대신 MS값을 (SS값을 n으로 나눈 값, 즉, variance값 혹은 Mean Square값) 구해서 보려고 하는데 이것은 같은 의미를 갖는다. 즉, 모든 SS값들에 n을 공토으로 나누어준 값을 저장하고 비교하려는 것이다.
> x.span <- seq(from = mean(x)-3*sd(x), + to = mean(x)+3*sd(x), + by = .1) > x.span [1] 35.0 35.1 35.2 35.3 35.4 35.5 35.6 35.7 35.8 [10] 35.9 36.0 36.1 36.2 36.3 36.4 36.5 36.6 36.7 [19] 36.8 36.9 37.0 37.1 37.2 37.3 37.4 37.5 37.6 [28] 37.7 37.8 37.9 38.0 38.1 38.2 38.3 38.4 38.5 [37] 38.6 38.7 38.8 38.9 39.0 39.1 39.2 39.3 39.4 [46] 39.5 39.6 39.7 39.8 39.9 40.0 40.1 40.2 40.3 [55] 40.4 40.5 40.6 40.7 40.8 40.9 41.0 41.1 41.2 [64] 41.3 41.4 41.5 41.6 41.7 41.8 41.9 42.0 42.1 [73] 42.2 42.3 42.4 42.5 42.6 42.7 42.8 42.9 43.0 [82] 43.1 43.2 43.3 43.4 43.5 43.6 43.7 43.8 43.9 [91] 44.0 44.1 44.2 44.3 44.4 44.5 44.6 44.7 44.8 [100] 44.9 45.0 45.1 45.2 45.3 45.4 45.5 45.6 45.7 [109] 45.8 45.9 46.0 46.1 46.2 46.3 46.4 46.5 46.6 [118] 46.7 46.8 46.9 47.0 47.1 47.2 47.3 47.4 47.5 [127] 47.6 47.7 47.8 47.9 48.0 48.1 48.2 48.3 48.4 [136] 48.5 48.6 48.7 48.8 48.9 49.0 49.1 49.2 49.3 [145] 49.4 49.5 49.6 49.7 49.8 49.9 50.0 50.1 50.2 [154] 50.3 50.4 50.5 50.6 50.7 50.8 50.9 51.0 51.1 [163] 51.2 51.3 51.4 51.5 51.6 51.7 51.8 51.9 52.0 [172] 52.1 52.2 52.3 52.4 52.5 52.6 52.7 52.8 52.9 [181] 53.0 53.1 53.2 53.3 53.4 53.5 53.6 53.7 53.8 [190] 53.9 54.0 54.1 54.2 54.3 54.4 54.5 54.6 54.7 [199] 54.8 54.9 55.0 55.1 55.2 55.3 55.4 55.5 55.6 [208] 55.7 55.8 55.9 56.0 56.1 56.2 56.3 56.4 56.5 [217] 56.6 56.7 56.8 56.9 57.0 57.1 57.2 57.3 57.4 [226] 57.5 57.6 57.7 57.8 57.9 58.0 58.1 58.2 58.3 [235] 58.4 58.5 58.6 58.7 58.8 58.9 59.0 59.1 59.2 [244] 59.3 59.4 59.5 59.6 59.7 59.8 59.9 60.0 60.1 [253] 60.2 60.3 60.4 60.5 60.6 60.7 60.8 60.9 61.0 [262] 61.1 61.2 61.3 61.4 61.5 61.6 61.7 61.8 61.9 [271] 62.0 62.1 62.2 62.3 62.4 62.5 62.6 62.7 62.8 [280] 62.9 63.0 63.1 63.2 63.3 63.4 63.5 63.6 63.7 [289] 63.8 63.9 64.0 64.1 64.2 64.3 64.4 64.5 64.6 [298] 64.7 64.8 64.9 65.0
x-mean(x) = residual = error
sum(residual^2) = SS (sum of square)
SS/n = variance, mean square (ms, MS)
이 residual을 구하기 위해서 쓰는 mean(x)의 대체값들을 (v값들) x.span에 모아 놓은 것이다.
이 값을 출력해보았는데 35.1 에서 시작하여 65에서 끝나며, 0.1씩 증가한다.
>
> residuals <- function(x, v) {
+ return(x - v)
+ }
>
> # sum of square residual 값을
> # 구하는 펑션
> ssr <- function(x, v) {
+ residuals <- (x - v)
+ return(sum(residuals^2))
+ }
>
> # mean square residual 값을
> # 구하는 펑션 (mean square
> # residual = variance)
> msr <- function(x, v) {
+ residuals <- (x - v)
+ return((mean(residuals^2)))
+ }
>
- 이 후 쓸 function들. (x-v) = residual이라고 부르니까 이 residual을 모으는 function
- function ssr = x 집합과 v값 (x.span의 한 숫자)를 인수를 주었을 때 구할 수 있는 Sum of Square값들 (실제로는 사용하지 않는다. 대신 msr 펑션으로 MS값을 구한다).
- function msr = 의 ssr을 n으로 나누어 구한 mean square residual을 (분산) 구하는 function
> ssrs <- c() # sum of square residuals
> msrs <- c() # mean square residuals = variance
> vs <- c() # the value of v in (x - v)
>
> # x.span의 값들을 v값으로 삼아 sum(x-x.span)^2 처럼 구하면
> # SS값을 구한 것이 된다. 우리가 배운 SS값은 x.span의 값으로
> # 샘플의 평균을 사용했을 때의 residual 값이다. x.span은
> # 샘플의 평균을 중심으로 여러가지 값을 사용하는 것을 가정한다.
>
> for (i in x.span) {
+ res.x <- residuals(x,i)
+ msr.x <- msr(x,i)
+ msrs <- append(msrs, msr.x)
+ vs <- append(vs, i)
+ }
> # 아래 plot은 SS값들이나 (두번째는) MS값들을 v값이 변화함에
> # 따라서 (x.span의 범위에 따라서 변화) 어떻게 변화하는지
> # 구한 것
>
> plot(vs, msrs)
>
comment
- x.span의 처음값인 35.1을 넣어서 (x-v)를 구한 후
- msr 펑션으로 mean square residual 값을 구한다.
- 그리고 이 값을 어딘가에 (msrs) 저장한 후
- 그 다음 value인 35.2값을 v 대신에 넣어 다시
- msr값을 구하여 위의 msrs에 추가하여 저장한다.
- 이것을 x.span의 모든값에 걸쳐 진행하면
- msrs에는 x.span의 모든 값을 대응하여 구한
- msr값들이 저장된다. 이 msr값 중에서 최소값을 찾고
- 이 최소값을 구할 때 쓴 v값을 (x.span 중 하나의 값) 찾고자 한다.
- 이를 위해서 for 를 이용한 loop문을 쓴다.
> # v값이 x.span에 따라서 변화하여 대입되었을 때의 > # MS값들을 (msr 펑션으로 구한 mean square값) > # 모아 놓은 값이 msrs > msrs [1] 249.975 246.985 244.015 241.065 238.135 [6] 235.225 232.335 229.465 226.615 223.785 [11] 220.975 218.185 215.415 212.665 209.935 [16] 207.225 204.535 201.865 199.215 196.585 [21] 193.975 191.385 188.815 186.265 183.735 [26] 181.225 178.735 176.265 173.815 171.385 [31] 168.975 166.585 164.215 161.865 159.535 [36] 157.225 154.935 152.665 150.415 148.185 [41] 145.975 143.785 141.615 139.465 137.335 [46] 135.225 133.135 131.065 129.015 126.985 [51] 124.975 122.985 121.015 119.065 117.135 [56] 115.225 113.335 111.465 109.615 107.785 [61] 105.975 104.185 102.415 100.665 98.935 [66] 97.225 95.535 93.865 92.215 90.585 [71] 88.975 87.385 85.815 84.265 82.735 [76] 81.225 79.735 78.265 76.815 75.385 [81] 73.975 72.585 71.215 69.865 68.535 [86] 67.225 65.935 64.665 63.415 62.185 [91] 60.975 59.785 58.615 57.465 56.335 [96] 55.225 54.135 53.065 52.015 50.985 [101] 49.975 48.985 48.015 47.065 46.135 [106] 45.225 44.335 43.465 42.615 41.785 [111] 40.975 40.185 39.415 38.665 37.935 [116] 37.225 36.535 35.865 35.215 34.585 [121] 33.975 33.385 32.815 32.265 31.735 [126] 31.225 30.735 30.265 29.815 29.385 [131] 28.975 28.585 28.215 27.865 27.535 [136] 27.225 26.935 26.665 26.415 26.185 [141] 25.975 25.785 25.615 25.465 25.335 [146] 25.225 25.135 25.065 25.015 24.985 [151] 24.975 24.985 25.015 25.065 25.135 [156] 25.225 25.335 25.465 25.615 25.785 [161] 25.975 26.185 26.415 26.665 26.935 [166] 27.225 27.535 27.865 28.215 28.585 [171] 28.975 29.385 29.815 30.265 30.735 [176] 31.225 31.735 32.265 32.815 33.385 [181] 33.975 34.585 35.215 35.865 36.535 [186] 37.225 37.935 38.665 39.415 40.185 [191] 40.975 41.785 42.615 43.465 44.335 [196] 45.225 46.135 47.065 48.015 48.985 [201] 49.975 50.985 52.015 53.065 54.135 [206] 55.225 56.335 57.465 58.615 59.785 [211] 60.975 62.185 63.415 64.665 65.935 [216] 67.225 68.535 69.865 71.215 72.585 [221] 73.975 75.385 76.815 78.265 79.735 [226] 81.225 82.735 84.265 85.815 87.385 [231] 88.975 90.585 92.215 93.865 95.535 [236] 97.225 98.935 100.665 102.415 104.185 [241] 105.975 107.785 109.615 111.465 113.335 [246] 115.225 117.135 119.065 121.015 122.985 [251] 124.975 126.985 129.015 131.065 133.135 [256] 135.225 137.335 139.465 141.615 143.785 [261] 145.975 148.185 150.415 152.665 154.935 [266] 157.225 159.535 161.865 164.215 166.585 [271] 168.975 171.385 173.815 176.265 178.735 [276] 181.225 183.735 186.265 188.815 191.385 [281] 193.975 196.585 199.215 201.865 204.535 [286] 207.225 209.935 212.665 215.415 218.185 [291] 220.975 223.785 226.615 229.465 232.335 [296] 235.225 238.135 241.065 244.015 246.985 [301] 249.975 >
comment
- msrs값에 저장된 msr값들 (mean square residual값들) 중에서
- 가장 작은 값을 찾아서 그 값을 구하도록 한 v값을 찾고자 한다.
- msrs값을 눈으로 살펴보기에는 너무 힘드므로 . . . .
> # 아래는 위에서 계산한 msr 값들을 저장한 msrs값들 중에서 최소값이 되는 것을 찾은 > # 것. 우리는 이것이 샘플의 평균임을 안다. > min(msrs) [1] 24.975 > # 최소값일 때의 위치 (msrs에서 몇번째인지) > min.pos.msrs <- which(msrs == min(msrs)) > min.pos.msrs [1] 151 > # msr 최소값이 구해졌을 때 사용된 v값 > vs[min.pos.msrs] [1] 50 > > plot(vs, msrs, cex=1, lwd=1, lty=3) > abline(v=vs[min.pos.msrs]) > text(x=50, y=150, "msr gets minimal value, when v = 50" ) > > > > > > >
comment
- msrs 의 min(msrs)값을 찾는다 24.975
- 이 최소값이 어느 위치에 있는지 (몇번째 자리에 있는지) 찾는다 which(msrs == min(msrs))
- 이 위치가 151이다
- 이 151번째 사용된 (최소값인 msr값 = mean(x-v)^2)을 결과한) vs값을 찾는다 (vs[151])
- 이 값을 출력하니 50 이고 이 값은 x의 평균값이다.
- 이를 plot으로 출력한 것
다음으로 MS값을 구하는 식인
- y = sum( (x-v)^2 ) / n 값을
v에 대해서v를 가지고 y를 미분하면 (dy/dv)- 변화하는 v값마다의 기울기 값을 구할 수 있는데
- 이 값이 0이 되는 지점의 v값이 무엇인지를 구하면 위의 R코드에서 구한 값을 구할 수 있게 된다. 아래는 그 과정이다.
- 그래프에서 가장 작은 기울기값을 갖는 v 값을 구한다고 (derivatives) 가정하고 이해를 하면 수학적으로 이해할 수 있다. 1)
\begin{eqnarray*}
\dfrac{dy}{dv} = \dfrac{\text{d}}{\text{dv}} \dfrac{\sum{(x-v)^2}}{n} & = & \dfrac {\sum{2(x-v)*(-1)}}{n} \\
& = & \dfrac{\sum{-2(x-v)}}{n} \\
& = & -\dfrac{2}{n} \sum{(x-v)} \\
\end{eqnarray*}
위의 식이 0이 (기울기가 0이 되는 부분) 될 때의 v 값을 찾아야 하므로
\begin{eqnarray*} -\dfrac{2}{n} \sum{(x-v)} & = & 0 \\ \sum{(x-v)} & = & 0 \\ \sum{x} - n*v & = & 0 \\ n*v & = & \sum{x} \\ v & = & \dfrac {\sum{x}}{n} \\ \end{eqnarray*}
즉, 기울기가 0이 될때의 v값은 x집합의 평균값일 때이다.
R에서 msr이 최소값이 되는 v값 효율적으로 찾기
위의 미분을 이용한 방법을 써서 MS값이 최소값이 되는 순간을 R에서 구현해 볼 수 있다. 아래 아웃풋의 코멘트 부분을 읽을 것.
why n-1 gradient explanation
# the above no gradient
gradient <- function(x, v){
residuals = x - v
# y = (sum(x-v)^2)/n 혹은 (mean(x-v)^2)
# 의 식이 ms값을 구하는 식인데
# 이를 v에 대해서 미분하면 chain rule을 써야 한다
# 자세한 것은 http://commres.net/wiki/estimated_standard_deviation
# 문서 중에 미분 부분 참조
# dy/dv = ( 2(x-v)*-1 ) / n chain rule
# dy/dv = -2 (x-v) / n = -2 (mean(residual))
dx = -2 * mean(residuals)
# return(list("ds" = dx))
return(dx)
} # function returns ds value
zx <- (x-mean(x))/sd(x)
# pick one random v in (x-v)
v <- rnorm(1)
# Train the model with scaled features
learning.rate = 1e-1
msrs <- c()
vs <- c()
nlen <- 75
for (epoch in 1:nlen) {
residual <- residuals(zx, v)
msr.x <- msr(zx, v)
msrs <- append(msrs, msr.x)
grad <- gradient(zx, v)
step.v <- grad * learning.rate #
v <- v - step.v # 그 다음 v값
vs <- append(vs, v) # v값 저장
}
tail(msrs)
tail(vs)
plot(vs, msrs, type='b')
# scaled
vs # 변화하는 v 값들의 집합
vs.orig <- (vs*sd(x))+mean(x)
vs.orig
# 마지막 v값이 최소값에 근접한 값
v
v.orig <- (v*sd(x))+mean(x)
v.orig
plot(vs.orig, msrs, type='b')
output
> # the above no gradient
>
> gradient <- function(x, v){
+ residuals = x - v
+ # y = (sum(x-v)^2)/n 혹은 (mean(x-v)^2)
+ # 의 식이 ms값을 구하는 식인데
+ # 이를 v에 대해서 미분하면 chain rule을 써야 한다
+ # 자세한 것은 http://commres.net/wiki/estimated_standard_deviation
+ # 문서 중에 미분 부분 참조
+ # dy/dv = ( 2(x-v)*-1 ) / n chain rule
+ # dy/dv = -2 (x-v) / n = -2 (mean(residual))
+ dx = -2 * mean(residuals)
+ # return(list("ds" = dx))
+ return(dx)
+ } # function returns ds value
>
>
comment
이 R script의 목적은 v값이 최소값이 되는 지점을 자동적으로 찾아보려는 것이다. 이것을 위해서 우선 v값으로 사용할 첫 점수를 랜덤하게 구한 후 (아래 그래프에서 빨간색 지점), 자동적으로 그 다음 v 점수를 찾고 (녹색지점), 그 다음 v 점수를 찾고 (황금색 지점), . . . 이런 과정을 계속하면서 각 v 점수에서의 msr값을 구해서 이에 해당하는 v값을 찾아 보려고 한다. 빨간색, 녹색, 황금색, . . . 이를 자동적으로 구하기 위해서 두가지 방법을 사용하는데 그 것이
- gradient function과
- learning_rate 값이다.
gradient 펑션은 dy/dv 의 연쇄 미분식인 (chain rules) -2(x-v) / n = -2 mean(res) 값을 구하는 것이다. 이렇게 구한 값에 learning_rate값을 곱한후, 이것을 먼저 사용한 v값에서 (빨간색 지점) 빼 주어 다음 v값으로 (녹색지점) 사용하려고 한다. 이 녹색지점에서의 v값을 사용했을 때의 gradient값을 구한 후 다시 이값에 learning_rate인 0.1을 곱하여 그다음 스텝의 값을 얻고, 이 값을 바로 전의 v값에서 빼 준 값을 그 다음 v값으로 사용한다. 이렇게 구하는 v값들은 0.1씩 곱해주는 효과때문에 오른 쪽으로 옮겨가는 지점이 “점진적으로 줄어들게 되고” 이 지점이 msr의 최소값이 되는 지점으로 가게 된다.
클릭하면 큰 이미지로 볼 수 있음
> zx <- (x-mean(x))/sd(x) > # pick one random v in (x-v) > v <- rnorm(1) > >
comment
- 랜덤하게 v값을 찾음
v ← rnorm(1) - 원래는 mean(x)값 근처의 값을 랜덤하게 골라야 하므로
v ← rnorm(1, mean(x), sd(x))와 같이 써야 하지만 (현재의 경우,rnorm(1, 50, 5)가 된다)- 그렇게 하질 않고 x 집합의 원소들을 표준점수화 한 후
zx ← (x-mean(x))/sd(x)(이렇게 표준점수화 하면 x 변인의 평균이 0, 표준편차가 1 이 되는 집합으로 변한다) v ← rnorm(1, 0, 1)로 구한다. 뒤의 인자인 0, 1 은 default이모로 생략.- 이렇게 하는 이유는 혹시 나중에 다른 x집합에 똑같은 작업을 하더라도 그 집합의 평균과 표준편차를 사용하지 않고
- 단순히
rnorm(1)을 이용해서 찾으려고 하는 것이다.
> # Train the model with scaled features > learning.rate = 1e-1 >
comment
- 이 0.1은 gradient function으로 구한 해당 v값에 대한 y 미분값을 (기울기 값을) 구한 후, 여기에 곱하기 위해서 지정한다.
> msrs <- c()
> vs <- c()
>
> nlen <- 75
> for (epoch in 1:nlen) {
+ residual <- residuals(zx, v)
+ msr.x <- msr(zx, v)
+ msrs <- append(msrs, msr.x)
+
+ grad <- gradient(zx, v)
+ step.v <- grad * learning.rate #
+ v <- v - step.v # 그 다음 v값
+ vs <- append(vs, v) # v값 저장
+ }
>
> tail(msrs)
[1] 0.999 0.999 0.999 0.999 0.999 0.999
> tail(vs)
[1] 6.415945e-08 5.132756e-08 4.106205e-08 3.284964e-08 2.627971e-08 2.102377e-08
>
> plot(vs, msrs)
>
>
comment
- nlen는 그냥 자의적으로 지정한다. 여기서는 75로 했다.
- for 문에서 처음 v값은 위에서 랜덤으로 구한 값이다 (v).
- 이 v값으로 gradient 펑션의 아웃풋 값을 구하고 (-2(residual)) =
grad ← gradient(x, v) - 이 값에 learning_rate값을 곱한 값을 구하여
step.v ← grad * learaning_rate - 이 값을 원래 v값에서 빼준 후에
- 다시 (for 문에서 반복하는 동안) v값으로 바꾼 후
v ← v - step.v - for문을 nlen번 만큼 반복한다.
- 이 과정에서 저장한
- msr값들과
- vs값들의 마지막 6개를 살펴본다.
- v 값으로 x집합의 평균값을 사용하는 것이 최소 msr값이 된다는 것이 맞다면
- v 값은 0이 될것이다 (왜냐하면 x집합을 zx로 바꿨기 때문에, 즉 평균이 0이고 sd값이 1일 집합으로 바꿨기 때문에)
- 아래 그래프의 각 포인트는 v값의 이동을 나타내는데 grad*learning_rate의 영향으로 점진적으로 하가하여 최소값으로 도달한다.
> # scaled > vs # 변화하는 v 값들의 집합 [1] 3.119260e-01 2.495408e-01 1.996326e-01 1.597061e-01 1.277649e-01 1.022119e-01 8.176952e-02 [8] 6.541561e-02 5.233249e-02 4.186599e-02 3.349279e-02 2.679424e-02 2.143539e-02 1.714831e-02 [15] 1.371865e-02 1.097492e-02 8.779935e-03 7.023948e-03 5.619158e-03 4.495327e-03 3.596261e-03 [22] 2.877009e-03 2.301607e-03 1.841286e-03 1.473029e-03 1.178423e-03 9.427383e-04 7.541907e-04 [29] 6.033525e-04 4.826820e-04 3.861456e-04 3.089165e-04 2.471332e-04 1.977066e-04 1.581652e-04 [36] 1.265322e-04 1.012258e-04 8.098061e-05 6.478449e-05 5.182759e-05 4.146207e-05 3.316966e-05 [43] 2.653573e-05 2.122858e-05 1.698286e-05 1.358629e-05 1.086903e-05 8.695226e-06 6.956181e-06 [50] 5.564945e-06 4.451956e-06 3.561565e-06 2.849252e-06 2.279401e-06 1.823521e-06 1.458817e-06 [57] 1.167054e-06 9.336428e-07 7.469143e-07 5.975314e-07 4.780251e-07 3.824201e-07 3.059361e-07 [64] 2.447489e-07 1.957991e-07 1.566393e-07 1.253114e-07 1.002491e-07 8.019931e-08 6.415945e-08 [71] 5.132756e-08 4.106205e-08 3.284964e-08 2.627971e-08 2.102377e-08 > vs.orig <- (vs*sd(x))+mean(x) > vs.orig [1] 51.55963 51.24770 50.99816 50.79853 50.63882 50.51106 50.40885 50.32708 50.26166 50.20933 [11] 50.16746 50.13397 50.10718 50.08574 50.06859 50.05487 50.04390 50.03512 50.02810 50.02248 [21] 50.01798 50.01439 50.01151 50.00921 50.00737 50.00589 50.00471 50.00377 50.00302 50.00241 [31] 50.00193 50.00154 50.00124 50.00099 50.00079 50.00063 50.00051 50.00040 50.00032 50.00026 [41] 50.00021 50.00017 50.00013 50.00011 50.00008 50.00007 50.00005 50.00004 50.00003 50.00003 [51] 50.00002 50.00002 50.00001 50.00001 50.00001 50.00001 50.00001 50.00000 50.00000 50.00000 [61] 50.00000 50.00000 50.00000 50.00000 50.00000 50.00000 50.00000 50.00000 50.00000 50.00000 [71] 50.00000 50.00000 50.00000 50.00000 50.00000 > > # 마지막 v값이 최소값에 근접한 값 > v [1] 2.102377e-08 > v.orig <- (v*sd(x))+mean(x) > v.orig [1] 50 > > plot(vs.orig, msrs, type='b') > > >
comment
만약에 처음에 구한 랜덤 v값이 평균의 오른 쪽에있었더라면, 아래 그림과 같이 평균에 접근했을 것이다.
그렇다면 왜 n-2 혹은 n-(1/2)가 아니고 n-1인가? 이를 수학적인 증명을 통해서 살펴보면 다음 장과 같다.
수학적 증명
우선,
\begin{eqnarray*} Var[X] & = & E[(X-\mu)^{2}] \\ & = & E[(X^{2} - 2 X \mu + \mu^{2})] \\ & = & E[X^{2}] - 2 \mu E[X] + E[\mu^2] \\ & = & E[X^{2}] - 2 \mu E[X] + E[\mu^{2}], \;\; \text{because}\; E[X] = \mu \text{, } \; E[\mu^2] = \mu^2, \\ & = & E[X^{2}] - 2 \mu^{2} + \mu^{2} \\ & = & E[X^{2}] - \mu^{2} \end{eqnarray*}
이므로
\begin{align} E\left[X^2\right] & = Var\left[X\right] + \mu^2 \nonumber \\ & = \sigma^{2} + \mu^2 \\ \end{align}
마찬가지로
\begin{align}
Var \left[ \overline{X}\right] & = E \left[\overline{X}^2 \right] - \left[E(\overline{X})\right]^2 \nonumber \\
& = E\left[\overline{X}^{2}\right] - \mu^{2} \nonumber
\end{align}
따라서
\begin{align}
E\left[\overline{X}^{2}\right] & = Var\left[\overline{X}\right] + \mu^2 \nonumber \\
& = \frac {\sigma^{2}} {n} + \mu^{2}
\end{align}
참고로 위에서 $Var\left[\overline{X}\right] = \dfrac {\sigma^{2}} {n} $ 에 해당하는 설명은 mean and variance of the sample mean 문서를 볼 것.
참고로 Expected value (기대값)와 Variance (분산)의 연산에 과한 법칙으로는 2)
X,Y are Independent variables.
\begin{align*} E[aX] & = a E[X] \\ E[X+Y] & = E[X] + E[Y] \\ Var[aX] & = a^{\tiny{2}} Var[X] \\ Var[X+Y] & = Var[X] + Var[Y] \end{align*}
우리가 알고자 하는 것은 아래의 식이 population의 parameter인 $\sigma^{2}$ 의 값과 같은가이다.
\begin{align*}
E[s^{2}] & = E \left[\frac{\displaystyle\sum_{i=1}^{n}(X_{i}-\overline{X})^{2}}{n-1} \right] \qquad
\cdot\cdot\cdot\cdot\cdot\cdot\cdot\cdot\cdot\cdot\cdot\cdot \;\; (a) \\
& = \sigma^{2}
\end{align*}
위의 식에서 일부만을 추출해서 먼저 보자.
\begin{align*} E \left[\sum{(X_{i}-\overline{X})^{2}} \right] & = E \left[\sum(X_{i}^{2}- 2 X_{i} \overline{X} + \overline{X}^{2})\right] \\ & = E \left[ \sum{X_{i}^2} - \sum{2X_{i} \overline{X}} + \sum {\overline{X^2}} \right] \\ & = E \left[ \sum{X_{i}^2} - 2 \overline{X} \sum{X_{i}} + \sum{\overline{X^2}} \right] \\ & = E \left[ \sum{X_{i}^2} - 2 \overline{X} \sum{X_{i}} + n \overline{X^2} \right] \\ & = E \left[ \sum{X_{i}^2} - 2 \overline{X} \cdot (n \overline{X}) + n \overline {X^2} \right] \\ & = E \left[ \sum{X_{i}^2} - n \overline{X}^2 \right] \\ & = \sum {E\left(X_{i}^2\right)} - E\left(n\overline{X}^2\right) \\ & = \sum {E\left(X_{i}^2\right)} - n E\left(\overline{X}^2\right) \;\;\; \dots\dots\dots\dots\dots (3) \end{align*}
한 편, 위의 $(1), (2)$에서
\begin{align*} E\left[X_{i}^{2}\right] & = \sigma^{2} + \mu^2 \;\;\; \dots\dots\dots\dots\dots (1) \\ E\left[\overline{X}^{2}\right] & = \dfrac {\sigma^{2}} {n} + \mu^{2} \;\;\; \dots\dots\dots\dots\dots (2) \end{align*}
위의 $(1), (2)$를 $(3)$에 대입해보면
\begin{align*} E \left[\sum{(X_{i}-\overline{X})^{2}} \right] & = \sum{E\left(X_{i}^{2}\right)} - n E\left(\overline{X}^{2}\right) \\ & = \sum{\left(\sigma^{2} + \mu^{2}\right)} - n \left(\dfrac{\sigma^2}{n} + \mu^2\right) \\ & = n\sigma^{2} + n\mu^{2} - \sigma^{2} - n\mu^{2} \\ & = \left(n-1\right) \sigma^{2} \end{align*}
위는 식 (a)의 일부이므로 이를 온전한 식에 대입해보면,
\begin{eqnarray*}
E \left[\sum{(X_{i}-\overline{X})^{2}} \right] & = & (n-1) \sigma^{2} \\
\end{eqnarray*}
\begin{eqnarray*} E[s^{2}] & = & E \left[ \frac{\displaystyle\sum_{i=1}^{n}(X_{i}-\overline{X})^{2}}{n-1} \right] \\ & = & \dfrac{1}{n-1} E \left[\sum{(X_{i}-\overline{X})^{2}} \right] \\ & = & \dfrac{1}{n-1} (n-1) \sigma^{2} \\ & = & \sigma^{2} \end{eqnarray*}
그러므로, n-1로 나눠 준 샘플분산의 (sample's variance) 기대값은
\begin{eqnarray*}
E(s^2) = \sigma^{2}
\end{eqnarray*}
만약에 우리가 population의 variance를 구하듯이 n을 이용한다고 하면,
\begin{eqnarray*} E[s^{2}] & = & E \left[ \frac{\displaystyle\sum_{i=1}^{n}(X_{i}-\overline{X})^{2}} {n} \right], \;\;\; \text{note that we use n instead of n-1} \\ & = & \dfrac{1}{n} E \left[\sum{(X_{i}-\overline{X})^{2}} \right] \\ & = & \dfrac{1}{n} (n-1) \sigma^{2} \\ & = & \left(\dfrac{n-1}{n}\right) \sigma^{2} \\ \end{eqnarray*}
즉, 원래 $\sigma^2$ 값보다 조금 작은 값을 갖게 될 것이다 (이를 biased result라고 한다). 따라서 샘플을 취한 후에 모집단의 분산을 추정할 때에는 n 대신에 n-1을 사용하는 것이 맞다. 그렇다면 모집단의 분산을 구할 때는 n으로 (N으로) 나누어 주면 된다고 생각된다. 그러나 일반적으로 모집단의 분산을 구할 때에도 N-1로 나누어 구하게 된다. 이유는 모집단의 경우에 N이 충분히 큰 경우인데 이 때에는 N으로 나누어 주나, N-1로 나누어주나 큰 차이가 없기 때문이다. 따라서, R에서 분산을 구하는 var(x)에는 x의 성격에 상관없이 SS를 n-1로 나누어 분산을 구하게 된다.
tags
estimated_standard_deviation.1759211676.txt.gz · Last modified: by hkimscil