统计学中的参数估计思想与参数估计方法python实现（以2010-2014上证指数收益率的均值区间估计为例）

参数估计是统计学中“推断统计” 的一部分，推断统计是统计学中最重要的部分之一，它源于统计学中一个最基本的目的，也是这个概率论与数理统计这一学科的数理统计部分存在的唯一目的：用样本推断总体情况。用样本推断总体情况。用样本推断总体情况。参数估计分为点估计和区间估计，高考后学生利用对答案，凭感觉，靠记忆等多种已知信息下，可以告诉父母自己大致能考个749分，而一些逻辑更严谨的学生会给出一个上下限的范围，比如语文作文的文章读起来就是满分作文，但写的时候字体不算出彩，那么就会保守的估计，说判卷老师眼神好又心情好，总分能拿个749，但要遇上了不对付的人，就只好得个745分了。给出的范围就是区间估计，确定的分数就是点估计。

点估计有它的局限性，点估计又叫做定值估计，以样本指标的数值作为作为总体指标的统计量。实际操作中，又往往以样本指标的实际值直接作为总体估计未知参数的估计值。像是要估计厂家生产钨丝灯泡的平均使用寿命，市场上控制除寿命外其他变量后的灯泡中，使用寿命高的更得到市场消费者的青睐。从经验上看，老师傅会猜测，注意我的用词，猜测。大概这批灯泡在3500H左右，或者会说，“充其量能用到个4000H，这灯泡就了不得了。”这批灯泡被买走后，等到最后用坏才能得知每个灯泡真实的使用寿命。老师傅有自己的生产经验，根据多年的客户反馈和自己对产品的把关，灯泡基本都在3500H后寿终正寝了，可偏受各种后期因素的影响，有的灯泡因为这些因素（使用环境，天气情况，生产材料）寿命出现了变化。显然3500H的估计值只能覆盖一部分的灯泡，而[3500,4000]在结果上往往会代表了更多灯泡的实际情况。但也有人较真，说我觉得灯泡的寿命一定在[0,∞)之间，这么说确实是对的，但没什么参考价值，因为这个结果人尽皆知。上面的老师傅提到的“这灯泡就了不得了”，其实是对灯泡的一种认可，也可以理解为，老师傅认为寿命达到4000H及以上的灯泡战胜了绝大多数的同类灯泡，这样的灯泡少见，假如生产100个灯泡的话，这样的灯能出来1个就烧高香了，那么换言之，老师傅因为对自己的材料把关自信，技术自信，以上的各种先决条件让他认为灯泡保底寿命3500H，里面最优质的灯泡可以用够4000H，那么老师傅就有99%的把握认为100个灯99个都在他给的范围内，剩下的那个在他的认知里已经超纲了。现在大家更普遍接受的是区间估计的方法，可以给出一个范围，还能给出一个估计范围的把握（置信水平）。

统计学靠的不是猜。而是科学的论证。

接下来的全部内容，也是当前各领域接受最多的估计方法，区间估计的全部实现过程。

# 上证指数的参数估计
import numpy as np
import pandas as pd
from scipy import stats
from matplotlib import pyplot as plt
from jqdatasdk import *
auth('ID','password')  # 出于保密目的，这里不填写自己的账号密码了

数据获取的方法，有很多的渠道，我这里使用聚宽数据接入所需数据：
预备了上证指数SH，深圳指数SZ只是写着玩的，有需要又有兴趣测试的话，可以用下面的这个SZ数据。

# 聚宽中获取上证指数的指数信息，一天一更新
SH_JQ = get_bars("000001.XSHG", 1213, unit='1d',
         fields=['date','open','high','low','close'],
         include_now=False, end_dt="2015-01-01", fq_ref_date=None,df=True)
SZ_JQ = get_bars("399001.XSHE", 1213, unit='1d',
         fields=['date','open','high','low','close'],
         include_now=False, end_dt="2015-01-01", fq_ref_date=None,df=True)

# 根据收盘价信息计算每日的收益率
SH_JQ['rateOfReturn'] = (SH_JQ.close-SH_JQ.close.shift())/SH_JQ.close.shift()
SZ_JQ['rateOfReturn'] = (SZ_JQ.close-SH_JQ.close.shift())/SZ_JQ.close.shift()

目标是对上证指数的收益率的参数估计，首先需要了解到数据的分布情况

# TRD_JQ.rateOfReturn.hist()
fig,ax = plt.subplots()
# 背景部分
    # 背景网格线
plt.grid()
# 主图部分
    # 频率分布直方图
n, bins, patches = ax.hist(SH_JQ.rateOfReturn[1:],density=True)
    # 区间内的正态分布曲线
plt.plot(np.arange(-0.06,0.06,0.002),stats.norm.pdf(np.arange(-0.06,0.06,0.002),SH_JQ.rateOfReturn[1:].mean(),SH_JQ.rateOfReturn[1:].std()),'--')
# 其他元素
ax.set_ylabel('Probability density')
ax.set_title(r'Histogram of TRD_JQ.rateOfReturn: $\mu={}$, $\sigma={}$'.format('%.2f' % SH_JQ.rateOfReturn[1:].mean(),'%.2f' % SH_JQ.rateOfReturn[1:].std()))
plt.show()

输出结果如下所示：

从可视化结果来看，近似服从正态分布吧。

非参数检验方法下，还可以使用K-S检验验证数据是否服从某类已知理论分布，或者极少数用于比较两个数据集分布的情况。

因为KS检验的局限性，经常使用的拟合优度检验和Kolmogorov-Smirnov检验的检验功效较低，在许多计算机软件的Kolmogorov-Smirnov检验无论是大小样本都用大样本近似的公式，很不精准，一般使用Shapiro-Wilk检验和Lilliefor检验。
Kolmogorov-Smirnov检验只能检验是否一个样本来自于一个已知样本，而Lilliefor检验可以检验是否来自未知总体。

stats.kstest(SH_JQ.rateOfReturn[1:],"norm",args=(SH_JQ.rateOfReturn[1:].mean(),SH_JQ.rateOfReturn[1:].std()))

KstestResult(statistic=0.054550249175892906, pvalue=0.0014160668616530904)

结果返回的stats表达了KS检验统计量，pvalue表达p值。KS检验统计量的值域为[0,1],值越大表达其分布间存在显著差异。实际上，尤其在机器学习的模型中，利用KS检验的思想引申至混淆矩阵下计算出的TPR和FPR差值的最大值，进而用于评价模型区分效果的好坏。举例说明的话，就是判别设定不同阈值时对应不同分类器下形成的TPR曲线与FPR曲线在y轴方向上的最大差值。如图所示：

每条红线代表不同的阈值，而每个阈值下设定的左右区间又可以看作是每个分类器产生的不同样本分类结果。假定1代表正例而0代表反例，{0，1}的二分类样本真实标签依概率从小到大从左至右排序。

• M1为第一条红线设定阈值下的分类器,此分类器预测结果为{0,0,1,1,1,1,1,1,1,1},真实结果为{0,0,0,0,0,1,1,1,1,1},预测正确的有样本中的观测1,2,6,7,8,9,10,共7个。其中真实为正例"1"的结果全部判断正确，反例0的五个结果中，判断正确2个，错误三个。混淆矩阵描述为TP = 5 , TN = 0 , FP = 3 , FN = 2。综上，计算得到TPR = 1 , TNR = 0 , FPR = 0.6 , FNR = 0.4
• 同理M2得：TPR = 1 , TNR = 0 , FPR = 0 , FNR = 1
• 同理M3得：TPR = 0.6 , TNR = 0.4 , FPR = 0 , FNR = 1
• 同理M4得：TPR = 0.2 , TNR = 0.8 , FPR = 0 , FNR = 1

TPR = [1,1,0.6,0.2]
FPR = [0.6,0,0,0]
# cluster_name = ['M1','M2','M3','M4']
cluster_name = [1,2,3,4]
# 主图部分
plt.plot(cluster_name,TPR,label = 'TPR')
plt.plot(cluster_name,FPR,label = 'FPR')
for i in range(len(TPR)):
    plt.plot([cluster_name[i],cluster_name[i]],[TPR[i],FPR[i]],color = 'green',lw =1.4,ls ='-.',alpha = 0.6)
    plt.scatter(cluster_name[i],TPR[i],cmap='o',color = 'r')
    plt.scatter(cluster_name[i],FPR[i],cmap='o',color = 'r')
# 其他元素
plt.legend()
plt.show()

四个分类器的对应输出结果如下所示：

可见，完全分类正确的M2分类器的区分效果最好，这也侧面解释了KS的值域[0,1]，且KS统计量越大，分类结果越好，也表示结果差异越显著。
上面的KS检验只是一个小插曲，KS的非参数检验方法，当然有它的局限性，还请按需使用。

# 对收益率的样本均值进行区间估计
    # alpha置信水平/自由度/样本均值/样本标准差
stats.t.interval(0.95,len(SH_JQ.rateOfReturn[1:])-1,SH_JQ.rateOfReturn[1:].mean(),stats.sem(SH_JQ.rateOfReturn[1:]))

(-0.0006112384423833683, 0.0007318937978484808)

即置信度在0.95时上证指数的收益率的置信区间为(-0.0006, 0.0007)