金融量化 - scipy 教程(01)

前篇已經大體介紹了NumPy，接下來讓咱們看看SciPy能作些什麼。NumPy替咱們搞定了向量和矩陣的相關操做，基本上算是一個高級的科學計算器。SciPy基於NumPy提供了更爲豐富和高級的功能擴展，在統計、優化、插值、數值積分、時頻轉換等方面提供了大量的可用函數，基本覆蓋了基礎科學計算相關的問題。

在量化分析中，運用最普遍的是統計和優化的相關技術，本篇重點介紹SciPy中的統計和優化模塊，其餘模塊在隨後系列文章中用到時再作詳述。

本篇會涉及到一些矩陣代數，如若感受不適，可考慮跳過第三部分或者在理解時簡單採用一維的標量代替高維的向量。

首先仍是導入相關的模塊，咱們使用的是SciPy裏面的統計和優化部分：

import numpy as np
import scipy.stats as stats
import scipy.optimize as opt

2、統計部分

2.1 生成隨機數

咱們從生成隨機數開始，這樣方便後面的介紹。生成n個隨機數可用rv_continuous.rvs(size=n)或rv_discrete.rvs(size=n)，其中rv_continuous表示連續型的隨機分佈，如均勻分佈（uniform）、正態分佈（norm）、貝塔分佈（beta）等；rv_discrete表示離散型的隨機分佈，如伯努利分佈（bernoulli）、幾何分佈（geom）、泊松分佈（poisson）等。咱們生成10個[0, 1]區間上的隨機數和10個服從參數a=4，b=2的貝塔分佈隨機數：

rv_unif = stats.uniform.rvs(size=10)
rv_unif
rv_beta = stats.beta.rvs(size=10, a=4, b=2)
rv_beta

在每一個隨機分佈的生成函數裏，都內置了默認的參數，如均勻分佈的上下界默認是0和1。但是一旦須要修改這些參數，每次生成隨機都要敲這麼老長一串有點麻煩，能不能簡單點？SciPy裏頭有一個Freezing的功能，能夠提供簡便版本的命令。SciPy.stats支持定義出某個具體的分佈的對象，咱們能夠作以下的定義，讓beta直接指代具體參數a=4和b=2的貝塔分佈。爲讓結果具備可比性，這裏指定了隨機數的生成種子，由NumPy提供。

np.random.seed(seed=2015)
rv_beta = stats.beta.rvs(size=10, a=4, b=2)

np.random.seed(seed=2015)
beta = stats.beta(a=4, b=2)

2.2 假設檢驗

好了，如今咱們生成一組數據，並查看相關的統計量（相關分佈的參數能夠在這裏查到）：

norm_dist = stats.norm(loc=0.5, scale=2)
n = 200
dat = norm_dist.rvs(size=n)
print("mean of data is: " + str(np.mean(dat)))
print("median of data is: " + str(np.median(dat)))
print("standard deviation of data is: " + str(np.std(d

假設這個數據是咱們獲取到的實際的某些數據，如股票日漲跌幅，咱們對數據進行簡單的分析。最簡單的是檢驗這一組數據是否服從假設的分佈，如正態分佈。這個問題是典型的單樣本假設檢驗問題，最爲常見的解決方案是採用K-S檢驗（ Kolmogorov-Smirnov test）。單樣本K-S檢驗的原假設是給定的數據來自和原假設分佈相同的分佈，在SciPy中提供了kstest函數，參數分別是數據、擬檢驗的分佈名稱和對應的參數：

mu = np.mean(dat)
sigma = np.std(dat)
stat_val, p_val = stats.kstest(dat, 'norm', (mu, sigma))
print('KS-statistic D = %6.3f p-value = %6.4f' % (stat_val

假設檢驗的p-value值很大（在原假設下，p-value是服從[0, 1]區間上的均勻分佈的隨機變量，可參考 http://en.wikipedia.org/wiki/P-value ），所以咱們接受原假設，即該數據經過了正態性的檢驗。在正態性的前提下，咱們可進一步檢驗這組數據的均值是否是0。典型的方法是t檢驗（t-test），其中單樣本的t檢驗函數爲ttest_1samp：

stat_val, p_val = stats.ttest_1samp(dat, 0)
print('One-sample t-statistic D = %6.3f, p-value = %6.

咱們看到p-value<0.05，即給定顯著性水平0.05的前提下，咱們應拒絕原假設：數據的均值爲0。咱們再生成一組數據，嘗試一下雙樣本的t檢驗（ttest_ind）：

norm_dist2 = stats.norm(loc=-0.2, scale=1.2)
dat2 = norm_dist2.rvs(size=100)
stat_val, p_val = stats.ttest_ind(dat, dat2, equal_var=False)
print('Two-sample t-statistic D = %6.3f, p-value = %6.4f'

注意，這裏咱們生成的第二組數據樣本大小、方差和第一組均不相等，在運用t檢驗時須要使用Welch's t-test，即指定ttest_ind中的equal_var=False。咱們一樣獲得了比較小的p-value$，在顯著性水平0.05的前提下拒絕原假設，即認爲兩組數據均值不等。

stats還提供其餘大量的假設檢驗函數，如bartlett和levene用於檢驗方差是否相等；anderson_ksamp用於進行Anderson-Darling的K-樣本檢驗等。

2.3 其餘函數

有時須要知道某數值在一個分佈中的分位，或者給定了一個分佈，求某分位上的數值。這能夠經過cdf和ppf函數完成：

g_dist = stats.gamma(a=2)
print("quantiles of 2, 4 and 5:")
print(g_dist.cdf([2, 4, 5]))
print("Values of 25%, 50% and 90%:")
print(g_dist.pdf([0.25, 0.5, 0.95]))

對於一個給定的分佈，能夠用moment很方便的查看分佈的矩信息，例如咱們查看N(0,1)的六階原點矩：

stats.norm.moment(6, loc=0, scale=1)

describe函數提供對數據集的統計描述分析，包括數據樣本大小，極值，均值，方差，偏度和峯度：

norm_dist = stats.norm(loc=0, scale=1.8)
dat = norm_dist.rvs(size=100)
info = stats.describe(dat)
print("Data size is: " + str(info[0]))
print("Minimum value is: " + str(info[1][0]))
print("Maximum value is: " + str(info[1][1]))
print("Arithmetic mean is: " + str(info[2]))
print("Unbiased variance is: " + str(info[3]))
print("Biased skewness is: " + str(info[4]))
print("Biased kurtosis is: " + str(info[5]))

當咱們知道一組數據服從某些分佈的時候，能夠調用fit函數來獲得對應分佈參數的極大似然估計（MLE, maximum-likelihood estimation）。如下代碼示例了假設數據服從正態分佈，用極大似然估計分佈參數：

norm_dist = stats.norm(loc=0, scale=1.8)
dat = norm_dist.rvs(size=100)
mu, sigma = stats.norm.fit(dat)
print("MLE of data mean:" + str(mu))
print("MLE of data standard deviation:" + str(sigma))

pearsonr和spearmanr能夠計算Pearson和Spearman相關係數，這兩個相關係數度量了兩組數據的相互線性關聯程度：

norm_dist = stats.norm()
dat1 = norm_dist.rvs(size=100)
exp_dist = stats.expon()
dat2 = exp_dist.rvs(size=100)
cor, pval = stats.pearsonr(dat1, dat2)
print("Pearson correlation coefficient: " + str(cor))
cor, pval = stats.pearsonr(dat1, dat2)
print("Spearman's rank correlation coefficient: " + str(cor))

其中的p-value表示原假設（兩組數據不相關）下，相關係數的顯著性。最後，在分析金融數據中使用頻繁的線性迴歸在SciPy中也有提供，咱們來看一個例子：

x = stats.chi2.rvs(3, size=50)
y = 2.5 + 1.2 * x + stats.norm.rvs(size=50, loc=0, scale=1.5)
slope, intercept, r_value, p_value, std_err = stats.linregress(x, y)
print("Slope of fitted model is:" , slope)
print("Intercept of fitted model is:", intercept)
print("R-squared:", r_value**2)

在前面的連接中，能夠查到大部分stat中的函數，本節權做簡單介紹，挖掘更多功能的最好方法仍是直接讀原始的文檔。另外，StatsModels（ http://statsmodels.sourceforge.net ）模塊提供了更爲專業，更多的統計相關函數。若在SciPy沒有知足需求，能夠採用StatsModels。