逐步判别f2.frm

<VB数理统计实用算法>书中的算法源程序

Attribute VB_Creatable = False
Attribute VB_PredeclaredId = True
Attribute VB_Exposed = False
'逐步判别
Option Explicit
Dim I As Integer, J As Integer, K As Integer
Dim sngH As Single, sngC As Single

Private Sub Form_Load()
    lblC(0).Visible = False: lblTV(0).Visible = False
    lblCh(0).Visible = False: lblTr(0).Visible = False
    lbltY(0).Visible = False
    Q(1) = 0 '用样本频率代替先验概率，如果 Q(1) = 1,则认为先验概率相等
End Sub

'计算
Private Sub cmdCalculate_Click()
    Dim G1 As Integer, M1 As Integer, H As Integer
    Dim step As Integer, Imin As Integer, Imax As Integer
    Dim R As Integer, R12 As Integer
    Dim D As Double, D1 As Double, eps As Double
    Dim Umax As Double, Umin As Double, Ui As Double
    Dim U As Double, X2 As Double
    Dim F12 As Double, F1 As Double, F2 As Double
    If cmdCalculate.Caption = "继续" Then
        Unload Me: frmContinue.Visible = True
    End If
    lblC(0).Visible = True: lblTV(0).Visible = True
    lblCh(0).Visible = True: lblTr(0).Visible = True
    lbltY(0).Visible = True
    F1 = Val(txtIn.Text): F2 = Val(txtOut.Text) '取得引入F和剔出F
    If F1 = F2 = 0 Then MsgBox "F1 = F2 =0，全部变量都引进，不剔除"
    If F1 < F2 Then
        MsgBox "F1不能小于F2！", vbOKOnly + vbExclamation, "参数错误"
        Exit Sub
    End If
    G1 = G + 1: M1 = M + 1
    For I = 1 To G
        Ng(I) = 0
    Next I
    For I = 1 To G1
        For J = 1 To M
            Mx(I, J) = 0
        Next J
    Next I
    For K = 1 To N
        H = XY(K, M1)
        If (H - 1) * (H - G) > 0 Then
            MsgBox ("原分类错，修改后重作")
            End
        End If
        Ng(H) = Ng(H) + 1                       'Ng为各种分类的样本数
        For I = 1 To M
            Mx(H, I) = XY(K, I) + Mx(H, I)
        Next I
    Next K
    For I = 1 To M
        D = 0
        For H = 1 To G
            D = D + Mx(H, I)
            Mx(H, I) = Mx(H, I) / Ng(H)         '平均值
        Next H
        Mx(G1, I) = D / N                       '总平均值
    Next I
    For I = 1 To M
        For J = 1 To M
            T(I, J) = 0
            W(I, J) = 0
        Next J
    Next I
    For K = 1 To N
        For I = 1 To M
            D = XY(K, I) - Mx(G1, I): X(I) = D
            For J = 1 To I
                T(I, J) = D * X(J) + T(I, J)    '总离差矩阵
            Next J
        Next I
    Next K
    For H = 1 To G
        D = Ng(H)
        For I = 1 To M
            D1 = Mx(H, I) - Mx(G1, I)
            X(I) = D * D1
            For J = 1 To I
                W(I, J) = D1 * X(J) + W(I, J)   '组间离差 B，先置于 W 的下三角
            Next J
        Next I
    Next H
'计算组内离差矩阵 W = T - B，并把 W、T 补成方阵
    For I = 1 To M
        For J = 1 To I
            W(I, J) = T(I, J) - W(I, J)
            W(J, I) = W(I, J)
            T(J, I) = T(I, J)
        Next J
    Next I
    D1 = Q(1)
    For H = 1 To G
        If D1 = 0 Then
            D = Ng(H) / N                       '用样本频率代替先验概率
        Else
            If D1 = 1 Then                      '先验概率相等
                D = 1 / G
            Else
                D = Q(H)
            End If
        End If
        LnQ(H) = Log(D) / Log(2.718282)
    Next H
'U记逐步计算中之Wilks量，U = |W| / |T|
'X记录选出的量，X(i) = 1 表示 Xi 已选中
L2: L = 0: step = 0: eps = 0.000001: U = 1
    For I = 1 To M
        X(I) = 0
    Next I
Lstep: step = step + 1
    Umax = eps: Umin = 1000
    For I = 1 To M
        If X(I) = 0 Then
            If T(I, I) >= eps Then
                Ui = W(I, I) / T(I, I)
                If Ui < eps Then
                    Ui = eps
                    If W(I, I) = 0 Then GoTo Lend
                End If
            End If
            If Ui < Umin Then
                Umin = Ui                   '从未选量中找出 Umin
                Imin = I
            End If
        Else
            Ui = T(I, I) / W(I, I)          '计算各变量的离差比 Ui
            If Ui > Umax Then
                Umax = Ui                   '从已选量中找出 Umax
                Imax = I
            End If
        End If
    Next I
'开始检验
    F12 = (1 - Umax) * (N - L - G + 1) / (Umax * (G - 1))
    If F12 < F2 Then
        L = L - 1: R = Imax
        R12 = -R: X(R) = 0
    Else
        F12 = (1 - Umin) * (N - L - G) / (Umin * (G - 1))
        If F12 <= (F1 + eps) Then
            GoTo L3
        Else
            L = L + 1: R = Imin: R12 = R
            X(R) = 1
        End If
    End If
    U = (W(R, R) / T(R, R)) * U             '逐步计算 U = |W| / |T|
'X2是U的近似值
    X2 = -(N - 1 - (L + G) / 2) * Log(U) / Log(2.718282)
    D = 1 / T(R, R)
    T(R, R) = 1
    D1 = 1 / W(R, R)
    W(R, R) = 1
    For J = 1 To M
        T(R, J) = T(R, J) * D
        W(R, J) = W(R, J) * D1
    Next J
    For I = 1 To M
        If I <> R Then
            D = T(I, R)
            T(I, R) = 0
            D1 = W(I, R)
            W(I, R) = 0
            For J = 1 To M
                T(I, J) = T(I, J) - D * T(R, J)
                W(I, J) = W(I, J) - D1 * W(R, J)
            Next J
        End If
    Next I
    CalCoeD2 N, L, X, W, Mx, C, D2gh, Ng, Fhg
    Class XY, X, C, LnQ, Pr, Dm, New_Sort
    GoTo Lstep
L3: If L = 0 Then GoTo Lend
    CalCoeD2 N, L, X, W, Mx, C, D2gh, Ng, Fhg
    Class XY, X, C, LnQ, Pr, Dm, New_Sort
'利用不显示的标签数组保存计算结果
    For I = 1 To N
        Load lblC1(I): Load lblTV1(I): Load lblTy1(I)
        Load lblCh1(I): Load lblTr1(I)
        sngC = Val(I): lblC1(I).Caption = sngC          '样本号
        sngC = New_Sort(I, 1): lblTV1(I).Caption = sngC '原有分类
        sngC = New_Sort(I, 2): lblTy1(I).Caption = sngC '判别分类
        If lblTV1(I).Caption = lblTy1(I).Caption Then lblCh1(I) = " " _
            Else lblCh1(I) = "*"                        '分类改变
'后验概率，取3位有效数字
        sngC = New_Sort(I, 3): lblTr1(I).Caption = Int(sngC * 1000 + 0.5) / 1000
    Next I
'利用标签数组显示判别分类结果
    sngH = lblC(0).Height                       '标签数组元素的高度
    For I = 1 To N                              '置放标签数组
        Load lblC(I): Load lblTV(I): Load lbltY(I)
        Load lblCh(I): Load lblTr(I)
        lblC(I).Move lblC(0).Left, lblC(0).Top + I * sngH
        lblTV(I).Move lblTV(0).Left, lblTV(0).Top + I * sngH
        lbltY(I).Move lbltY(0).Left, lbltY(0).Top + I * sngH
        lblCh(I).Move lblCh(0).Left, lblCh(0).Top + I * sngH
        lblTr(I).Move lblTr(0).Left, lblTr(0).Top + I * sngH
        lblC(I).Visible = True: lblTV(I).Visible = True
        lbltY(I).Visible = True: lblCh(I).Visible = True
        lblTr(I).Visible = True
        lblC(I) = lblC1(I): lblTV(I) = lblTV1(I): lbltY(I) = lblTy1(I)
        lblCh(I) = lblCh1(I): lblTr(I) = lblTr1(I)
    Next I
    cmdCalculate.Caption = "继续"
    Exit Sub
Lend:
    MsgBox "判别分类数L为0或组内离差矩阵W(I,I)=0"
End Sub

'在右侧控制判别结论的垂直滚动条
Private Sub VScroll1_Change()
    Dim V As Integer
    On Error Resume Next
    V = VScroll1.Value
    lblC(0) = "样本号": lblTV(0) = "原有分类": lbltY(0) = "判别分类"
    lblCh(0) = "分类改变": lblTr(0) = "后验概率"
    For I = 1 To intRow
        If I + V > intRow Then GoTo 10
        lblC(I) = lblC1(I + V)                    '编号
        lblTV(I) = lblTV1(I + V)                  '原有分类
        lbltY(I) = lblTy1(I + V)                  '判别分类
        lblCh(I) = lblCh1(I + V)                  '分类改变
        lblTr(I) = lblTr1(I + V)                  '后验概率
        GoTo 20
10:     lblC(I) = "": lblTV(I) = "": lbltY(I) = ""
        lblCh(I) = "": lblTr(I) = ""
20: Next I
End Sub

'结束
Private Sub cmdExit_Click()
    Unload Me
    End
End Sub