主题：【第十届原创】利用VBA进行能力验证Z值的自动计算

利用VBA进行能力验证Z值的自动计算

在能力验证时，大多数采用Z值来评价检测结果。Z值的计算公式为：Z=(x-X)/σ，式中：x为实验室检测结果，X为指定值，σ为能力评定标准差。一般能力验证采用稳健（Robust）统计技术确定指定值和能力评定标准差，即采用稳健统计的中位值作为指定值，尺度化中位绝对差（MADe）或标准化四分位距（NIQR）作为能力评定标准差。评价结果一般为：│Z│≤2为满意结果，2＜│Z│＜3为有问题结果（可疑值），│Z│≥3为不满意结果（离群值）。
Z值的计算方法，在CNAS-GL02:2014《能力验证结果的统计处理和能力评价指南》和CNAS-GL40《能力验证的选择核查与利用指南》中均有介绍，但相对来说，Z值的计算还是比较繁琐。本文介绍一种利用VBA自动计算Z值的方法，可以在输入实验室代码和检测结果后，达到一键计算Z值的效果。
1前言
该部分对稳健统计方法作一背景介绍。
1.1采用稳健统计方法的理由
由于经典统计方法对离群值敏感，因此通常优先采用对离群值相对不敏感的稳健统计方法。中位值、尺度化中位绝对差（MADe）和标准化四分位距（nIQR）均是简易稳健统计量。算法A通过迭代方法转化原始数据，为近似正态分布提供均值和标准偏差的替代计算方法，这种方法在预期离群值比例低于20%的情况下非常有用。
1.2对总体平均值和标准偏差的简单估计方法
1.2.1中位值
中位值是对称分布总体平均值的一种简单估计，该方法对离群值不敏感。
假设参加者提交的p个数据按递增顺序表示为：x₁,x₂,…x_i,…x_p，则中位值med(x)为：
med(x)=x_{(p+1)/2}              当p为奇数时
med(x)=[x_{p/2}+x_{1+p/2}]/2      当p为偶数时
1.2.2尺度化中位绝对差（MADe）
MADe是正态分布数据的总体标准偏差的估计值，MADe计算方法对较高比例（50%）的离群值不敏感。当p个数据递增排列并计算出med(x)后，计算p个数据中每个数据与中位值的绝对差di(i=1到p)，再计算绝对差的中位值，将得到的中位值乘以1.483即可得到MADe。
d_i=│x_i-med(x)│
MADe(x)=1.483med(d)
1.2.3标准化四分位距（nIQR）
nIQR是一种类似于MADe的稳健统计方法，该方法相对简单并使用广泛。可将参加者结果递增排列，计算第75百分位和第25百分位参加者结果的差值，然后乘以系数0.7413即可得到nIQR。
1.2.4算法A
应用此法计算可得到总体平均值和标准差的稳健值。
p个数据按递增顺序表示为：x₁,x₂,…x_i,…x_p。这些数据的稳健平均值和稳健标准差记为x*和s*。先计算p个数据的中位值作为初始稳健平均值（x*），计算其绝对中位差作为初始稳健标准差（s*）。
x* =medx_i
s*=1.483×med│x_i-x*│
根据以下步骤更新x*和s*的值：
δ=1.5s*
对于每个x_i来说，若
x_i<x*-δ，则x_i*=x*-δ
x_i>x*+δ，则x_i*=x*+δ
否则，x_i保持不变。
然后计算x*和s*的新的取值：



稳健估计值x*和s*可由迭代计算得出，例如用新取值数据更新x*和s*，直至过程收敛。当稳健平均值和稳健标准差的第三位有效数字在连续两次迭代中不再变化时，即可认为过程是收敛的。
2 EXCEL计算过程
以《三种特定过敏原免疫球蛋白E抗体（d1，f1和e3）的浓度》之27个实验室报告的d1数据为例（来自于GB/T 28043-2011《利用实验室间比对进行能力验证的统计方法》表2）。
2.1在EXCEL表适当位置（A9:B35）分别输入实验室代码和d1数据；
2.2将实验室代码及d1数据复制到D9:E35区域，并以列E为主要关键字对D9:E35区域进行升序排序；计算均值和标准差，以及稳健平均值（中位值）和稳健标准差的初始值（绝对中位差）；
2.3将E9:E35数据复制到F9:F35区域，计算δ，x*-δ，x*+δ的值，将超出截止值范围的值，用截止值替换。即小于x*-δ的用x*-δ替换，大于x*+δ的用x*+δ替换。计算新的均值和标准差，以及稳健平均值（与均值相同）和稳健标准差（标准差*1.134），此即为第1次迭代计算；
2.4重复上述迭代过程，将第1次迭代计算的结果复制到下一列，计算δ，x*-δ，x*+δ的值，将超出截止值范围的值，用截止值替换。计算新的均值和标准差，以及稳健平均值（与均值相同）和稳健标准差（标准差*1.134）。当稳健平均值和稳健标准差的第三位有效数字在连续两次迭代中不再变化时，终止迭代。用公式Z=(x-X)/σ计算Z值，其中x为各实验室检测结果，X为稳健平均值，σ为稳健标准差。
2.5得到各实验室检测结果Z值后，利用EXCEL的插入图表功能，绘制各实验室的Z值图，并对图表进行适当修饰，图表中Z值从低到高依次排列，各实验室检测结果的偏离情况一目了然。
3 VBA编程
3.1编程思路：根据上述EXCEL计算过程，进行VBA编程。
3.2优化过程：对迭代计算采用循环控制，第3次迭代后增加判断语句，当稳健平均值和稳健标准差的第三位有效数字在连续两次迭代中不再变化时，终止迭代。
3.3自动制作图表：利用VBA的强大功能，对检测结果Z值的图表制作进行编程，基本达到预期目标，大大节约图表制作的工作量，对图表制作予以规范化和格式化。制作出的图表效果如下：

4 程序代码
以下代码在WinXP Professional（版本号5.1.2600 SP3）+EXCEL2003（版本号11.0）测试通过。为方便大家理解和探讨，在此公布所有程序代码。
4.1总体架构
程序分为三部分设计，分别是“计算Z值”、“自动绘图”和“清除数据”。使用时先点击“清除数据”按钮，可以清除上一次使用时遗留下来的数据；当输入实验室代码和检测结果后，点击“计算Z值”按钮即可计算出各实验室检测结果的Z值；然后点击“自动绘图”按钮，即可绘制出检测结果Z值的图表，供制作能力验证报告材料等场合使用。
4.2 计算Z值的代码
Private Sub CommandButton1_Click()
Dim i, j, k, R As Integer
  R =Range("B65536").End(xlUp).Row
    If R <20 Then
      MsgBox "数据太少，请核查"
      Exit Sub
    End If
  '将实验室代码及检测结果复制到指定区域并进行递增排序
Range(Cells(9, 1), Cells(R, 2)).Copy Destination:=Range(Cells(9, 4),Cells(R, 5))
Range(Cells(9, 4), Cells(R, 5)).Select
Selection.Sort Key1:=Range("E9"), Order1:=xlAscending,Header:=xlGuess, _
OrderCustom:=1, MatchCase:=False, Orientation:=xlTopToBottom, _
  SortMethod:=xlPinYin,DataOption1:=xlSortNormal
  '计算均值和标准差，以及稳健平均值（中位值）和稳健标准差的初始值（绝对中位差）
  Cells(1, 5) =Application.Average(Range(Cells(9, 5), Cells(R, 5)))
  Cells(2, 5) =Application.StDev(Range(Cells(9, 5), Cells(R, 5)))
  Cells(3, 5) ="0"
  Cells(4, 5) =Application.Median(Range(Cells(9, 5), Cells(R, 5)))
  Cells(8, 14)= "第0次x-x*"
    For i = 9To R
      IfCells(i, 5) <> 0 Then Cells(i, 14) = Abs(Cells(i, 5) - Cells(4, 5))
      Next i
    Cells(5, 5)= 1.483 * Application.Median(Range(Cells(9, 14), Cells(R, 14)))
'进行迭代计算，并设定满足条件后退出迭代计算
  For j = 1 To8
    If j >=3 Then
      If _
    Abs(Cells(4, j + 4) / Application.Power(10,Application.RoundDown(Log(Cells(4, j + 4)) / Log(10#), 0)) _
      -Cells(4, j + 3) / Application.Power(10, Application.RoundDown(Log(Cells(4, j +3)) / Log(10#), 0))) < 0.01 And _
    Abs(Cells(4, j + 3) / Application.Power(10,Application.RoundDown(Log(Cells(4, j + 3)) / Log(10#), 0)) _
      -Cells(4, j + 2) / Application.Power(10, Application.RoundDown(Log(Cells(4, j +2)) / Log(10#), 0))) < 0.01 And _
      Abs(Cells(5,j + 4) / Application.Power(10, Application.RoundDown(Log(Cells(5, j + 4)) /Log(10#), 0)) _
      -Cells(5, j + 3) / Application.Power(10, Application.RoundDown(Log(Cells(5, j +3)) / Log(10#), 0))) < 0.01 And _
    Abs(Cells(5, j + 3) / Application.Power(10,Application.RoundDown(Log(Cells(5, j + 3)) / Log(10#), 0)) _
      -Cells(5, j + 2) / Application.Power(10, Application.RoundDown(Log(Cells(5, j +2)) / Log(10#), 0))) < 0.01 Then
      Exit For
      Else
      End If
    End If
    Range(Cells(9,j + 4), Cells(R, j + 4)).Copy Destination:=Range(Cells(9, j + 5), Cells(R, j +5))
      Cells(6,j + 5) = 1.5 * Cells(5, j + 4)
      Cells(7,j + 5) = Cells(4, j + 4) - Cells(6, j + 5)
      Cells(8,j + 5) = Cells(4, j + 4) + Cells(6, j + 5)
      For i = 9 To R
        IfCells(i, 5) <> 0 And Cells(i, j + 5) < Cells(7, j + 5) Then Cells(i, j+ 5) = Cells(7, j + 5)
        IfCells(i, 5) <> 0 And Cells(i, j + 5) > Cells(8, j + 5) Then Cells(i, j+ 5) = Cells(8, j + 5)
        Next i
      Cells(1,j + 5) = Application.Average(Range(Cells(9, j + 5), Cells(R, j + 5)))
      Cells(2,j + 5) = Application.StDev(Range(Cells(9, j + 5), Cells(R, j + 5)))
      Cells(3,j + 5) = j
      Cells(4,j + 5) = Application.Average(Range(Cells(9, j + 5), Cells(R, j + 5)))
      Cells(5,j + 5) = 1.134 * Application.StDev(Range(Cells(9, j + 5), Cells(R, j + 5)))
      '计算Z值
      Cells(8,16) = "Z值＝(x-x*)/s*"
    Range(Cells(9, 4), Cells(R, 4)).Copy Destination:=Range(Cells(9, 15),Cells(R, 15))
      For i = 9 To R
      Cells(i, 16) = (Cells(i, 5) - Cells(4,Range("IV1").End(xlToLeft).Column)) / Cells(5,Range("IV1").End(xlToLeft).Column)
      Next i
    Next j
    '将Z值复制到检测结果右侧
  Range(Cells(9, 15), Cells(R, 16)).Select
  Selection.Sort Key1:=Range("O9"), Order1:=xlAscending,Header:=xlGuess, OrderCustom:=1, MatchCase:=False, _
  Orientation:=xlTopToBottom, SortMethod:=xlPinYin,DataOption1:=xlSortNormal
  Range(Cells(9, 16), Cells(R, 16)).Copy Destination:=Range(Cells(9, 3),Cells(R, 3))
  Selection.Sort Key1:=Range("P9"), Order1:=xlAscending,Header:=xlGuess, OrderCustom:=1, MatchCase:=False, _
  Orientation:=xlTopToBottom, SortMethod:=xlPinYin,DataOption1:=xlSortNormal
End Sub
4.3 清除数据的代码
Private Sub CommandButton2_Click()
  '清除图表和单元格数据
ActiveSheet.ChartObjects.Delete
  Cells.Clear
  '在部分单元格输入项目名称
  Cells(8, 1) ="实验室代码"
  Cells(8, 2) ="检测结果"
  Cells(8, 3) ="Z值"
  Cells(1, 4) ="平均值"
  Cells(2, 4) ="标准差"
  Cells(3, 4) ="迭代步骤"
  Cells(4, 4) ="新的x*"
  Cells(5, 4) ="新的s*"
  Cells(6, 4) ="δ=1.5s*"
  Cells(7, 4) ="x*-δ"
  Cells(8, 4) ="x*+δ"
  '设置单元格字体、字号、对齐方式
Cells.Font.Name = "宋体"
Cells.Font.Name = "Arial Narrow"
Cells.Font.Size = 12
Cells.HorizontalAlignment = xlCenter
  Cells(9,1).Select
End Sub
4.4自动绘图的代码
Private Sub CommandButton3_Click()
  Dim mychartAs ChartObject
  Dim R AsInteger
  R =Range("B65536").End(xlUp).Row
  With Sheet1
    '在指定区域绘制图表，并设置各类参数
  .ChartObjects.Delete
    Set mychart= .ChartObjects.Add(120, 40, 600, 300)
    Withmychart.Chart
      .SetSourceData Source:=Range(Cells(9, 15),Cells(R, 16)), PlotBy:=xlColumns
    .HasLegend = False
    .Axes(xlValue).TickLabels.Font.Name = "arial narrow"    '设置Y轴字体
    .Axes(xlCategory).MajorTickMark = xlNone    '无X轴主要刻度线
    .Axes(xlCategory).TickLabelPosition = xlTickLabelPositionNone    '不显示刻度线标签
    .SeriesCollection(1).ApplyDataLabels AutoText:=False, LegendKey:= _
      False,ShowSeriesName:=False, ShowCategoryName:=True, ShowValue:=False, _
    ShowPercentage:=False, ShowBubbleSize:=False
    .SeriesCollection(1).DataLabels.AutoScaleFont = True
    .SeriesCollection(1).DataLabels.Orientation = 90
    .SeriesCollection(1).DataLabels.Font.Name = " Arial Narrow "
    .SeriesCollection(1).DataLabels.Font.FontStyle = "Bold"
      .SeriesCollection(1).DataLabels.Font.Size= 12
      .PlotArea.Interior.ColorIndex = 0
      End With
    End With
End Sub
5存在问题
5.1 考虑到不同检测项目的结果数值可能相差很大，因此未对计算过程的数值进行有效数字的修约，最终Z值与理论结果可能存在±0.01左右的偏差。
5.2 图表区域设置得较小，如果检测结果较多导致图表柱形图间距偏小的话，可以在图表绘制后将图表区适当拉长，还可以通过修改代码中的字号适当缩小数据标志的大小，以使图表更加美观。
5.3限于水平，未将数值轴主要网格线等予以进一步美化。同时由于本人VBA功底薄弱，完全是从零基础起步，通过网络查找资料摸索着编程。文中肯定存在着诸多不当之处，希望与众多版友共同探讨提高。