Differences

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

--- z-test_and_t-test_simulation_in_r [2024/09/12 11:51] – created hkimscil
+++ z-test_and_t-test_simulation_in_r [2024/09/13 10:13] (current) – hkimscil
@@ Line 1: / Line 1: @@
 <code>
-n.ajstu <- 100000
+n.ajstu <- 100000 # 모집단 100000
-mean.ajstu <- 110
+mean.ajstu <- 110 # 모집단 평균 110이라 가정
-sd.ajstu <- 10
+sd.ajstu <- 10 # 표준편차 10이라 가정
 set.seed(1024)
+# rnorm2 펑션을 이용해서 모집단 만들기
+rnorm2 <- function(n,mean,sd) { mean+sd*scale(rnorm(n)) }
 ajstu <- rnorm2(n.ajstu, mean=mean.ajstu, sd=sd.ajstu)
@@ Line 10: / Line 12: @@
 sd(ajstu)
 var(ajstu)
 iter <- 10000 # # of sampling
-n.4 <- 4
+# 이 후 n.# 은 샘플사이즈를 말함
+n.4 <- 4   # 샘플사이즈 4
 means4 <- rep (NA, iter)
+# sample 평션을 이용해서 4개 샘플을 ajoust에서
+# 취한 후, 평균을 내서 이를 means[i]에 기록한다
+# 이것을 iteration 숫자만큼 한다
 for(i in 1:iter){
   means4[i] = mean(sample(ajstu, n.4))
@@ Line 54: / Line 61: @@
   means2500[i] = mean(sample(ajstu, n.2500))
 }
-# n.3600 <- 3600
-# means3600 <- rep (NA, iter)
-# for(i in 1:iter){
-#  means3600[i] = mean(sample(ajstu, n.3600))
-# }
 h4 <- hist(means4)
@@ Line 68: / Line 69: @@
 h1600 <- hist(means1600)
 h2500 <- hist(means2500)
-# h3600 <- hist(means3600)
 plot(h4, ylim=c(0,3000), col="red")
@@ Line 76: / Line 76: @@
 plot(h900, add = T, col="yellow")
 # plot(h2500, add = T, col="brown")
-# plot(h3600, add = T, col="orange")
-sss <- c(4,25,100,400,900,1600,2500,3600) # sss sample sizes
+# 아래는 종합적인 출력을 위해서 좀 복잡하게
-ses <- rep (NA, length(sss)) # std errors
+# 코딩을 하였지만, 요지는 특정 샘플사이즈에서
+# standard deviation of sample means 즉
+# standard error를 구하는 과정
+sss <- c(4,25,100,400,900,1600,2500,3600) # 8 종류의 샘플사이즈
+ses <- rep (NA, length(sss)) # 각 샘플 사이즈에 대응하는 standard error 값 저장을 위한 공간
 for(i in 1:length(sss)){
   ses[i] = sqrt(var(ajstu)/sss[i])
@@ Line 85: / Line 88: @@
 data.frame(ses)
-se.1 <- ses
+se.1 <- ses # one standard deviation
-se.2 <- 2 * ses
+se.2 <- 2 * ses  # two standard deviation
 # alt.2 <- qnorm(.975)
 # se.2 <- alt.2 * ses
-se.2
+# 2 standard deviation of sample means를 (se) 썼을 때
+# 아랫쪽, 위쪽 경계
 lower.s2 <- mean(ajstu)-se.2
 upper.s2 <- mean(ajstu)+se.2
-sample <- rnorm2(100, 112.41, 10)
+# 이제 머리가 좋아지는 약을 먹은 하나의 샘플을
+# 가정하고 이 샘플의 평균이  113라고 가정하면
+smean <- 113
-ses <- round(ses,3)
+# z-test를 위한 계산에서 분자는 difference라고 했고
-smean <- mean(sample)
+# 이를 구한 것
 diff <- smean - mean(ajstu)
 diff
-# z-test
+# z-test diff를 stanard error 값으로 나눈 것
+# 이 diff를 샘플사이즈 별로 (8종류) 구한 ses로
+# 나누어 본다
 z.sc <- (diff)/ses
 z.sc
-# z.crit <- 2
+# 위의 z.sc 점수와 강사가 이야기한 2점과 비교
-z.crit <- qnorm(.975)
+# 이 2점을 (standard error 2개) z critical score라고
-z.crit
+# 부른다
+z.crit <- 2
+# z.crit <- qnorm(.975)
+# z.crit
+# 위에서 구한 모든 것들을 종합해서 출력
 data.frame(cbind(sss, smean, mean(ajstu), diff, ses,
                  lower.s2, upper.s2,
                  z.sc, z.crit, z.sc>z.crit))
-# t-test
-ses <- round(ses,3)
+# t-test 의 상황에서는 z critical value가 2였던 것에
-smean <- smean
+# 비해서 샘플사이즈에 따라서 critical value가 달라지게
+# 된다. 이를 구한 것이 아래의 qt평션을 이용한 t.crit
 diff <- smean - mean(ajstu)
 dfs <- sss - 1
 t.sc <- diff/ses
-perc <- round(pt(t.sc, dfs, lower.tail = FALSE),9)
+# 아래는 샘플사이즈 별로 각각 구한
+# t calculated value 들에 해당하는
+# t distribution에서의 퍼센트
+perc <- round(pt(t.sc, dfs, lower.tail = FALSE), 5)
+# 아래는 degrees of freedom이 dfs 이고
+# 97.5 퍼센트일 때 봐야 할 (구해야 할)
+# t 값을 구한 것. 이것을 샘플사이즈 별로 한 것
+# 4, 25, . . . .
 t.crit <- qt(.975, dfs, lower.tail = TRUE)
+# 모두 모아서 출력
 data.frame(cbind(sss, smean, mean(ajstu), diff, ses,
                  t.sc, dfs, perc, t.crit, t.sc>t.crit))
@@ Line 127: / Line 151: @@
 qnorm(.975)
+# 하나만 따로 생각해 봄
+# 평균 113, 표준편차가 10인 100개로 이루어진
+# 샘플을 구하여 sa 에 지정
+sa <- rnorm2(100, 113, 10)
+# t test 수행
 ?t.test
 t.test(sample, mu=mean(ajstu))
+# two sample t test 의 경우
 sa <- rnorm2(100, 110, 10)
-sb <- rnorm2(100, 115, 10)
+sb <- rnorm2(100, 113, 10)
 sa.mean <- mean(sa)
@@ Line 146: / Line 177: @@
 ssa <- var(sa)*dfa
 ssb <- var(sb)*dfb
+ssa
+ssb
+ssa+ssb
+dfa+dfb
 var.pooled <- (ssa+ssb)/(dfa+dfb)
 vp <- var.pooled
+vp
 se <- sqrt((vp/na)+(vp/nb))
-diff <- sb.mean - sa.mean
+diff <- sa.mean - sb.mean
 t.cal <- diff/se
 se
 diff
 t.cal
+t.crit <- qt(.975, 198)
+t.crit
+pt(t.cal, 198)*2
+t.test(sa, sb)
-t.test(s1, s2)
 </code>