Aim: เพื่อจะเข้าใจที่มาของ confidence interval อย่างถ่องแท้ จึงขอทุ่มเทเวลาเพื่อทบทวน principle ของสถิติและความน่าจะเป็นเสียหน่อย post นี้จะ nerd มาก จริงๆทำไว้ให้ตัวเองกลับมาอ่านเองในภายหลังแล้วเข้าใจ

ขอย้อนกลับจากปลายทาง เรารู้ว่าเมื่อเราทำการศึกษาหนึ่งๆขึ้นมาโดยการสุ่มกลุ่มตัวอย่าง เราไม่สามารถมั่นใจได้ 100% ว่าคำตอบของเราตรงกับความจริงแท้ในธรรมชาติ เราจึงคำนวณช่วงตัวเลขขึ้นมาช่วงหนึ่ง เพื่อจะบอกว่าเรามั่นใจ 95% ว่า่ความจริงแท้ต้องอยู่ในช่วงนี้แหละ ช่วงนั้นคือ 95%CI ซึ่ง = Mean ± 1.96 * sd/Sqrt(n)

คำถามคาใจก็คือ ทำไม ??? สูตรถึงเป็นเช่นนี้ แล้วมันใช้ได้จริงๆเหรอ เพื่อจะตอบคำถามนี้ จึงต้องย้อนกลับไปสู่หลักพื้นฐานที่สุด

    1. ค่า z กับ Standard normal distribution คืออะไร?
      normal distribution ที่ standard คือชุดข้อมูลที่มี mean = 0 และ sd = 1 ที่เป็นเลขนี้ก็เพราะมันจำง่ายดีซึ่งกราฟจะประมาณนี้ (แกน y เป็นความถี่/ความหนาแน่น ของข้อมูลที่มีค่าตามแกน x ภาษา stat เรียก density function)และเราสามารถระบุข้อมูลแต่ละตัวในกราฟได้ด้วย +/- จำนวนเท่าของ sd ดังภาพซึ่งจริงๆแล้ว จำนวนเท่าที่ว่า ก็คือค่า z นั่นเอง
      z1
    2. ข้อดีของกราฟมาตรฐานนี้คือ ในแต่ละค่า z เราสามารถเทียบตารางได้ว่าค่า z นี้อยู่ percentile ที่เท่าไรในชุดข้อมูล สังเกตหัวตารางจะเขียนว่า 1-tail z ซึ่งแปลผลตรงไปตรงมาแบบที่ผมบอก แต่ในทางแพทย์เราจะใช้แบบ 2-tail มากกว่า ซึ่งก็แปลผลเทียบกันไม่ยาก (z = 1.96 จะตรงกับ 97.5% และ 95% พอดี และเป็นเลขที่เรานิยมใช้)
    3. สมมติเรามีข้อมูลความสูงของประชากรมนุษย์ ซึ่งเป็น normal distribution mean=170 sd=10 แล้วเราสงสัยว่า ผู้ชายสูง 185 อยู่ percentile ที่เท่าไรในประชากรนี้ สิ่งที่เราต้องทำคือการคำนวณ z score ของความสูง 185 = เท่าไร ซึ่ง สูตรก็คือ

z = (x-mean)/sd

 ในกรณีตัวอย่างนี้จะพบว่า z = (185-170)/10 = 1.5 ซึ่งเราสามารถเทียบได้เป็น percentile ที่ 0.933 (คำถามนี้ต้องใช้ 1-tail z ส่วนวิธีเทียบโดยทั่วไปใช้ z table ส่วนผมใช้ function pnorm() ใน R)

 4. ในทางกลับกัน เราอาจจะอยากรู้ว่าความสูง mean +/- เท่าไรจึงจะครอบคลุม 95% ของประชากรทั้งหมด เราก็สามารถคำนวณได้ไม่ยากโดยกลับสมการนิดหน่อย ดังนี้

  x = mean ± z×sd

เมื่อ z ในที่นี้คือ 2 tail z ที่ cover 95% ของ standard curve ซึ่งก็คือ 1.96

==========ตัดจบครึ่งนึงก่อน สรุปใจความดังนี้===========
1. standard normal distribution มี mean=0 sd=1 และค่า z ก็คือจำนวนเท่าของ sd ที่ห่างออกจาก mean นั่นเอง
2. ถ้ามีชุดข้อมูลที่เป็น normal distribution อยู่ เราสามารถเทียบค่า x ใดใดในชุดข้อมูลของเราเป็นค่า z ได้ด้วย z = (x-mean)/sd และกลับกัน x = mean ± z×sd
3. ค่า z ที่ได้สามารถนำไปคิดหา percentile ได้ แล้วแต่ว่าเป็น 1-tail หรือ 2-tail

ดังนั้นเงื่อนไขสำคัญของการใช้ z ก็คือความเป็น normal distribution ของชุดข้อมูลนั่นเองครับ
ทีนี้เรามีสิ่งที่คลับคล้ายจะเป็น 95% CI แล้ว แต่ต้องบอกไว้ไม่ให้สับสนว่ามันยังไม่ใชj 95%CI นะครับ

ตอนต่อไปจะมาต่อกันที่ว่าทำไมค่าจริงส่วนใหญ่ในธรรมชาติจึงเป็น normal distribution และทำไมเราถึงนำหลักการของ z นี้มาคิดหา 95%CI ได้ในกรณีส่วนใหญ่ *แม้ข้อมูลจะไม่เป็น normal distribution ก็ตาม*

———————–R Code—————————
rnorm ใช้สร้าง set ข้อมูลที่เป็น normal distribution
dnorm ใช้สร้าง density function ก็คือแกน y ของ curve นั่นเอง เมื่อกำหนดค่า x ต่างๆ
pnorm ใช้คำนวณว่าค่า x (หรือ z) cover กี่ percentile
qnorm เป็น inverse ของ pnorm คือเราป้อน 1-tail percentile (เช่น 0.95) แล้วจะได้ค่า x (หรือ z) กลับมา

# means comment
# <- means 'assign value' ~ =
#plot1 ::
x <- seq(from=-4,to=4,by=0.001)
y <- dnorm(x,0,1)
plot(y=y,x=x,type="l",main="Standard Normal Distribution (mean=0 sd=1)", xlab="Z score = amount of SD away from mean")
  lines(y=c(0,0.3),x=c(-2,-2),lwd=2,lty=3,col="#6666BB")
  text(y=0.32,x=-2,"-2SD",col="#6666BB")
  lines(y=c(0,0.3),x=c(1,1),lwd=2,lty=3,col="#33AA33")
  text(y=0.32,x=1.2,"+1SD",col="#33AA33")

#plot2 ::
x <- seq(from=-4,to=4,by=0.001)
y <- dnorm(x,0,1)
plot(y=y,x=x,type="l",main="1 tail Z",xlab="Z score = amount of SD away from mean")
  seg1 <- seq(from=-4, to=1.96, by=0.01)
  polygon(c(-4,seg1,1.96),c(0,dnorm(seg1),0),col="#FFAA44")
  text(y=0.25,x=2.5,"Z at 1.96 \n cover 97.5%")

#plot3 ::
x <- seq(from=-4,to=4,by=0.001)
y <- dnorm(x,0,1)
plot(y=y,x=x,type="l",main="2 tail Z",xlab="Z score = amount of SD away from mean")
  seg1 <- seq(from=-1.96, to=1.96, by=0.01)
  polygon(c(-1.96,seg1,1.96),c(0,dnorm(seg1),0),col="#44AAFF")
  text(y=0.25,x=2.5,"Z at 1.96 \n cover 95%")

#plot4 ::
hist(Man,50,xlab="Height",col="#DDDD88")
  lines(y=dnorm(140:200,170,10)*2000,x=140:200)
  lines(y=c(0,80),x=c(185,185),lwd=3,lty=3,col="red")
  text(x=192,y=60,"Height 185\nis at 93 percentile")

#plot5 ::
  Man <- rnorm(2000,mean=170,sd=10)
  hist(Man,50,xlab="Height",col="#DDDD88",ylim=c(0,100))
  lines(y=dnorm(140:200,170,10)*2000,x=140:200)
  lines(y=c(0,80),x=c(qnorm(0.975,170,10),qnorm(0.975,170,10)),lwd=3,lty=3,col="red")
  lines(y=c(0,80),x=c(qnorm(0.025,170,10),qnorm(0.025,170,10)),lwd=3,lty=3,col="red")
  text(x=170,y=100,"Heigth 150.4-189.6 (170+/-19.6) covers 95% of population")

Advertisements