Fisher最优分割法（附python实现）

1 最优分割法简介

最有分割法是对有序样品的一种聚类方法。当样品是按顺序排列，在分类中不允许打破样品的顺序。即 ，对 个有序样品进行分割，就可能有 种划分方法，这每一种分法称为一种分割。在所有的这些分割中，找到一种分割法，这种分割法使得各段内样品之间的差异最小，而各段之间的差异最大。这种对 个样品分段并使组内离差平方和最小的分割方法，称为最优分割法。段内数值变化可用变差或者极差来表示，比如样品段 变差：

  ，

表示样品段内样本间的差异情况，越小表示段内各样品之间差异性越小；因此，把 称为段直径。

2 Fisher最优分割法简介

Fisher最优分割法是用离差平方和来表示同类样本之间的差异程度，通过简便的计算步骤和作图，确定最优分类数，使同类样本间的差异最小，各类别样本间的差异最大，并用F检验法检验最优分类数的合理性。

设某一类 包含的样品有 ，记为 ，该类的均值向量 为   ， ， 表示这一类的直径，常用的直径有 ，当是单变量时，也可定义直径为­ ， 是中位数。

3 Fisher最优分割递推公式

用 表示将n 个有序样品分为 类的一种分法，定义损失函数为 ，当 和 固定时，越小表示各类的离差平方和越小，因此要寻找一种分法 ，使分类损失函数 达到最小，记为 。

递推公式：

             

公式含义：如果要找到n 个样品分为 个类的最优分割，应建立在将 个样品分为k-1 类的最优分割的基础上。

4 Fisher最优分割计算流程

若分类数 已知，求分类法 使它在损失函数意义下达最小，其求法如下：首先找分点 使递推公式达极小，即 。于是得第k类 。然后找 ，使它满足 得到第 类。类似的方法一次可以得到所有类 ,这就是所求的最优解，即： 。

总之，为了求最优解，主要是计算 和 }

5 Fisher最优分割-Python完整代码：

import numpy as np

def G(X, i, j):

    #tex:

    #$$ \overline{x_G}=\frac{1}{j-i+1}\sum^{j}_{t=i}X_t $$

    xg = 0

    for t in range(i-1, j):

        xg += X[t]

    xg /= (j - i + 1)

    return xg

def D(X, i, j):

    #tex:

    #$$ D(i,j)=\sum^{j}_{t=i}(X_t-\overline{x_G})'(X_t-\overline{x_G}) $$

    d = 0

    xg = G(X, i, j)

    for t in range(i-1, j):

        d += (X[t] - xg) ** 2

    return d

def Lq(X, N, K):

    #tex:

    #$$ L(b(n,k))=\sum^k_{i=1}D(i_t,i_{t+1}-1) $$

    #$$ \begin{equation} \left\{ \begin{array}{lr} L(p(n,2))=1 \\ L(p(n,k))=\underset{k \le j \le n}{min} D(1,j-1)+D(n,j) \end{array} \right. \end{equation} $$

    l = np.zeros((N+1, K+1))

    for n in range(1, N+1):

        l[n, 1] = D(X, 1, n)

    for k in range(2, K+1):

        for n in range(k, N+1):

            l[n, k] = min([l[j-1, k-1] + D(X, j, n) for j in range(k, n+1)])

    return l

def B(Lp):

    N = Lp.shape[0] - 1

    K = Lp.shape[1] - 1

    Bl = np.zeros(K+1)

    for k in range(1, K):

        Bl[k] = Lp[N, k] / Lp[N, k+1]

    return Bl

def split_class(X, K, l=None):

    if l is None:

        l = Lq(X, len(X), K * 2 + 1)

        Bl = B(l)

        for i in range(1, K):

            if Bl[i]-1 <= 0.05:

                K = i+1

                l = Lq(X, len(X), K * 2 + 1)

                break

    N = len(X)

    if K == 1:

        return []

    else:

        dist = np.array([abs(l[t-1, K-1] + D(X, t, N) - l[N, K]) for t in range(K, N+1)])

        t = list(range(K, N+1))[int(np.argmin(dist))]

        r = split_class(X[:t-1], K-1, l)

        r.append(t - 1)

        return r

直接调用split_class函数即可

6 Fisher最优分割-确定分类数K

分类个数的确定如果能从实际问题中首先确定当然最好，如果不能则可以从 随k的变化趋势图中找到拐点处，作为确定k的依据。当曲线拐点很平缓时，可选择的k较多，这时需要用其他方法来确定，如均方比法和特征根法。

稳定分类数可以这样定义: 当分类数小于此值时,  较大, 类内样本间的差别仍然显著, 需要继续再细分, 否则分类选取时, 代表性和信息量不能保证; 当分类数大于此值时, 类内样本间的差别已经不显著, 再细分的意义不大。稳定分类数代表了 值由较大值向较小值变化的转折点。按分类选取的原则在每一类中选取一个代表点, 则在优化布点时, 分类数应不小于稳定分类数, 才能保证优选点的代表性。