บทที่ 5
Probability and
Distributions
5.1 Randon
Sampling
ในยุคต้นๆของทฤษฎีความน่าจะเป็น เรื่องที่ศึกษามักเป็นเรื่องเกี่ยวกับเกมส์และการพนัน ดังนั้นตัวอย่างสุ่มจึงมักเกี่ยวข้องกับไพ่ที่หยิบขึ้นมาจากสำรับ หรือการหยิบลูกบอลจากกล่อง เป็นต้นในR เราสามารถจำลองสถานการณ์ (simulation) เหล่านี้ได้ด้วยฟังก์ชันsample() ซึ่งมีรูปแบบดังนี้
sample(x, size, replace = FALSE, prob = NULL)
เมื่อx คือเวคเตอร์ข้อมูลที่ประกอบด้วยค่าของข้อมูลที่ เป็นไปได้ทั้งหมด,size คือจำนวนตัวอย่างสุ่มที่ต้องการ, replace=เป็น argument ที่ใช้เพื่อกำหนดว่าการสุ่มตัวอย่างที่กำลังจะทำนั้นเป็ นการสุ่มตัวอย่างแบบใส่คืน (with replacement) หรือไม่ใส่คืน (without replacement) และprob= เป็ น argument ที่ใช้เพื่อกำหนดค่าความน่าจะเป็นของข้อมู ลแต่ละค่าในเวคเตอร์ x หากเราต้องการสุ่มตัวเลข 5 ตัว จากตัวเลขชุดหนึ่งที่มีค่าตั้งแต่ 1 ถึง 40 สามารถทำได้ดังนี้
> sample(1:40, 5)
[1] 20 13 23 18 14
argument แรก1:40 หมายถึงเวกเตอร์ของตัวเลข 1 ถึง 40 จำนวน 40 ค่า ซึ่งถือเป็นประชากรของข้อมูล ส่วนargument ที่สองคือ 5 หมายถึงขนาดของตัวอย่างสุ่มที่เราต้องการ ในความเป็นจริงเราสามารถใช้ sample(40,5)ก็ได้ซึ่งจะให้ผลเหมือนกัน เนื่องจากตัวเลข 40 หมายถึงช่วงของเลขจำนวนเต็มที่มีค่าตั้งแต่ 1 ถึง 40 ให้สังเกตว่า default ของฟังก์ชันsample() จะเป็ นการสุ่มแบบไม่ใส่คืน ดังนั้นหากต้องการให้เป็ นการสุ่มแบบใส่คืนสามารถเปลี่ยนได้โดยเปลี่ยน argumentreplace=TRUE เช่น การจำลองข้อมูลของผลของการโยนเหรียญ 1เหรียญ 10 ครั้ง
> sample(c("H", "T"), 10, replace=TRUE)
[1] "T" "T" "T" "H" "H" "H" "T" "H" "H" "T"
การจำลองข้อมูลในตัวอย่างข้างต้นเป็นกรณีที่ ความน่าจะเป็นในการเกิดแต่ละผลลัพธ์เท่ากัน แต่เราสามารถจำลองข้อมูลที่ความน่าจะเป็นของการเกิดแต่ละผลลัพธ์ไม่เท่ากันได้ด้วย เช่น การจำลองข้อมูลผลการผ่าตัด ซึ่งผลของการผ่าตัดมีสองแบบคือ สำเร็จ หรือล้มเหลว โดยที่ความน่าจะเป็ นที่การผ่าตัดจะประสบความสำเร็จคือ 90 เปอร์เซ็นต์ ทำได้ดังนี้
> sample(c("succ", "fail"), 10, replace=TRUE, prob = c(0.90, 0.10))
[1] "succ" "succ" "succ" "succ" "succ" "succ" "succ" "succ" "succ" "succ"
จากที่เราทราบดีอยูแล้วว่า เราสามารถอธิบายลักษณะของตัวเลขสุ่มได้ด้วยการแจกแจง (distribution) ซึ่งการแจกแจงก็คือฟังก์ชันที่ ใช้ในการหาค่าความน่าจะเป็นเมื่อตัวแปรสุ่มมีค่าใดค่าหนึ่ง หรือมีค่าอยู่ ในช่วงใดๆ เช่นP(a≤ 67 X≤ b) เป่นต้นซึ่งในกรณีของตัวแปรสุ่มแบบต่อเนื่องฟังก์ชันการแจกแจงของตัวแปรสุ่มจะอยู่ในรูปของฟังก์ชัน
ความหนาแน่น (density function) แต่ถ้าตัวแปรสุ่มเป็นแบบไม่ต่อเนื่องฟังก์ชันจะอยู่ในรูปของฟังก์ชันการแจกแจง(distribution function) เช่นP(X= k) เป็ นต้น ซึ่งR มีคำสั่งที่ใช้ในการสุ่มค่าสังเกตหรือข้อมูลที่มีการแจกแจงแบบต่างๆได้ ต่อไปนี้เป็นตัวอย่างของการใช้ R ในการสุ่มตัวอย่างจากการแจกแจงแบบต่างๆ
Uniform.รูปแบบของคำสั่งที่ใช้ในการสุ่มตัวอย่างจากการแจกแจงแบบ uniform คือ
runif(n, min=0, max=1)
ตัวอย่างเช่น หากเวลาในการรอสัญญานไฟจราจรมีการแจกแจงแบบเอกรูป (uniform distribution) ในช่วง [0,2]
> runif(1,0,2) # time at light
[1] 0.7325975
> runif(5,0,2) # time at 5 lights
[1] 0.5658379 0.1512934 1.3920678 1.6725483 0.9159859
> runif(5) # 5 random number in [0,1]
[1] 0.3441802 0.2687535 0.5044450 0.8939212 0.7787676
ในยุคต้นๆของทฤษฎีความน่าจะเป็น เรื่องที่ศึกษามักเป็นเรื่องเกี่ยวกับเกมส์และการพนัน ดังนั้นตัวอย่างสุ่มจึงมักเกี่ยวข้องกับไพ่ที่หยิบขึ้นมาจากสำรับ หรือการหยิบลูกบอลจากกล่อง เป็นต้นในR เราสามารถจำลองสถานการณ์ (simulation) เหล่านี้ได้ด้วยฟังก์ชันsample() ซึ่งมีรูปแบบดังนี้
sample(x, size, replace = FALSE, prob = NULL)
เมื่อx คือเวคเตอร์ข้อมูลที่ประกอบด้วยค่าของข้อมูลที่ เป็นไปได้ทั้งหมด,size คือจำนวนตัวอย่างสุ่มที่ต้องการ, replace=เป็น argument ที่ใช้เพื่อกำหนดว่าการสุ่มตัวอย่างที่กำลังจะทำนั้นเป็ นการสุ่มตัวอย่างแบบใส่คืน (with replacement) หรือไม่ใส่คืน (without replacement) และprob= เป็ น argument ที่ใช้เพื่อกำหนดค่าความน่าจะเป็นของข้อมู ลแต่ละค่าในเวคเตอร์ x หากเราต้องการสุ่มตัวเลข 5 ตัว จากตัวเลขชุดหนึ่งที่มีค่าตั้งแต่ 1 ถึง 40 สามารถทำได้ดังนี้
> sample(1:40, 5)
[1] 20 13 23 18 14
argument แรก1:40 หมายถึงเวกเตอร์ของตัวเลข 1 ถึง 40 จำนวน 40 ค่า ซึ่งถือเป็นประชากรของข้อมูล ส่วนargument ที่สองคือ 5 หมายถึงขนาดของตัวอย่างสุ่มที่เราต้องการ ในความเป็นจริงเราสามารถใช้ sample(40,5)ก็ได้ซึ่งจะให้ผลเหมือนกัน เนื่องจากตัวเลข 40 หมายถึงช่วงของเลขจำนวนเต็มที่มีค่าตั้งแต่ 1 ถึง 40 ให้สังเกตว่า default ของฟังก์ชันsample() จะเป็ นการสุ่มแบบไม่ใส่คืน ดังนั้นหากต้องการให้เป็ นการสุ่มแบบใส่คืนสามารถเปลี่ยนได้โดยเปลี่ยน argumentreplace=TRUE เช่น การจำลองข้อมูลของผลของการโยนเหรียญ 1เหรียญ 10 ครั้ง
> sample(c("H", "T"), 10, replace=TRUE)
[1] "T" "T" "T" "H" "H" "H" "T" "H" "H" "T"
การจำลองข้อมูลในตัวอย่างข้างต้นเป็นกรณีที่ ความน่าจะเป็นในการเกิดแต่ละผลลัพธ์เท่ากัน แต่เราสามารถจำลองข้อมูลที่ความน่าจะเป็นของการเกิดแต่ละผลลัพธ์ไม่เท่ากันได้ด้วย เช่น การจำลองข้อมูลผลการผ่าตัด ซึ่งผลของการผ่าตัดมีสองแบบคือ สำเร็จ หรือล้มเหลว โดยที่ความน่าจะเป็ นที่การผ่าตัดจะประสบความสำเร็จคือ 90 เปอร์เซ็นต์ ทำได้ดังนี้
> sample(c("succ", "fail"), 10, replace=TRUE, prob = c(0.90, 0.10))
[1] "succ" "succ" "succ" "succ" "succ" "succ" "succ" "succ" "succ" "succ"
จากที่เราทราบดีอยูแล้วว่า เราสามารถอธิบายลักษณะของตัวเลขสุ่มได้ด้วยการแจกแจง (distribution) ซึ่งการแจกแจงก็คือฟังก์ชันที่ ใช้ในการหาค่าความน่าจะเป็นเมื่อตัวแปรสุ่มมีค่าใดค่าหนึ่ง หรือมีค่าอยู่ ในช่วงใดๆ เช่นP(a≤ 67 X≤ b) เป่นต้นซึ่งในกรณีของตัวแปรสุ่มแบบต่อเนื่องฟังก์ชันการแจกแจงของตัวแปรสุ่มจะอยู่ในรูปของฟังก์ชัน
ความหนาแน่น (density function) แต่ถ้าตัวแปรสุ่มเป็นแบบไม่ต่อเนื่องฟังก์ชันจะอยู่ในรูปของฟังก์ชันการแจกแจง(distribution function) เช่นP(X= k) เป็ นต้น ซึ่งR มีคำสั่งที่ใช้ในการสุ่มค่าสังเกตหรือข้อมูลที่มีการแจกแจงแบบต่างๆได้ ต่อไปนี้เป็นตัวอย่างของการใช้ R ในการสุ่มตัวอย่างจากการแจกแจงแบบต่างๆ
Uniform.รูปแบบของคำสั่งที่ใช้ในการสุ่มตัวอย่างจากการแจกแจงแบบ uniform คือ
runif(n, min=0, max=1)
ตัวอย่างเช่น หากเวลาในการรอสัญญานไฟจราจรมีการแจกแจงแบบเอกรูป (uniform distribution) ในช่วง [0,2]
> runif(1,0,2) # time at light
[1] 0.7325975
> runif(5,0,2) # time at 5 lights
[1] 0.5658379 0.1512934 1.3920678 1.6725483 0.9159859
> runif(5) # 5 random number in [0,1]
[1] 0.3441802 0.2687535 0.5044450 0.8939212 0.7787676
เพื่อแสดงให้เห็นถึงลักษณะการแจกแจงแบบ
uniform ที่มีค่าต่ำสุดเท่ากับ 0 และค่าสูงสุดเท่ากับ 1 ทำได้ดังนี้
> x = runif(100)
> hist(x,prob=TRUE, col=gray(.9), main="Uniform on [0,1]")
> curve(dunif(x,0,1),add=TRUE)
Normal.การแจกแจงแบบปกตินั้นมีพารามิเตอร์ที่เกี่ยวข้องอยู่ 2 ตัว คือµ และσ เช่น IQ มีการแจกแจงปกติด้วยแจงแบบปกติโดยมีค่าเฉลี่ ยเท่ากับ 280 และค่าเบี่ ยงเบนมาตรฐานเท่ากับ 10 เราสามารถทำให้ค่าของตัวแปรสุ่มปกตินี้เป็ นค่ามาตรฐานZ= (X− µ)/σได้ ซึ่งค่ามาตรฐานนี้จะมีการแจกแจงแบบปกติมาตรฐานที่มีค่าเฉลี่ยเท่ากับ 0 และค่าเบี่ยงเบนมาตรฐานเท่ากับ
rnorm(n, mean=0, sd=1)
ตัวอย่างต่อไปนี้แสดงการสุ่มตัวอย่างจากการแจกแจงแบบปกติ
> rnorm(1,100,16) # an IQ score
[1] 122.7393
> rnorm(1,mean=280,sd=10) # how long for a baby (10 days early)
[1] 276.9401
เพื่อให้เห็นถึงลักษณะของการแจกแจงแบบปกติที่มีค่าเฉลี่ย 0 และค่าเบี่ยงเบนมาตรฐานเท่ากับ 1 ทำได้ดังนี้
> x = rnorm(100)
> hist(x,prob=TRUE, col= "blanchedalmond", main="Normal(0,1)")
> curve(dnorm(x),add=TRUE)
ได้กราฟดังรูปที่ 5.2
รูปที่ 5.2: Normal(0,1)
Binomial.ตัวแปรสุ่มทวินามเป็ นตัวแปรสุ่มแบบไม่ต่อเนื่อง โดยที่ค่าของตัวแปรสุ่มทวินามคือจำนวนของความสำเร็จที่เกิดขึ้นในการทดลองที่ เป็นอิสระกัน n ครั้ง เราสามารถสุ่มตัวอย่างจากการแจกแจงทวินามได้โดยใช้คำสั่ง rbinom(n, size, prob) ในการสร้างตัวเลขสุ่ มทวินามทำได้ดังนี้
> n = 10; p=.5
> rbinom(1,n,p) # 6 successes in 10 trials
[1] 4
> rbinom(5,n,p) # 5 binomial numbers
[1] 5 5 4 5 6
อย่างไรก็ตามเมื่อ n มีขนาดใหญ่ตัวแปรสุ่มทวินามจะมีลักษณะเข้าใกล้ตัวแปรสุ่มปกติ ซึ่งเป็ นผลมาจากทฤษฎี ขีดจำกัดกลาง (Central limit theorem) โดยจะได้ว่า
จะมีการแจกแจงประมาณได้ด้วยการแจกแจงแบบปกติ
รูปที่ 5.3 แสดงการแจกแจงของตัวแปรสุ่มทวินามที่ n= 5,15,50 เมื่อp= 0.25 ซึ่งจากกราฟจะเห็นได้ว่าเมื่อ
n มีขนาดเพิ่มขึ้นรูปร่างของการแจกแจงทวินามจะใกล้เคียงกับลักษณะของการแจกแจงแบบปกตินั่นคือมี
ลักษณะเป็นรูประฆังคว่ำ สมมาตรซึ่งกราฟเหล่านี้ได้จากคำสั่ง
> x = runif(100)
> hist(x,prob=TRUE, col=gray(.9), main="Uniform on [0,1]")
> curve(dunif(x,0,1),add=TRUE)
Normal.การแจกแจงแบบปกตินั้นมีพารามิเตอร์ที่เกี่ยวข้องอยู่ 2 ตัว คือµ และσ เช่น IQ มีการแจกแจงปกติด้วยแจงแบบปกติโดยมีค่าเฉลี่ ยเท่ากับ 280 และค่าเบี่ ยงเบนมาตรฐานเท่ากับ 10 เราสามารถทำให้ค่าของตัวแปรสุ่มปกตินี้เป็ นค่ามาตรฐานZ= (X− µ)/σได้ ซึ่งค่ามาตรฐานนี้จะมีการแจกแจงแบบปกติมาตรฐานที่มีค่าเฉลี่ยเท่ากับ 0 และค่าเบี่ยงเบนมาตรฐานเท่ากับ
rnorm(n, mean=0, sd=1)
ตัวอย่างต่อไปนี้แสดงการสุ่มตัวอย่างจากการแจกแจงแบบปกติ
> rnorm(1,100,16) # an IQ score
[1] 122.7393
> rnorm(1,mean=280,sd=10) # how long for a baby (10 days early)
[1] 276.9401
เพื่อให้เห็นถึงลักษณะของการแจกแจงแบบปกติที่มีค่าเฉลี่ย 0 และค่าเบี่ยงเบนมาตรฐานเท่ากับ 1 ทำได้ดังนี้
> x = rnorm(100)
> hist(x,prob=TRUE, col= "blanchedalmond", main="Normal(0,1)")
> curve(dnorm(x),add=TRUE)
ได้กราฟดังรูปที่ 5.2
รูปที่ 5.2: Normal(0,1)
Binomial.ตัวแปรสุ่มทวินามเป็ นตัวแปรสุ่มแบบไม่ต่อเนื่อง โดยที่ค่าของตัวแปรสุ่มทวินามคือจำนวนของความสำเร็จที่เกิดขึ้นในการทดลองที่ เป็นอิสระกัน n ครั้ง เราสามารถสุ่มตัวอย่างจากการแจกแจงทวินามได้โดยใช้คำสั่ง rbinom(n, size, prob) ในการสร้างตัวเลขสุ่ มทวินามทำได้ดังนี้
> n = 10; p=.5
> rbinom(1,n,p) # 6 successes in 10 trials
[1] 4
> rbinom(5,n,p) # 5 binomial numbers
[1] 5 5 4 5 6
อย่างไรก็ตามเมื่อ n มีขนาดใหญ่ตัวแปรสุ่มทวินามจะมีลักษณะเข้าใกล้ตัวแปรสุ่มปกติ ซึ่งเป็ นผลมาจากทฤษฎี ขีดจำกัดกลาง (Central limit theorem) โดยจะได้ว่า
จะมีการแจกแจงประมาณได้ด้วยการแจกแจงแบบปกติ
รูปที่ 5.3 แสดงการแจกแจงของตัวแปรสุ่มทวินามที่ n= 5,15,50 เมื่อp= 0.25 ซึ่งจากกราฟจะเห็นได้ว่าเมื่อ
n มีขนาดเพิ่มขึ้นรูปร่างของการแจกแจงทวินามจะใกล้เคียงกับลักษณะของการแจกแจงแบบปกตินั่นคือมี
ลักษณะเป็นรูประฆังคว่ำ สมมาตรซึ่งกราฟเหล่านี้ได้จากคำสั่ง
par(mfrow=c(1,3))
n=5;p=.25
x=rbinom(100,n,p)
hist(x,prob=TRUE,col="lightblue", main="n=5, p=.25")
xvals=0:n
points(xvals,dbinom(xvals,n,p),type="h",lwd=2)
points(xvals,dbinom(xvals,n,p),type="p",lwd=2)
#++++++++++++++++++++++++++++++++++++
n=15;p=.25
x=rbinom(100,n,p)
hist(x,prob=TRUE,col="lightgreen", main="n=15, p=.25")
xvals=0:n
points(xvals,dbinom(xvals,n,p),type="h",lwd=2)
points(xvals,dbinom(xvals,n,p),type="p",lwd=2)
#++++++++++++++++++++++++++++++++++++
n=50;p=.25
x=rbinom(100,n,p)
hist(x,prob=TRUE,col="brown1", main="n=50, p=.25")
xvals=0:n
points(xvals,dbinom(xvals,n,p),type="h",lwd=2)
points(xvals,dbinom(xvals,n,p),type="p",lwd=2)
Exponential.การแจกแจงแบบ exponentail เป็นการแจกแจงที่สำคัญอีกการแจกแจงหนึ่ง ซึ่งมักใช้ในการอธิบาย
ถึงลักษณะการแจกแจงของอายุการใช้งานของอุปกรณ์ไฟฟ้ า เช่น ถ้าค่าเฉลี่ยของอายุการใช้งานของหลอดไฟเท่าากับ2500 ชั่วโมง เราอาจจะคิดได้ว่าอายุการใช้งานนี้มีการแจกแจงแบบ exponentail ที่มีค่าเฉลี่ยเท่ากับ 2500 พารามิเตอร์สำหรับการแจกแจงนี้คือrate = 1/mean รูปแบบคำสั่งสำหรับการสุ่มตัวอย่างที่มีการแจกแจงแบบ exponentailคือ
rexp(n, rate = 1)
ตัวอย่างต่อไปนี้เป็นการสุ่มตัวอย่างข้อมูลที่มีการแจกแจงแบบ exponentail ที่มี rate เท่ากับ 1/2500 และได้กราฟของการแจกแจงของข้อมูลชุดนี้ดังรูปที่ 5.4
n=5;p=.25
x=rbinom(100,n,p)
hist(x,prob=TRUE,col="lightblue", main="n=5, p=.25")
xvals=0:n
points(xvals,dbinom(xvals,n,p),type="h",lwd=2)
points(xvals,dbinom(xvals,n,p),type="p",lwd=2)
#++++++++++++++++++++++++++++++++++++
n=15;p=.25
x=rbinom(100,n,p)
hist(x,prob=TRUE,col="lightgreen", main="n=15, p=.25")
xvals=0:n
points(xvals,dbinom(xvals,n,p),type="h",lwd=2)
points(xvals,dbinom(xvals,n,p),type="p",lwd=2)
#++++++++++++++++++++++++++++++++++++
n=50;p=.25
x=rbinom(100,n,p)
hist(x,prob=TRUE,col="brown1", main="n=50, p=.25")
xvals=0:n
points(xvals,dbinom(xvals,n,p),type="h",lwd=2)
points(xvals,dbinom(xvals,n,p),type="p",lwd=2)
Exponential.การแจกแจงแบบ exponentail เป็นการแจกแจงที่สำคัญอีกการแจกแจงหนึ่ง ซึ่งมักใช้ในการอธิบาย
ถึงลักษณะการแจกแจงของอายุการใช้งานของอุปกรณ์ไฟฟ้ า เช่น ถ้าค่าเฉลี่ยของอายุการใช้งานของหลอดไฟเท่าากับ2500 ชั่วโมง เราอาจจะคิดได้ว่าอายุการใช้งานนี้มีการแจกแจงแบบ exponentail ที่มีค่าเฉลี่ยเท่ากับ 2500 พารามิเตอร์สำหรับการแจกแจงนี้คือrate = 1/mean รูปแบบคำสั่งสำหรับการสุ่มตัวอย่างที่มีการแจกแจงแบบ exponentailคือ
rexp(n, rate = 1)
ตัวอย่างต่อไปนี้เป็นการสุ่มตัวอย่างข้อมูลที่มีการแจกแจงแบบ exponentail ที่มี rate เท่ากับ 1/2500 และได้กราฟของการแจกแจงของข้อมูลชุดนี้ดังรูปที่ 5.4
> x =
rexp(100,1/2500)
> hist(x, prob=TRUE, col="yellow3", main="Exponentail mean=2500")
> curve(dexp(x,1/2500),add=TRUE)
รูปที่ 5.3: Random binomial data with the theoretical distribution
รูปที่ 5.4: Random exponential data with the theoretical density
> hist(x, prob=TRUE, col="yellow3", main="Exponentail mean=2500")
> curve(dexp(x,1/2500),add=TRUE)
รูปที่ 5.3: Random binomial data with the theoretical distribution
รูปที่ 5.4: Random exponential data with the theoretical density
5.2 Probability
Calculation and Combinatorics
ในกรณีที่การสุ่มตัวอย่างเป็ นแบบไม่ใส่คืน เช่นsample(1:40), 5 เราทราบว่าความน่าจะเป็ นของการที่จะได้ค่าใดค่าหนึ่งเป็นค่าแรกคือ1/40 ค่าที่ สองเป็ น1/39 ตามลำดับ ดังนั้นในตัวอย่างนี้ความน่าจะเป็ นที่ จะได้ตัวอย่างสุ่มชุดใดชุดหนึ่งจะเท่ากับ1/(40× 39× 38××37× 36) หรือในR เราสามารถใช้ฟังก์ชันprod() ในการคำนวณหาค่าความน่าจะเป็นนี้ได้ดังนี้
> 1/prod(40:36)
[1] 1.266449e-08
อย่างไรก็ตามจะสังเกตเห็นได้ว่าค่าความน่าจะเป็นที่จะได้หมายเลขชุดใดๆนั้นเกิดขึ้นตามลำดับที่กำหนด (nPr) ถ้าในกรณี ที่ ลำดับของหมายเลขไม่ สำคัญหมายความว่าตัวเลขชุดเดียวกันถึงแม้ว่าจะมีลำดับการเรียงที่แตกต่างกันก็ ยังถือว่าเป็นตัวเลขชุดเดียวกัน (nCr) เราสามารถหาความน่าจะเป็นได้ดังนี้
> prod(5:1)/prod(40:36)
[1] 1.519738e-06
หรือใช้ฟังก์ชันchoose() ใน R คำนวณคานี้ได้ดังนี้
> 1/choose(40,5)
[1] 1.519738e-06
ในกรณีที่การสุ่มตัวอย่างเป็ นแบบไม่ใส่คืน เช่นsample(1:40), 5 เราทราบว่าความน่าจะเป็ นของการที่จะได้ค่าใดค่าหนึ่งเป็นค่าแรกคือ1/40 ค่าที่ สองเป็ น1/39 ตามลำดับ ดังนั้นในตัวอย่างนี้ความน่าจะเป็ นที่ จะได้ตัวอย่างสุ่มชุดใดชุดหนึ่งจะเท่ากับ1/(40× 39× 38××37× 36) หรือในR เราสามารถใช้ฟังก์ชันprod() ในการคำนวณหาค่าความน่าจะเป็นนี้ได้ดังนี้
> 1/prod(40:36)
[1] 1.266449e-08
อย่างไรก็ตามจะสังเกตเห็นได้ว่าค่าความน่าจะเป็นที่จะได้หมายเลขชุดใดๆนั้นเกิดขึ้นตามลำดับที่กำหนด (nPr) ถ้าในกรณี ที่ ลำดับของหมายเลขไม่ สำคัญหมายความว่าตัวเลขชุดเดียวกันถึงแม้ว่าจะมีลำดับการเรียงที่แตกต่างกันก็ ยังถือว่าเป็นตัวเลขชุดเดียวกัน (nCr) เราสามารถหาความน่าจะเป็นได้ดังนี้
> prod(5:1)/prod(40:36)
[1] 1.519738e-06
หรือใช้ฟังก์ชันchoose() ใน R คำนวณคานี้ได้ดังนี้
> 1/choose(40,5)
[1] 1.519738e-06
5.3 Buil-In
Distributions inR
ถูกสร้างเป็นฟังก์ชันไว้ในR เรียบร้อยแล้ว เราจึงสามารถใช้ฟังก์ชันเหล่านี้ในการหาค่าความน่าจะเป็นหรือค่าวิกฤตสำหรับการทดสอบสมมติฐานได้
โดยใช้แทนตารางสถิติแบบดั้งเดิม ในที่นี้เราจะพิจารณาการแจกแจงแบบปกติ และการแจกแจงแบบทวินาม
อย่างไรก็ตามสำหรับการแจกแจงแบบอื่นๆก็มีรูปแบบเดียวกัน
ต่อไปนี้เป็ นสิ่งที่เราสามารถหาได้จากการแจกแจงทางสถิติ
ต่อไปนี้เป็ นสิ่งที่เราสามารถหาได้จากการแจกแจงทางสถิติ
• Density หรือ point probability (ขึ้นอยู่กับเป็ นการแจกแจงไม่ต่อเนื่องหรือต่อเนื่อง)
• Cumulated probability, distribution function
• Quantiles
• Pseudo-random numbers
• Cumulated probability, distribution function
• Quantiles
• Pseudo-random numbers
Rมีฟังก์ชันสำหรับหาค่าเหล่านี้ได้ เช่น สำหรับการแจกแจงปกติจะมีฟังก์ชัน dnorm(), pnorm(), qnorm()และ rnorm()ตามลำดับ (density,probability,quantile,
และrandom)
5.3.1 Densities
Density สำหรับการแจกแจงแบบต่อเนื่องคือค่าของความน่าจะเป็ นที่ ตัวแปรสุ่มจะมีค่าเข้าใกล้ค่า x
ซึ่ งความน่าจะเป็นที่ตัวแปรสุ่มจะมีค่าอยู่ในช่วงที่กำหนดนี้คือพื้นที่ใต้กราฟที่สอดคล้องกับการแจกแจงของตัวแปรสุ่มนี้ สำหรับการแจกแจงแบบไม่ต่อเนื่อง คำว่า density หมายถึงความน่าจะเป็ นแบบจุด ซึ่งคือความน่าจะเป็นที่ตัวแปรสุ่มจะมีค่าเท่ากับ x density function เป็นการแจกแจงปกติ สามารถทำได้ดังนี้
> x = seq(-4,4,0.1)
> plot(x,dnorm(x), type="l")
ไดกราฟดังรูปที่ 5.5 ซึ่ง argument ค่า -4 ถึง 4 และ 0.1 หมายถึงเราสร้างให้เวกเตอร์x ประกอบด้วยข้อมูลที่มีค่าตั้งแต่ -4 ถึง 4 โดยที่แต่ละค่าห่างกัน 0.1 นั่นคือเวกเตอร์ x ประกอบด้วย (−4.0, −3.9, −3.8,..., 3.8, 3.9, 4.0) ในฟังก์ชันplot() เราใช้ argumenttype = "l" เพื่อให้มีการสร้างกราฟเป็นอีกวิธีหนึ่ งที่ สามารถสร้างกราฟในแบบเดียวกันนี้ได้คือใช้ฟงก์ชันcurve() ดังนี้
รูปที่ 5.5: Density of normal distribution
> curve(dnorm(x), from=-4, to=4)
คำสั่งนี้จะสะดวกมากกว่าในการสร้างกราฟรูปแบบนี้ แต่ ฟังก์ชันcurve() ต้องการค่าของy ในรูปฟังก์ชันของx
สำหรับการแจกแจงแบบไม่ต่อเนื่ อง ค่าของตัวแปรสุ่มจะเป็ นค่ามีจุดตัดที่ ชัดเจน จึงนิยมที่ จะสร้างกราฟแบบเส้นตรงตั้งฉากกับแกนx เช่นเราต้องการกราฟของการแจกแจงทวินามที่มีค่าn= 50 และp= 0.33 ทำได้ดังนี้
> x = 0:50
> plot(x,dbinom(x,size=50,prob=0.33), type="h")
ได้กราฟดังรูปที่ 5.6 ในฟังก์ชันdnorm() นี้เราใช้ argument 3 arguments โดย argument แรกคือx ซึ่งเป็น
ค่าของตัวแปรสุ่ม argument ที่สองคือsize=50 เป็นการกำหนดขนาดของการทดลอง และ argument ที่สามคือ
prob=0.33เป็ นการกำหนดความน่าจะเป็ นของการเกิดความสำเร็จในการทดลองแต่ละครั้ง ส่วน argument type
= "h"ในฟังก์ชัน plot()เป็ นการสร้างกราฟเส้นตรงตั้งฉากกับแกนx
5.3.2 Cumulative Distribution Functions
ฟังก์ชันการแจกแจงความน่าจะเป็ นสะสม (cumulative distribution function) หมายถึงค่าความน่าจะเป็นที่ตัวแปรสุ่มมีค่าเท่ากับหรือน้อยกว่าค่าที่กำหนด ฟังก์ชันในR ที่ใช้ในการคำนวณค่าความน่าจะเป็ นสะสม
คือฟังก์ชันที่ขึ้นต้นด้วย"p" เช่นpnorm และpbinom เป็ นต้น
หากเราทราบว่าตัวแปรสุ่มตัวหนึ่งมีการแจกแจงแบบปกติ มีค่าเฉลี่ยเท่ากับ 132 และค่าเบี่ยงเบนมาตรฐานเท่ากับ 13 และอยากทราบว่าความน่าจะเป็ นที่ตัวแปรสุ่มนี้จะมีค่ามากกว่า 160 จะเท่ากับเท่าไร ทำได้ดังนี้
5.3.1 Densities
Density สำหรับการแจกแจงแบบต่อเนื่องคือค่าของความน่าจะเป็ นที่ ตัวแปรสุ่มจะมีค่าเข้าใกล้ค่า x
ซึ่ งความน่าจะเป็นที่ตัวแปรสุ่มจะมีค่าอยู่ในช่วงที่กำหนดนี้คือพื้นที่ใต้กราฟที่สอดคล้องกับการแจกแจงของตัวแปรสุ่มนี้ สำหรับการแจกแจงแบบไม่ต่อเนื่อง คำว่า density หมายถึงความน่าจะเป็ นแบบจุด ซึ่งคือความน่าจะเป็นที่ตัวแปรสุ่มจะมีค่าเท่ากับ x density function เป็นการแจกแจงปกติ สามารถทำได้ดังนี้
> x = seq(-4,4,0.1)
> plot(x,dnorm(x), type="l")
ไดกราฟดังรูปที่ 5.5 ซึ่ง argument ค่า -4 ถึง 4 และ 0.1 หมายถึงเราสร้างให้เวกเตอร์x ประกอบด้วยข้อมูลที่มีค่าตั้งแต่ -4 ถึง 4 โดยที่แต่ละค่าห่างกัน 0.1 นั่นคือเวกเตอร์ x ประกอบด้วย (−4.0, −3.9, −3.8,..., 3.8, 3.9, 4.0) ในฟังก์ชันplot() เราใช้ argumenttype = "l" เพื่อให้มีการสร้างกราฟเป็นอีกวิธีหนึ่ งที่ สามารถสร้างกราฟในแบบเดียวกันนี้ได้คือใช้ฟงก์ชันcurve() ดังนี้
รูปที่ 5.5: Density of normal distribution
> curve(dnorm(x), from=-4, to=4)
คำสั่งนี้จะสะดวกมากกว่าในการสร้างกราฟรูปแบบนี้ แต่ ฟังก์ชันcurve() ต้องการค่าของy ในรูปฟังก์ชันของx
สำหรับการแจกแจงแบบไม่ต่อเนื่ อง ค่าของตัวแปรสุ่มจะเป็ นค่ามีจุดตัดที่ ชัดเจน จึงนิยมที่ จะสร้างกราฟแบบเส้นตรงตั้งฉากกับแกนx เช่นเราต้องการกราฟของการแจกแจงทวินามที่มีค่าn= 50 และp= 0.33 ทำได้ดังนี้
> x = 0:50
> plot(x,dbinom(x,size=50,prob=0.33), type="h")
ได้กราฟดังรูปที่ 5.6 ในฟังก์ชันdnorm() นี้เราใช้ argument 3 arguments โดย argument แรกคือx ซึ่งเป็น
ค่าของตัวแปรสุ่ม argument ที่สองคือsize=50 เป็นการกำหนดขนาดของการทดลอง และ argument ที่สามคือ
prob=0.33เป็ นการกำหนดความน่าจะเป็ นของการเกิดความสำเร็จในการทดลองแต่ละครั้ง ส่วน argument type
= "h"ในฟังก์ชัน plot()เป็ นการสร้างกราฟเส้นตรงตั้งฉากกับแกนx
5.3.2 Cumulative Distribution Functions
ฟังก์ชันการแจกแจงความน่าจะเป็ นสะสม (cumulative distribution function) หมายถึงค่าความน่าจะเป็นที่ตัวแปรสุ่มมีค่าเท่ากับหรือน้อยกว่าค่าที่กำหนด ฟังก์ชันในR ที่ใช้ในการคำนวณค่าความน่าจะเป็ นสะสม
คือฟังก์ชันที่ขึ้นต้นด้วย"p" เช่นpnorm และpbinom เป็ นต้น
หากเราทราบว่าตัวแปรสุ่มตัวหนึ่งมีการแจกแจงแบบปกติ มีค่าเฉลี่ยเท่ากับ 132 และค่าเบี่ยงเบนมาตรฐานเท่ากับ 13 และอยากทราบว่าความน่าจะเป็ นที่ตัวแปรสุ่มนี้จะมีค่ามากกว่า 160 จะเท่ากับเท่าไร ทำได้ดังนี้
> 1 -
pnorm(160, mean = 132, sd=13)
[1] 0.01562612
[1] 0.01562612
รูปที่ 5.6: Point probabilities in binom(50,0.33)
จากค่าความน่าจะเป็
นที่ได้หมายความว่ามีเพียง 1.5 เปอร์เซ็นต์
ของประชากรเท่านั้นที่มีค่ามากกว่า 160 ซึ่งฟังก์ชัน pnorm()จะให้ค่าความน่าจะเป็ นที่ตัวแปรสุ่มจะมีค่าน้อยกว่าหรือเท่ากับค่าที่กำหนดไว้ใน
argument แรกในฟังก์ชันด้วยค่าเฉลี่ยและค่าเบี่ยงเบนมาตรฐานที่กำหนด โดยปกติเราจะนำการหาค่าความน่าจะเป็ นนี้ไปใช้ในการทดสอบสมมุ ติฐานทางสถิติต่างๆ
เช่น ใน Simple sign
test:มีคนไข้ 20 คน แต่ละคนได้รับทรีทเมนต 2 ชนิด และผู้ทดลองได้ถามคนไข้แต่ละคนว่าชอบทรีทเมนต์ใดมากกว่ากันปรากฏว่ามีคนไข้ 16 คน ที่ชอบทรีทเมนต์ A มากกว่า จากข้อมูลที่ได้เราต้องการทดสอบว่าทรีทเมนต์
A ดีกว่าทรีทเมนต์ B อย่างมีนัยสำคัญทางสถิติหรือไม่
ซึ่งถ้าทรีทเมนต์ทั้งสองดีพอๆกัน เราย่อมคาดหมายว่าจำนวนคนไข้ชอบทรีทเมนต์ A จะมีการแจกแจงแบบทวินามด้วยค่าp= 0.5 และn= 20 ดังนั้นความน่าจะเป็ นที่ จำนวนคนไข้ที่ ชอบทรีทเมนต์ A จะน้อยกว่าหรือเท่ากับ 16 คน คือ
> pbinom(16, size=20,prob=0.5)
[1] 0.9987116
เราต้องการค่าความน่าจะเป็ นที่X≥16 ถ้าเรานำค่าความน่าจะเป็ นที่ ได้นี้ไปลบออกจาก 1 เราจะไม่ได้ค่าความน่าจะเป็ น ที่ถูกต้อง เราควรทำดังนี้
> 1-pbinom(15,size=20,prob=0.5)
[1] 0.005908966
หากต้องการหาค่าความน่าจะเป็ นของการทดสอบสมมติฐานแบบสองทาง เช่น ความน่าจะเป็ นที่ มีจำนวนคนไข้ที่ ชอบทรีทเมนต์ A มากกว่าหรือเท่ากับ 16 คน หรือน้อยกว่าเท่ากับ 4 คน
> 1-pbinom(15,size=20,prob=0.5) + pbinom(4,size=20,prob=0.5)
[1] 0.01181793
ซึ่งเท่ากับสองเท่าของความน่าจะเป็ นทางเดียว
5.3.3 Quantiles
Quantile function เป็ นส่วนกลับของฟังก์ชันความน่าจะเป็ นสะสม นั่นคือp-quantile คือค่าที่มีคุณสมบัติที่ทำให้ค่าความน่าจะเป็ นที่ ตัวแปรสุ่มจะมีค่าเท่ากับหรือน้อยกว่าค่าๆนี้เท่ากับp เช่น ตามนิยามแล้วค่ามัธยฐานจะเป็ นค่าควอนไทล์ที่ 50%
โดยทั่วไปตารางของการแจกแจงทางสถิติมักจะให้ค่าที่อยู่ในรูปของควอนไทล์ สำหรับชุดของค่าความน่าจะเป็นที่กำหนด ตารางสถิติจะแสดงขอบเขตที่ค่าสถิติทดสอบต้องมีค่ามากกว่าค่านี้เพื่อการปฏิเสธสมมติฐานหลักที่ระดับนัยสำคัญที่ กำหนด ในทางทฤษฎีมักใช้ค่าควอนไทล์ในการคำนวณหาช่วงความเชื่ อมั่นของพารามิเตอร์ เช่น ต้องการหาช่วงความเชื่ อมั่น 95%ของµ จาก
x+ σ/√n× N0.025≤ µ≤¯x+ σ/√n× N0.975
> pbinom(16, size=20,prob=0.5)
[1] 0.9987116
เราต้องการค่าความน่าจะเป็ นที่X≥16 ถ้าเรานำค่าความน่าจะเป็ นที่ ได้นี้ไปลบออกจาก 1 เราจะไม่ได้ค่าความน่าจะเป็ น ที่ถูกต้อง เราควรทำดังนี้
> 1-pbinom(15,size=20,prob=0.5)
[1] 0.005908966
หากต้องการหาค่าความน่าจะเป็ นของการทดสอบสมมติฐานแบบสองทาง เช่น ความน่าจะเป็ นที่ มีจำนวนคนไข้ที่ ชอบทรีทเมนต์ A มากกว่าหรือเท่ากับ 16 คน หรือน้อยกว่าเท่ากับ 4 คน
> 1-pbinom(15,size=20,prob=0.5) + pbinom(4,size=20,prob=0.5)
[1] 0.01181793
ซึ่งเท่ากับสองเท่าของความน่าจะเป็ นทางเดียว
5.3.3 Quantiles
Quantile function เป็ นส่วนกลับของฟังก์ชันความน่าจะเป็ นสะสม นั่นคือp-quantile คือค่าที่มีคุณสมบัติที่ทำให้ค่าความน่าจะเป็ นที่ ตัวแปรสุ่มจะมีค่าเท่ากับหรือน้อยกว่าค่าๆนี้เท่ากับp เช่น ตามนิยามแล้วค่ามัธยฐานจะเป็ นค่าควอนไทล์ที่ 50%
โดยทั่วไปตารางของการแจกแจงทางสถิติมักจะให้ค่าที่อยู่ในรูปของควอนไทล์ สำหรับชุดของค่าความน่าจะเป็นที่กำหนด ตารางสถิติจะแสดงขอบเขตที่ค่าสถิติทดสอบต้องมีค่ามากกว่าค่านี้เพื่อการปฏิเสธสมมติฐานหลักที่ระดับนัยสำคัญที่ กำหนด ในทางทฤษฎีมักใช้ค่าควอนไทล์ในการคำนวณหาช่วงความเชื่ อมั่นของพารามิเตอร์ เช่น ต้องการหาช่วงความเชื่ อมั่น 95%ของµ จาก
x+ σ/√n× N0.025≤ µ≤¯x+ σ/√n× N0.975
เมื่อN0.025 คือควอนไทล์ที่2.5% ในการแจกแจงปกติ
ถ้าσ= 12 และn= 5 และเราหาค่าเฉลี่ย
x= 83 เราสามารถหาช่วงความเชื่ อมั่นได้ดังนี้
> xbar = 83
> sigma = 12
> n = 5
> sem = sigma / sqrt(n)
> sem
[1] 5.366563
> xbar + sem * qnorm(0.025)
[1] 72.48173
> xbar + sem * qnorm(0.975)
[1] 93.51827
ดังนั้นช่วงความเชื่อมั่น95% ของµ คือ 72.48 ถึง 93.52
โดยปกติ ค่าควอนไทล์มักเขียนอยู่ ในรูปΦ−1 (0.975) เมื่อΦ คือสัญลักษณ์มาตรฐานของฟังก์ชันการแจกแจง
ความน่ าจะเป็ นสะสมของการแจกแจงแบบปกติ
x= 83 เราสามารถหาช่วงความเชื่ อมั่นได้ดังนี้
> xbar = 83
> sigma = 12
> n = 5
> sem = sigma / sqrt(n)
> sem
[1] 5.366563
> xbar + sem * qnorm(0.025)
[1] 72.48173
> xbar + sem * qnorm(0.975)
[1] 93.51827
ดังนั้นช่วงความเชื่อมั่น95% ของµ คือ 72.48 ถึง 93.52
โดยปกติ ค่าควอนไทล์มักเขียนอยู่ ในรูปΦ−1 (0.975) เมื่อΦ คือสัญลักษณ์มาตรฐานของฟังก์ชันการแจกแจง
ความน่ าจะเป็ นสะสมของการแจกแจงแบบปกติ
5.4 Simulations
หากเราสามารถจำลองข้อมูลสุ่มชนิดต่างๆได้จะทำให้เราสามารถทำการทดลองเพื่อตอบคำถามเกี่ยวกับการแจกแจง ของข้อมูลได้เนื่องจากมีฟังก์ชันมากมายที่สามารถใช้ในการสร้างตัวเลขสุ่มได้ เราจะสามารถทราบถึงลักษณะการแจกแจงของค่าตัวเลขสุ่มเหล่านี้ได้โดยพิจารณาจาก histogram หรือวิธีการอื่นๆ ในที่นี้เราต้องการสร้างตัวเลขสุ่มชนิดใหม่ขึ้นมาและตรวจสอบว่ าการแจกแจงของตัวเลขสุ่ มนี้เป็ นการแจกแจงแบบใด
5.4.1 The Central Limit Theorem
เราเริ่มต้นด้วยตัวอย่างของ central limit theorem (CLT) ซึ่งเป็ นทฤษฎีที่สำคัญมากในทางสถิติ ซึ่งทฤษฎีนี้ได้กล่าวไว้ว่า หากXi ถูกสุ่มมาอย่างเป็ นอิสระกันจากประชากรที่มีค่าเฉลี่ยµ และค่าเบี่ยงเบนมาตรฐานσ ที่ทราบค่าแล้วค่ามาตรฐานของค่าเฉลี่ยของX คือ
จะมีการแจกแจงเข้าใกล้การแจกแจงแบบปกติที่มีค่าเฉลี่ยเท่ากับ 0 และค่าเบี่ยงเบนมาตรฐานเท่ากับ 1 นั่นคือหากn มีขนากใหญ่พอแล้ว ค่าเฉลี่ยของตัวอย่างจะมีการแจกแจงประมาณได้ด้วยการแจกแจงแบบปกติ โดยมีค่าเฉลี่ยเท่ากับ
จะมีการแจกแจงประมาณได้ด้วยN(0,1) และเราจะใช้R ในการตรวจสอบทฤษฎีนี้ดังต่ อไปนี้
> n = 10
> p=.25
> s=rbinom(100,n,p)
> x = (s - n*p)/sqrt(n*p*(1-p))
> hist(x,prob=TRUE,col="gold2")
> curve(dnorm(x), col="blue", add=TRUE)
จากคำสั่งในข้างต้นจะได้กราฟดังรูปที่ 5.7 ซึ่งจะเห็นได้ว่าการแจกแจงของข้อมูลนั้นใกล้เคียงกับN(0,1)
รูปที่ 5.7: Scaled binomial data is approximatelyN(0,1)
5.4.2 For Loops
หากเราสามารถจำลองข้อมูลสุ่มชนิดต่างๆได้จะทำให้เราสามารถทำการทดลองเพื่อตอบคำถามเกี่ยวกับการแจกแจง ของข้อมูลได้เนื่องจากมีฟังก์ชันมากมายที่สามารถใช้ในการสร้างตัวเลขสุ่มได้ เราจะสามารถทราบถึงลักษณะการแจกแจงของค่าตัวเลขสุ่มเหล่านี้ได้โดยพิจารณาจาก histogram หรือวิธีการอื่นๆ ในที่นี้เราต้องการสร้างตัวเลขสุ่มชนิดใหม่ขึ้นมาและตรวจสอบว่ าการแจกแจงของตัวเลขสุ่ มนี้เป็ นการแจกแจงแบบใด
5.4.1 The Central Limit Theorem
เราเริ่มต้นด้วยตัวอย่างของ central limit theorem (CLT) ซึ่งเป็ นทฤษฎีที่สำคัญมากในทางสถิติ ซึ่งทฤษฎีนี้ได้กล่าวไว้ว่า หากXi ถูกสุ่มมาอย่างเป็ นอิสระกันจากประชากรที่มีค่าเฉลี่ยµ และค่าเบี่ยงเบนมาตรฐานσ ที่ทราบค่าแล้วค่ามาตรฐานของค่าเฉลี่ยของX คือ
จะมีการแจกแจงเข้าใกล้การแจกแจงแบบปกติที่มีค่าเฉลี่ยเท่ากับ 0 และค่าเบี่ยงเบนมาตรฐานเท่ากับ 1 นั่นคือหากn มีขนากใหญ่พอแล้ว ค่าเฉลี่ยของตัวอย่างจะมีการแจกแจงประมาณได้ด้วยการแจกแจงแบบปกติ โดยมีค่าเฉลี่ยเท่ากับ
จะมีการแจกแจงประมาณได้ด้วยN(0,1) และเราจะใช้R ในการตรวจสอบทฤษฎีนี้ดังต่ อไปนี้
> n = 10
> p=.25
> s=rbinom(100,n,p)
> x = (s - n*p)/sqrt(n*p*(1-p))
> hist(x,prob=TRUE,col="gold2")
> curve(dnorm(x), col="blue", add=TRUE)
จากคำสั่งในข้างต้นจะได้กราฟดังรูปที่ 5.7 ซึ่งจะเห็นได้ว่าการแจกแจงของข้อมูลนั้นใกล้เคียงกับN(0,1)
รูปที่ 5.7: Scaled binomial data is approximatelyN(0,1)
5.4.2 For Loops
โดยทั่วไปกระบวนการในการสร้างตัวเลขสุ่มจำนวน
100 ค่า อาจทำได้ไม่ง่ายนัก ดังนั้นเราอาจจำเป็นต้องที่จะสร้าง
ตัวเลขสุ่มนี้ที่ละค่า
ซึ่งทำได้โดยการใช้ "for loop" ตัวอย่างต่อไปนี้แสดงการสร้างตัวเลขสุ่มจำนวน 100 ค่า ที่ละค่า
> results =
numeric(0)
# a place to
store the results
> for (i in
1:100){
# the for loop
+ s =
rbinom(1,n,p)
# just 1 this
time
+ results[i] =
(s - n*p)/sqrt(n*p*(1-p))
# store the
answer
+ }
>
hist(results)
จากคำสั่งเหล่านี้เราได้สร้างตัวแปรresults ไว้สำหรับเก็บข้อมูลที่ เราสร้างขึ้น เราสามารถดูลักษณะของข้อมูลที่ เรา
สร้างขึ้นได้โดยใช้คำสั่ง hist(results) ซึ่งได้ดังรูปที่ 5.8
CLT with normal
data.CLT สามารถใช้กับตัวแปรสุ่มปกติได้ด้วย พิจารณาตัวอย่างต่อไปนี้ ให้Xi เป็นตัวเลขสุ่ม
ปกติที่มีค่าเฉลี่ยµ= 5 และσ= 5
>
results=c()
> mu=0;sigma=1
> for(i in 1:200) {
+ x=rnorm(100,mu,sigma)
> mu=0;sigma=1
> for(i in 1:200) {
+ x=rnorm(100,mu,sigma)
รูปที่ 5.8: Simulation of CLT with binomial data
+ results[i]=(mean(x) - mu)/(sigma/sqrt(100))
+ }
>
hist(results,prob=TRUE,ylim=c(0,.4), col="green")
>
curve(dnorm(x), add=TRUE)
จากคำสั่งข้างบนได้กราฟดังรูปที่ 5.9 โดย histogram ที่ได้ชี้ให้เห็นว่าข้อมูลมีการแจกแจงประมาณไดด้วยการแจกแจง
ปกติ
Histogram of
results
รูปที่ 5.9: Simulation of CLT with normal data
5.4.3 Normal Probability Plots
กราฟอีกชนิดหนึ่งที่
สามารถนำมาใชในการพิจารณาวาขอมู ลมี การแจกแจงแบบปกติ หรือไม ไดแกnormal probability
plotซึ่งกราฟชนิดนี้ใชไดดีกวา histogram ในการตรวจสอบลักษณะการแจกแจงปกติของขอมูล
แนวคิดของการสราง
กราฟรูปแบบนี้มาจากการพลอตควอนไทลของขอมูลของเรากับควอนไทลของการแจกแจงแบบปกติ
ซึ่งควอนไทลที่qคือคาของขอมู ลที่ ขอมู ลจำนวนq×100% ของขอมู ลทั้งหมดที่ มีคานอยกวาคานี้ ดังนั้นควอนไทล 0.25 คือ ควอไทลที่1, ควอนไทล 0.5 คือควอไทลที่ 2 หรือ median และควอนไทล 0.75 คือควอไทลที่ 3 ควอนไทลของการแจกแจงทางทฤษฎีก็คลายคลึงกันเพียงแต แทนที่ จะเป
นจำนวนของขอมูลที่ มีคานอยกวา ก็กลายเป นพื้นที่ ที่อยู ทางดานซายของคาที่กำหนด
เชน median จะแบงพื้นที่ใตโคงความหนาแนน (density curve) ออกเป นสองสวนเทาๆกัน
การแปลความหมายของ
normal probability plot ทำไดโดยพิจารณาการกระจายของจุดบนกราฟ
หากการกระจายมีแนวโนมเป นเสนตรงแลวแสดงวาขอมูลมีการแจกแจงประมาณไดดวยการแจกแจงแบบปกติ
หากการกระจายเป นเสนโคงแสดงว าการแจกแจงนั้นมีหางยาวหรือหางสั้นได
ในR มีฟงกชันที่ใชในการสราง normal probability plot ไดแกqqnorm() หรือฟงกชันทั่วไปคือqqplot()
และqqline() ที่ ใชในการลากเสนอางอิง
ตัวอยางตอไปนี้แสดงการสราง
normal probability plot สำหรับขอมูลที่มีการแจกแจงแบบตางๆ
>
par(mfrow=c(1,4))
> x =
rnorm(100);qqnorm(x, main="Normal(0,1)");qqline(x)
> x = rnorm(100,10,15);qqnorm(x, main="Normal(10,15)");qqline(x)
> x = rexp(100,1/10);qqnorm(x, main="Exponential mu = 10");qqline(x)
> x = runif(100,0,1);qqnorm(x, main="Uniform(0,1)");qqline(x)
> x = rnorm(100,10,15);qqnorm(x, main="Normal(10,15)");qqline(x)
> x = rexp(100,1/10);qqnorm(x, main="Exponential mu = 10");qqline(x)
> x = runif(100,0,1);qqnorm(x, main="Uniform(0,1)");qqline(x)
ไดกราฟดังรูปที่ 5.10 จะเห็นไดวากราฟ 2 รูปแรกมีลักษณะของการกระจายคลายเสนตรง ในขณะที่กราฟสองรูปหลัง
ไม
รูปที่ 5.10: Some normal plots
รูปที่ 5.10: Some normal plots
ไม่มีความคิดเห็น:
แสดงความคิดเห็น