期末考 2022

• 把整張期末考卷該計算的部分用 R (或 Python) 寫一遍，老狗學點新把戲，一樂。

1 [Q1]

• 1. 請寫出 Binomial Distribution, X的 機率質量函數 (pmf), \(f_X(x)\)。設其參數為 (n=10, p=0.5)，計算 \(P(X≤2)\)=?
• 2. 請寫出 Poisson Distribution, Y 的 機率質量函數 (pmf) \(f_Y (y)\)。設其參數(平均數)為 (λ=5)，計算 Var(Y)=?
• 3. 請證明，當 λ=n⋅p， 而 $ n → ∞$ 時，\(f_X (s)→f_Y (s)\)
• 4. 令 n=100, p=0.1, λ=10， 取 s=10, 計算 \(f_X(s)\) 與 \(f_Y(s)\)，以及其相對誤差 \(|f_X(s)-f_Y(s)|/(f_Y(s))\)，以其驗證 (3) 之正確性。
• 5. Gamma 函數 \(Γ(s)\)是俗稱的「廣義」階乘函數，其定義為： \(Γ(s)=\int_0^∞ x^{s-1} e^{-x} dx\)，請計算 \(Γ(0.5)= ?\) 以及 \(Γ(5.5)= ?\)
• 6. 根據「中央極限定理」 (Central Limit Theorem , CLT)，取 n 個獨立且分布相同 (iid) 的隨機變數 \(X_1,X_2, …. X_n\) 其共同的平均數及變異數設定為 \(E(X_i )=μ_X , Var(X_i)=σ_X^2\), 把這些iid的隨機變數作如何處理就可以令其機率分布趨近於標準型的常態分佈 \(Norm(μ=0,σ=1)？\)

1.1 使用 R

quit
#Q1.1 

q101= pbinom(2, 10, .5)

#Q1.4

n= 100
p= .1
λ= 10
s= 10
dbin= dbinom(s,n,p)
dpoi= dpois(s,λ)
err=  abs(dbin-dpoi)/dpoi

q104= c(dbin, dpoi, err)

#Q1.5

q105= c(gamma(.5), gamma(5.5))

# ans ---
q101

[1] 0.0546875

q104

[1] 0.13186535 0.12511004 0.05399496

q105

[1]  1.772454 52.342778

1.2 使用 Python

reticulate::repl_python()
import math
import numpy as np
import scipy.stats as st


#Q1.1 

#q101= pbinom(2, 10, .5)
q101= st.binom(10, .5).cdf(2)

#Q1.4

n= 100
p= .1
λ= 10
s= 10

#dbin= dbinom(s,n,p)
#dpoi= dpois(s,λ)

dbin= st.binom(n,p).pmf(s)
dpoi= st.poisson(λ).pmf(s)

err=  abs(dbin-dpoi)/dpoi

q104= (dbin, dpoi, err)

#Q1.5

q105= (math.gamma(.5), math.gamma(5.5))

# ans ---
q101

0.0546875

q104

(0.13186534682448794, 0.12511003572113372, 0.053994957833851165)

q105

(1.7724538509055159, 52.34277778455352)

2 [Q2]

n 個獨立且分布相同的隨機變數 \(X_1,X_2, …, X_n\)， 其共同的期望值為 μ, 標準差為 σ, 統計量 \(\bar{X} , S^2\) 分別代表它們的【樣本平均數】以及【樣本變異數】。 數學公式定義如下：

• \(\bar{X}= {1\over n} \sum_{i=1}^n X_i\) ,
• \(S^2={1\over {n-1}} \sum_{i=1}^n {(X_i-\bar{X})}^2\)

又令統計量 T 代表 樣本平均 \(\bar{X}\) 與母體平均 μ 之 標準化距離，數學公式定義如下：

• \(T={\bar{X}-μ\over S/\sqrt{n}}\)

令 n = 100, μ = 50, σ = 10

請回答以下問題：

• 7. \(E(\bar{X})= ?\)
• 8. \(Var(\bar{X})= ?\)
• 9. \(E(S^2)= ?\)
• 10. \(\bar{X} 的 機率密度函數 (pdf), f_{\bar{X}}(x)= ?\)
• 11. \(T 的 機率密度函數 (pdf), f_T (t)= ?\)
• 12. \(P(48 ≤ \bar{X} ≤ 52)= ?\)
• 13. \(P(-2 ≤ T≤ 2)=?\)
• 14. \(若 P(-∞ < \bar{X} ≤ x)= 95% 那麼 x = ?\)
• 15. \(若 P(t≤ T <∞)= 95% 那麼 t = ?\)

2.1 使用 R

quit
n= 100
μ= 50  
σ= 10

# Q2.12
q212= pnorm(52, mean= μ, sd= σ/sqrt(n)) -
  pnorm(48, mean= μ, sd= σ/sqrt(n))

# Q2.13
q213= pt(2, df=n-1) - pt(-2, df=n-1)

# Q2.14
q214= qnorm(.95, mean= μ, sd= σ/sqrt(n))

# Q2.15
q215= qt(.95, df=n-1, lower.tail=FALSE)

q212

[1] 0.9544997

q213

[1] 0.9517603

q214

[1] 51.64485

q215

[1] -1.660391

2.2 使用 Python

reticulate::repl_python()

n= 100
μ= 50  
σ= 10

# 2.12
# q212= pnorm(52, mean= μ, sd= σ/sqrt(n)) -
#  pnorm(48, mean= μ, sd= σ/sqrt(n))

q212= st.norm(μ, σ/n**.5).cdf(52) - st.norm(μ, σ/n**.5).cdf(48)

# Q2.13
# q213= pt(2, df=n-1) - pt(-2, df=n-1)
q213= st.t(df=n-1).cdf(2) - st.t(df=n-1).cdf(-2)

# Q2.14
# q214= qnorm(.95, mean= μ, sd= σ/sqrt(n))
q214= st.norm(μ, σ/n**.5).ppf(.95)

# Q2.15
# q215= qt(.95, df=n-1, lower.tail=FALSE)
q215= st.t(df=n-1).isf(.95)

q212

0.9544997361036416

q213

0.951760306627367

q214

51.64485362695147

q215

-1.6603911559963895

3 [Q3]

本校某社團根據現今社會氛圍，認為大學生「不反對婚前性行為」的比例p很可能超過50%，為了支持這個想法，社長隨機挑選校內n = 100位學生做問卷調查，詢問同學們是否「不反對婚前性行為？」。 回答【是】的學生人數令為隨機變數X。 令p之估測值( estimate)為 \(\hat{p}=X/n\)，假設每位學生皆是獨立作答，不互相影響。請回答以下問題：

• 16. 試寫出估測值 \(\hat{p}\) 之期望值\(E(\hat{p})\) 以及標準差 \(std(\hat{p})\) 之理論值 = ?
      
• 17. 這次調查得到結果為 55份回答【是】， 45份回答【否】， 請問 \(\hat{p}\) = ?
• 18. 根據這份調查結果，請寫出 p 之 95% 之信心區間。
• 19. 以這份調查結果，能否以顯著水準(significance level) α = 5% 依假設檢定之原理來做成接受 \(H_1: “p > 50\%”\) 的結論，請寫出檢定過程。並說明所獲結論的理由。
• 20. 如果把問卷份數 從 n = 100 增加到 n = 400，重作檢定，在調查結果的比例維持不變的情形之下，結論是否有不同？ 請說明你的結論以及推論的過程。

3.1 使用 R

quit
n= 100
X=  55
α= .05

# Q3.17
q317= pH= X/n

# Q3.18

#se= 1.96*(pH*(1-pH)/n)^.5
se=  -qnorm(α/2) * (pH*(1-pH)/n)^.5
q318= pH + c(-se, se)

q317

[1] 0.55

q318

[1] 0.452493 0.647507

#
# 參考的解答
#

# Exact binomial test

q1= binom.test(X, n, p= .5)
q318_1= q1$conf.int
q318_1

[1] 0.4472802 0.6496798
attr(,"conf.level")
[1] 0.95

# proportions test with continuity correction

q2= prop.test(X, n,  p= .5) #, correct=FALSE)
q318_2= q2$conf.int
q318_2

[1] 0.4475426 0.6485719
attr(,"conf.level")
[1] 0.95

4 [Q3.19, .20]

\(H_1: "p > .5"\)

\(H_0: "p \leq .5"\)

Set α = 5%

\(H_0: "p \leq .5"\)

\(T: "{{\hat{p}-p_0} \over \sqrt{p_0(1-p_0)\over n}} > Q^{-1}(α)"\) ?

or

\(T: "pValue \equiv Q({{\hat{p}-p_0} \over \sqrt{p_0(1-p_0)\over n}}) < α"\) ?

where,

\(z_α = Q^{-1}(α)\)

or

\(Q(z_α) = α\)

The decision Rule:

If the test (T) is True, then we can Reject \(H_0\), and Accept \(H_1\)

else,
we can NOT Reject \(H_0\), can NEITHER Accept \(H_1\)

---

4.1 使用 R

# Q3.19

p0= .5
z_α= qnorm(α, lower.tail = FALSE)

theStats= (pH-p0)/(p0*(1-p0)/n)^.5
theTest= (  theStats> z_α )

theStats

[1] 1

theTest

[1] FALSE

if (theTest==TRUE) { theDecision= "Reject H0"
  } else {theDecision= "Can NOT Reject H0"}

q319= theDecision
q319

[1] "Can NOT Reject H0"

# Q3.20

nLarger= 400
theStats= (pH-p0)/(p0*(1-p0)/nLarger)^.5
theTest= (  theStats> z_α )

theStats

[1] 2

theTest

[1] TRUE

if (theTest==TRUE) { theDecision= "Reject H0"
  } else {theDecision= "Can NOT Reject H0"}


q320= theDecision
q320

[1] "Reject H0"

#
# 參考的解答，來自 R ，仍需進一步解讀
#


q2= prop.test( pH*n, n,       p= .5, alternative = 'greater')
q4= prop.test( pH*nLarger, nLarger, p= .5, alternative = 'greater')


q2

    1-sample proportions test with continuity correction

data:  pH * n out of n, null probability 0.5
X-squared = 0.81, df = 1, p-value = 0.1841
alternative hypothesis: true p is greater than 0.5
95 percent confidence interval:
 0.462983 1.000000
sample estimates:
   p 
0.55 

q4

    1-sample proportions test with continuity correction

data:  pH * nLarger out of nLarger, null probability 0.5
X-squared = 3.8025, df = 1, p-value = 0.02559
alternative hypothesis: true p is greater than 0.5
95 percent confidence interval:
 0.5076335 1.0000000
sample estimates:
   p 
0.55 

# 進一步解讀
# 運用 pValue

theStats= q2$statistic
thePV=    q2$p.value

theTest= (  thePV < α )

if (theTest==TRUE) { theDecision= "Reject H0"
  } else {theDecision= "Can NOT Reject H0"}

q320_1= theDecision

theStats

X-squared 
     0.81 

thePV

[1] 0.1840601

q320_1

[1] "Can NOT Reject H0"

# 進一步解讀
# 運用 pValue

theStats= q4$statistic
thePV=    q4$p.value

theTest= (  thePV < α )

if (theTest==TRUE) { theDecision= "Reject H0"
  } else {theDecision= "Can NOT Reject H0"}

q320_2= theDecision

theStats

X-squared 
   3.8025 

thePV

[1] 0.02558806

q320_2

[1] "Reject H0"

4.2 使用 Python

reticulate::repl_python()
n= 100
X=  55
α= .05

# Q3.17
q317= pH= X/n

# Q3.18

#se= np.array([-1.96, 1.96])*(pH*(1-pH)/n)**.5
ci= pH + st.norm.ppf([α/2, 1-α/2]) * (pH*(1-pH)/n)**.5
q318= ci
 
q317

0.55

q318

array([0.45249302, 0.64750698])

from statsmodels.stats.proportion import proportions_ztest

pH= .55
n=  100

zstat1, pValue1= proportions_ztest(
  count= pH*n, 
  nobs= n, 
  value= 0.5, 
  alternative= 'larger')
  
theTest1= (pValue1 < α)

q319= theDecision1= "Reject H0" if theTest1  else "can NOT Reject H0"
q319

'can NOT Reject H0'

nLarger= 400

zstat2, pValue2= proportions_ztest(
  count= pH*nLarger, 
  nobs= nLarger, 
  value= 0.5, 
  alternative= 'larger')


theTest2= (pValue2 < α)

q320= theDecision2= "Reject H0" if theTest2  else "can NOT Reject H0"

q320

'Reject H0'

5 [Ref]

5.1 手稿1

手稿1

5.2 手稿2

手稿2

5.3 手稿3

手稿3