当我们实验测得一组数据的时候，总有一些数据是有偏差的，我们就需要提出这些有偏差的奇异值。其中一种方法是格拉布斯
原理
在做测量不确定度的评定时，对于测量结果进行数据处理之前，往往要进行异常值的剔除工作。超出在规定条件下预期的误差叫做异常值。产生异常值的原因一般是由于疏忽、失误或突然发生的不该发生的原因造成的，如读错、记错、仪器示值突然跳动、突然震动、操作失误等。所以必须在计算测量结果及不确定度评定中要考虑异常值的判别和剔除。
异常值的判别方法也叫异常值检验法，即：判断异常值的统计检验法。其方法有很多种，例如格拉布斯法、狄克逊法、偏度-峰度法、拉依达法、奈尔法等等。每种方法都有其适用范围和优缺点。每种统计检验法都会犯错误1和错误2。但是有人做过统计，在所有方法中，格拉布斯法犯这两种错误的概率最小，所以本文介绍如何使用格拉布斯法来剔除异常值，其判别步骤如下：
1、假设现在有一组测量数据为：例如测量10次(n＝10)，获得以下数据：8.2、5.4、14.0、7.3、4.7、9.0、6.5、10.1、7.7、6.0。
2、排列数据：将上述测量数据按从小到大的顺序排列，得到4.7、5.4、6.0、6.5、7.3、7.7、8.2、9.0、10.1、14.0。可以肯定，可疑值不是最小值就是最大值。
3、计算平均值x-和标准差s：x-＝7.89；标准差s＝2.704。计算时，必须将所有10个数据全部包含在内。
4、计算偏离值：平均值与最小值之差为7.89－4.7＝3.19；最大值与平均值之差为14.0－7.89＝6.11。
5、确定一个可疑值：比较起来，最大值与平均值之差6.11大于平均值与最小值之差3.19，因此认为最大值14.0是可疑值。
6、计算Gi值：Gi＝(xi－x- )/s；其中i是可疑值的排列序号——10号；因此G10＝( x10－x- )/s＝(14.0－7.89)/2.704＝2.260。由于 x10－x-是残差，而s是标准差，因而可认为G10是残差与标准差的比值。下面要把计算值Gi与格拉布斯表给出的临界值GP(n)比较，如果计算的Gi值大于表中的临界值GP(n)，则能判断该测量数据是异常值，可以剔除。但是要提醒，临界值GP(n)与两个参数有关：检出水平α (与置信概率P有关)和测量次数n (与自由度f有关)。
7、定检出水平α：如果要求严格，检出水平α可以定得小一些，例如定α＝0.01，那么置信概率P＝1－α＝0.99；如果要求不严格，α可以定得大一些，例如定α＝0.10，即P＝0.90；通常定α＝0.05，P＝0.95。
8、查格拉布斯表获得临界值：根据选定的P值(此处为0.95)和测量次数n(此处为10)，查格拉布斯表，横竖相交得临界值G95(10)＝2.176。
9、比较计算值Gi和临界值G95(10)：Gi＝2.260，G95(10)＝2.176，Gi＞G95(10)。
10、判断是否为异常值：因为Gi＞G95(10)，可以判断测量值14.0为异常值，将它从10个测量数据中剔除。
11、余下数据考虑：剩余的9个数据再按以上步骤计算，如果计算的Gi＞G95(9)，仍然是异常值，剔除；如果Gi＜G95(9)，不是异常值，则不剔除。本例余下的9个数据中没有异常值。
实现
def Average_Number(List):
return sum(List) / len(List)
def Variance(List):
temp_List = []
ave = Average_Number(List)
for num in List:
temp_List.append(((num - ave) ** 2) / len(List))
return sum(temp_List)
def Standard_Deviation(List):
return Variance(List) ** (1 / 2)
def Gi(number, lst):
if len(lst) <= 17:
average = Average_Number(lst)
SD = Standard_Deviation(lst)
gi = (abs(number - average)) / SD
return gi
else:
exit(0)
def Grubbs_Table(n):
if n == 3:
return [1.148, 1.153, 1.155, 1.155, 1.155]
elif n == 4:
return [1.425, 1.463, 1.481, 1.492, 1.496]
elif n == 5:
return [1.602, 1.672, 1.715, 1.749, 1.764]
elif n == 6:
return [1.729, 1.822, 1.887, 1.944, 1.973]
elif n == 7:
return [1.828, 1.938, 2.020, 2.097, 2.139]
elif n == 8:
return [1.909, 2.032, 2.126, 2.220, 2.274]
elif n == 9:
return [1.977, 2.110, 2.215, 2.323, 2.387]
elif n == 10:
return [2.036, 2.176, 2.290, 2.410, 2.482]
elif n == 11:
return [2.088, 2.234, 2.355, 2.485, 2.564]
elif n == 12:
return [2.134, 2.285, 2.412, 2.550, 2.636]
elif n == 13:
return [2.175, 2.331, 2.462, 2.607, 2.699]
elif n == 14:
return [2.213, 2.371, 2.507, 2.659, 2.755]
elif n == 15:
return [2.247, 2.409, 2.549, 2.705, 2.806]
elif n == 16:
return [2.279, 2.443, 2.585, 2.747, 2.852]
elif n == 17:
return [2.309, 2.475, 2.620, 2.785, 2.894]
def Outliers_Remover(lst, confidence_level=1):
length = len(lst)
lst = lst
mark = []
for i in range(0, length):
n = Gi(lst[i], lst)
if n > Grubbs_Table(len(lst))[confidence_level - 1]:
mark.append(i)
for j in mark:
lst.pop(j)
return lst
def main():
lst = [1, 2, 3, 4, 10000]
lst_out = Outliers_Remover(lst, confidence_level=5)
print(lst_out)
return 0
if __name__ == "__main__":
main()置信等级为1，则是相对应格拉布斯临界值表的90.00%的置信概率
置信等级为2，则是相对应格拉布斯临界值表的95.00%的置信概率
置信等级为3，则是相对应格拉布斯临界值表的97.50%的置信概率
置信等级为4，则是相对应格拉布斯临界值表的99.00%的置信概率
置信等级为5，则是相对应格拉布斯临界值表的99.50%的置信概率
