Differences

This shows you the differences between two versions of the page.

--- estimated_standard_deviation [2025/09/05 20:21] – [R에서 msr이 최소값이 되는 v값 효율적으로 찾기] hkimscil
+++ estimated_standard_deviation [2025/09/07 06:24] (current) – [수학적 증명] hkimscil
@@ Line 28: / Line 28: @@
 ====== 실험적, R에서 시뮬레이션으로 이해 ======
+아래 output의 코멘트를 읽을 것
 <code>
 rm(list=ls())
@@ Line 137: / Line 138: @@
 </WRAP>
 <WRAP half column>
-SS = sum(x-mean(x))^2 인데, mean(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 = 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을 공토으로 나누어준 값을 저장하고 비교하려는 것이다.
@@ Line 417: / Line 418: @@
 \end{eqnarray*}
+즉, 기울기가 0이 될때의 v값은 x집합의 평균값일 때이다.
 ===== R에서 msr이 최소값이 되는 v값 효율적으로 찾기 =====
-위에 따르면, 우리가  찾는 v 값은 샘플의 평균값이 ($\frac {\sum{x}}{n}$) 된다. 따라서, 평균값으로 MS값을 구하게 되면 언제나 가장 작은  값을 취하게 되는 결과를 갖는다. 이렇게 작은 값을 갖는 현상을 보정하려고 n 대신에 n보다 조금 작은 숫자인 n-1을 가지고 SS 부분을 나누어 준다.
+위의 미분을 이용한 방법을 써서 MS값이 최소값이 되는 순간을 R에서 구현해 볼 수 있다. 아래 아웃풋의 코멘트 부분을 읽을 것.
-위의 미분을 이용한 방법을 R에서 구현해 볼 수 있다. 이 때 위의 미분식을 사용한다.
 <del>[[why n-1 gradient explanation]]
 </del>
@@ Line 465: / Line 465: @@
 tail(vs)
-plot(msrs)
+plot(vs, msrs, type='b')
-plot(vs)
 # scaled
@@ Line 477: / Line 476: @@
 v.orig <- (v*sd(x))+mean(x)
 v.orig
+plot(vs.orig, msrs, type='b')
 </code>
@@ Line 528: / Line 529: @@
   * 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값들은 0.1씩 곱해주는 효과때문에 오른 쪽으로 옮겨가는 지점이 "<fc #ff0000>**점진적으로 줄어들게 되고**</fc>" 이 지점이 msr의 최소값이 되는 지점으로 가게 된다.
+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씩 곱해주는 효과때문에 오른 쪽으로 옮겨가는 지점이 "<fc #ff0000>**점진적으로 줄어들게 되고**</fc>" 이 지점이 msr의 최소값이 되는 지점으로 가게 된다.
-{{:pasted:20250905-190437.png}}
+{{:pasted:20250905-202627.png}}
+클릭하면 큰 이미지로 볼 수 있음
 </WRAP>
 </WRAP>
@@ Line 564: / Line 564: @@
   * 원래는 mean(x)값 근처의 값을 랜덤하게 골라야 하므로
   * ''v <- rnorm(1, mean(x), sd(x))'' 와 같이 써야 하지만 (현재의 경우, '' rnorm(1, 50, 5) ''가 된다)
-  * 그렇게 하질 않고 x 집합의 원소들을 표준점수화 한 후 '' zx <- (x-mean(x))/sd(x) ''
+  * 그렇게 하질 않고 x 집합의 원소들을 표준점수화 한 후 '' zx <- (x-mean(x))/sd(x) '' (이렇게 표준점수화 하면 x 변인의 평균이 0, 표준편차가 1 이 되는 집합으로 변한다)
-  * (이렇게 표준점수화 하면 x 변인의 평균이 0, 표준편차가 1 이 되는 집합으로 변한다)
   * '' v <- rnorm(1, 0, 1) ''로 구한다. 뒤의 인자인 0, 1 은 default이모로 생략.
   * 이렇게 하는 이유는 혹시 나중에 다른 x집합에 똑같은 작업을 하더라도 그 집합의 평균과 표준편차를 사용하지 않고
@@ Line 583: / Line 582: @@
 <WRAP half column>
 comment
-  * 이 0.1은 gradient function으로 구한 해단 v값에 대한 y 기울기 값을 (미분값)
+  * 이 0.1은 gradient function으로 구한 해당 v값에 대한 y 미분값을 (기울기 값을) 구한 후, 여기에 곱하기 위해서 지정한다.
-  * 구한 후, 여기에 곱하기 위해서 지정한다.
 </WRAP>
@@ Line 613: / Line 611: @@
 [1] 6.415945e-08 5.132756e-08 4.106205e-08 3.284964e-08 2.627971e-08 2.102377e-08
 >
-> plot(msrs)
+> plot(vs, msrs)
-> plot(vs)
+>
 >
 </code>
@@ Line 620: / Line 618: @@
 <WRAP half column>
 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의 영향으로 점진적으로 하가하여 최소값으로 도달한다.
+{{:pasted:20250905-214631.png}}
 </WRAP>
 </WRAP>
@@ Line 657: / Line 670: @@
 > v.orig
 [1] 50
+>
+> plot(vs.orig, msrs, type='b')
 >
 >
@@ Line 664: / Line 679: @@
 <WRAP half column>
 comment
+{{:pasted:20250905-231742.png}}
+만약에 처음에 구한 랜덤 v값이 평균의 오른 쪽에있었더라면, 아래 그림과 같이 평균에 접근했을 것이다.
+{{:pasted:20250905-231513.png}}
 </WRAP>
 </WRAP>
@@ Line 787: / Line 806: @@
 \end{eqnarray*}
-즉, 원래 $\sigma^2$ 값보다 조금 작은 값을 갖게 될 것이다 (이를 biased result라고 한다).
+즉, 원래 $\sigma^2$ 값보다 조금 작은 값을 갖게 될 것이다 (이를 biased result라고 한다). 따라서 샘플을 취한 후에 모집단의 분산을 추정할 때에는 n 대신에 n-1을 사용하는 것이 맞다. 그렇다면 모집단의 분산을 구할 때는 n으로 (N으로) 나누어 주면 된다고 생각된다. 그러나 일반적으로 모집단의 분산을 구할 때에도 N-1로 나누어 구하게 된다. 이유는 모집단의 경우에 N이 충분히 큰 경우인데 이 때에는 N으로 나누어 주나, N-1로 나누어주나 큰 차이가 없기 때문이다. 따라서, R에서 분산을 구하는 var(x)에는 x의 성격에 상관없이 SS를 n-1로 나누어 분산을 구하게 된다.
-====== R 에서 SS값이 최소인 v값을 찾기 ======
-<code>
-#library(ggplot2)
-#library(ggpmisc)
-rm(list=ls())
-rnorm2 <- function(n,mean,sd){
-  mean+sd*scale(rnorm(n))
-}
-# set.seed(191)
-nx <- 20
-mx <- 50
-sdx <- mx * 0.15
-x <- rnorm2(nx, mx, sdx)
-mean(x)
-sd(x)
-length(x)
-hist(x)
-x.span <- seq(from = mean(x)-6*sd(x),
-              to = mean(x)+6*sd(x),
-              by = 0.1)
-residuals <- function(x, v) {
-  return(x - v)
-}
-ssr <- function(x, v) {
-  residuals <- (x - v)
-  return(sum(residuals^2))
-}
-msr <- function(x, v) {
-  residuals <- (x - v)
-#  return((sum(residuals^2))/(length(x)-1))
-  return((mean(residuals^2)))
-}
-srs <- c() # sum of residuals
-ssrs <- c() # sum of square residuals
-msrs <- c() # mean square residuals = variance
-vs <- c() # the value of v in (x - v)
-for (i in x.span) {
-  res.x <- residuals(x,i)
-  srs.x <- sum(res.x)
-  ssr.x <- ssr(x,i)
-  msr.x <- msr(x,i)
-  srs <- append(srs, srs.x)
-  ssrs <- append(ssrs, ssr.x)
-  msrs <- append(msrs, msr.x)
-  vs <- append(vs, i)
-}
-plot(srs)
-plot(msrs)
-plot(srs)
-min(msrs)
-min.pos.msrs <- which(msrs == min(msrs))
-min.pos.msrs
-print(vs[min.pos.msrs])
-plot(vs, msrs)
-plot(vs, srs)
-# the above no gradient
-# mse 값으로 계산 rather than sse
-# 후자는 값이 너무 커짐
-gradient <- function(x, v){
-  residuals = x - v
-  dx = -2 * mean(residuals)
-  return(list("ds" = dx))
-} # function returns ds value
-residuals <- function(x, v) {
-  return(x - v)
-}
-ssr <- function(x, v) {
-  residuals <- (x - v)
-  return(sum(residuals^2))
-}
-msr <- function(x, v) {
-  residuals <- (x - v)
-  return((sum(residuals^2))/(length(x)-1))
-#  return(mean(residuals^2))
-}
-# pick one random v in (x-v)
-v <- rnorm(1)
-# Train the model with scaled features
-learning.rate = 1e-1
-ssrs <- c()
-msrs <- c()
-mres <- c()
-vs <- c()
-# Record Loss for each epoch:
-zx <- (x-mean(x))/sd(x)
-nlen <- 100
-for (epoch in 1:nlen) {
-  residual <- residuals(zx, v)
-  ssr.x <- ssr(zx, v)
-  msr.x <- msr(zx, v)
-  ssrs <- append(ssrs, ssr.x)
-  msrs <- append(msrs, msr.x)
-  grad <- gradient(zx, v)
-  step.v <- grad$ds * learning.rate
-  v <- v - step.v
-  vs <- append(vs, v)
-}
-tail(srs)
-tail(msrs)
-tail(ssrs)
-tail(vs)
-plot(srs)
-plot(msrs)
-plot(ssrs)
-plot(vs)
-# scaled
-v
-v.orig <- (v*sd(x))+mean(x)
-v.orig
-</code>
-====== output ======
-<code>
-> #library(ggplot2)
-> #library(ggpmisc)
->
-> rm(list=ls())
-> rnorm2 <- function(n,mean,sd){
-+   mean+sd*scale(rnorm(n))
-+ }
->
-> # set.seed(191)
-> nx <- 20
-> mx <- 50
-> sdx <- mx * 0.15
-> x <- rnorm2(nx, mx, sdx)
->
-> mean(x)
-[1] 50
-> sd(x)
-[1] 7.5
-> length(x)
-[1] 20
-> hist(x)
->
-> x.span <- seq(from = mean(x)-6*sd(x),
-+               to = mean(x)+6*sd(x),
-+               by = 0.1)
->
-> residuals <- function(x, v) {
-+   return(x - v)
-+ }
->
-> ssr <- function(x, v) {
-+   residuals <- (x - v)
-+   return(sum(residuals^2))
-+ }
->
-> msr <- function(x, v) {
-+   residuals <- (x - v)
-+ #  return((sum(residuals^2))/(length(x)-1))
-+   return((mean(residuals^2)))
-+ }
->
-> srs <- c() # sum of residuals
-> ssrs <- c() # sum of square residuals
-> msrs <- c() # mean square residuals = variance
-> vs <- c() # the value of v in (x - v)
->
-> for (i in x.span) {
-+   res.x <- residuals(x,i)
-+   srs.x <- sum(res.x)
-+   ssr.x <- ssr(x,i)
-+   msr.x <- msr(x,i)
-+   srs <- append(srs, srs.x)
-+   ssrs <- append(ssrs, ssr.x)
-+   msrs <- append(msrs, msr.x)
-+   vs <- append(vs, i)
-+ }
-> plot(srs)
-> plot(msrs)
-> plot(srs)
->
-> min(msrs)
-[1] 53.4375
-> min.pos.msrs <- which(msrs == min(msrs))
-> min.pos.msrs
-[1] 451
-> print(vs[min.pos.msrs])
-[1] 50
->
-> plot(vs, msrs)
-> plot(vs, srs)
->
->
->
-> # the above no gradient
-> # mse 값으로 계산 rather than sse
-> # 후자는 값이 너무 커짐
->
-> gradient <- function(x, v){
-+   residuals = x - v
-+   dx = -2 * mean(residuals)
-+   return(list("ds" = dx))
-+ } # function returns ds value
->
-> residuals <- function(x, v) {
-+   return(x - v)
-+ }
->
-> ssr <- function(x, v) {
-+   residuals <- (x - v)
-+   return(sum(residuals^2))
-+ }
->
-> msr <- function(x, v) {
-+   residuals <- (x - v)
-+   return((sum(residuals^2))/(length(x)-1))
-+ #  return(mean(residuals^2))
-+ }
->
-> # pick one random v in (x-v)
-> v <- rnorm(1)
-> # Train the model with scaled features
-> learning.rate = 1e-1
->
-> ssrs <- c()
-> msrs <- c()
-> mres <- c()
-> vs <- c()
-> # Record Loss for each epoch:
-> zx <- (x-mean(x))/sd(x)
->
-> nlen <- 100
-> for (epoch in 1:nlen) {
-+   residual <- residuals(zx, v)
-+   ssr.x <- ssr(zx, v)
-+   msr.x <- msr(zx, v)
-+   ssrs <- append(ssrs, ssr.x)
-+   msrs <- append(msrs, msr.x)
-+
-+   grad <- gradient(zx, v)
-+
-+   step.v <- grad$ds * learning.rate
-+   v <- v - step.v
-+   vs <- append(vs, v)
-+ }
->
-> tail(srs)
-[1] -890 -892 -894 -896 -898 -900
-> tail(msrs)
-[1] 1 1 1 1 1 1
-> tail(ssrs)
-[1] 19 19 19 19 19 19
-> tail(vs)
-[1] 2.936258e-11 2.349006e-11 1.879204e-11 1.503363e-11 1.202690e-11 9.621523e-12
->
-> plot(srs)
-> plot(msrs)
-> plot(ssrs)
-> plot(vs)
-> # scaled
-> v
-[1] 9.621523e-12
-> v.orig <- (v*sd(x))+mean(x)
-> v.orig
-[1] 50
->
->
-</code>
 ====== tags ======
 {{tag>"research methods" "조사방법론" "표준편차" "n-1" "자유도" "degrees of freedom" "n-1" "표준오차"}}