NumPy庫已經提供了一些基本的統計函式，如求期望、方差、中位數、最大值和最小值等。示例程式碼：

import numpy as np

#構建一個1000個隨機變數的陣列
x = np.random.randn(1000)

#對陣列元素的值進行統計
mean = x.mean()
std = x.std()
var = x.var()

print(mean,std,var)

執行結果：

(0.02877273942510088, 0.97623362287515114, 0.95303208643194282)

mean是期望值，std是標準差，var是方差，使用numpy.array物件已有的方法獲得統計指標快速有效，而SciPy庫則提供了更高階的統計工具，它的Stats模組包含了多種概率分佈的隨機變數（隨機變數是指概率論中的概念，不是Python中的變數），其中隨機變數又分為連續和離散兩種。所有的連續隨機變數都是rv_continuous的派生類的物件，而所有的離散隨機變數都是rv_discrete的派生類的物件。

4.5.1 連續概率分佈

SciPy的stats模組提供了大約80種連續隨機變數和10多種離散分佈變數，這些分佈都依賴於numpy.random函式。可以通過如下語句獲得stats模組中所有的連續隨機變數，示例程式碼：

from scipy import stats
[k for k, v in stats.__dict__.items() if isinstance(v, stats.rv_continuous)]

執行結果：

['ksone', 'kstwobign', 'norm', 'alpha', 'anglit', 'arcsine', 'beta', 'betaprime', 'bradford', 'burr', 'burr12', 'fisk', 'cauchy', 'chi', 'chi2', 'cosine', 'dgamma', 'dweibull', 'expon', 'exponnorm', 'exponweib', 'exponpow', 'fatiguelife', 'foldcauchy', 'f', 'foldnorm', 'frechet_r', 'weibull_min', 'frechet_l', 'weibull_max', 'genlogistic', 'genpareto', 'genexpon', 'genextreme', 'gamma', 'erlang', 'gengamma', 'genhalflogistic', 'gompertz', 'gumbel_r', 'gumbel_l', 'halfcauchy', 'halflogistic', 'halfnorm', 'hypsecant', 'gausshyper', 'invgamma', 'invgauss', 'invweibull', 'johnsonsb', 'johnsonsu', 'laplace', 'levy', 'levy_l', 'levy_stable', 'logistic', 'loggamma', 'loglaplace', 'lognorm', 'gilbrat', 'maxwell', 'mielke', 'kappa4', 'kappa3', 'nakagami', 'ncx2', 'ncf', 't', 'nct', 'pareto', 'lomax', 'pearson3', 'powerlaw', 'powerlognorm', 'powernorm', 'rdist', 'rayleigh', 'reciprocal', 'rice', 'recipinvgauss', 'semicircular', 'skewnorm', 'trapz', 'triang', 'truncexpon', 'truncnorm', 'tukeylambda', 'uniform', 'vonmises', 'vonmises_line', 'wald', 'wrapcauchy', 'gennorm', 'halfgennorm']
 

連續隨機變數物件主要使用如下方法，下表所示：

下面以標準正態分佈（函式表示f(x)=(1/√2π)exp(-x^2/2)）為例，簡單介紹隨機變數的用法。示例程式碼：

from scipy import stats
# 設定正態分佈引數，其中loc是期望值引數，scale是標準差引數
X = stats.norm(loc=1.0, scale=2.0)

# 計算隨機變數的期望值和方差
print(X.stats())

執行結果：

(array(1.0), array(4.0))

以上程式碼說明，norm可以像函式一樣呼叫，通過loc和scale引數可以指定隨機變數的偏移和縮放參數。對於正態分佈的隨機變數來說，這兩個引數相當於指定其期望值和標準差，標準差是方差的算術平方根。X的stats()方法，可以計算隨機變數X分佈的特徵值，如期望值和方差。
此外，通過呼叫隨機變數X的rvs()方法，可以得到包含一萬次隨機取樣值的陣列x，然後呼叫NumPy的mean()和var()計算此陣列的均值和方差，其結果符合隨機變數X的特性，示例程式碼：

#對隨機變數取10000個值
x = X.rvs(size=10000)
print(np.mean(x), np.var(x))

執行結果

(1.0287787687588861, 3.9944276709242805)

使用fit()方法對隨機取樣序列x進行擬合，它返回的是與隨機取樣值最吻合的隨機變數引數，示例程式碼：

#輸出隨機序列的期望值和標準差
print(stats.norm.fit(x))

執行結果：

(1.0287787687588861, 1.998606432223283)

在下面的例子中，計算取樣值x的直方圖統計以及累計分佈，並與隨機變數的概率密度函式和累積分佈函式進行比較。示例程式碼：

pdf, t = np.histogram(x, bins=100, normed=True)
t = (t[:-1]+t[1:])*0.5
cdf = np.cumsum(pdf) * (t[1] - t[0])
p_error = pdf - X.pdf(t)
c_error = cdf - X.cdf(t)
print("max pdf error: {}, max cdf error: {}".format(np.abs(p_error).max(), np.abs(c_error).max()))

執行結果：

max pdf error: 0.0208405611169, max cdf error: 0.0126874590568

通過繪圖的方式檢視概率密度函式求得的理論值（theory value）和直方圖統計值（statistic value）的結果，可以看出二者是一致的，示例程式碼：

import pylab as pl
pl.plot(t, pdf, color="green", label = "statistic value")
pl.plot(t, X.pdf(t), color="yellow", label ="theory value")
pl.legend(loc = "best")
pl.show()

繪圖如下：

也可以用同樣的方式顯示隨機變數X的累積分佈和陣列pdf的累加結果，示例程式碼：

import pylab as pl
pl.plot(t, cdf, color="green", label = "statistic value")
pl.plot(t, X.cdf(t), color="yellow", label ="theory value")
pl.legend(loc = "best")
pl.show()

繪圖如下：

4.5.2 離散概率分佈

投擲有六個面的骰子時，只能獲得1到6的整數，因此所得到的概率分佈是離散的。我們以值域離散的分佈稱為離散概率分佈，包括泊松分佈、二項分佈、幾何分佈等。通常使用概率質量函式（PMF）描述其分佈情況，如幾何分佈函式PMF=(1-p)^(k-1)p。
在stats模組中所有描述離散分佈的隨機變數都從rv_discrete類繼承，也可以直接用rv_discrete類自定義離散概率分佈。假設有一個不均勻的骰子，其各點出現的概率不相等，我們用如下程式碼定義其分佈，示例程式碼：

# 陣列x儲存骰子的所有可能值，陣列p儲存每個值出現的概率
x = range(1, 7)
p = (0.4, 0.2, 0.1, 0.1, 0.1, 0.1)

# 建立表示這個骰子的隨機變數dice，呼叫其rvs()方法投擲此骰子20次，獲得符合概率p的隨機數
dice = stats.rv_discrete(values=(x, p))
print(dice.rvs(size=20))

執行結果：

array([3, 6, 4, 5, 5, 2, 1, 3, 3, 1, 1, 3, 1, 5, 1, 3, 4, 1, 2, 2])

除了自定義離散概率分佈，我們也可以利用stats模組裡的函式定義各種分佈。下面以生成幾何分佈為例，其函式是geom()，示例程式碼：

import numpy as np
from scipy.stats import geom

# 設定幾何分佈的引數
p = 0.5
dist = geom(p)

# 設定樣本區間  
x = np.linspace(0, 5, 1000)  

# 得到幾何分佈的 PMF 和CDF  
pmf = dist.pmf(x) 
cdf = dist.cdf(x)  

# 生成500個隨機數  
sample = dist.rvs(500)

4.5.3 描述與檢驗函式

SciPy中有超過60種統計函式。stats模組包括了諸如kstest 和normaltest等樣本測試函式，用來檢測樣本是否服從某種分佈。在使用這些工具前，要對資料有較好的理解，否則可能會誤讀它們的結果。樣本分佈檢驗為例，示例程式碼：

import numpy as np 
from scipy import stats 

# 生成包括100個服從正態分佈的隨機數樣本
sample = np.random.randn(100) 

# 用normaltest檢驗原假設
out = stats.normaltest(sample) 
print('normaltest output') 
print('Z-score = ' + str(out[0])) 
print('P-value = ' + str(out[1])) 

# kstest 是檢驗擬合度的Kolmogorov-Smirnov檢驗，這裡針對正態分佈進行檢驗
# D是KS統計量的值，越接近0越好
out = stats.kstest(sample, 'norm') 
print('\nkstest output for the Normal distribution') 
print('D = ' + str(out[0])) 
print('P-value = ' + str(out[1])) 

# 類似地可以針對其他分佈進行檢驗，例如Wald分佈
out = stats.kstest(sample, 'wald') 
print('\nkstest output for the Wald distribution') 
print('D = ' + str(out[0])) 
print('P-value = ' + str(out[1]))

SciPy的stats模組中還提供了一些描述函式，如幾何平均（gmean）、偏度（skew
）、樣本頻數（itemfreq）等。示例程式碼：

import numpy as np 
from scipy import stats 

# 生成包括100個服從正態分佈的隨機數樣本
sample = np.random.randn(100) 

# 調和平均數，樣本值須大於0 
out = stats.hmean(sample[sample > 0]) 
print('Harmonic mean = ' + str(out)) 

# 計算-1到1之間樣本的均值
out = stats.tmean(sample, limits=(-1, 1)) 
print('\nTrimmed mean = ' + str(out)) 

# 計算樣本偏度
out = stats.skew(sample) 
print('\nSkewness = ' + str(out)) 

# 函式describe可以一次給出樣本的多種描述統計結果
out = stats.describe(sample) 
print('\nSize = ' + str(out[0])) 
print('Min = ' + str(out[1][0])) 
print('Max = ' + str(out[1][1])) 
print('Mean = ' + str(out[2])) 
print('Variance = ' + str(out[3])) 
print('Skewness = ' + str(out[4])) 
print('Kurtosis = ' + str(out[5]))