返回目录：excel表格制作
第七章 代码实战
第一节 算法
7.1.1 什么是算法
算法是计算机程序的灵魂，是实现某种特定工作的“方法”。假设以日常做菜为例，我们总是按照一定的步骤（例如按照菜谱）将油盐酱醋等材料进行“加工”。其中，这种步骤（菜谱）就是做这道菜的“算法”。可见，尽管油盐酱醋等材料是一样的，但是由于做法（算法）的不同而可以呈现不同的风味。
7.1.2 计算机中对程序的描述
数据+算法=程序。在计算机中，数据和算法是不可分割的两部分，根据数据的结构不同，算法也曾现千差万别。一个算法往往只能适应于一种数据结构，或者几种数据结构，当数据结构发生变化的时候，算法需要发生相应的变化。例如，可以用3个介于0~255之间的数字来表示一种RGB颜色。现在如果数据增加到了4个，用4个数字表示颜色，就需要用到CMYK系统了。在Excel中，定义了某列或者某行作为数据的一个索引，现在如果需要改变索引的位置，往往需要重新构建VBA代码。
第二节 数理统计基础知识
7.2.1 名词解释
数理统计：数理统计是以概率论为基础，研究社会和自然界中大量随机现象数量变化基本规律的一种方法。
事件：就是发生的某一件事情。例如，“火车以300公里的时速前进”，“水烧到100℃沸腾”等都是事件。事件分为：绝对事件和概率事件，绝对事件是指，在某种情况下绝对会发生的事件，或者在某种情况下绝对不可能发生的事情，例如“在一个标准大气压下，纯水在100℃时达到沸腾”就是一个绝对事件。概率事件是指，在一定情况下，某种事件可能发生或者可能不发生，例如，“抛一枚硬币，正面向上”。
概率：概率事件发生的可能性，叫做概率事件的“概率”。例如上述例子中，抛硬币的结果，可能是正面向上，也有可能是背面向上，”正面向上“这件事情发生的概率为1/2。
试验：为了观察事件而进行的动作，例如，掷一个骰子。每掷一下就是一个“试验”。
观察：某一次试验的具体结果。
总体：所研究的对象所有观察结果。
样本：从总体中抽取出来的，用于了解总体情况的一部分观察。
样本空间：即样本中的观察的次数。
统计量：样本中的观察值经过一定的计算得到的，可以反映总体在某些方面具体情况的一个数值。
估计：样本的一个统计量，用于对总体的某方面特征进行评估。
无偏估计：该统计量在样本空间无限大时与总体的特征相互一致。
描述性统计：用于描述样本所在总体的某方面特征的统计过程。
统计推断：用于对不同样品所在总体方面的某种特征进行推断的过程。
描述性统计
描述性统计用于对样品所在总体的某一方面特性进行描述，描述性统计一般牵涉到对总体的分布、平均数、方差、峰度、偏度等方面的描述。
7.2.2 描述性统计
描述性统计是描述总体的某一种特征的统计过程。在描述性统计中，非常重要的一个概念是“分布”。
观察值与其取值时的概率有着一定的联系，这种联系在数理统计中以分布的形式表达出来。下面举一些日常生活中常见的例子说明。
（1）离散变量的概率分布
掷一个骰子，正面向上为1到6的概率均为1/6，那么这个出现的点数和概率之间的关系可以用下表表示：
点数123456
概率1/61/61/61/61/61/6
两个骰子，正面向上的点数之和与概率之间也有一定的关系，这种关系如下：
点数23456789101112
概率1/362/363/364/365/366/365/364/363/362/361/36
以上举的例子，事件A和其发生的概率之间的这种一一对应的关系就叫做分布列（分布率）
以上的事件A可以成为离散型变量事件，即这个观察值可以一个一个地“列举“出来，下面我们再看两个例子，看一看连续型的观察值事件。
（2）连续变量的概率分布
电灯泡的使用寿命t（小时）与其发生的概率之间的关系。
第一，这种关系，正确表述应为”电灯泡至少工作t小时与概率间的关系“。由于电灯泡的工作时间是一个连续的变量，即在一定分范围内，可以取任意的数值。所以电灯炮恰好工作t小时，这样的事件在日常生活中是不会遇到的，在理论上，一个连续变量在某一点的概率也为0.
第二，我们将这种电灯泡至少工作T小时发生的概率表述为下面这样的表达式：
即电灯泡至少工作T小时的概率等于一个函数f(x)在负无穷大到T上的积分。把函数f(x)成为分布函数，简称分布。事实证明，这个例子中，f(x)具有下列表达式的形式：
第二个例子是数理统计的一个最为基础的公式，即平常实验中我们所得到的实验数据，在一般情况下都是服从正态分布的，正态分布的分布率为：
若一个变量服从正太分布，正态分布函数公式中，μ为总体的平均数，σ为总体的方差。可以简记为N(μ, σ )
（2）分布函数的特征
分布函数具有以下特征：
分布函数总是大于或等于零。
分布函数在负无穷大到正无穷大上的积分等于1
（3）平均数
用于描述总体中观察值的“均值”。在数理统计中，“平均数”一词包含我们平时所说的“算数平均数”，但不仅仅指代“算数平均数”。平均数还有如下的形式：
算术平均数：即我们平时所说的平均数，等于一组数的和除以个数。算数平均数用M或者M1表示。
平方平均数：一组数中每一个数的平方和除以个数再开方。平方平均数用M2表示。
调和平均数：个数除以一组数中每一个数的倒数之和，调和平均数用H表示。
众数：一组数中出现次数最多或者频率最高的数，众数一般用Mod表示。
中位数：一组数按顺序排列，位于中间的数。当这组数个数为奇数时，即为中间位置的数，当个数为偶数时，等于位于中间位置的两个数的算数平均数。
需要指出的是，当变量服从正态分布时，其算术平均数=众数=中位数。
（4）方差和标准差
标准差用于描述一组数的离散程度。其定义是，一组数与其平均数差的平方平均数。方差定义是：标准差的平方。
在总体中，总体的标准差等于
通常，用样本的方差来估计总体的方差，样本的方程公式为：
（5）偏度、峰度
偏度用于描述一个分部与正态分布相偏的程度，峰度用于描述一个分部中心聚集的程度。这两个统计量的公式这里不再给出。
（6）求一组数的众数VBA程序
现在给出一组数字，求出这一组数字的众数。根据上面的定义，众数就是出现次数最多的数字。那么这个程序要怎么写？
最简单的想法，先让众数等于第一个数，把数组中全部等于这个数的数字给“数”一遍；然后如果后面的数字的出现个数大于前面的数，就把众数置换为第二个数字。
Function M(Rng As Range)
Dim N As Integer, i As Integer, j As Integer, Tmp As Double, T1 As Integer, T2 As Integer, mt As Double
N = Rng.Columns.Count * Rng.Rows.Count
T2 = 0
mt = Rng(1)
For i = 1 To N
Tmp = Rng(i)
T1 = 0
For j = 1 To N
If Tmp = Rng(j) Then T1 = T1 + 1
Next j
If T1 > T2 Then T2 = T1: mt = Rng(i)
Next i
M = mt
End Function
但是，这个程序有着明显的一个漏洞，即如果后面的数字的个数等于前面的数字，那么众数并不会改为后面的数字，而是始终保持前面的数字。
读者可以想一想如果遇到有相同个数的数字时如何处理？
（7）求一个数组的中位数的VBA程序
求一个数组的中位数，首先需要对这个数组进行排序。下面我们讲解排序的一个经典的算法：冒泡算法。
冒泡算法：根据“轻者向上，重者向下”的原理，完成对数据的排序。现在不妨假设一个简单的序列3、5、1、9、6、7。共6个数字。冒泡算法可以有两种算法，一种是固定比较，一种是移动比较。我们先说一下固定比较的算法：
第一轮，第一个数字与后面所有的数字进行比较，如果第一个数字大于后面的某一个数字，就把第一个数字与后面该位置的数字做交换。于是，第一个数字为3，第3个数字为1,3大于1，将3和1做交换，变成了一个新的序列：
1、5、3、9、6、7
第二轮，用第二个数字与后面的数字进行比较，如果大于则交换位置，于是5和3进行交换，序列变成：1、3、5、9、6、7
……
如此进行五轮交换，最后得到排序的结果1、3、5、6、7、9.
移动比较算法：
第一轮，用第一个数字与第二个数字进行比较，第二个数字与第三个数字进行比较，……若大于则交换，于是变成了：
3、1、5、6、7、9
第二轮，再用第一个数字与第二个数字进行比较，第二个数字与第三个数字进行比较……
……
如此进行5轮比较，得到排序结果。
排序以后，如果数组个数为奇数，则直接取（数组个数+1）/2位置上的数即可。如果数组个数为偶数，则取 （数组个数/2位置上的数+数组个数/2+1位置上的数）/2.
Function Med(Rng As Range)
Dim A() As Double, N As Integer, i As Integer, j As Integer, Temp As Double, R
N = 1
For Each R In Rng
ReDim Preserve A(N)
A(N) = R.Value
N = N + 1
Next R
N = N – 1
For i = 1 To N – 1
For j = i + 1 To N
If A(i) > A(j) Then Temp = A(i): A(i) = A(j): A(j) = Temp
Next j
Next i
If N Mod 2 = 1 Then
Med = A((N + 1) / 2)
Else
Med = (A(N / 2) + A(N / 2 + 1)) / 2
End If
End Function
7.2.3 假设检验
假设检验是用于统计分析中分析某一个命题是否成立的方法。
假设检验中有也有一些专业术语：
原假设：即原命题，即需要统计分析的命题。记作H0。
备择假设：即当原命题被推翻的时候选择的假设，通常为原假设的否定形式。
概率p：即通过统计分析得到了原假设成立的概率。
显著水平Alpha：统计分析并不能完全证实或者证伪原假设，通常情况下，我们按照“小概率事件在一次实验中出现的可能性为0”的假说对原假设进行接受或者否定。例如，天气预报中，如果明天下雨的概率小于0.0001，就会预报说不会下雨。同样的道理，我们设定一个比较小的概率值（通常为0.05或者0.01），当原假设的概率小于这个概率时，就推翻原假设，接受备择假设，若大于这个概率值，则接受原假设。这个人为设定的概率值就是“显著水平”。
显著、不显著：当原假设的概率值小于显著水平时，就说达到“显著”，否则为“不显著”。
第一类错误：也叫“拒真”错误，当统计结果根据p值拒绝了实际上成立的原假设时，为“第一类错误”，第一类错误也可以叫做“假阳性”错误。
第二类错误：也叫“纳伪”错误，即当统计结果根据p值接受了实际上不成立的原假设时，为“第二类错误”，第二类错误也可以叫做“假阴性”错误。
假设检验就是假设某命题成立，然后根据假设计算一个统计量，进而计算出这种假设成立的概率值，如果这个概率值小于某一个具体的数值（一般情况下为5%或者1%）就推翻原假设，接受备择假设。
（1）单样本总体均数的T检验
已知总体的均数，现抽样得到一个样本，判断样本是否来源于总体。这一类问题为单样本总体均数的T检验。下面举例说明。
在实际工作中，已知某机床正常工作时，加工出的零件直径平均数为10.5cm。现用该机床加工了10个零件，测得直径的数值分别为10.2、10.6、10.6、10.4、10.6、10.7、10.7、10.7、10.4、10.7。试问该机床的工作状态是否正常。
在单样本总体均数T检验中，我们需要构建一个统计量T，这个统计量T的公式如下：
根据公式可以编辑VBA代码如下：
Function TSingle(Rng As Range, u As Double) As Double
Dim A() As Double, Mean As Double, S As Double, SS As Double, df As Integer, N As Integer, i As Integer
N = 1
Dim R
For Each R In Rng
ReDim Preserve A(N)
A(N) = R.Value
N = N + 1
Next R
N = N – 1
df = N – 1
For i = 1 To N
S = S + A(i)
SS = SS + A(i) * A(i)
Next i
SS = SS - S * S / N
Mean = S / N
TSingle = (Mean - u) / (Sqr(SS / (df * N)))
End Function
统计过程，先做原假设为，该机器运转正常，即现在加工的零件直径与平均数无显著性差异。
根据上面的叙述和公式，以及VBA程序，可以计算出这个统计量T的数值。
计算结果，T值为1.108，然后再计算这个数值在T分布上的积分，即从无穷大积分到-1.108，以及1.108积分到正无穷。等于1减去从-1.108到1.108上T分布的积分（自由度为10-1=9），算出概率为0.2967,又0.2967大于0.05，所以接受原假设，即该机器运转正常。
（2）方差其次的两样本总体均数T检验
当一个样本扩展到两个样本的时候，统计量变得相当复杂。这时候应当做两种假设：两个样本所在的总体方差相等（其次）以及两样本所在总体的方程不等。
当两样本所在总体方差其次时，两样本总体均数T检验用来比较两样本所在总体均数是否相等，下面举例说明。
在技能大赛中，甲乙两名操作工进行零件加工的比赛，比赛分5轮进行。在规定时间内，甲完成的符合规范要求的件数为3件、4件、3件、5件、4件；乙完成的符合规范要求的零件数为2件、4件、4件、5件、2件。问甲乙两人在操作技术上时候存在一定的优劣差别。
在两样本总体均数T检验时，方差其次情况下，也是构建一个统计量T，T的公式为：
根据公式可以编写VBA代码如下：
Function TDEqual(rng1 As Range, rng2 As Range) As Double
Dim R
Dim S1 As Double, N1 As Integer, M1 As Double, SS1 As Double
Dim S2 As Double, N2 As Integer,M2 as Double,SS2 As Double
N1 = 1
For Each R In rng1
S1 = S1 + R
SS1 = SS1 + R * R
N1 = N1 + 1
Next R
N1 = N1 – 1
M1 = S1 / N1
SS1 = SS1 - S1 * S1 / N1
N2 = 1
For Each R In rng2
S2 = S2 + R
SS2 = SS2 + R * R
N2 = N2 + 1
Next R
N2 = N2 – 1
M2 = S2 / N2
SS2 = SS2 - S2 * S2 / N2
Dim SS As Double
SS =(SS1 + SS2) * (1 / N1 + 1 / N2) / (N1 + N2 - 2)
TDEqual = (M1 - M2) / Sqr(SS)
End Function
统计过程，原假设，甲乙两人技能无明显差别，即单位时间内甲完成的零件个数和乙完成的零件个数无显著差异。
根据我们编写的VBA代码进行计算，T值为0.566，计算出p值等于0.2，大于0.05，所以，甲乙两人完成零件的个数无显著差异。
（3）方差不齐次的两样本T检验
当两个样本所在总体方差不等的时候，需要用到方差不其次的样本T检验。方差不齐次时，统计量T的公式为：
根据公式可以编写VBA代码如下
Function TDDif(Rng1 As Range, Rng2 As Range)
Dim S1 As Double, N1 As Integer, SS1 As Double, M1 As Double
Dim S2 As Double, N2 As Integer, SS2 As Double, M2 As Double
Dim R
N1 = 1
For Each R In Rng1
S1 = S1 + R
SS1 = SS1 + R * R
N1 = N1 + 1
Next R
N1 = N1 – 1
N2 = 1
For Each R In Rng2
S2 = S2 + R
SS2 = SS2 + R * R
N2 = N2 + 1
Next R
N2 = N2 – 1
M1 = S1 / N1
M2 = S2 / N2
SS1 = SS1 - S1 * S1.N1
SS2 = SS2 - S2 * S2 / N2
TDDif = (M1 - M2) / Sqr(SS1 / (N1 * N1 - N1) + SS2 / (N2 * N2 - N2))
End Function
（4）方差其次性检验
方差其次性检验用来检验两个样本方差是否其次，用两个样本的均方差中较大的一个除以较小的一个，即可以得到统计量F，根据统计量F的数值在F分布中求得积分，就得到两组方差相等的概率值。
方差其次性检验的公式如下：
根据公式编写VBA代码
Function FValue(Rng1 as Range,Rng2 as Range)
Dim S1 As Double,SS1 As Double,N1 As Integer
Dim S2 As Double,SS2 As Double,N2 As Integer
Dim R
For Each R In Rng1
S1=S1+R
SS1=SS1+R*R
N1=N1+1
Next R
For Each R In Rng2
S2=S2+R
SS2=SS2+R*R
N2=N2+1
Next R
SS1=SS1-S1*S1/N1
SS2=SS2-S2*S2/N2
Fvalue=(SS1/(N1-1))/(SS2/(N2-1))
If FValue<1 Then FValue=1/Fvalue
End Function
（5）配对样本T检验
配对样本T检验用于统计配对样本之间的总体均数是否相等。在实验设计时，设计配对样本可以尽量的减少实验误差，消除个体差异。下距离说明。
为了验证某种饲料对小鼠的血红蛋白含量的作用，现将小鼠用普通饲料饲喂一周，同一只小鼠再用饲料A饲喂一周，采集得到不同饲料饲喂时，每只小鼠在不同饲料饲喂后血红蛋白含量（g/L）见下表，问饲料是否会增加小鼠血红蛋白含量？
表一
小鼠编号
1
2
3
4
5
普通饲料
102
127
110
105
120
A饲料
110
132
121
112
129
配对样本的T检验过程，是划归为单一样本的总体均数T检验问题来解决的。这个问题中，先将小鼠不同饲料饲喂情况下的血红蛋白的差值计算出来，然后再将问题简化为：一个样本的总体均数是否为0，统计量T的公式为：
根据公式可以编写VBA代码如下：
Function TParied(Rng1 As Range, Rng2 As Range)
Dim D() As Double
Dim S As Double, SS As Double, N As Integer, i As Integer
Dim R1, R2
N = 1
For Each R1 In Rng1
ReDim Preserve D(N)
D(N) = R1
N = N + 1
Next R1
N = N – 1
i = 1
For Each R2 In Rng2
D(i) = D(i) - R2
i = i + 1
Next R2
For i = 1 To N
S = S + D(i)
SS = SS + D(i) * D(i)
Next i
SS = SS - S * S / N
TParied = (S / N) / Sqr(SS / (N * N - N))
End Function
统计过程，假设饲料A对小鼠血红蛋白含量无显著影响，根据配对T检验公式，构造统计量，并计算出T值等于-8.0，然后计算自由度为4的T分布积分，得到概率值为p=0.002，小于0.05。所以推翻原来饲料对小鼠血红蛋白无显著影响的假设，饲料对小鼠的血红蛋白含量有明显的增加作用。
（6）单因素方差分析
单因素方差分析，用于3组或者3组以上的比较，当实验的组数达到3组或者3组以上时，就不再适合用T检验来做统计分析了，此时的统计量为F。统计量F的定义是，组间数据的均方除以组内数据的均方。
平方和和自由度的剖分：在不同的实验组之间，得到一个数据的变异程度。数据的变异程度我们用平方和SS来表示，SS的定义是，一组数的平方和等于这一组数与其平均数的差值的平方之和。而一组数的均方MS的定义是，均方等于平方和除以自由度。在不同的组之间，我们可以得到组间的平方和SSA。此时，我们把总的平方和剖分为组间的平方和和组内的平方和，即SST=SSA+SSE，将总的自由度剖分为组间自由度和组内自由度。关于总平方和、组间平方和、组内平方和、总自由度、组间自由度，组内自由度的公式如下：
假设现在进行N组实验，每组实验有M个重复，则：
对于总平方和和总自由度
对于组间平方和和组间自由度：
对于组内平方和和自由度
对于F的计算：
这里我们给出一个Sub过程，读者可以仿照这个过程写一个计算F的函数
Sub Anova()
Dim Rng As Range, i As Integer, j As Integer, M As Integer, N As Integer, SST As Double, SSA As Double, SSE As Double
Dim dfT As Integer, dfA As Integer, dfE As Integer, CS As Double, SA As Double
Set Rng = Application.InputBox("选择连续的区域", "以行分组", , , , , , 8)
M = Rng.Columns.Count
N = Rng.Rows.Count
dfT = N * M – 1
dfA = N – 1
dfE = dfT – dfA
For i = 1 To N
SA = 0
For j = 1 To M
SST = SST + Rng.Cells(i, j) * Rng.Cells(i, j)
CS = CS + Rng.Cells(i, j)
SA = SA + Rng.Cells(i, j)
Next j
SSA = SSA + SA * SA
Next i
CS = CS * CS / (N * M)
SST = SST – CS
SSA = SSA / M – CS
SSE = SST – SSA
Sheet2.Cells(2, 2) = SSA
Sheet2.Cells(2, 3) = dfA
Sheet2.Cells(2, 4) = SSA / dfA
Sheet2.Cells(2, 5) = (SSA / dfA) / (SSE / dfE)
Sheet2.Cells(3, 2) = SSE
Sheet2.Cells(3, 3) = dfE
Sheet2.Cells(3, 4) = SSE / dfE
Sheet2.Cells(4, 2) = SST
Sheet2.Cells(4, 3) = dfT
Sheet2.Cells(4, 4) = SST / dfT
Sheet2.Range("F2").Formula = "=F.DiST.RT(E2,C2,C3)"
Sheet2.Activate
End Sub
程序作用，在Sheet2中显示单因素方差分析（Anova）的结果。
（7）相关分析和回归分析
作为统计中一对非常重要的概念，相关和回归在统计工作中占有相当大的比例。相关分析揭示两个实验变量直接的“相关性”，回归分析则是通过某种代数表达式表示二者直接的“关系”。相关和回归用于对某种事物的推测，即预测(Forecast)某个实验的结果。
在实际的科学实验中，我们也可以通过某些实验结果来观察两个变量之间的关系，在一定范围内，某种溶液对光的吸收量总是随着溶液浓度的增加而增加。这种一个变量总随着另一个变量变化而变化的特性就叫做“相关性”。
相关系数：相关系数是一个介于-1~1之间的数，表示两个变量之间变化的关联程度，相关系数的绝对值越高，表面二者之间的“相关性”越大；反之，相关系数的绝对值越小，二者之间的“相关性”也就越低。科学工作中，经常用相关系数的平方来判断两个变量之间的相关性。
当相关系数为一个正数时，两个变量之间的作用为“正相关”，即一个变量的增加将导致另一个变量的增加；当相关系数为负数时，两个变量之间的相关性为“负相关”，即一个变量的增加将导致另一个变量的减少。
相关系数的计算：两个变量x，y之间的相关系数r的公式为：
这里必须要说的是，在Excel中，已经含有计算相关系数的公式了，但是不妨碍我们自己写一个计算公式。
Function Correlate(Rng1 As Range, Rng2 As Range)
Dim A() As Double, B() As Double
Dim S1 As Double, SS1 As Double, S2 As Double, SS2 As Double, SP As Double
Dim R, N As Integer
N = 1
For Each R In Rng1
ReDim Preserve A(N)
A(N) = R
N = N + 1
Next R
N = 1
For Each R In Rng2
ReDim Preserve B(N)
B(N) = R
N = N + 1
Next R
N = N – 1
For i = 1 To N
S1 = S1 + A(i)
SS1 = SS1 + A(i) * A(i)
S2 = S2 + B(i)
SS2 = SS2 + B(i) * B(i)
SP = SP + A(i) * B(i)
Next i
SS1 = SS1 - S1 * S1 / N
SS2 = SS2 - S2 * S2 / N
SP = SP - S1 * S2 / N
Correlate = SP / Sqr(SS1 * SS2)
End Function
“相关性”表示两个变量之间的关乎关系，如果两个变量之间具有一定的相关性，则其中一个变量可以用另一个变量的函数进行表示，这种一个变量用另一个变量或者另外变量表示的过程就叫做“回归”。
在回归问题中，需要表示的变量被称为“因变量”，而另外的变量成为“自变量”，因变量总是随着自变量的变化而变化。我们把只有一个自变量的回归成为“一元回归”，把自变量和因变量之间的线性关系的回归叫做“线性回归”。
回归的最简单模型便是一元线性回归问题，在一元线性回归中，假设因变量y随着自变量x的变化而曾现线性的变化。一元线性回归问题，无非是找到两个合适的系数，使得Y=bX+a这个算式计算出来的Y与实际观察值y之间的“差异性”最小。
对于一元线性回归，具有下列公式：
可以仿照这个公式写一个VBA程序，同时输出回归方程和相关系数的平方R的值。
Function Reg(Rng1 As Range, Rng2 As Range) As String
Dim X() As Double, Y() As Double
Dim S1 As Double, SS1 As Double, S2 As Double, SS2 As Double, SP As Double
Dim R, N As Integer
Dim b As Double, a As Double
N = 1
For Each R In Rng1
ReDim Preserve X(N)
X(N) = R
N = N + 1
Next R
N = 1
For Each R In Rng2
ReDim Preserve Y(N)
Y(N) = R
N = N + 1
Next R
N = N – 1
For i = 1 To N
S1 = S1 + X(i)
SS1 = SS1 + X(i) * X(i)
S2 = S2 + Y(i)
SS2 = SS2 + Y(i) * Y(i)
SP = SP + X(i) * Y(i)
Next i
SS1 = SS1 - S1 * S1 / N
SS2 = SS2 - S2 * S2 / N
SP = SP - S1 * S2 / N
R = SP * SP / (SS1 * SS2)
b = SP / SS1
a = (S2 - S1 * b) / N
Reg = "y=" & b & "x+" & a & ",R=" & R
End Function
第三节 微积分算法的VBA代码实现
7.3.1 数值估算法
在计算机发展历程中，微积分思想解决了估算一个数值方法。在未发明“无穷大”和“无穷小”时，人们对于一些数值的计算显得无能为力。例如，求10.022这样的数值，按照一般的算法，就需要计算一个4位数的乘法。
在微积分的数值估算法中，有一个著名的泰勒公式：
、 这个公式就像从天而降一样，对于复杂的计算来说简直就是福音。它告诉我们，当我们对于结果的要求不是那么精确时，可以省略掉后面的计算部分而只保留前面的结果。
以10.022这个计算为例，说一下在微积分估算数值中的算法。
（1）将问题抽象成一个函数，现在这个函数就是f(x)=x2
（2）取合适的X0值，X0这里可以取10，ΔX=0.02
（3）对于y=f(X0+ΔX),在ΔX非常小的时候，可以约等于
f(X0)+f’(X0)ΔX
（4）又f(x)=x2,f’(x)=2x所以，f(X0+ΔX)=f(X0)+2X0ΔX
（5）带入，可以得到10.022=f(10.02)=f(10)+2×10×0.02=102+0.4=100.4
用计算器计算一下，10.022=100.4004，两者相差仅有0.0004.
再举一个例子，现在要求计算1000.1的开3次方，即1000.11/3.在没有计算器的时候，这个开立方的计算将无比复杂。但是任然可以利用上面的思路，这里，f(x)=x1/3.f’(x)=1/3x-2/3.那么，取X0=1000，ΔX=0.1.
f(X0+ΔX)=f(X0)+f’(X0) ΔX=10001/3+1/3×1000-2/3×0.1
10001/3=10,1000-2/3=0.01,所以，1000.11/3=10+0.00033=10.00033。跟计算器计算出的数值几乎一样。
类似的，我们可以写一个VBA程序来计算例如像4.02的开平方这样的计算。
Function MySqr(A As Double) As Double
Dim X0 As Integer
Dim DX As Double
X0=Fix(A)
DX=A-X0
MySqr=Sqr(X0)+0.5*DX/Sqr(X0)
End Function
读者可能会问，这里还是利用了系统的Sqr函数，那么可不可以不利用系统的Sqr函数来计算一个数的开平方呢？答案是肯定的，这里我们就需要用到上一章所将到的函数的递归调用了。因为我们已经知道，0的开平方是0, 1的开平方是1，因为现在已经知道了初始的数值，不妨让函数递归时，每一次减去0.1，一直减少到整数部分为0或者为1。这个程序就留给读者完成了！
7.3.2 曲边梯形面积的计算
函数y=f(x)与x轴，以及x=a，x=b构成的图形的面积（曲边梯形）被定义为f(x)在a到b上的定积分。定积分的计算方法是，将x轴上a到b之间的距离N等分，每一个小份就可以近似看出一个矩形，然后分别求取每一小份的面积，再对每个小份面积求和。当N趋近于无穷大时，求和的面积就等于曲边梯形面积。在高等数学中，我们通过计算一个函数原函数的方法来求定积分的精确值，但是，计算机中完全没有这个必要，因为计算机并不知道函数的原函数是什么，更不会计算例如负无穷大到正无穷大这样的反常积分。计算机只会做的事情，就是像上面那样的“分割”和“求和”，至于说取极限，也只要利用上一章所将的过程，即判断两次计算结果的差的绝对值是否小于一个确切的数值。
例如：利用定积分求sin(x)在0到π/2上的定积分。
Function IntG(N As Integer)
Dim a As Double
Dim b As Double
Dim StepI As Double
Dim i As Integer
Dim y As Double
a = 0
b = Application.WorksheetFunction.Pi() / 2
StepI = (b - a) / N
For i = 0 To N
y = Sin(a + StepI * i) * StepI
IntG = IntG + y
Next i
End Function
Function ResT()
Dim N As Integer
N = 100
While Abs(IntG(N) - IntG(N + 10)) > 0.00001
N = N + 10
Wend
ResT = IntG(N)
End Function
这里的两个函数，一个用于“分割求和”，一个用于“求极限”。直接调用ResT就可以得到结果1.000008.
7.3.3 计算概率值
根据上述“分割求和”与“取极限”的思想，我们可以求出一个概率的数值。在之前已经谈到，一个时间发生的概率等于这个事件的某个统计量在其分布函数上的一个定积分。现在我们暂且不谈诸如t检验中用到的t分布这样的太过于复杂的表达式，先看一看标准正态分布情况下的概率值的计算。在之前的讨论中，我们当时给出的概率值均是“两尾概率”，所谓的两尾概率指的是，当两个样本（或者更多的样本）的总体均数“不相等”时的概率，“不相等”的概率应该等于第一个变量大于第二个变量的概率加上第一个变量小于第二个变量的概率。因此，“两尾”概率等于“一侧”概率的两倍。
现在我们来看一下标准正太分布的两尾概率求解。
标准正态分布的分布函数为:
Function InTG(X As Double, N As Integer) As Double
Dim a As Double
Dim b As Double
Dim StepI As Double
Dim i As Integer
Dim y As Double
a = -X
b = X
StepI = (b - a) / N
For i = 0 To N
y = Exp(-(a + StepI * i) ^ 2 / 2) * StepI / Sqr(2 * 3.14)
InTG = InTG + y
Next i
End Function
Function ResT(X As Double) As Double
Dim N As Integer
N = 100
While Abs(InTG(X, N) - InTG(X, N + 10)) > 0.00001
N = N + 10
Wend
ResT = InTG(X, N)
ResT=1-ResT
End Function
最后这个ResT=1-ResT，因为整个积分等于1，两尾的概率等于1减去中间的概率。
大家可以用这个来计算一下，当x取值为2.95时，两尾的概率小于0.5%，达到显著水平，这就是工厂上以平均数±3倍标准差作为判断正常运行标准的原因。
第四节 按照条件搜索相应的记录
7.4.1 按精确条件搜索记录
平时工作中，总会遇到按条件搜索记录的要求，例如，现在我们想在工作表1（Sheet1）中的第N列（编号为A~Z的列）中找到与一个“目标”完全相同的行并将其数据拷贝到同一个工作簿的工作表2（Sheet2）中。
Sub MyFind()
Dim Rng As Range,N As Integer
Dim RngColumn As Range,Aim As String,Colu As String
Dim i As Integer, j As Integer
Set Rng=Sheet1.UsedRange
N=Rng.Rows.Count
L1:
Colu=Inputbox(“输入需要查找的列A~Z”)
If Asc(Colu)<65 or Asc(Colu)>90 Then Goto L1
Colu=Colu & “:” & Colu
Set RngColumn=Sheet1.Range(Colu)
Aim=Inputbox(“输入查找项”)
j=1
For i=1 to N
If RngColumn(i)=Aim Then Sheet1.Rows(i).Copy Sheet2.Rows(j): j=j+1
Next i
End Sub
7.4.2 按照范围条件搜索
与精确条件不同，范围条件往往是一个范围，对于这些范围的条件，我们现在的搜索内容不再是某一个具体的值，而是需要拓展到一个比如大于、小于、不等于等的条件。
例如，如果我们想找出表1（Sheet1）中某一列的数大于或小于一个数字的所有的行并拷贝到工作表2（Sheet2）中。这里，我们需要进行两个操作，第一个操作，确定需要提取的目标是“大于”或者“小于”，第二，确定目标是什么。对于第一个操作，我们可以通过一个字符串来解决，即向计算机输入大于号“＞”或者小于号“＜”，对于第二个操作，仍然可以用Aim表示。
Sub MyFind()
Dim Rng As Range, N As Integer
Dim RngC As Range, Aim As Double, Colu As String
Dim Op As String, Cmd As String
Dim i As Integer, j As Integer
Set Rng=Sheet1.UsedRange
N=Rng.Rows.Count
L1:
Colu=Inputbox(“输入需要查找的列A~Z”)
If Asc(Colu)<65 or Asc(Colu)>90 Then Goto L1
Colu=Colu & “:” & Colu
Set RngC=Sheet1.Range(Colu)
Cmd=Inputbox(“输入条件，符号和数字之间留一个空格”)
Op=Split(Cmd,” “)(0)
Aim=Val(Split(Cmd,” “)(1))
j=1
If Op=”>” Then
For i=1 To N
If RngC(i)>Aim Then Sheet1.Rows(i).Copy Sheet2.Rows(j) : j=j+1
Next I
Else
For i=1 to N
If RngC(i)<Aim Then Sheet1.Rows(i).Copy Sheet2.Rows(j) : j=j+1
Next i
End If
End Sub
请大家先用笔或者什么的把这个程序的代码记录一下，在下一章代码的健壮性中我们将对这个程序进行讨论。
7.4.3 多个条件的搜索
在多条件搜索的情况下，搜索条件不再是一个，而是可以拓展到多个条件，条件与条件之间可以是“或者”、“并且”等关系。现在先看一下两个条件的搜索。例如，想搜索表1（Sheet1）中符合B列大于10并且C列小于6的所有的行，并复制到表2（Sheet2）。注意到，这里考虑到通用性，不一定是B列和C列，中间的连接符也不一定是“并且”也可以是“或者”。
Sub MyFind()
Dim Cond As String
Dim Cond1 As String, Op As String,Cond2 As String, i As Integer, j As Integer
j=1
Cond=Inputbox(“输入条件”)
Cond1=Split(Cond,”.”)(0)
Op=Split(Cond,”.”)(1)
Cond2=Split(Cond,”.”)(2)
Select Case Op
Case “And”
For i=1 to Sheet1.UsedRange.Rows.Count
If MyComp(i,Cond1) And MyComp(I,Cond2) Then Sheet1.Rows(i).Copy Sheet2.Rows(j) : j=j+1
Next i
Case “Or”
For i=1 to Sheet1.UsedRange.Rows.Count
If MyComp(i,Cond1) OR MyComp(I,Cond2) Then Sheet1.Rows(i).Copy Sheet2.Rows(j) : j=j+1
Next i
Case Else
End Select
End Sub
Function MyComp(RowI As Integer,Cond As String)As Boolean
Dim C1 As String,Col As String
Dim Op As String
Dim C2 As String,V As Double
Dim Rng As Range
C1=Split(Cond,” “)(0)
C2=Split(Cond,” “)(1)
C3=Split(Cond,” “)(2)
V=Val(C3)
Col=C1 & “:” & C1
Set Rng=Sheet1.Range(Col)
Select Case C2
Case”>”
If Rng(RowI)>V Then MyComp=True Else MyComp=False
Case “<”
If Rng(RowI)<V Then MyComp=True Else MyComp=False
Case “=”
If Rng(RowI)=V Then MyComp=True Else MyComp=False
Case “!=”
If Rng(RowI)<>V Then MyComp=True Else MyComp=False
Case “>=”
If Rng(RowI)>=V Then MyComp=True Else MyComp=False
Case “<=”
If Rng(RowI)<=V Then MyComp=True Else MyComp=False
Case Else
MyComp=False
End Select
End Function
请读者同样也将这个程序“打包”复制一下，我们在第8章将讨论这个程序中可能出现的问题。
7.4.4 模糊搜索
模糊搜索是一个非常大的难题，计算机程序的编写法则是：不怕繁琐，就怕模糊。当一个条件模糊不清的时候，我们不清除需要计算机做些什么，这样的问题是最为头疼的。遇到此类问题时，先把问题尽量地量化和精确化，然后再给出具体解决方案。
例如：找出B列中数字较大的行。这里的问题就是，如何定义“较大？”假如我们将B列的最大数和最小数之间的差定义为B列数据的范围，那么，这个数据范围的前10%算不算“较大？”还是前20%算作“较大？”当然，这里需要根据B列共有多少个数来定。假设我们现在让用户自己输入这个百分比。
Sub MyFind()
Dim Rng As Range, Max As Double, Min As Double, i As Integer, j as Integer
Dim Rt As Double, Psent As Double
Set Rng=Sheet1.Range(“B:B”)
Max=Rng(1)
Min=Rng(1)
j=1
For i=1 To Sheet1.UsedRange.Rows.Count
If Rng(i)> Max Then Max=Rng(i)
If Rng(i)< Min Then Min=Rng(i)
Next i
Rt=Max-Min
Psent=Inputbox(“输入百分比%”)
Psent=Psent/100
For i=1 to Sheet1.Usedrange.Rows.Count
If Rng(i)>Max-Psent*Rt Then Sheet1.Rows(i).Copy Sheet2.Rows(j) : j=j+1
Next i
End Sub
当然这个问题也可以稍微改一下，即我们搜索全部排名的前5%或者10%，请读者思考一下这个搜索排名前5%的程序应该如何实现。
第五节 几个简单统计过程的VBA实现
7.5.1 统计过程
计算出每个统计量，并列表，然后计算p值的过程成为统计过程。在统计过程中涉及到的内容有：
（1）对该统计过程的每一个统计量进行计算。
（2）计算p值。
我们可以看到，这里不再是只求出一个具体的数值，而是需要显示整个统计的过程，所以这里用Sub过程比较好，我们利用Excel现有的单元格来显示我们统计的结果。这一节中，将不再详细讨论过程是如何实现的，而是只给出过程的VBA程序，对于VBA程序的语句，大家可以自行理解。
7.5.2 两样本总体均数T检验（方差其次）的Sub过程
Sub TTestFor()
Dim Rng1 As Range, Rng2 As Range, S1 As Double, S2 As Double, SS1 As Double, SS2 As Double, df As Integer
Dim N1 As Integer, N2 As Integer, i, SS As Double, T As Double
Set Rng1 = Application.InputBox("选择第一个区域", "第一个区域", , , , , , 8)
Set Rng2 = Application.InputBox("选择第二个区域", "第二个区域", , , , , , 8)
For Each i In Rng1
S1 = S1 + i
SS1 = SS1 + i * i
N1 = N1 + 1
Next i
SS1 = SS1 - S1 * S1 / N1
For Each i In Rng2
S2 = S2 + i
SS2 = SS2 + i * i
N2 = N2 + 1
Next i
SS2 = SS2 - S2 * S2 / N2
SS = SS1 + SS2
SS = SS / (N1 + N2 - 2)
SS = SS * (1 / N1 + 1 / N2)
SS = Sqr(SS)
T = (S1 / N1 - S2 / N2) / SS
df = N1 + N2 – 2
Sheet1.Range("B2") = S1 / N1
Sheet1.Range("B3") = S2 / N2
Sheet1.Range("B4") = S1 / N1 - S2 / N2
Sheet1.Range("C2") = Sqr(SS1 / (N1 - 1))
Sheet1.Range("C3") = Sqr(SS2 / (N2 - 2))
Sheet1.Range("C4") = SS
Sheet1.Range("D2") = N1 – 1
Sheet1.Range("D3") = N2 – 1
Sheet1.Range("D4") = N1 + N2 – 2
Sheet1.Range("E4") = T
Sheet1.Range("F4").Formula = "=T.DIST.2T(ABS(E4),D4)"
Sheet1.Activate
End Sub
作用：将一个Sheet中的两个单元格区域的数字做T检验。检验结果放在Sheet1中A1到F4单元格区域。
7.5.3 相关系数矩阵的Sub过程
Sub Cor()
Dim Rng As Range
Set Rng = Application.InputBox("选择区域", "连续区域", , , , , , 8)
Dim i, j As Integer
For i = 1 To Rng.Columns.Count
For j = 1 To Rng.Columns.Count
Sheet2.Cells(i, j) = WorksheetFunction.Correl(Rng.Columns(i), Rng.Columns(j))
Next j
Next i
End Sub
7.5.4求矩阵的逆矩阵的Sub过程
Sub Solve()
Dim RngO As Range, Rng1 As Range, i As Integer, j As Integer, k As Integer, N As Integer, M As Integer
Set RngO = Application.InputBox("请选择区域", "选择方形区域", , , , , , 8)
N = RngO.Rows.Count
M = RngO.Columns.Count
If N <> M Then
MsgBox ("区域选择有误！")
Else
For k = 1 To N
Set Rng1 = RngO.Offset(RngO.Rows.Count + 1, 0)
Rng1.Cells(k, k) = 1 / RngO.Cells(k, k)
For j = 1 To N
If j <> k Then Rng1.Cells(k, j) = RngO.Cells(k, j) / RngO.Cells(k, k): Rng1.Cells(j, k) = -RngO.Cells(j, k) / RngO.Cells(k, k)
Next j
For i = 1 To N
For j = 1 To N
If i <> k And j <> k Then Rng1.Cells(i, j) = RngO.Cells(i, j) - RngO.Cells(i, k) * Rng1.Cells(k, j)
Next j
Next i
Set RngO = Rng1
Next k
End If
End Sub