描述性統(tǒng)計分析

描述性統(tǒng)計所提取的統(tǒng)計的信息稱為統(tǒng)計量，包括頻數(shù)與頻率，反映集中趨勢的均值、中位數(shù)、眾數(shù)和分位數(shù)，反映離散程度的極差、方差和標(biāo)準(zhǔn)差，反映分布形狀（相對于正態(tài)分布）的偏度和峰度。

變量分為類別變量和數(shù)值變量，類別變量往往被作為維度，數(shù)值變量往往被作為指標(biāo)。類別可以經(jīng)過特定的轉(zhuǎn)換轉(zhuǎn)換為數(shù)值，從而作為指標(biāo)，數(shù)值變量也可以經(jīng)過特定的分箱或轉(zhuǎn)換轉(zhuǎn)換為文本型變量，從而作為類別或維度。

頻數(shù)與頻率

最基本的統(tǒng)計量就是頻數(shù)與頻率，它們適用于類別變量。

頻數(shù)，指數(shù)據(jù)中類別變量每個不同取值出現(xiàn)的次數(shù)。

頻率，指每個類別變量的頻數(shù)與總次數(shù)的比值，通常采用百分?jǐn)?shù)表示。

下面我們以鳶尾花數(shù)據(jù)集為例說明這些概念，首先導(dǎo)包并讀取數(shù)據(jù)：

import numpy as np
import pandas as pd
import matplotlib.pyplot as plt
import seaborn as sns
from sklearn.datasets import load_iris

# 設(shè)置seaborn繪圖的樣式，并顯示中文
sns.set(style="darkgrid")
plt.rcParams["font.family"] = "SimHei"
plt.rcParams["axes.unicode_minus"] = False

iris = load_iris()
data = np.column_stack([iris.data, iris.target])
data = pd.DataFrame(data, columns=["sepal_length", "sepal_width", "petal_length", "petal_width", "type"])
data.sample(10)

注意：column_stack會先將一維數(shù)組轉(zhuǎn)換為2維列向量后，按列水平進(jìn)行拼接

numpy拼接小知識補充

np.column_stack([iris.data, iris.target])

等價于：

np.hstack([iris.data, iris.target.reshape(-1, 1)])

或：

np.concatenate([iris.data, iris.target.reshape(-1, 1)], axis=1)

或：

np.c_[iris.data, iris.target.reshape(-1, 1)]

而iris.target.reshape(-1, 1)也可以用新增軸來表示，等價于：

iris.target[:, np.newaxis]

np.newaxis的本質(zhì)等于None，可以直接用None替換，即：

iris.target[:, None]

個人覺得column_stack最方便，因為實現(xiàn)了將一維數(shù)組自動轉(zhuǎn)換為2維列向量。

下面計算鳶尾花數(shù)據(jù)中，每個類別出現(xiàn)的頻數(shù)：

frequency = data["type"].value_counts()
frequency

2.0    50
1.0    50
0.0    50
Name: type, dtype: int64

將頻數(shù)除以總數(shù)即表示每個類別出現(xiàn)的頻率，使用百分比表示：

percentage = frequency * 100 / len(data)
percentage

2.0    33.333333
1.0    33.333333
0.0    33.333333
Name: type, dtype: float64

反映趨中趨勢的幾個指標(biāo)

有均值、中位數(shù)、眾數(shù)和分位數(shù)。

均值、中位數(shù)和眾數(shù)

均值，即平均值，其為—組數(shù)據(jù)的總和除以數(shù)據(jù)的個數(shù)。

中位數(shù)，將一組數(shù)據(jù)升序排列，位于該組數(shù)據(jù)最中間位置的值，就是中位數(shù)。如果數(shù)據(jù)個數(shù)為偶數(shù)，則取中間兩個數(shù)值的均值。

眾數(shù)，一組數(shù)據(jù)中出現(xiàn)次數(shù)最多的值

從上圖可以看到在正態(tài)分布下，三者是相同的，在偏態(tài)分布下，三者會有所不同。

數(shù)值變量在正態(tài)分布時，可以使用均值與中值表示集中趨勢。在偏態(tài)分布下，均值容易受極端值的影響，所以一般使用中值表示集中趨勢。

類別變量通常使用眾數(shù)表示集中趨勢，但眾數(shù)在一組數(shù)據(jù)中可能不是唯一的。

舉個例子，要統(tǒng)計居民的總體收入水平，使用哪項指標(biāo)衡量更合適呢？首先收入屬于數(shù)值變量，可以使用均值與中位數(shù)表示集中趨勢。

但20%的人掌握著80%的人財富，居民收入是個嚴(yán)重右偏的分布，均值會受極端值的影響，所以使用中位數(shù)指標(biāo)更合適。

下面我們計算花萼長度的均值，中位數(shù)以及眾數(shù)：

# 計算花萼長度的均值。
mean = data["sepal_length"].mean()
# 計算花萼長度的中位數(shù)。
median = data["sepal_length"].median()
# 計算花萼長度的眾數(shù)。
mode = data["sepal_length"].mode()
print(mean, median)
# mode方法返回的是Series類型。
print(mode)

5.843333333333335 5.8
0    5.0
dtype: float64

也可以使用scipy的stats模塊來求一組數(shù)據(jù)的眾數(shù)。

from scipy import stats
stats.mode(data["sepal_length"])

ModeResult(mode=array([5.]), count=array([10]))

同時會返回該眾數(shù)出現(xiàn)的頻次。

看看分布：

# 繪制數(shù)據(jù)的分布（直方圖 + 密度圖）。
sns.distplot(data["sepal_length"])
# 繪制垂直線。
plt.axvline(mean, ls="-", color="r", label="均值")
plt.axvline(median, ls="-", color="g", label="中值")
plt.axvline(mode, ls="-", color="indigo", label="眾數(shù)")
plt.legend()

Serise的mode方法和stats.mode()方法的區(qū)別

Serise的mode方法的返回值類型是Serise。stats.mode()方法的返回值類型是 ModeResult

如果眾數(shù)的值不唯一，Series的mode()方法會顯示所有眾數(shù)，而stats.mode()方法只顯示其中一個，但同時能知道該眾數(shù)的個數(shù)。

s = pd.Series([1, 1, 1, 2, 2, 2, 3, 3, 3, 4, 5])

對于Series的mode方法：

s.mode()

0    1
1    2
2    3
dtype: int64

對于stats.mode()方法：

stats.mode(s)

ModeResult(mode=array([1], dtype=int64), count=array([3]))

分位數(shù)

分位數(shù)，通過n-1個分位將數(shù)據(jù)劃分為n個區(qū)間，使得每個區(qū)間的數(shù)值個數(shù)相等（或近似相等）。其中，n為分位數(shù)的數(shù)量．常用的分位數(shù)有四分位數(shù)與百分位數(shù)。

以四分位數(shù)為例，通過3個分位，將數(shù)據(jù)劃分為4個區(qū)間（百分位數(shù)可根據(jù)四分位數(shù)對比理解）。

• 第1個分位稱為1/4分位（下四分位）數(shù)據(jù)中1/4的數(shù)據(jù)小于該分位值。
• 第2個分位稱為2/4分位（中四分位）數(shù)據(jù)中2/4的數(shù)據(jù)小于該分位值。
• 第3個分位稱為3/4分位（上四分位）數(shù)據(jù)中3/4的數(shù)據(jù)小于該分位值。

使用Numpy中計算分位數(shù)：

x = [1, 3, 10, 15, 18, 20, 23, 40]
# quantile與percentile都可以計算分位數(shù)，不同的是，quantile方法，
# q（要計算的分位數(shù)）的取值范圍為[0, 1]，而percentile方法，q的取值范圍為[0, 100]。
print(np.quantile(x, q=[0.25, 0.5, 0.75]))
print(np.percentile(x, q=[25, 50, 75]))

[ 8.25 16.5  20.75]
[ 8.25 16.5  20.75]

使用Pandas中計算分位數(shù)：

x = [1, 3, 10, 15, 18, 20, 23, 40]
s = pd.Series(x)
s.describe()

count     8.000000
mean     16.250000
std      12.395276
min       1.000000
25%       8.250000
50%      16.500000
75%      20.750000
max      40.000000
dtype: float64

describe方法支持自定義分位位置：

s.describe(percentiles=[0.25, 0.9])

count     9.000000
mean     16.777778
std      11.702326
min       1.000000
25%      10.000000
50%      18.000000
90%      26.400000
max      40.000000
dtype: float64

分位數(shù)計算的原理

首先計算分位點所在的索引位置：

x = np.array([1, 3, 10, 15, 18, 20, 23, 40])
n = len(x)
# 計算四分位的索引（index）。
q1_index = (n - 1) * 0.25
q2_index = (n - 1) * 0.5
q3_index = (n - 1) * 0.75
print(q1_index, q2_index, q3_index)

1.75 3.5 5.25

索引位置不是整數(shù)時，使用最近位置的兩個整數(shù)，加權(quán)計算來得到四分位的位置，每個整數(shù)的權(quán)重為距離的反比。加權(quán)計算：

index = np.array([q1_index, q2_index, q3_index])
# 計算左邊元素的值。
left = np.floor(index).astype(np.int8)
# 計算右邊元素的值。
right = np.ceil(index).astype(np.int8)
# 獲取index的小數(shù)部分與整數(shù)部分。
weight, _ = np.modf(index)
# 根據(jù)左右兩邊的整數(shù)，加權(quán)計算四分位數(shù)的值。權(quán)重與距離成反比。
q = x[left] * (1 - weight) + x[right] * weight
print(q)

[ 8.25 16.5  20.75]

當(dāng)索引位置是整數(shù)時，計算過程可以簡化為：

x = np.array([1, 3, 10, 15, 18, 20, 21, 23, 40])
n = len(x)
# 計算四分位的索引（index）。
index = (np.array([0.25, 0.5, 0.75])*(n - 1)).astype(np.int8)
print(x[index])

[10, 18, 21]

反映離散程度的極差、方差和標(biāo)準(zhǔn)差

極差指一組數(shù)據(jù)中，最大值與最小值之差。

方差體現(xiàn)的是一組數(shù)據(jù)中，每個元素與均值偏離的大小。

σ 2 = 1 n ? 1 ∑ i = 1 n ( x i ? x ˉ ) 2 \Huge{\sigma^{2}=\frac{1}{n-1} \sum_{i=1}^{n}\left(x_{i}-\bar{x}\right)^{2}} σ2=n?11?∑i=1n?(xi??xˉ)2

• x i x_i xi?：數(shù)組中的每個元素。
• n：數(shù)組元素的個數(shù)。
• x ˉ \bar{x} xˉ：數(shù)組中所有元素的均值。

標(biāo)準(zhǔn)差為方差的開方。

σ = 1 n ? 1 ∑ i = 1 n ( x i ? x ˉ ) 2 \Huge{\sigma=\sqrt{\frac{1}{n-1} \sum_{i=1}^{n}\left(x_{i}-\bar{x}\right)^{2}}} σ=n?11?∑i=1n?(xi??xˉ)2 ?

關(guān)于三者，說明如下：

• 極差的計算非常簡單，但是極差沒有充分的利用數(shù)據(jù)信息。
• 方差（標(biāo)準(zhǔn)差）可以體現(xiàn)數(shù)據(jù)的分散性，方差（標(biāo)準(zhǔn)差）越大，數(shù)據(jù)越分散，方差（標(biāo)準(zhǔn)差）越小，數(shù)據(jù)越集中。
• 方差（標(biāo)準(zhǔn)差）也可以體現(xiàn)數(shù)娼的波動性（穩(wěn)定性）。方差（標(biāo)準(zhǔn)差）越大，數(shù)據(jù)波動性越大，方差（標(biāo)淮差）越小，數(shù)據(jù)波動性越小。
• 當(dāng)數(shù)據(jù)較大時，也可以使用n代替n-1

# 計算極差。
sub = data["sepal_length"].max() - data["sepal_length"].min()
# 計算方差。
var = data["sepal_length"].var()
# 計算標(biāo)準(zhǔn)差。
std = data["sepal_length"].std()
print(sub, var, std)

3.6000000000000005 0.6856935123042505 0.8280661279778629

花瓣長度的方差較大，花瓣寬度的方差較小，繪圖對比：

plt.figure(figsize=(15, 4))
plt.ylim(-0.5, 1.5)
plt.plot(data["petal_length"], np.zeros(len(data)), ls="",
         marker="o", ms=10, color="g", label="花瓣長度")
plt.plot(data["petal_width"], np.ones(len(data)), ls="",
         marker="o", ms=10, color="r", label="花瓣寬度")
plt.axvline(data["petal_length"].mean(), ls="--", color="g", label="花瓣長度均值")
plt.axvline(data["petal_width"].mean(), ls="--", color="r", label="花瓣寬度均值")
plt.legend()

反映分布形狀的偏度和峰度

偏度

偏度是統(tǒng)計數(shù)據(jù)分布偏斜方向和程度的度量，是統(tǒng)計數(shù)據(jù)分布非對稱程度的數(shù)字特征。

• 如果數(shù)據(jù)對稱分布（例如正態(tài)分布），則偏度為0。
• 如果數(shù)據(jù)左偏分布．則偏度小于O。
• 如果數(shù)據(jù)右偏分布．則偏度大于0。

# 構(gòu)造左偏分布數(shù)據(jù)。
t1 = np.random.randint(1, 11, size=100)
t2 = np.random.randint(11, 21, size=500)
t3 = np.concatenate([t1, t2])
left_skew = pd.Series(t3)
# 構(gòu)造右偏分布數(shù)據(jù)。
t1 = np.random.randint(1, 11, size=500)
t2 = np.random.randint(11, 21, size=100)
t3 = np.concatenate([t1, t2])
right_skew = pd.Series(t3)
# 計算偏度。
print(left_skew.skew(), right_skew.skew())
# 繪制核密度圖。
sns.kdeplot(left_skew, shade=True, label="左偏")
sns.kdeplot(right_skew, shade=True, label="右偏")
plt.legend()

-0.9911238058650503 0.7820903371872946

峰度

峰度是描述總體中所有取值分布形態(tài)陡緩程度的統(tǒng)計量。可以將峰度理解為數(shù)據(jù)分布的高矮程度。峰度的比較是相對于標(biāo)準(zhǔn)正態(tài)分布的。

• 對于標(biāo)準(zhǔn)正態(tài)分布，峰度為0．
• 如果峰度大于0，則密度圖高于標(biāo)準(zhǔn)正態(tài)分布。
  • 數(shù)據(jù)在分布上比標(biāo)準(zhǔn)正態(tài)分布密集，方差（標(biāo)淮羞）較小．
• 如果峰度小于0，則密度圖低于標(biāo)準(zhǔn)正態(tài)分布。
  • 說明數(shù)據(jù)在分布上比標(biāo)準(zhǔn)正態(tài)分布分散，方差（標(biāo)準(zhǔn)差）較大。

# 標(biāo)準(zhǔn)正態(tài)分布。
standard_normal = pd.Series(np.random.normal(0, 1, size=10000))
print("標(biāo)準(zhǔn)正態(tài)分布峰度：", standard_normal.kurt(), "標(biāo)準(zhǔn)差：", standard_normal.std())
print("花萼寬度峰度：", data["sepal_width"].kurt(), "標(biāo)準(zhǔn)差：", data["sepal_width"].std())
print("花瓣長度峰度：", data["petal_length"].kurt(), "標(biāo)準(zhǔn)差：", data["petal_length"].std())
sns.kdeplot(standard_normal, label="標(biāo)準(zhǔn)正態(tài)分布")
sns.kdeplot(data["sepal_width"], label="花萼寬度")
sns.kdeplot(data["petal_length"], label="花瓣長度")

標(biāo)準(zhǔn)正態(tài)分布峰度： 0.02338847301358893 標(biāo)準(zhǔn)差： 0.9980947521404823
花萼寬度峰度： 0.2907810623654279 標(biāo)準(zhǔn)差： 0.4335943113621737
花瓣長度峰度： -1.4019208006454036 標(biāo)準(zhǔn)差： 1.7644204199522617

繪制對數(shù)正態(tài)分布的圖像

如果一個分布取對數(shù)后為正態(tài)分布，則該分布稱為對數(shù)正態(tài)分布。

import random
import numpy as np

logdata = [random.lognormvariate(0, 0.6) for i in range(100000)]
data = np.random.normal(0, 0.6, size=100000)
sns.kdeplot(data, label="正態(tài)分布")
sns.kdeplot(np.log(logdata), label="對數(shù)正態(tài)分布取對數(shù)")
sns.kdeplot(np.exp(data), label="正態(tài)分布取冪")
sns.kdeplot(logdata, label="對數(shù)正態(tài)分布")
plt.xlim(-5, 5)

推斷統(tǒng)計分析

推斷統(tǒng)計的概念

總體，就是被研究的全部數(shù)據(jù)，總體中的某個數(shù)據(jù)，就是個體。從總體中抽取部分個體．就構(gòu)成了樣本，樣本中包含的個體數(shù)量，稱為樣本容量。

在實際的研究中，往往無法獲取全部數(shù)據(jù)，只能對總體進(jìn)行抽樣。推斷統(tǒng)計就是研究根據(jù)樣本數(shù)據(jù)去推斷總體數(shù)量特征的方法，它在對樣本數(shù)據(jù)進(jìn)行描述的基礎(chǔ)上，對統(tǒng)計總體的未知數(shù)量特征做出以概率形式表述的推斷，從而通過樣本統(tǒng)計量來估計總體參數(shù)。

推斷統(tǒng)計分析分為參數(shù)估計和假設(shè)檢驗。

參數(shù)估計

點估計與區(qū)間估計

點估計，就是使用樣本的統(tǒng)計量去代替總體參數(shù)。例如，我們要求鳶尾花的平均花瓣長度，就可以使用樣本的均值來估計總體的均值：

import numpy as np
import pandas as pd
import matplotlib.pyplot as plt
import seaborn as sns
from sklearn.datasets import load_iris

# 設(shè)置seaborn繪圖的樣式。
sns.set(style="darkgrid")
plt.rcParams["font.family"] = "SimHei"
plt.rcParams["axes.unicode_minus"] = False

iris = load_iris()
data = np.column_stack((iris.data, iris.target))
data = pd.DataFrame(data, 
        columns=["sepal_length", "sepal_width", "petal_length", "petal_width", "type"])
print(data["petal_length"].mean())

3.7586666666666693

點估計實現(xiàn)簡單，但是容易受到隨機抽樣的影響，無法保證結(jié)論的準(zhǔn)確性。

區(qū)間估計則根據(jù)樣本的統(tǒng)計量，計算出一個可能的區(qū)間與概率，表示總體的參數(shù)會有多少概率位于該區(qū)間中。區(qū)間估計指定的區(qū)間，我們稱為置信區(qū)間，而區(qū)間估計指定的概率，稱為置信度。例如，鳶尾花花瓣長度有70%的可能性在3.4cm-3.8cm之間，3.4cm- 3.8cm就是置信區(qū)間，而70%就是置信度。

總之，點估計是使用一個值來代替總體的參數(shù)值，能夠給出具體的估計值但缺乏準(zhǔn)確性。區(qū)間估計是使用的一個置信區(qū)間與置信度，表示總體參數(shù)有多少可能（置信度）會在該范圍（置信區(qū)間）內(nèi)，能夠給出合理的范圍和支持概率。

經(jīng)過抽樣，獲取—個樣本之后，該如何才能確定置信區(qū)間與置信度呢？區(qū)間估計的基石就是根據(jù)中心極限定理。

中心極限定理

定理內(nèi)容：如果總體（分布不重要）均值為μ，方差為 σ 2 \sigma^{2} σ2。我們進(jìn)行隨機抽樣，樣本容量為n，當(dāng)n增大時，則樣本均值逐漸趨近服從正態(tài)分布： X ˉ ～ N ( μ , σ 2 / n ) \bar{X}\sim N\left(\mu, \sigma^{2} / n\right) Xˉ～N(μ,σ2/n)

該定理說明了總體與樣本之間，在分布上的聯(lián)系。該定理說明在抽樣的樣本容量n足夠大時，進(jìn)行多次抽樣．則每次抽樣會得到—個均值，這些均值會圍繞在總體均值左右，呈正態(tài)分布。均值等于總體的均值，標(biāo)準(zhǔn)差等于總體標(biāo)準(zhǔn)差 σ \sigma σ除以 n \sqrt{n} n ?

程序模擬

下面模擬總體的均值為30，標(biāo)準(zhǔn)差為80，抽樣的樣本容量n為64，看看1000次抽樣的樣本均值是否構(gòu)成均值為30，標(biāo)準(zhǔn)差為10的正態(tài)分布：

data = np.random.normal(30, 80, 100000)
mean_arr = np.zeros(1000)
for i in range(1000):
    mean_arr[i] = np.random.choice(data, size=64, replace=False).mean()
print("均值：", np.mean(mean_arr))
print("標(biāo)準(zhǔn)差：", np.std(mean_arr), 80/np.sqrt(64))
print("偏度：", pd.Series(mean_arr).skew())
print("峰度：", pd.Series(mean_arr).kurt())
sns.distplot(mean_arr)

均值： 29.98515540795136
標(biāo)準(zhǔn)差： 10.387123594380887 10.0
偏度： -0.0499634202251343
峰度： 0.1875444656831493

從上述結(jié)果可以看到，樣本的均值和標(biāo)準(zhǔn)差接近于理論值。

正態(tài)分布的特性

在正態(tài)分布中，數(shù)據(jù)的分布比例如下：

• 以均值為中心，在1倍標(biāo)準(zhǔn)差內(nèi) ( x ˉ ? σ , x ˉ + σ ) (\bar{x}-\sigma, \bar{x}+\sigma) (xˉ?σ,xˉ+σ)，包含約68%的樣本數(shù)據(jù)。
• 以均值為中心，在2倍標(biāo)準(zhǔn)差內(nèi) ( x ˉ ? 2 σ , x ˉ + 2 σ ) (\bar{x}-2\sigma, \bar{x}+2\sigma) (xˉ?2σ,xˉ+2σ)，包含約95%的樣本數(shù)據(jù)。
• 以均值為中心，在3倍標(biāo)準(zhǔn)差內(nèi) ( x ˉ ? 3 σ , x ˉ + 3 σ ) (\bar{x}-3\sigma, \bar{x}+3\sigma) (xˉ?3σ,xˉ+3σ)，包含約99.7%的樣本數(shù)據(jù)。

可以用程序模擬一下：

scale = 50
x = np.random.normal(0, scale, size=10000000)
for times in np.arange(1, 4):
    y = x[(x > -times * scale) & (x < times * scale)]
    print(f"{times}倍標(biāo)準(zhǔn)差：{len(y)/len(x):.1%}")

1倍標(biāo)準(zhǔn)差：68.3%
2倍標(biāo)準(zhǔn)差：95.4%
3倍標(biāo)準(zhǔn)差：99.7%

二分搜索查找當(dāng)正態(tài)分布覆蓋99%樣本數(shù)據(jù)時，大概是多少倍標(biāo)準(zhǔn)差：

import numpy as np
scale = 1
x = np.random.normal(0, scale, size=100000000)
confidenceLevel = 99
arr = np.arange(2, 3, 0.00001)
iterations = 0

low, high = 0, len(arr) - 1
while low <= high:
    mid = (high + low) >> 1
    times = arr[mid]
    rate = len(x[(x > -times * scale) & (x < times * scale)]) * 100 / len(x)
    if rate >= confidenceLevel:
        high = mid - 1
    else:
        low = mid + 1
    iterations += 1
result = round(arr[low], 4)
print(result, "迭代次數(shù)：", iterations)

2.5761 迭代次數(shù)： 17

可以看到標(biāo)準(zhǔn)差大概是2.58倍時，正態(tài)分布覆蓋99%樣本數(shù)據(jù)。

可以通過scipy獲取準(zhǔn)確值：

from scipy.stats import norm

def calcTimes(confidenceLevel):
    alpha = 1 - confidenceLevel
    return norm.ppf(1 - alpha / 2)

calcTimes(0.99)

結(jié)果：

2.5758293035489004

執(zhí)行calcTimes(0.95)的結(jié)果是1.959963984540054，說明1.96倍標(biāo)準(zhǔn)差能覆蓋95%的樣本數(shù)據(jù)。

計算在正態(tài)分布情況下，指定倍數(shù)的標(biāo)準(zhǔn)差能覆蓋多大比例的樣本數(shù)據(jù)可以使用以下命令：

from scipy.stats import norm
# norm.cdf(x=?)計算正態(tài)分布概率（面積）：P(X<x) 

def calcArea(times):
    return (norm.cdf(x=times) - norm.cdf(x=-times))*100

calcArea(2.58)

結(jié)果：

99.01199684844586

置信度與置信區(qū)間

根據(jù)中心極限定理，如果多次抽樣（總體均值 μ \mu μ，方差為 σ 2 \sigma^2 σ2）,則樣本均值( x ˉ \bar{x} xˉ)構(gòu)成正態(tài)分布，滿足 X ˉ ～ N ( μ , σ 2 / n ) \bar{X}\sim N\left(\mu, \sigma^{2}/n\right) Xˉ～N(μ,σ2/n)。如果我們對總體進(jìn)行一次抽樣，則該樣本的均值有95%的概率會在 ( μ ? 2 σ , μ + 2 σ ) (\mu-2\sigma, \mu+2\sigma) (μ?2σ,μ+2σ)，僅有5%的概率會在 ( μ ? 2 σ , μ + 2 σ ) (\mu-2\sigma, \mu+2\sigma) (μ?2σ,μ+2σ)范圍外。根據(jù)小概率事件（很小的概率在一次抽樣中基本不會發(fā)生），如果抽樣的樣本均值在 ( μ ? 2 σ , μ + 2 σ ) (\mu-2\sigma, \mu+2\sigma) (μ?2σ,μ+2σ)之外，我們就可以認(rèn)為，本次抽樣來自總體的均值并非我們所期望的均值。

通常，我們以2倍標(biāo)準(zhǔn)差作為判斷依據(jù)，則以樣本均值為中心，正負(fù)2倍標(biāo)準(zhǔn)差構(gòu)成的區(qū)間，就是置信區(qū)間。而2倍標(biāo)準(zhǔn)差包含了95%的數(shù)據(jù)，因此此時置信度為95%。即我們有95%的信心認(rèn)為，總體均值會在置信區(qū)間之內(nèi)。

下面，我們模擬進(jìn)行一次抽樣，看看實際的總體均值是否會在置信區(qū)間之內(nèi)呢：

mean = np.random.randint(-10000, 10000)  #總體均值
std = 50  # 定義總體標(biāo)準(zhǔn)差
n = 50  # 樣本容量
data = np.random.normal(loc=mean, scale=std, size=10000)
sample = np.random.choice(data, size=n, replace=False)
sample_mean = sample.mean()
print("總體均值：", mean)
print("一次抽樣的樣本均值：", sample_mean)
plt.plot(mean, 0, marker="*", color="orange", ms=15, label="總體均值")
plt.plot(sample_mean, 0, marker="o", color="r", label="樣本均值")

se = std / np.sqrt(n)
min_ = sample_mean - 1.96 * se
max_ = sample_mean + 1.96 * se
print("置信區(qū)間（95%置信度）", (min_, max_))
plt.hlines(0, xmin=min_, xmax=max_, colors="b", label="置信區(qū)間")
plt.axvline(min_, ymin=0.4, ymax=0.6, color="r", ls="--", label="左邊界")
plt.axvline(max_, ymin=0.4, ymax=0.6, color="g", ls="--", label="右邊界")
plt.legend()

總體均值： 8819
一次抽樣的樣本均值： 8824.987475450054
置信區(qū)間（95%置信度） (8811.128182538798, 8838.84676836131)

還可以直接通過scipy計算：

from scipy import stats
stats.norm.interval(0.95, loc=sample.mean(), scale=std / np.sqrt(n))

結(jié)果同樣為(8811.128437206558, 8838.84651369355)

實例

經(jīng)過大量長期統(tǒng)計，A公司商品日均生產(chǎn)量為50，標(biāo)準(zhǔn)差為15。公司領(lǐng)導(dǎo)對最近半個月（15天）的生產(chǎn)量進(jìn)行突擊檢查，請問：

• 如果最近半個月日均產(chǎn)量為45．是否可能存在消極變數(shù)?例如，機器老化．員工怠工等。
• 如果最近半個月日均產(chǎn)量為59，是否可能存在積極改進(jìn)?例如，提高生產(chǎn)效率，加大人力投入等。

mean, std = 50, 15
n = 15
sigma = std / np.sqrt(n)
start, end = stats.norm.interval(0.95, loc=mean, scale=sigma)
print(f"樣本均值有95%的概率在{start:.1f}-{end:.1f}范圍內(nèi)")

樣本均值有95%的概率在42.4-57.591范圍內(nèi)

因此日均產(chǎn)量為45在合理范圍內(nèi)，日均產(chǎn)量為59很可能存在積極改進(jìn)。

總體標(biāo)準(zhǔn)差未知時的區(qū)間估計

上述區(qū)間估計的方法，適用于已知總體標(biāo)準(zhǔn)差且樣本量大于30的情況下。但正常情況下，我們對總體的標(biāo)準(zhǔn)差都是未知的，有時抽樣樣本量也會少于30。

不過，進(jìn)行多次抽樣，樣本的均值服從t分布，當(dāng)每次抽樣樣本量n增大時（達(dá)到30以上），t分布將逐漸接近于正態(tài)分布。

使用 t 分布的好處是不需要知道總體的標(biāo)準(zhǔn)差，可以直接根據(jù)樣本的標(biāo)準(zhǔn)差進(jìn)行計算。

于是我們可以根據(jù)鳶尾花花瓣長度的樣本數(shù)據(jù)，去估算在95%置信度下全球鳶尾花花瓣長度的均值。

首先導(dǎo)包并加載數(shù)據(jù)：

import pandas as pd
from sklearn.datasets import load_iris
from scipy import stats

iris = load_iris()
data = pd.DataFrame(
    iris.data,
    columns=["sepal_length", "sepal_width", "petal_length", "petal_width"])

計算：

confidenceLevel = 0.95
mean = data.petal_length.mean()
start, end = stats.t.interval(alpha=confidenceLevel, df=len(data) - 1, loc=mean, scale=stats.sem(data.petal_length))
print(f"樣本均值是{mean:.2f}，在{confidenceLevel:.0%}置信度下，全球鳶尾花花瓣長度的均值在{start:.2f}-{end:.2f}cm的范圍")

結(jié)果：樣本均值是3.76，在95%置信度下，全球鳶尾花花瓣長度的均值在3.47-4.04cm的范圍

假設(shè)檢驗

示例

某車間用一臺包裝機包裝葡萄糖，袋裝糖的凈重是—個隨機變量．它服從正態(tài)分布。當(dāng)機器正常時，其均值為0.5kg，標(biāo)準(zhǔn)差為0.015kg。某日開工后為檢驗包裝機是否正常，隨機地抽取它所包裝的糖9袋，稱得凈重為(kg)

0.497, 0.506, 0.518, 0.524, 0.498, 0.511, 0.520, 0.515, 0.512

問：機器是否正常？

按照前面區(qū)間估計的做法，我們可以根據(jù)總體的特征計算出樣本均值的置信區(qū)間，然后判斷當(dāng)前樣本的均值是否存在于這個范圍之內(nèi)：

import numpy as np
from scipy import stats

# 總體的均值和標(biāo)準(zhǔn)差
mean, std = 0.5, 0.015
# 樣本
a = np.array([0.497, 0.506, 0.518, 0.524, 0.498, 0.511, 0.520, 0.515, 0.512])
left, right = stats.norm.interval(0.95, loc=mean, scale=std / np.sqrt(len(a)))
print(f"95%置信度下，樣本均值的置信區(qū)間：({left:0.3f},{right:0.3f})")
print("樣本均值：", a.mean())

結(jié)果：

95%置信度下，樣本均值的置信區(qū)間：(0.490,0.510)
樣本均值： 0.5112222222222224

顯然樣本均值不在根據(jù)總體估算的置信區(qū)間之內(nèi)，所以在95%置信度下，機器是不正常的。

另一種方案就是根據(jù)樣本去估計總體，由于樣本量為9小于30，所以我們使用t分布來估計總體的均值：

import numpy as np
from scipy import stats

# 總體的均值和標(biāo)準(zhǔn)差
mean, std = 0.5, 0.015
# 樣本
a = np.array([0.497, 0.506, 0.518, 0.524, 0.498, 0.511, 0.520, 0.515, 0.512])
left, right = stats.t.interval(0.95, df=len(a) - 1, loc=a.mean(), scale=stats.sem(a))
print(f"95%置信度下，總體均值的置信區(qū)間：({left:0.3f},{right:0.3f})")
print("總體實際均值：", mean)

結(jié)果：

95%置信度下，總體均值的置信區(qū)間：(0.504,0.518)
總體實際均值： 0.5

顯然95%置信度下總體實際均值也不在總體均值的置信區(qū)間之內(nèi)，所以我們可以認(rèn)為機器是不正常的。

假設(shè)檢驗的概念

在假設(shè)檢驗的思路下，則是先給出一個原假設(shè)，然后根據(jù)反證法看是否出現(xiàn)小概率事件，如果沒有出現(xiàn)小概率事件的情況下則接受原假設(shè)，若出現(xiàn)小概率事件情況下則拒絕原假設(shè)接受備擇假設(shè)。

對于上面這個問題的思路就是，先作原假設(shè)：機器是正常的，備擇假設(shè)：機器是不正常的。結(jié)果樣本均值不在置信區(qū)間內(nèi)，出現(xiàn)了小概率事件，所以拒絕原假設(shè)接受備擇假設(shè)：機器是不正常的。

假設(shè)檢驗，也稱為顯著性檢驗，是通過樣本的統(tǒng)計量，來判斷與總體參數(shù)之間的差異是否顯著。原假設(shè)也稱為零假設(shè)或 H 0 H_0 H0?，備則假設(shè)也稱為對立假設(shè)或 H 1 H_1 H1?。

根據(jù)樣本信息進(jìn)行分析判斷，是選擇接受（維持）原假設(shè)還是拒絕原假設(shè)（接受備擇假設(shè)）。

在假設(shè)檢驗中，我們認(rèn)為，小概率事件在—次試驗中是不會發(fā)生的。—旦小慨率事件發(fā)生，則我們就有理由拒絕原假設(shè)。

假設(shè)檢驗遵循"疑罪從無"的原則，接受原假設(shè)，并不代表原假設(shè)—定是正確的，只是沒有充分的證據(jù)，去證明原假設(shè)是錯誤的，因此，只能維持原假設(shè)。

P-Value與顯著性水平

P-Value是一個概率值，表示支持原假設(shè)的概率。在原假設(shè)為等值假設(shè)時，P-Value也表示樣本統(tǒng)計量與總體參數(shù)無差異的概率。

顯著性水平則是我們預(yù)先設(shè)定的—個閾值，使用 α \alpha α表示，通常 α \alpha α的取值為0.05（1-a為置信度）。當(dāng)P-Value的值大于 α \alpha α時，支持原假設(shè)，否則，拒絕原假設(shè)。

假設(shè)檢驗的步驟如下:

1. 設(shè)置原假設(shè)與備擇假設(shè)。
2. 設(shè)置顯著性水平 α \alpha α（通常選擇 α \alpha α= 0.05）。
3. 根據(jù)問題選擇假設(shè)檢驗的方式。
4. 計算統(tǒng)計量，并通過統(tǒng)計量獲取 P P P值。
5. 根據(jù) P P P值與 α \alpha α值，決定接受原假設(shè)還是備擇假設(shè)。

Z檢驗與t檢驗

Z檢驗用來判斷樣本均值是否與總體均值具有顯著性差異。Z檢驗要求總體呈正態(tài)分布、總體方差已知，否則無法使用Z檢驗，另外要求樣本容量最好大于30，否則建議使用t檢驗。

Z統(tǒng)計量的計算方式如下：

Z = x ˉ ? μ 0 S x ˉ = x ˉ ? μ 0 σ / n \Huge{Z=\frac{\bar{x}-\mu_{0}}{S_{\bar{x}}}=\frac{\bar{x}-\mu_{0}}{\sigma / \sqrt{n}}} Z=Sxˉ?xˉ?μ0??=σ/n ?xˉ?μ0??

• x ˉ \bar{x} xˉ：樣本均值。
• μ 0 \mu_{0} μ0?：待檢驗的總體均值（假設(shè)的總體均值）。
• S x ˉ S_{\bar{x}} Sxˉ?：樣本均值分布的標(biāo)準(zhǔn)差（標(biāo)準(zhǔn)誤差）。
• σ \sigma σ：總體的標(biāo)準(zhǔn)差。
• n n n：樣本容量。

通過假設(shè)檢驗來求解前面的示例：

1. 設(shè)置原假設(shè)與備擇假設(shè)。
   • 原假設(shè)： μ = μ 0 = 0.5 k g \mu=\mu_0 = 0.5kg μ=μ0?=0.5kg（機器運作正常）
   • 備擇假設(shè)： μ ≠ μ 0 ≠ 0.5 k g \mu \neq \mu_0 \neq 0.5 kg μ?=μ0??=0.5kg（機器運作不正常）
2. 設(shè)置顯著性水平。
   • 設(shè)置 α \alpha α= 0.05
3. 根據(jù)問題選擇假設(shè)檢驗的方式。
   • 袋裝糖的凈重呈正態(tài)分布，總體標(biāo)準(zhǔn)差已知，選擇 Z Z Z檢驗。
4. 計算統(tǒng)計量，并通過統(tǒng)計量獲取 P P P值。
   • 根據(jù)選擇的假設(shè)檢驗，進(jìn)行計算。

import numpy as np
from scipy import stats

# 總體的均值和標(biāo)準(zhǔn)差
mean, std = 0.5, 0.015
# 樣本
a = np.array([0.497, 0.506, 0.518, 0.524, 0.498, 0.511, 0.520, 0.515, 0.512])
# 樣本均值
sample_mean = a.mean()
sigma = std / np.sqrt(len(a))

Z = (sample_mean - mean) / sigma
print("統(tǒng)計量Z:", Z)
P = 2 * stats.norm.sf(abs(Z))
print("P-Value值：", P)

結(jié)果：

統(tǒng)計量Z: 2.244444444444471
P-Value值： 0.02480381963225589

5. 根據(jù)P值與α值，決定接受原假設(shè)還是備擇假設(shè)。

   • 由結(jié)果可知， P P P值即支持原假設(shè)的概率小于顯著性水平α的值0.05，故我們拒絕原假設(shè)。接受備擇假設(shè)，即我們認(rèn)為機器運作是不正常的。

這就是假設(shè)檢驗的步驟。

t檢驗則適用于總體呈正態(tài)分布且方差未知的情況。

隨著樣本容量的增大（樣本量達(dá)到30以上時），t分布逐漸接近于正態(tài)分布。此時，t檢驗也就近似應(yīng)用于Z檢驗。

t統(tǒng)計量的計算方式如下：

t = x ˉ ? μ 0 S x ˉ = x ˉ ? μ 0 S / n \Huge{t=\frac{\bar{x}-\mu_{0}}{S_{\bar{x}}}=\frac{\bar{x}-\mu_{0}}{S / \sqrt{n}}} t=Sxˉ?xˉ?μ0??=S/n ?xˉ?μ0??

• x ˉ \bar{x} xˉ：樣本均值。
• μ 0 \mu_{0} μ0?：待檢驗的總體均值（假設(shè)的總體均值）。
• S x ˉ S_{\bar{x}} Sxˉ?：樣本均值的標(biāo)準(zhǔn)差（標(biāo)準(zhǔn)誤差）。
• S S S：樣本的標(biāo)準(zhǔn)差。
• n n n：樣本容量。

鳶尾花花瓣長度的均值為3.5cm，根據(jù)假設(shè)檢驗的步驟判斷是否正確：

1. 設(shè)置原假設(shè)與備擇假設(shè)。
   • 原假設(shè)：總體均值 μ = μ 0 = 3.5 c m \mu=\mu_0 = 3.5cm μ=μ0?=3.5cm（該說法正確）
   • 備擇假設(shè)：總體均值 μ = μ 0 ≠ 3.5 c m \mu = \mu_0 \neq 3.5cm μ=μ0??=3.5cm（該說法不正確）
2. 設(shè)置顯著性水平。
   • 設(shè)置 α \alpha α= 0.05
3. 根據(jù)問題選擇假設(shè)檢驗的方式。
   • 鳶尾花呈正態(tài)分布，總體標(biāo)準(zhǔn)差未知，選擇 t t t檢驗。
4. 計算統(tǒng)計量，并通過統(tǒng)計量獲取 P P P值。

import pandas as pd
from sklearn.datasets import load_iris
import numpy as np

iris = load_iris()
data = pd.DataFrame(
    iris.data,
    columns=["sepal_length", "sepal_width", "petal_length", "petal_width"])
mean = data["petal_length"].mean()
std = data["petal_length"].std()
print("樣本均值：", mean, "樣本標(biāo)準(zhǔn)差:", std)

t = (mean - 3.5) / (std / np.sqrt(len(data)))
print("t統(tǒng)計量:", t)

# 計算p值
# df:自由度，即變量可以自由取值的個數(shù)
P = 2 * stats.t.sf(abs(t), df=len(data) - 1)
print("P-Value值：", P)

# 還可以通過scipy提供的相關(guān)方法來進(jìn)行t檢驗的計算，無需自行計算。
t, p_twoTail = stats.ttest_1samp(data['petal_length'], 3.5)
print(f"t統(tǒng)計量:{t}, P-Value值：{p_twoTail}")

樣本均值： 3.7586666666666693 樣本標(biāo)準(zhǔn)差: 1.7644204199522617
t統(tǒng)計量: 1.7954942587239626
P-Value值： 0.07460161706985045

還可以直接通過scipy提供的相關(guān)方法來進(jìn)行t檢驗的計算：

t, p_twoTail = stats.ttest_1samp(data['petal_length'], 3.5)
print(f"t統(tǒng)計量:{t}, P-Value值：{p_twoTail}")

結(jié)果：

t統(tǒng)計量:1.79549425872394, P-Value值：0.07460161706985409

P P P值即支持原假設(shè)的概率大于顯著性水平α的值0.05，故我們沒有充足的理由拒絕原假設(shè)，則接受原假設(shè)，即我們認(rèn)為鳶尾花花瓣長度的均值為3.5cm。

雙邊檢驗與單邊檢驗

• 原假設(shè)： μ = μ 0 = 3.5 c m \mu=\mu_0 = 3.5cm μ=μ0?=3.5cm
• 備擇假設(shè)： μ ≠ μ 0 ≠ 3.5 c m \mu \neq \mu_0 \neq 3.5cm μ?=μ0??=3.5cm

對于上面的等值假設(shè)，檢驗的是總體均值( μ \mu μ)與假設(shè)均值( μ 0 \mu_0 μ0?)是否相等，當(dāng) μ ≠ μ 0 \mu≠\mu_0 μ?=μ0?時．總體均值可以大于假設(shè)均值，也可以小于假設(shè)均值，像這樣的檢驗稱為雙邊假設(shè)檢驗（雙邊檢驗）。非等值假設(shè)時，總體參數(shù)大于或者小于假設(shè)參數(shù)值．像這樣的檢驗稱為單邊假設(shè)檢驗(單邊檢驗)。例如，我們僅關(guān)注技術(shù)改進(jìn)后是否比以前有所提高，而不是關(guān)注是否與以前不同。

以均值為例，設(shè)總體均值為 μ \mu μ，假設(shè)均值為 μ 0 \mu_0 μ0?

如果設(shè)立：

• 原假設(shè): μ ≤ μ 0 \mu \leq \mu_{0} μ≤μ0?
• 備擇假設(shè): μ > μ 0 \mu > \mu_{0} μ>μ0?
  則稱這樣的假設(shè)為右邊假設(shè)檢驗(右邊檢驗)。

如果設(shè)立：

• 原假設(shè): μ ≥ μ 0 \mu \geq \mu_{0} μ≥μ0?
• 備擇假設(shè): μ < μ 0 \mu < \mu_{0} μ<μ0?
  則稱這樣的假設(shè)為左邊假設(shè)檢驗(左邊檢驗)。

說明：在單邊檢驗中，原假設(shè)為維持現(xiàn)狀，備則假設(shè)為改變現(xiàn)狀。

我們可以根據(jù)計算P值的方式來記憶是左邊假設(shè)檢驗還是右邊假設(shè)檢驗，

例如 μ ≤ μ 0 \mu \leq \mu_{0} μ≤μ0?等價于 μ ? μ 0 ≤ 0 \mu - \mu_{0} \leq 0 μ?μ0?≤0，那么應(yīng)該使統(tǒng)計量越小P值越大，就應(yīng)該計算統(tǒng)計量右邊圍成的面積。

而 μ ≥ μ 0 \mu \geq \mu_{0} μ≥μ0?等價于 μ ? μ 0 ≥ 0 \mu - \mu_{0} \geq 0 μ?μ0?≥0則應(yīng)該計算左邊圍成的面積。

所以對于單邊檢驗，左邊和右邊是指應(yīng)該計算統(tǒng)計量左邊還是右邊圍成的面積。

而對于等值假設(shè)，即雙邊檢驗，應(yīng)該使統(tǒng)計量越接近0，P值越大，應(yīng)該計算統(tǒng)計量對應(yīng)左右兩側(cè)位置與兩邊圍成的面積。

例如，判斷鳶尾花的平均花瓣長度等于3.5cm的說法是否正確，計算出統(tǒng)計量為1.8后，再計算左右兩邊圍成的面積即為P值：

右邊假設(shè)檢驗

判斷鳶尾花的平均花瓣長度不超過3.5cm的說法是否正確

1. 設(shè)置原假設(shè)與備擇假設(shè)。
   • 原假設(shè)： μ ≤ μ 0 \mu \leq \mu_0 μ≤μ0?
   • 備擇假設(shè)： μ > μ 0 \mu > \mu_0 μ>μ0?
2. 設(shè)置顯著性水平。
   • 設(shè)置 α \alpha α= 0.05
3. 根據(jù)問題選擇假設(shè)檢驗的方式。
   • 鳶尾花呈正態(tài)分布，總體標(biāo)準(zhǔn)差未知，選擇 t t t檢驗。
4. 計算統(tǒng)計量，并通過統(tǒng)計量獲取 P P P值。

計算下方紅色部分（右邊）的面積：

print("t統(tǒng)計量:", t)
P = stats.t.sf(t, df=len(data) - 1)
print("P-Value值：", P)

t統(tǒng)計量: 1.79549425872394
P-Value值： 0.037300808534927045

5. 根據(jù)P值與α值，決定接受原假設(shè)還是備擇假設(shè)。

   • P < α P<\alpha P<α，因此拒絕原假設(shè)，我們認(rèn)為鳶尾花的平均花瓣長度超過3.5cm。

左邊假設(shè)檢驗

假如判斷鳶尾花的平均花瓣長度不小于3.5cm的說法是否正確

1. 設(shè)置原假設(shè)與備擇假設(shè)。
   • 原假設(shè)： μ ≥ μ 0 \mu \geq \mu_0 μ≥μ0?
   • 備擇假設(shè)： μ < μ 0 \mu < \mu_0 μ<μ0?
2. 設(shè)置顯著性水平。
   • 設(shè)置 α \alpha α= 0.05
3. 根據(jù)問題選擇假設(shè)檢驗的方式。
   • 鳶尾花呈正態(tài)分布，總體標(biāo)準(zhǔn)差未知，選擇 t t t檢驗。
4. 計算統(tǒng)計量，并通過統(tǒng)計量獲取 P P P值。

計算下方紅色部分（左邊）的面積：

print(t)
P = stats.t.cdf(t, df=len(data) - 1)
print("P-Value值:", P)

1.79549425872394
P-Value值: 0.962699191465073

5. 根據(jù)P值與𝛼值，決定接受原假設(shè)還是備擇假設(shè)。

   • P > α P>\alpha P>α，因此維持原假設(shè)，我們認(rèn)為鳶尾花的平均花瓣確實長度不小于3.5cm

示例

某公司要求，平均日投訴量均值不得超過1%?，F(xiàn)檢查—個部門的服務(wù)情況。在該部門維護(hù)的—個500人客戶群中，近7天的投訴量分別為5, 6, 8, 4, 4, 7, 0。請問該部門是否達(dá)標(biāo)？

原假設(shè)平均日投訴量均值小于等于1%，是個右邊假設(shè)檢驗，總體標(biāo)準(zhǔn)差未知，樣本量小于30，選擇 t t t檢驗：

data = np.array([5, 6, 8, 4, 4, 7, 0]) / 500 * 100
print(data)
# 假設(shè) 平均日投訴量 <=1%
mean = data.mean()
std = data.std()
print("樣本均值：", mean, "樣本標(biāo)準(zhǔn)差:", std)
t = (mean - 1) / (std / np.sqrt(len(data)))
print("t統(tǒng)計量:", t)
P = stats.t.sf(t, df=len(data) - 1)
print("P-Value值：", P)

[1.  1.2 1.6 0.8 0.8 1.4 0. ]
樣本均值： 0.9714285714285715 樣本標(biāo)準(zhǔn)差: 0.4831867007225075
t統(tǒng)計量: -0.1564465546936854
P-Value值： 0.5595938210714403

支持原假設(shè)的概率超過顯著性水平，故接受原假設(shè)，即平均日投訴量均值沒有超過1%，該部門達(dá)標(biāo)。