Copying "Introduction to Statistical Analysis" in Python Vol. 8 - Chapter 4 "First Statistical Tests" ① Tests for Population Mean, Population Variance, and Population Proportion

Chapter 4 "First Statistical Tests"

Author of the book: Dr. Sadao Ishimura

This article covers the "Introduction to Statistical Analysis" Chapter 4 "First Statistical Tests" Python copying activity .

This is a copying series where I calmly convert the book's figures, tables, and calculations into Python.
In this article, among the statistical test themes in Chapter 4, I will tackle tests for population mean, population variance, and population proportion.
I will continue to make use of ChatGPT!

Now, let's open the book and set off on a journey of statistical analysis 🚀

Introduction

This blog series introduces the "joy of statistical analysis" gained through copying the Python code from the book "Introduction to Statistical Analysis" (Tokyo Tosho, referred to as the "Text").

The introduction of the book and citation notation are posted in the linked article.

Chapter 4 First Statistical Tests

This article covers the following sections of Chapter 4.

The data used in the article cites the data published in the text itself.
For items with a small amount of data, I register the data in the code, and for items with a large amount of data, I convert them into CSV files and read them.

I will import the libraries used in Chapter 4.

### インポート

# 数値計算
import math                                                # python標準ライブラリ
import numpy as np
import pandas as pd

# 統計
import scipy.stats as stats
import pingouin as pg
from statsmodels.stats.proportion import proportions_ztest # 母比率のz検定
from statsmodels.stats.weightstats import ttest_ind        # 2標本のt検定

# 描画
import matplotlib.pyplot as plt
import seaborn as sns
plt.rcParams['font.family'] = 'Meiryo'

Introduction

As I begin learning statistical tests, ChatGPT sent me some encouragement.

When the persuasiveness of a presentation really goes up, it's not because of a "hunch," but because there is "evidence in the data."
Statistical testing is a tool that provides that "evidence".
It gives numerical answers to doubts like "Is there really a difference?" or "Isn't it just a coincidence?"
Let's cultivate the judgment to not get lost in data together from here!

"Statistical testing" is also a difficult concept, but let's do our best to tackle it.
First, please enjoy a "gentle introduction" to statistical testing by ChatGPT.

🎯 What is statistical testing? - First, the idea of "questioning the hypothesis"
For things you think "this is probably the case" (= hypothesis),
asking "is it okay as it is?" and making a judgment using data.
That is the concept of statistical testing.

🍰 A gentle introduction to each test

🔸 Testing Population Means
When you want to verify whether the mean value differs from a standard, such as asking, 'Can I say the average satisfaction for this product is 80 points?'

🔸 Testing Population Variances
When you want to evaluate stability or the magnitude of variation, such as asking, 'Does the variation in quality exceed the standard?'

🔸 Testing Population Proportions
When you want to compare changes in proportions, such as asking, 'Has the response rate for this campaign increased compared to the last one?'

🧾 In closing the intro
As shown, testing is a tool for 'using data to re-examine whether things are fine as they are.'
What is the truth of the matter? Let's think about it carefully together.

The procedure for statistical testing is roughly as follows.
For details, please be sure to read 📘 Text Section 4.1!

1. [Planning Phase] Things to decide before research/study

• Selection of test method
  → Choose a test that fits the target you want to compare and the content of your hypothesis
  (e.g., test for population mean, test for population variance, test for population proportion, etc.)
  *Once the test method is decided, the test statistic is also determined

• Setting the hypothesis

  • Null hypothesis $${H_0}$$: The hypothesis you want to overturn (e.g., 'there is no difference')

  • Alternative hypothesis $${H_1}$$: The hypothesis you want to support (e.g., 'there is a difference')

• Setting the significance level
  → The threshold for judgment (e.g., 5% = 0.05) is decided in advance

2. [Execution Phase] Things to do after obtaining data

• Calculation of test statistic
  → Based on the obtained data, calculate the test statistic under the null hypothesis

• Judgment of the null hypothesis
  → Compare the test statistic with the rejection region, or compare the $${p}$$-value with the significance level to determine 'whether the null hypothesis is rejected'

Section 4.2 Testing Population Means

[Summary of Testing Population Means]

[Hypothesis]
Set the hypothesis for the population mean $${\mu_0}$$ that you want to overturn as follows.

The hypothesis you want to support = alternative hypothesis, is set as one of the following three patterns.

◆ ◆ ◆

■ Testing the population mean of a normal population: Formula for test statistic p.146

I will borrow the formula for the test statistic in population mean testing from the text.

Given sample mean $$\bar{x}$$, sample variance $${s^2}$$, sample size $${N}$$, and null hypothesis population mean $$\mu_0$$, the test statistic for population mean $$\mu$$ using sample $${{x_1, x_2, \ldots, x_N}}$$ is

follows a $${t}$$-distribution with $${N-1}$$ degrees of freedom.

I will define a "test function for the population mean of a normal population" according to the formula.
I will use scipy.stats for probability calculations.

### 正規母集団の母平均の検定関数（母分散未知） p.146

# 正規母集団の母平均のt検定統計量の算出関数
def t_stat_pop_mean(x_list, mu0):
    N = len(x_list)
    x_bar = sum(x_list) / N
    s2 = sum([(x - x_bar)**2 for x in x_list]) / (N - 1)
    return (x_bar - mu0) / math.sqrt(s2 / N)

# 正規母集団の母平均のt検定関数 ※確率計算はscipy.stats利用
def pop_mean_ttest(x_list, mu0, alpha=0.05, alternative='two-sided'):
    N = len(x_list)
    t_val = t_stat_pop_mean(x_list, mu0)
    match alternative:
        case 'two-sided':
            c_val = stats.t.ppf(q=1 - alpha/2, df=N - 1)
            c_val = -c_val, c_val
            p_val = stats.t.sf(x=abs(t_val), df=N - 1) * 2
        case 'less':
            c_val = stats.t.ppf(q=alpha, df=N - 1)
            p_val = stats.t.cdf(x=t_val, df=N - 1)
        case 'greater':
            c_val = stats.t.ppf(q=1 - alpha, df=N - 1)
            p_val = stats.t.sf(x=t_val, df=N - 1)
    return {'t_value': t_val, 'c_value': c_val, 'alpha': alpha, 'p_value': p_val}

For the test data $${[14.1, 14.0, 14.1, 14.0]}$$, I will perform a two-tailed test for the population mean with the null hypothesis $${H_0}$$ "population mean $$\mu=\mu_0=14$$" and a significance level of 5% ($${\alpha=0.05}$$).

# テスト
x_list = [14.1, 14.0, 14.1, 14.0]
pop_mean_ttest(x_list, mu0=14, alpha=0.05)

【Execution Result】
"It cannot be said to be significant at the 5% significance level, and the null hypothesis cannot be rejected."

I will briefly explain the output content.

• t_value: Realized value of test statistic $${T}$$

• c_value: Rejection limit value
  Threshold for test statistic $${T}$$ at which the null hypothesis can be rejected
  When the realized value of $${T}$$ is outside the rejection limit, the null hypothesis can be rejected

• alpha: Significance level $$\alpha$$

• p_value: $${p}$$-value
  Probability that the test statistic becomes the realized value under the null hypothesis
  When the $${p}$$-value is less than or equal to the significance level, the null hypothesis can be rejected

◆ ◆ ◆

■ Population mean test (two-tailed test) example p.147
In the example on p.147, I will check whether the diameter data measured by randomly sampling from recent manufactured products meets specifications using a population mean test.

The hypotheses are as follows:
・Null hypothesis $${H_0}$$ "population mean $$\mu=15.4$$"
・Alternative hypothesis $${H_1}$$ "population mean $$\mu \neq 15.4$$"
　⇒Two-tailed test: Because it is NG if it is either larger or smaller than the specification

① Perform test with custom function
First, I will set the data and the null hypothesis $$\mu_0=15.4$$.

### 母平均の検定（両側検定） 部品の直径 p.147

# データの設定
samples = [15.5, 15.7, 15.4, 15.4, 15.6, 15.4, 15.6, 15.5, 15.4]  # 標本
mu0 = 15.4  # 帰無仮説の母平均μ₀

Next, we will perform a population mean test with a significance level of 5% using a two-tailed test.

# 関数利用
result = pop_mean_ttest(samples, mu0=mu0)
result

【Execution Result】
The test statistic T=2.683 is greater than the rejection limit (upper threshold) of 2.306 (it falls into the rejection region according to the text), so the null hypothesis is rejected at the 5% significance level.

In other words, the diameter of the recent manufactured products is different from the specifications!

② Visualization of the test
Let's visualize the t-distribution to see how it is rejected.

### 図示 p.148

## 設定と準備
# 統計関連の値
N = len(samples)                  # 標本サイズ
alpha = result['alpha']           # 有意水準
lower, upper = result['c_value']  # 棄却限界値
t_val = result['t_value']         # 検定統計量
t_dist = stats.t(df=N - 1)        # t分布
# グラフ描画用設定
color = 'tab:blue'                            # 基本の色
x_min, x_max = -5, 5                          # x軸の最小値、最大値
x_val = np.linspace(x_min, x_max, 1001)       # 確率密度関数用のx軸の値
x_val_lower = np.linspace(x_min, lower, 101)  # 下側棄却域用のx軸の値
x_val_upper = np.linspace(upper, x_max, 101)  # 上側棄却域用のx軸の値

## 描画処理
# 描画領域の設定
plt.figure(figsize=(7, 3))
# t分布の確率密度関数の描画
plt.plot(x_val, t_dist.pdf(x_val), color=color, label=f'自由度({N}-1)の$t$分布')
# 下側棄却域の塗りつぶし
plt.fill_between(x_val_lower, 0, t_dist.pdf(x_val_lower),
                 color=color, alpha=0.3, label=f'有意水準{alpha}の棄却域')
# 上側棄却域の塗りつぶし
plt.fill_between(x_val_upper, 0, t_dist.pdf(x_val_upper),
                 color=color, alpha=0.3)
# 検定統計量の垂直線の描画
plt.axvline(t_val, color='tab:red', ls='--', label=f'検定統計量={t_val:.3f}')
# 修飾
plt.xticks([lower, 0, upper])
plt.legend();

【Execution Result】

The bell-shaped curve is the probability density function of the t-distribution with degrees of freedom N-1=8.
The blue regions at both ends are the rejection regions for a 5% significance level.
Since it is a two-tailed test, there is a 2.5% rejection region at each end.
The red dotted line is the test statistic T.

Since T falls into the right (upper) rejection region, the null hypothesis is rejected.

③ Performing the test with Python libraries
Usually, we should utilize well-established libraries.
Here, we will use scipy.stats and pingouin.
Both pass the sample data and the null hypothesis population mean μ0 as arguments.

🖲️scipy.stats

# scipy.stats利用
stats.ttest_1samp(samples, popmean=mu0)

【Execution Result】
We obtained the test statistic T, the p-value, and the degrees of freedom of the t-distribution.

🖲️pingouin

# pingouin利用
pg.ttest(x=samples, y=mu0)

【Execution Result】
The main outputs are T: test statistic T, p-val: p-value, and CI95%: 95% confidence interval.

◆ ◆ ◆

■ Population Mean Test (One-tailed test) Example p.150
This is the potato harvest data from the example on p.150.
We will use a population mean test to investigate whether the harvest has "increased" (become larger) compared to before due to the new fertilizer.

The hypotheses are as follows:
・Null hypothesis H0: "Population mean μ=41.4"
・Alternative hypothesis H1: "Population mean μ > 41.4"
⇒ One-tailed test (greater): Interested in whether it has become larger than before

① Performing the test with a custom function
First, we set the data and the null hypothesis μ0=41.4.

### 母平均の検定（片側検定） ジャガイモの収穫量 p.150

# データの設定
samples = [42.9, 43.7, 41.2, 40.8, 41.3, 44.2]  # 標本
mu0 = 41.4  # 帰無仮説の母平均μ₀

Next, we will perform a population mean test with a significance level of 5% using a one-tailed test (greater).

# 関数利用
result = pop_mean_ttest(samples, mu0=mu0, alternative='greater')
result

【Execution Result】
The test statistic T=1.615 is smaller than the rejection limit (upper threshold) of 2.015 and falls outside the rejection region, so the null hypothesis cannot be rejected at the 5% significance level.

In other words, we cannot say that the fertilizer is effective.

② Visualization of the test
Let's visualize the $${t}$$ distribution and see how it cannot be rejected.

### 図示 p.151

## 設定と準備
# 統計関連の値
N = len(samples)                  # 標本サイズ
alpha = result['alpha']           # 有意水準
upper = result['c_value']         # 棄却限界値
t_val = result['t_value']         # 検定統計量
t_dist = stats.t(df=N - 1)        # t分布
# グラフ描画用設定
color = 'tab:blue'                            # 基本の色
x_min, x_max = -5, 5                          # x軸の最小値、最大値
x_val = np.linspace(x_min, x_max, 1001)       # 確率密度関数用のx軸の値
x_val_upper = np.linspace(upper, x_max, 101)  # 上側棄却域用のx軸の値

## 描画処理
# 描画領域の設定
plt.figure(figsize=(7, 3))
# t分布の確率密度関数の描画
plt.plot(x_val, t_dist.pdf(x_val), color=color, label=f'自由度({N}-1)の$t$分布')
# 上側棄却域の塗りつぶし
plt.fill_between(x_val_upper, 0, t_dist.pdf(x_val_upper),
                 color=color, alpha=0.3, label=f'有意水準{alpha}の棄却域')
# 検定統計量の垂直線の描画
plt.axvline(t_val, color='tab:red', ls='--', label=f'検定統計量={t_val:.3f}')
# 修飾
plt.xticks([0, upper])
plt.legend();

【Execution Result】
Since it is an upper one-tailed test, the blue rejection region is in the upper $${5\%}$$.

The bell-shaped curve is the probability density function of the $${t}$$ distribution with degrees of freedom $${N-1=5}$$.
The blue region on the upper side is the rejection region at the $${5\%}$$ significance level.
The red dotted line is the test statistic $${T}$$.

Since $${T}$$ is outside the rejection region, the null hypothesis cannot be rejected.

③ Executing the test with Python libraries
We use scipy.stats and pingouin.
Both pass the sample data, the null hypothesis population mean $${\mu_0}$$, and the one-tailed test (greater) as arguments.

🖲️scipy.stats

# scipy.stats利用
stats.ttest_1samp(samples, popmean=mu0, alternative='greater')

【Execution Result】
We obtained the test statistic $${T}$$, the $${p}$$ value, and the degrees of freedom of the $${t}$$ distribution.

🖲️pingouin

# pingouin利用
pg.ttest(x=samples, y=mu0, alternative='greater')

【Execution Result】
The main outputs are T: test statistic $${T}$$, p-val: $${p}$$ value, and CI95%: 95% confidence interval.

Understanding Check: Population Mean Test p.152

This is a population mean test regarding product strength.
With a $${5\%}$$ significance level and a two-tailed test, we will solve it calmly with code.

Set the sample data and the null hypothesis population mean $${\mu_0}$$.

### 母平均の検定（片側検定） ザイルの破断強度 p.152

# データの設定
samples = [4480, 4510, 4570, 4360, 4240, 4520,
           4260, 4650, 4380, 4130, 4530, 4290]  # 標本
mu0 = 4500  # 帰無仮説の母平均μ₀

Execute the test using a custom function.

# 関数利用
result = pop_mean_ttest(samples, mu0=mu0)
result

【Execution Result】
The null hypothesis could not be rejected.

Execute the test using scipy.stats.

# scipy.stats利用
stats.ttest_1samp(samples, popmean=mu0)

【Execution Result】

Execute the test using pingouin.

# pingouin利用
pg.ttest(x=samples, y=mu0)

【Execution Result】

Let's visualize the $${t}$$ distribution and see how it cannot be rejected.

### 図示 p.153

## 設定と準備
# 統計関連の値
N = len(samples)                  # 標本サイズ
alpha = result['alpha']           # 有意水準
lower, upper = result['c_value']  # 棄却限界値
t_val = result['t_value']         # 検定統計量
t_dist = stats.t(df=N - 1)        # t分布
# グラフ描画用設定
color = 'tab:blue'                            # 基本の色
x_min, x_max = -5, 5                          # x軸の最小値、最大値
x_val = np.linspace(x_min, x_max, 1001)       # 確率密度関数用のx軸の値
x_val_lower = np.linspace(x_min, lower, 101)  # 下側棄却域用のx軸の値
x_val_upper = np.linspace(upper, x_max, 101)  # 上側棄却域用のx軸の値

## 描画処理
# 描画領域の設定
plt.figure(figsize=(7, 3))
# t分布の確率密度関数の描画
plt.plot(x_val, t_dist.pdf(x_val), color=color, label=f'自由度({N}-1)の$t$分布')
# 下側棄却域の塗りつぶし
plt.fill_between(x_val_lower, 0, t_dist.pdf(x_val_lower),
                 color=color, alpha=0.3, label=f'有意水準{alpha}の棄却域')
# 上側棄却域の塗りつぶし
plt.fill_between(x_val_upper, 0, t_dist.pdf(x_val_upper),
                 color=color, alpha=0.3)
# 検定統計量の垂直線の描画
plt.axvline(t_val, color='tab:red', ls='--', label=f'検定統計量={t_val:.3f}')
# 修飾
plt.xticks([lower, 0, upper])
plt.legend();

[Execution Results]

Section 4.3 Test for Population Variance

[Summary of Test for Population Variance]

[Hypothesis]
Set the population variance hypothesis $${\sigma_0^2}$$ that you wish to overturn as follows.

The hypothesis you wish to support, i.e., the alternative hypothesis, is set as one of the following three patterns.

◆ ◆ ◆

■ Test for population variance of a normal population: Formula for test statistic p.154
I will borrow the formula for the test statistic in the test for population variance from the text.

Given sample variance $${s^2}$$, sample size $${N}$$, and null hypothesis population variance $${\sigma^2_0}$$, the test statistic for population variance $${\mu}$$ using sample $${{x_1, x_2, \ldots, x_N}}$$ is

which follows a $${\chi^2}$$ distribution with $${N-1}$$ degrees of freedom.

I will define the "test function for population variance of a normal population" according to the formula.
I will use scipy.stats for probability calculations.

### 正規母集団の母分散の検定関数 p.154 ※両側検定のp値が自信ない

# 正規母集団の母分散のχ²検定統計量の算出関数
def chi2_stat_pop_variance(x_list, sigma20):
    N = len(x_list)
    x_bar = sum(x_list) / N
    s2 = sum([(x - x_bar)**2 for x in x_list]) / (N - 1)
    return ((N - 1) * s2) / sigma20

# 正規母集団の母分散のχ²検定関数 ※確率計算はscipy.stats利用
def pop_variance_chi2test(x_list, sigma20, alpha=0.05, alternative='two-sided'):
    N = len(x_list)
    chi2_val = chi2_stat_pop_variance(x_list, sigma20)
    match alternative:
        case 'two-sided':
            c_val_lower = stats.chi2.ppf(q=alpha/2, df=N - 1)
            c_val_upper = stats.chi2.ppf(q=1 - alpha/2, df=N - 1)
            c_val = c_val_lower, c_val_upper
            half_val = stats.chi2.ppf(q=0.5, df=N - 1)
            if chi2_val <= half_val:
                p_val = stats.chi2.cdf(x=chi2_val, df=N - 1) * 2
            else:
                p_val = stats.chi2.sf(x=chi2_val, df=N - 1) * 2
        case 'less':
            c_val = stats.chi2.ppf(q=alpha, df=N - 1)
            p_val = stats.chi2.cdf(x=chi2_val, df=N - 1)
        case 'greater':
            c_val = stats.chi2.ppf(q=1 - alpha, df=N - 1)
            p_val = stats.chi2.sf(x=chi2_val, df=N - 1)
    return {'chi2_value': chi2_val, 'c_value': c_val, 'alpha': alpha,
            'p_value': p_val}

For test data $${[1, 2, 1, 2, 3]}$$, I will perform a two-tailed test for population variance with null hypothesis $${H_0}$$ "population variance $${\sigma^2=\sigma^2_0=6}$$" and a significance level of 5% ($${\alpha=0.05}$$).

# テスト
x_list = [1, 2, 1, 2, 3]
pop_variance_chi2test(x_list, sigma20=6, alpha=0.05)

[Execution Results]
"It is significant at the 5% significance level, and the null hypothesis can be rejected."

I will briefly explain the output content.

• chi2_value: Realized value of the test statistic $${T}$$

• c_value: Rejection limit value
  Threshold of the test statistic $${T}$$ at which the null hypothesis can be rejected
  When the realized value of $${T}$$ is outside the rejection limit, the null hypothesis can be rejected

• alpha: Significance level $${\alpha}$$

• p_value: $${p}$$ value
  The probability that the test statistic becomes the realized value under the null hypothesis
  When the $${p}$$ value is less than or equal to (or less than) the significance level, the null hypothesis can be rejected

◆ ◆ ◆

■ Test of population variance (one-sided test) Example p.155
In the example on p.155, we use a test of population variance to investigate whether the variation in the diameter of manufactured products has improved (decreased) due to a new manufacturing method.

The hypotheses are as follows:
・Null hypothesis $${H_0}$$: "Population variance $${\sigma^2 = 0.4^2}$$"
・Alternative hypothesis $${H_1}$$: "Population variance $${\sigma^2 < 0.4^2}$$"
　⇒One-sided test (less): Interested in whether the variation has decreased

① Execute test with custom function
Set the data and the null hypothesis $${\sigma^2_0=0.4^2}$$, then execute the test.

### 母分散の検定（片側検定）  リベットの直径 p.155

# 関数利用
samples = [35.2, 34.5, 34.9, 35.2, 34.8, 35.1, 34.9, 35.2, 34.9, 34.8]
result = pop_variance_chi2test(samples, sigma20=0.4**2, alpha=0.05,
                               alternative='less')
result

【Execution Results】
The test statistic $${T=2.906}$$ is smaller than the rejection limit (lower threshold) $${3.325}$$ (according to the text, it falls into the rejection region), so the null hypothesis is rejected at the $${5\%}$$ significance level.

In other words, the variation has decreased due to the new manufacturing method!

② Visualization of the test
Let's visualize the $${chi^2}$$ distribution to see how it is rejected.

### 図示 p.153

## 設定と準備
# 統計関連の値
N = len(samples)                  # 標本サイズ
alpha = result['alpha']           # 有意水準
lower = result['c_value']         # 棄却限界値
chi2_val = result['chi2_value']   # 検定統計量
chi2_dist = stats.chi2(df=N - 1)  # カイ二乗分布
# グラフ描画用設定
color = 'tab:blue'                            # 基本の色
x_min, x_max = 0, 25                          # x軸の最小値、最大値
x_val = np.linspace(x_min, x_max, 1001)       # 確率密度関数用のx軸の値
x_val_lower = np.linspace(x_min, lower, 101)  # 下側棄却域用のx軸の値

## 描画処理
# 描画領域の設定
plt.figure(figsize=(7, 3))
# カイ二乗分布の確率密度関数の描画
plt.plot(x_val, chi2_dist.pdf(x_val), color=color, 
         label=f'自由度({N}-1)の$\chi^2$分布')
# 下側棄却域の塗りつぶし
plt.fill_between(x_val_lower, 0, chi2_dist.pdf(x_val_lower),
                 color=color, alpha=0.3, label=f'有意水準{alpha}の棄却域')
# 検定統計量の垂直線の描画
plt.axvline(chi2_val, color='tab:red', ls='--', label=f'検定統計量={chi2_val:.3f}')
# 修飾
plt.xticks([lower, 0])
plt.legend();

【Execution Results】

The curve skewed to the left is the probability density function of the $${\chi^2}$$ distribution with degrees of freedom $${N-1=9}$$.
The blue region on the lower side is the rejection region at the $${5\%}$$ significance level.
The red dotted line is the test statistic $${T}$$.

Since $${T}$$ is in the rejection region, the null hypothesis is rejected.

Section 4.4 Test of Population Proportion

【Summary of Test of Population Proportion】

【Hypothesis】
Set the hypothesis $${p_0}$$ for the population proportion you want to overturn as follows:

The hypothesis you want to support, the alternative hypothesis, is set as one of the following three patterns.

◆ ◆ ◆

■ Testing population proportion of a binomial population: Formula for test statistic p.158
I will borrow the formula for the test statistic in population proportion testing from the text.

Given sample proportion $${m/N}$$, sample size $${N}$$, and null hypothesis population proportion $${p_0}$$, the test statistic for population proportion $${p}$$ is

follows a standard normal distribution.

I will define a "test function for the population proportion of a binomial population" according to the formula.
I will use scipy.stats for probability calculations.

### 二項母集団の母比率の検定関数 p.158

# 二項母集団の母分散のz検定統計量の算出関数
def z_stat_pop_ratio(m, N, p0):
    p_hat = m / N
    return (p_hat - p0) / math.sqrt(p0 * (1 - p0) / N)

# 二項母集団の母比率のz検定関数 ※確率計算はscipy.stats利用
def pop_ratio_ztest(m, N, p0, alpha=0.05, alternative='two-sided'):
    z_val = z_stat_pop_ratio(m, N, p0)
    std_norm_dist = stats.norm(loc=0, scale=1)
    match alternative:
        case 'two-sided':
            c_val = std_norm_dist.ppf(q=1 - alpha/2)
            c_val = -c_val, c_val
            p_val = std_norm_dist.sf(x=abs(z_val)) * 2
        case 'less':
            c_val = std_norm_dist.ppf(q=alpha)
            p_val = std_norm_dist.cdf(x=z_val)
        case 'greater':
            c_val = std_norm_dist.ppf(q=1 - alpha)
            p_val = std_norm_dist.sf(x=z_val)
    return {'z_value': z_val, 'c_value': c_val, 'alpha': alpha, 'p_value': p_val}

For a sample size of $${N=346}$$ and a count of data in the category of interest of $${m=29}$$, I will perform a population proportion test with the null hypothesis $${H_0}$$ "population proportion $${p=p_0=0.12}$$", a significance level of 5% ($${\alpha=0.05}$$), and a one-tailed test (less).

# テスト
pop_ratio_ztest(m=29, N=346, p0=0.12, alternative='less')

[Execution Result]
"It is significant at the 5% significance level, and the null hypothesis can be rejected."

I will briefly explain the output content.

• z_value: Realized value of the test statistic $${T}$$

• c_value: Rejection limit value
  The threshold for the test statistic $${T}$$ at which the null hypothesis can be rejected
  When the realized value of $${T}$$ is outside the rejection limit, the null hypothesis can be rejected

• alpha: Significance level $${\alpha}$$

• p_value: $${p}$$ value
  The probability that the test statistic will be the realized value under the null hypothesis
  When the $${p}$$ value is less than or equal to the significance level, the null hypothesis can be rejected

💡 Reference information: Points to note when using Python libraries
I searched for libraries that could be used for testing population means and found statsmodels, so I will try it out.

### statsmodelsで答え合わせしようと思ったが、異なる値が算出される
z_val, p_val = proportions_ztest(count=29, nobs=346, value=0.12,
                                 alternative='smaller')
print('z_value:', z_val, ', p_value:', p_val)

[Execution Result]
A different result was output compared to the custom function based on the formula in the text.

After researching various things, I assumed the following state.

• proportions_ztest is likely intended primarily for "two-sample tests for the difference in population proportions"

• It incorporates the pooled sample proportion into the calculation

• The one-sample population proportion test, like the one in this article, also uses the pooled sample proportion calculation.

I tested the hypothesis above.

### 二項母集団の母比率の検定関数
# statsmodelsと同じロジック（プールした標本比率の公式を使っていると思われる）

# 二項母集団の母分散のz検定統計量の算出関数2
def z_stat_pop_ratio2(m, N, p0):
    p_hat = m / N
    return (p_hat - p0) / math.sqrt(p_hat * (1 - p_hat) / N) # ★相違点

# 二項母集団の母比率のz検定関数2 ※確率計算はscipy.stats利用
def pop_ratio_ztest2(m, N, p0, alpha=0.05, alternative='two-sided'):
    z_val = z_stat_pop_ratio2(m, N, p0)
    std_norm_dist = stats.norm(loc=0, scale=1)
    match alternative:
        case 'two-sided':
            c_val = std_norm_dist.ppf(q=1 - alpha/2)
            c_val = -c_val, c_val
            p_val = std_norm_dist.sf(x=abs(z_val)) * 2
        case 'less':
            c_val = std_norm_dist.ppf(q=alpha)
            p_val = std_norm_dist.cdf(x=z_val)
        case 'greater':
            c_val = std_norm_dist.ppf(q=1 - alpha)
            p_val = std_norm_dist.sf(x=z_val)
    return {'z_value': z_val, 'c_value': c_val, 'alpha': alpha, 'p_value': p_val}

# テスト
pop_ratio_ztest2(m=29, N=346, p0=0.12, alternative='less')

[Execution Results]
The hypothesis is highly likely to be correct.

When trying a one-sample population proportion test with statsmodels, keep in mind that there is a difference in the calculation method!

◆ ◆ ◆

■ Population Proportion Test (One-tailed test) Example p.159
We will examine the effect of reducing the retention rate using the new teaching method from the example on p.159 using a population proportion test.

The hypotheses are as follows:
・Null hypothesis $${H_0}$$ "Population proportion $${p = 0.12}$$"
・Alternative hypothesis $${H_1}$$ "Population proportion $${p < 0.12}$$"
　⇒One-tailed test (less): Interested in whether it has decreased

Set the data and the null hypothesis $${σ^2_0=0.4^2}$$, and execute the test using a custom function.

### 母比率の検定（片側検定） 新指導法の効果 p.159

# 関数利用
result =pop_ratio_ztest(m=29, N=346, p0=0.12, alternative='less')
result

[Execution Results]
The test statistic $${T=-2.071}$$ is smaller than the rejection limit (lower threshold) $${-1.645}$$ (it falls into the rejection region according to the text), so the null hypothesis is rejected at the 5% significance level.

In other words, the retention rate has decreased due to the new teaching method!

Let's visualize the standard normal distribution to see how it is rejected.

### 図示 p.160

## 設定と準備
# 統計関連の値
alpha = result['alpha']           # 有意水準
lower = result['c_value']         # 棄却限界値
z_val = result['z_value']         # 検定統計量
std_norm_dist = stats.norm(loc=0, scale=1)  # 標準正規分布
# グラフ描画用設定
color = 'tab:blue'                            # 基本の色
x_min, x_max = -3, 3                          # x軸の最小値、最大値
x_val = np.linspace(x_min, x_max, 1001)       # 確率密度関数用のx軸の値
x_val_lower = np.linspace(x_min, lower, 101)  # 下側棄却域用のx軸の値

## 描画処理
# 描画領域の設定
plt.figure(figsize=(7, 3))
# 標準正規分布の確率密度関数の描画
plt.plot(x_val, std_norm_dist.pdf(x_val), color=color, label=f'標準正規分布')
# 下側棄却域の塗りつぶし
plt.fill_between(x_val_lower, 0, std_norm_dist.pdf(x_val_lower),
                 color=color, alpha=0.3, label=f'有意水準{alpha}の棄却域')
# 検定統計量の垂直線の描画
plt.axvline(z_val, color='tab:red', ls='--', label=f'検定統計量={z_val:.3f}')
# 修飾
plt.xticks([lower, 0])
plt.legend();

[Execution Results]

The bell-shaped curve is the probability density function of the standard normal distribution.
The blue area at the bottom is the rejection region for the 5% significance level.
The red dotted line is the test statistic $${T}$$.

Since $${T}$$ is in the rejection region, the null hypothesis is rejected.

Let's have ChatGPT wrap up the end of the article!

📘 A word from ChatGPT:

That is all for this copying session.

Series Articles

Next Article

Previous Article

Table of Contents

Blog Introduction

I am writing seven series of articles on note.
Please come and take a look!

1. Relaxed Statistics

This is a blog that roughly explores probability and statistics using the Statistical Test Grade 2 problem collection as a guide.
Feel free to treat it like casual conversation. Please come and take a look.
It corresponds to the Statistical Test Grade 2 Official Problem Collection CBT Edition.
Sample code for Python and EXCEL is also available.

2. Experiment! Fun Bayesian Modeling 1 & 2 with PyMC Ver. 5

I will draw and analyze the Bayesian models used in the psychological research from the books "Fun Bayesian Modeling" and "Fun Bayesian Modeling 2" using PyMC Ver. 5.
Starting with these books, many Bayesian models are written in R + Stan.
I strive to explore the possibilities of PyMC and make Bayesian modeling easy to practice.
Since these are familiar and easy-to-visualize themes, please run them with PyMC and let's enjoy it together!

3. Experiment! Iwanami Data Science 1 Bayesian Modeling with PyMC Ver. 5

I will draw and analyze the Bayesian models by four Bayesians from the book "Experiment! Iwanami Data Science Vol. 1" using PyMC Ver. 5.
This book is a great resource for learning the basics of Bayesian programming.
I feel like I've become friends with Bayesian methods by happily running PyMC models.
Everyone, please run them with PyMC and let's play and learn together!

4. Fun Copying: Bayesian, Python, etc.

I will blog about the results of my "book copying activities" for Bayesian, Python, and others.
I am mainly working on translations into Python.
I hope this serves as sample code for fellow learners who are also copying code 🍀

5. Introduction to Time Series Analysis for Psychology with R and Stan, using Python and PyMC Ver. 5

I will practice the time series analysis from the book "Introduction to Time Series Analysis for Psychology with R and Stan" using Python and PyMC Ver. 5.
This book is packed with time series analysis themes!
I realized the depth of time series analysis.
I will enjoy learning time series analysis with my favorite language, Python.

6. Writing about Data Science-like things

I write columns on statistics, data analysis, AI, machine learning, and Python on an irregular basis.
There are many articles related to statistics and data science books.
Series on "Statistics," "Python," "Math and Python," and "R" have been created.

7. Python Machine Learning Programming Practice Log

I have written articles about my various thoughts while studying the book 'Python Machine Learning Programming: PyTorch & scikit-learn Edition'.
This book is a textbook for scikit-learn and PyTorch.
Please feel free to give it a try if you like.

Thank you very much for reading until the end.