문제.

우리는 모집단의 (population) 평균값을 알고 있다면 샘플의 분산값을 다음과 같이 구할 수 있다. 그리고, 이것을 모집단의 분산값으로 추정할 수 있다.

\begin{eqnarray*} \widehat{\sigma^2} = \dfrac {\displaystyle\sum_{i=1}^{n}{(X_{i}-\mu)}} {n} \end{eqnarray*}

그러나, 현재 우리가 가지고 있는 것은 샘플 밖에 없다. 즉, 모집단의 평균은 알지 못하는 상태이기에 모집단 분산을 추정하는 계산에 사용할 수 없다. 따라서 샘플의 평균을 사용한다. 그런데, 샘플의 평균을 사용할 때는 분모에 N 대신에 n-1을 사용해야 한다. 왜 n-1을 사용하는것이 모집단의 분산값 추정에 도움이 되는가가 문제이다.

\begin{eqnarray*} \widehat{\sigma^2} \neq \frac {\displaystyle\sum_{i=1}^{n}{(X_{i}-\overline{X})}} {n} \end{eqnarray*}

\begin{eqnarray*} \widehat{\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 (개인의 점수를 평균으로 추측했을 때 틀린 만큼의 에러값) 하였다. 이렇게 평균을 가지고 개인점수를 예측하는 것이 가장 작은 오차를 갖는 방법이라고 하였다 (개인점수 중에서 평균이 제일 많이 나오므로, 평균으로 개인점수를 예측하면 제일 들 틀린다). 따라서 어느 한 집합에서 (샘플에서) 개인점수에서 평균을 빼고 이를 제곱하여 모두 더한 값은 (평균이 아닌) 다른 값으로 구한 값보다 항상 작게 된다 (최소값을 갖는다).

이것이 의미하는 것은 샘플의 평균으로 분산을 구하면 분자에 해당하는 SS값이 최소값이 되어 버리므로 예측하는 분산의 숫자가 항상 너무 작은 숫자가 되어 버린다.

이를 그림으로도 설명할 수 있다. 아래에서 녹색의 세로선은 모집단의 평균값이고, 붉은색의 세로선은 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
>

붉은 선 s1 평균 (10)
초록 선 모집단 평균 12 (혹은 그냥 12 가치라고 해도 됨)
세개 검은 선은 각각 6, 9, 15

붉은 선에서 세개 검은 선까지의 길이는 4, 1, 5가 되고 이를 더하면 10
초록 선에서 세개 검은 선까지의 길이는 6, 3, 3이 되고 이를 더하면 12

아래 output의 코멘트를 읽을 것

다음으로 MS값을 구하는 식인

• y = sum( (x-v)^2 ) / n 값을
• v에 대해서 v를 가지고 y를 미분하면 (dy/dv)
• 변화하는 v값마다의 기울기 값을 구할 수 있는데
• 이 값이 0이 되는 지점의 v값이 무엇인지를 구하면 위의 R코드에서 구한 값을 구할 수 있게 된다. 아래는 그 과정이다.
• 그래프에서 가장 작은 기울기값을 갖는 v 값을 구한다고 (derivatives) 가정하고 이해를 하면 수학적으로 이해할 수 있다. ¹⁾

\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집합의 평균값일 때이다.

위의 미분을 이용한 방법을 써서 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)
[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)
> 
> 

> # 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')
> 
> 
> 

그렇다면 왜 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 (분산)의 연산에 과한 법칙으로는 ²⁾

우리가 알고자 하는 것은 아래의 식이 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)$에서

위의 $(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로 나누어 분산을 구하게 된다.