【Python】一个简陋的基于混淆矩阵计算遥感分类精度（总体精度、Kappa系数、用户精度、生产者精度等）的代码

本文链接：https://blog.csdn.net/leaze932822995/article/details/123180360

1 简介

有一篇文章返修了，由于文章的重点不在分类所以我之前就只写了个Kappa系数上去，没想到审稿人居然要求我提供其他参数ಥ_ಥ可是我只大概存了个各类型的分类。。。虽然后来从垃圾堆里翻了下数据，但也只能得到一个混淆矩阵。。。
说实话分类我做得也不多，作为纯ArcGIS党以前只知道有软件可以跑各种精度或GEE用现成代码跑。。这一块更是纯纯地没有经验。。
本着靠自己不靠别人的原则，这次我就想写一个基于混淆矩阵计算总体精度（Overall Accuracy，OA）、Kappa系数（Kappa Coefficient）、错分误差（Commission Error，CE）、漏分误差（Omission Error，OE）、生产者精度（Producer‘s Accuracy，PA）、用户精度（User’s Accuracy，UA）的代码。
在此之前，我先找呼院士咨询相关公式，它给了我一篇徐小洋的文章。所以这次代码里的公式都是基于这篇文章。

2 案例和代码说明

为了避免一些不必要的麻烦。。我这里拿出我现场随便瞎写的案例来做示例。
如果需要从excel还是哪里读表，可以自行检索相关的读取方法。
首先，我的遥感分类场景里有8类地物，这里敲进去一个矩阵（实际上是8个list）。

a = [200, 2  , 20 , 1  , 9  , 000, 000, 1  ]
b = [6  , 50 , 10 , 2  , 20 , 000, 4  , 6  ]
c = [12 , 6  , 300, 5  , 5  , 000, 000, 11 ]
d = [1  , 1  , 3  , 50 , 9  , 000, 000, 18 ]
e = [4  , 2  , 6  , 2  , 300, 7  , 1  , 22 ]
f = [000, 1  , 2  , 1  , 30 , 80 , 2  , 9  ]
g = [1  , 5  , 3  , 2  , 10 , 10 , 24 , 4  ]
h = [1  , 3  , 15 , 5  , 30 , 5  , 2  , 300]

到时候如果类别个数不一样可以直接改。
在这里插入图片描述
然后我再写一行把它合成一个二维矩阵。

Confusion_Matrix = [
    a, b, c, d, e, f, g, h
]

上面就是需要结合实际情况改的内容了。。也就是说整个代码只有上面的要动。
我怕有的人不懂改，那么我基于上述8类的场景把它改成一个4类的场景。注意：我只删去了后4类。

a = [200, 2  , 20 , 1  ]
b = [6  , 50 , 10 , 2  ]
c = [12 , 6  , 300, 5  ]
d = [1  , 1  , 3  , 50 ]

Confusion_Matrix = [
    a, b, c, d
]

好的，接下来我写函数部分，函数部分不用改直接抄就行了。

import copy

def OA(Confusion_Matrix):
    CM = Confusion_Matrix
    xia_list = []
    for i in range(len(CM)):
        xia_list.append(sum(CM[i]))
    xia = sum(xia_list)
    shang_list = []
    for i in range(len(CM)):
        shang_list.append(CM[i][i])
    shang = sum(shang_list)
    OA = shang / xia
    return OA

def KP(Overall_Accuracy, Confusion_Matrix):
    CM = Confusion_Matrix
    p0 = Overall_Accuracy
    pe_xia_list = []
    for i in range(len(CM)):
        pe_xia_list.append(sum(CM[i]))
    pe_xia = sum(pe_xia_list) ** 2
    CMT = list(map(list, zip(*CM))) # 转置
    pe_shang_list = []
    for i in range(len(CM)):
        pe_shang_list.append(sum(CM[i]) * sum(CMT[i]))
    pe_shang = sum(pe_shang_list)
    pe = pe_shang / pe_xia
    KP = (p0 - pe) / (1 - pe)
    return KP

def CE(Class, Confusion_Matrix):
    CM = Confusion_Matrix
    Class_list = CM[Class]
    Commission_list = copy.deepcopy(Class_list)
    del Commission_list[Class]
    CE = sum(Commission_list) / sum(Class_list)
    return CE

def OE(Class, Confusion_Matrix):
    CM = Confusion_Matrix
    CMT = list(map(list, zip(*CM))) # 转置
    Class_list = CMT[Class]
    Ommission_list = copy.deepcopy(Class_list)
    del Ommission_list[Class]
    OE = sum(Ommission_list) / sum(Class_list)
    return OE

从上到下依次是总体精度（Overall Accuracy，OA）、Kappa系数（Kappa Coefficient）、错分误差（Commission Error，CE）、漏分误差（Omission Error，OE），而生产者精度（Producer‘s Accuracy，PA）和用户精度（User’s Accuracy，UA）基于CE和OE就可以直接得出，所以就不写复杂的函数了（懒）。
像OA这个函数，参数位把矩阵Confusion_Matrix扔进去就行了；
KP这个函数，需要填你计算好的OA和Confusion_Matrix；
CE和OE这两个函数是基于某一地类的，因此需要填入一个地类的序号，比如说第一类就填0，第二类就填1；
UA和PA是基于CE和OE得出的。。。所以我就没写，直接看输出部分吧。
最后写一个输出把计算的这些精度指标print出来。

# 总体精度
Overall_Accuracy = OA(Confusion_Matrix)
print(Overall_Accuracy)
# 卡帕系数
Kappa = KP(Overall_Accuracy, Confusion_Matrix)
print(Kappa)

# 针对某一类
# 错分误差
Commission_Error = CE(0, Confusion_Matrix) # 第一个参数位写第几类，0是第一类，1是第二类
print(Commission_Error)
# 漏分误差
Omission_Error = OE(0, Confusion_Matrix) # 第一个参数位写第几类，0是第一类，1是第二类
print(Omission_Error)
# 用户精度
Users_Accuracy = 1 - Commission_Error
print(Users_Accuracy)
# 生产者精度/制图精度
Producers_Accuracy = 1 - Omission_Error
print(Producers_Accuracy)

下面是8地类案例print出来的结果。一共有6个指标。
在这里插入图片描述

3 完整代码

a = [200, 2  , 20 , 1  , 9  , 000, 000, 1  ]
b = [6  , 50 , 10 , 2  , 20 , 000, 4  , 6  ]
c = [12 , 6  , 300, 5  , 5  , 000, 000, 11 ]
d = [1  , 1  , 3  , 50 , 9  , 000, 000, 18 ]
e = [4  , 2  , 6  , 2  , 300, 7  , 1  , 22 ]
f = [000, 1  , 2  , 1  , 30 , 80 , 2  , 9  ]
g = [1  , 5  , 3  , 2  , 10 , 10 , 24 , 4  ]
h = [1  , 3  , 15 , 5  , 30 , 5  , 2  , 300]

Confusion_Matrix = [
    a, b, c, d, e, f, g, h
]

import copy

def OA(Confusion_Matrix):
    CM = Confusion_Matrix
    xia_list = []
    for i in range(len(CM)):
        xia_list.append(sum(CM[i]))
    xia = sum(xia_list)
    shang_list = []
    for i in range(len(CM)):
        shang_list.append(CM[i][i])
    shang = sum(shang_list)
    OA = shang / xia
    return OA

def KP(Overall_Accuracy, Confusion_Matrix):
    CM = Confusion_Matrix
    p0 = Overall_Accuracy
    pe_xia_list = []
    for i in range(len(CM)):
        pe_xia_list.append(sum(CM[i]))
    pe_xia = sum(pe_xia_list) ** 2
    CMT = list(map(list, zip(*CM))) # 转置
    pe_shang_list = []
    for i in range(len(CM)):
        pe_shang_list.append(sum(CM[i]) * sum(CMT[i]))
    pe_shang = sum(pe_shang_list)
    pe = pe_shang / pe_xia
    KP = (p0 - pe) / (1 - pe)
    return KP

def CE(Class, Confusion_Matrix):
    CM = Confusion_Matrix
    Class_list = CM[Class]
    Commission_list = copy.deepcopy(Class_list)
    del Commission_list[Class]
    CE = sum(Commission_list) / sum(Class_list)
    return CE

def OE(Class, Confusion_Matrix):
    CM = Confusion_Matrix
    CMT = list(map(list, zip(*CM))) # 转置
    Class_list = CMT[Class]
    Ommission_list = copy.deepcopy(Class_list)
    del Ommission_list[Class]
    OE = sum(Ommission_list) / sum(Class_list)
    return OE

# 总体精度
Overall_Accuracy = OA(Confusion_Matrix)
print(Overall_Accuracy)
# 卡帕系数
Kappa = KP(Overall_Accuracy, Confusion_Matrix)
print(Kappa)

# 针对某一类
# 错分误差
Commission_Error = CE(0, Confusion_Matrix) # 第一个参数位写第几类，0是第一类，1是第二类
print(Commission_Error)
# 漏分误差
Omission_Error = OE(0, Confusion_Matrix) # 第一个参数位写第几类，0是第一类，1是第二类
print(Omission_Error)
# 用户精度
Users_Accuracy = 1 - Commission_Error
print(Users_Accuracy)
# 生产者精度/制图精度
Producers_Accuracy = 1 - Omission_Error
print(Producers_Accuracy)