Python实现Fleiss Kappa一致性分析,并计算Z值和p值等相关统计量

news2024/11/15 14:01:14

参考资料

Fleiss Kappa的定义

Fleiss Kappa的原论文因为要付费才能阅读,我这里就不放链接了 

Fleiss' kappa - Wikipediahttps://en.wikipedia.org/wiki/Fleiss%27_kappa

Fleiss Kappa相关统计量 Z值,p值,95%置信区间

属性一致性分析 的 kappa 统计量的方法和公式 - Minitab请选择您所选的方法或公式。https://support.minitab.com/zh-cn/minitab/20/help-and-how-to/quality-and-process-improvement/measurement-system-analysis/how-to/attribute-agreement-analysis/attribute-agreement-analysis/methods-and-formulas/kappa-statistics/#testing-significance-of-fleiss-kappa-unknown-standard刘伟. 属性值测量系统的相关性研究[D].南京理工大学,2010.https://kns.cnki.net/kcms2/article/abstract?v=30M9_8jKqe0eeyUFGF2Hn7wI0LSs-Y73fT2Qnj25ZztRFGn-k7w38HpuuprNL7uDobRQegHTkIu0rsYpWZAWQXSLhgCo5vf_lwfoH8muLibc8UPwv_ZzuvwiTEcRGaXb&uniplatform=NZKPT&language=CHSFleiss' kappa in SPSS Statistics | Laerd Statisticshttps://statistics.laerd.com/spss-tutorials/fleiss-kappa-in-spss-statistics.php

第一个和第二个链接提供了计算公式,第三个链接是SPSS的计算Fleiss Kappa的图文教程,没有涉及公式推理。

前情提要

网上关于Fleiss Kappa的资料比较多,但是相关统计量,如Z值和P值的资料比较少。维基百科只提供了Fleiss Kappa的计算方法,没有提及相关统计量的计算。

期间找了很多工具,如SPSS,Minitab和SPSSAU。SPSS和Minitab是付费软件,SPSSAU是免费的在线分析平台。但SPSSAU只允许5万条及以下的数据进行分析,而我的数据量超过了5万条。

找遍全网好像没有人公开计算Fleiss Kappa相关统计量的代码,而可以计算这些统计量的软件要么是付费的,要么是存在数据量的限制,我就想能不能自己实现一个计算工具。

Fleiss Kappa原理

关于Fleiss Kappa的原理部分,我觉得还是看维基百科比较好,他写的比较清楚,还提供了数据示例。下面是关于Fleiss Kappa的计算代码。

import numpy as np
def fleiss_kappa(data: np.array):
    """
    Calculates Fleiss' kappa coefficient for inter-rater agreement.

    Args:
        data: numpy array of shape (subjects, categories), where each element represents
              the number of raters who assigned a particular category to a subject.

    Returns:
        kappa: Fleiss' kappa coefficient.
    """
    subjects, categories = data.shape
    n_rater = np.sum(data[0])

    p_j = np.sum(data, axis=0) / (n_rater * subjects)
    P_e_bar = np.sum(p_j ** 2)

    P_i = (np.sum(data ** 2, axis=1) - n_rater) / (n_rater * (n_rater - 1))
    P_bar = np.mean(P_i)

    K = (P_bar - P_e_bar) / (1 - P_e_bar)

subjects是样本数量,相当于维基百科中的N

categories是类别数量,相当于维基百科中的k

n_rater是投票者的数量,相当于维基百科中的n

data是输入的矩阵,第i行j列的元素是维基百科中的n_{ij}

这里我提供维基百科的公式截图,可以对照看一下

 

相关统计量的计算

下图是Minitab提供的公式,链接我放在了博客开头

 

z值算出来后,p-value和95%置信区间就水到渠成了。

在代码中我使用tmp进行了替换,简化了表达式 

    tmp = (1 - P_e_bar) ** 2
    var = 2 * (tmp - np.sum(p_j * (1 - p_j) * (1 - 2 * p_j))) / (tmp * subjects * n_rater * (n_rater - 1))
    
    # standard error
    SE = np.sqrt(var) 

    Z = K / SE

    p_value = 2 * (1 - norm.cdf(np.abs(Z)))

    ci_bound = 1.96 * SE / subjects
    lower_ci_bound = K - ci_bound
    upper_ci_bound = K + ci_bound

下面是完整代码

import numpy as np
from scipy.stats import norm

def fleiss_kappa(data: np.array):
    """
    Calculates Fleiss' kappa coefficient for inter-rater agreement.

    Args:
        data: numpy array of shape (subjects, categories), where each element represents
              the number of raters who assigned a particular category to a subject.

    Returns:
        kappa: Fleiss' kappa coefficient.
    """
    subjects, categories = data.shape
    n_rater = np.sum(data[0])

    p_j = np.sum(data, axis=0) / (n_rater * subjects)
    P_e_bar = np.sum(p_j ** 2)

    P_i = (np.sum(data ** 2, axis=1) - n_rater) / (n_rater * (n_rater - 1))
    P_bar = np.mean(P_i)

    K = (P_bar - P_e_bar) / (1 - P_e_bar)

    tmp = (1 - P_e_bar) ** 2
    var = 2 * (tmp - np.sum(p_j * (1 - p_j) * (1 - 2 * p_j))) / (tmp * subjects * n_rater * (n_rater - 1))
    
    # standard error
    SE = np.sqrt(var) 

    Z = K / SE

    p_value = 2 * (1 - norm.cdf(np.abs(Z)))

    ci_bound = 1.96 * SE / subjects
    lower_ci_bound = K - ci_bound
    upper_ci_bound = K + ci_bound

    print("Fleiss Kappa: {:.3f}".format(K))
    print("Standard Error: {:.3f}".format(SE))
    print("Z: {:.3f}".format(Z))
    print("p-value: {:.3f}".format(p_value))
    print("Lower 95% CI Bound: {:.3f}".format(lower_ci_bound))
    print("Upper 95% CI Bound: {:.3f}".format(upper_ci_bound))
    print()

这个函数只能处理格式形如维基百科示例的数据,对于其他格式的数据,需要相关的转换函数。

这里提供了两个转换函数,和对应的测试数据

def transform(*raters):
    """
    Transforms the ratings of multiple raters into the required data format for Fleiss' Kappa calculation.

    Args:
        *raters: Multiple raters' ratings. Each rater's ratings should be a list or array of annotations.

    Returns:
        data: numpy array of shape (subjects, categories), where each element represents the number of raters
              who assigned a particular category to a subject.

    """
    assert all(len(rater) == len(raters[0]) for rater in raters), "Lengths of raters are not consistent."
    
    subjects = len(raters[0])
    categories = max(max(rater) for rater in raters) + 1
    data = np.zeros((subjects, categories))

    for i in range(subjects):
        for rater in raters:
            data[i, rater[i]] += 1
    
    return data

def tranform2(weighted):
    """
    Transforms weighted data into the required data format for Fleiss' Kappa calculation.

    Args:
        weighted: List of weighted ratings. Each row represents [rater_0_category, rater_1_category, ..., rater_n_category, weight].

    Returns:
        data: numpy array of shape (subjects, categories), where each element represents the number of raters
              who assigned a particular category to a subject.

    """
    n_rater = len(weighted[0]) - 1
    raters = [[] for _ in range(n_rater)]
    for i in range(len(weighted)):
        for j in range(len(raters)):
            raters[j] = raters[j] + [weighted[i][j] for _ in range(weighted[i][n_rater])]
    
    data = transform(*raters)
    
    return data

def test():
    # Example data provided by wikipedia https://en.wikipedia.org/wiki/Fleiss_kappa
    data = np.array([
        [0, 0, 0, 0, 14],
        [0, 2, 6, 4, 2],
        [0, 0, 3, 5, 6],
        [0, 3, 9, 2, 0],
        [2, 2, 8, 1, 1],
        [7, 7, 0, 0, 0],
        [3, 2, 6, 3, 0],
        [2, 5, 3, 2, 2],
        [6, 5, 2, 1, 0],
        [0, 2, 2, 3, 7]
    ])

    fleiss_kappa(data)

    # need transform
    rater1 = [1, 2, 2, 1, 2, 2, 1, 1, 3, 1, 2, 2]
    rater2 = [1, 2, 1, 2, 1, 2, 3, 2, 3, 2, 3, 1]
    rater3 = [1, 2, 2, 1, 3, 3, 3, 2, 1, 2, 3, 1]

    data = transform(rater1, rater2, rater3)
    fleiss_kappa(data)

    # The first row indicates that both rater 1 and 2 rated as category 0, this case occurs 8 times.
    # need transform2
    weighted_data = [
        [0, 0, 8],
        [0, 1, 2],
        [0, 2, 0],
        [1, 0, 0],
        [1, 1, 17],
        [1, 2, 3],
        [2, 0, 0],
        [2, 1, 5],
        [2, 2, 15]
    ]
    data = tranform2(weighted_data)
    fleiss_kappa(data)

test()

代码准确性

对于测试的3个数据,我同时使用了SPSSAU和SPSS进行统计分析,它们的结果和我代码计算结果一致。

对于我自身的42万条数据,使用SPSS和我代码计算的结果一致。

代码和测试数据已上传到github,欢迎下载和打星

Fleiss-Kappa/ at main · Lucienxhh/Fleiss-Kappa · GitHubicon-default.png?t=N4P3https://github.com/Lucienxhh/Fleiss-Kappa

本文来自互联网用户投稿,该文观点仅代表作者本人,不代表本站立场。本站仅提供信息存储空间服务,不拥有所有权,不承担相关法律责任。如若转载,请注明出处:http://www.coloradmin.cn/o/623908.html

如若内容造成侵权/违法违规/事实不符,请联系多彩编程网进行投诉反馈,一经查实,立即删除!

相关文章

如何使用切片辅助超推理 SAHI 技术对 YOLOv8 进行推理过程和代码实现

前面章节已经详细描述了 切片辅助超推理(SAHI )技术原理介绍 引入SAHI,这是一种专为小物体检测而设计的尖端流水线。SAHI 利用切片辅助推理和微调技术的力量,彻底改变了检测对象的方式。SAHI 物体检测的与众不同之处在于它与任何物体检测器的无缝集成,无需进行繁琐的微调…

Nmap安装

Nmap 文章目录 Nmap简介下载安装Zenmapnmap 配置环境变量检查是否安装成功界面 简介 Nmap是一款非常强大实用的工具,可用于检测网络上的存活主机,检测目标主机的开放端口,检测端口上相应服务上网版本,主机操作系统等信息&#xf…

基于计量学角度对传感器的灵敏度的理解和举例

基于计量学角度对传感器的灵敏度的理解和举例 灵敏度指标是考察传感器特性的主要指标之一,是对传感器设计、生产和选型过程中非常重要的参数。本文将基于计量学知识对灵敏度进行举例介绍。 一、灵敏度定义 灵敏度是传感器、测量装置或仪器的响应的变化除以对应的…

设计模式(八):结构型之装饰器模式

设计模式系列文章 设计模式(一):创建型之单例模式 设计模式(二、三):创建型之工厂方法和抽象工厂模式 设计模式(四):创建型之原型模式 设计模式(五):创建型之建造者模式 设计模式(六):结构型之代理模式 设计模式…

中国计算机学会CCF推荐的国际会议(图像处理方向)

CCF推荐的国际会议(医学图像处理方向) 1 介绍2 最新目录3 投了会议可以再投期刊吗?4 个人感想 1 介绍 CCF根据论文的质量和影响力,对国际期刊和国际会议进行了评估和分类,以便研究者在选择发表论文或参与学术交流时有…

亚马逊美国站 隐形墨水笔CPC认证 儿童玩具

隐形墨水笔 一种采用透明墨水并集成紫外线灯的笔。 隐形笔是指书写后的字迹在正常光线下不能看见(隐形)的笔,只有在紫外灯的照射下才能看见。一般用来标记物品或写不想轻易让人看见的语句。 市面上的隐形笔分为水性隐形墨水笔、油性隐形墨水…

【CSS 04】Zoro 外边距 外边距合并 内边距 内容高度与宽度 框模型 轮廓

CSS 说在前面外边距 margin外边距合并 margin_collapse内边距 padding高度与宽度 dimension框(盒子)模型 boxmodel轮廓 outline 说在前面 最近发现一个有趣的事情,就是CSDN会把我写在【】中的Zoro当做文章主要技术关键词,尽管我在…

科技政策 | 国家、广东省、深圳市的制造业单项冠军企业(产品)遴选申报指南来啦!

原创 | 文 BFT机器人 近日深圳市工业和信息化局开始了2023年技术创新项目扶持计划制造业单项冠军奖励项目申报。 那么制造业单项冠军,国家、广东省、深圳市各级申请需要什么条件呢?有那些注意事项呢?今天,一文带大家读懂。 01 申…

设计模式(九):结构型之桥接模式

设计模式系列文章 设计模式(一):创建型之单例模式 设计模式(二、三):创建型之工厂方法和抽象工厂模式 设计模式(四):创建型之原型模式 设计模式(五):创建型之建造者模式 设计模式(六):结构型之代理模式 设计模式…

赛效:如何在线给图片降噪

1:在电脑上搜索“改图鸭”网站,登录账号后在“图片修复”菜单里点击“图片降噪”。 2:点击降噪工具页面中间的,将图片添加到该页面。 3:图片添加上去后会自动降噪处理,对比原图和效果图,如果对效…

【Vue】预渲染之prerender-spa-plugin解析,方便搜索引擎的抓取

prerender-spa-plugin解析 项目背景:对于那些需要推广,希望能在百度搜索时排名靠前的网站而言,使用单页面应用的无法被抓取背景,VUE项目想SEO优化,但vue是单页面应用,不利于搜索引擎的抓取实现过程&#xf…

pwn入门(0.0)

单字节:byte 双字节:词或字word 四字节:双词Dword,doubleword 八字节:四词Qword,quadword 16字节:double quadword RBX:存储地址 RCX:计数器循环 RDX:整除取余 R…

onnx模型转TensorRT模型时出错

把onnx模型转TensorRT模型出错 错误一:WARNING: onnx2trt_utils.cpp:366: Your ONNX model has been generated with INT64 weights, while TensorRT does not natively support INT64. Attempting to cast down to INT32. 解决方法 1、安装onnx-simplifier pip…

ts 装饰器

使用装饰器前,需要把 tsconfig.json 中 experimentalDecorators 设置为 true学习了小满的B站课程:https://www.bilibili.com/video/BV1wR4y1377K?p24 前言 ts中有几种装饰器类型: 类装饰器 ClassDecorator方法装饰器 MethodDecorator参数装…

NFTScan | 05.29~06.04 NFT 市场热点汇总

欢迎来到由 NFT 基础设施 NFTScan 出品的 NFT 生态热点事件每周汇总。 周期:2023.05.29 ~ 2023.06.04 NFT Hot News:NFT 热点资讯 01/ 数据:NFT 巨鲸 Pranksy 以 52.17 枚 ETH 抛售 25 枚 Doodles 5 月 29 日,据 NFT Whale Aler…

Java + lua

luaj 主要特征 luaj 用法示例 luaj 实现原理 查找并调用指定的 Java 方法 从 Java 方法获取返回值 将 Lua function 作为参数传递给 Java 方法 在某些业务场景下,我们可能会遇到 lua 中要调用 java 代码情况,当然这个用 JNI 肯定是可以做到的&…

三十六、数学知识——组合数(递推法 + 预处理法 + 卢卡斯定理 + 分解质因数求解组合数 + 卡特兰数)

组合数算法主要内容 一、基本思路1、组合数基本概念2、递推法——询问次数多 a b 值较小 模处理(%mod)3、预处理阶乘方法——询问次数较多 a b 值很大 模处理(%mod)4、卢卡斯定理——询问次数较少 (a b 值很大&am…

泰克AFG31051信号发生器产品参数

AFG31000系列任意波函数发生器 概述 验证连接 DUT 后输出波形 InstaView? 技术用在任意波函数发生器上可直接查看连接 DUT 后的实时波形,无需使用示波器或其他设备,节省测试时间并避免因阻抗不匹配导致的实验错误。 高保真度信号与高级模式 在连续模式…

软考高级架构师笔记-6计算机系统性能评价信息系统基础知识

目录 1. 前言 & 考情分析2. 系统配置与性能评价1. 性能指标2. 性能指标3. 阿姆达尔解决方案3.信息系统基础知识1.信息系统的分类2.信息系统的生命周期3.信息系统战略规划3.常见系统介绍1.客户关系管理CRM2.供应链管理SCM3.企业应用集成EAI4.结语1. 前言 & 考情分析 前…

VR全景创业是否真的赚钱?项目真的靠谱吗?

说到创业,也许你心里会觉得这类项目对于技术和资金等要求都是比较高的,先别急着反驳,这是大多数人的心理。我们选择某个创业项目时,都需要从创业成本和盈利利润来做预算,现在告诉你有一个VR全景创业项目,几…