从0开始学习R语言--Day14--贝叶斯统计与结构方程模型

news2025/6/5 10:37:24

贝叶斯统计

在很多时候,我们经常会看到在统计分析中出现很多反直觉的结论,比如假如有一种病,人群中的患病率为1%,患者真患病时,检测结果为阳性的概率是99%,如果没有,则检测结果为阳性的概率是5%,但是如果让我们计算一个人检测结果为阳性时真正患病的概率为多少,根据概率公式P(患病|检测为阳性) = [P(检测为阳性|患病)× P(患病)]/P(检测为阳性) 我们可以计算,而P(检测为阳性)可以通过全概率公式P(检测为阳性)=P(检测为阳性|患病)P(患病) + P(检测为阳性|没患病)P(没患病)得出,计算得到的结果为16.67%,可能有人会疑惑,为什么患真病的前提下,检测为阳性时99%,反过来就只有16.67%,其实会有这样的错句,一是因为没抓住基础的概率,即本身该病的患病率就是很低的1%,二是在人群中没病但检测为阳性的人(5%×99%≈4.95%)比有病且检测为阳性的人(1%×99%≈0.99%)多很多。

这种先有对数据样本的初始判断,观察数据在给定假设下的表现,通过后续检验后得出新判断的过程,就是贝叶斯统计,其中,初始的判断叫做先验概率(Prior),即人群中的患病率1%;在某个假设成立的情况下,观察到的表现叫做似然值(Likehood),即患者有病时检测为阳性的概率和没病时检测为阳性的概率,这里是否有病就是假设;在看到计算结果后得出的新判断叫做后验概率(Posterior)。

下面我们依然是给出一个例子来解释:

# 加载必要的库
library(tidyverse)
library(ggplot2)

# 定义参数
p_disease <- 0.01       # 患病先验概率(1%)
p_true_positive <- 0.99 # 真阳性率(99%)
p_false_positive <- 0.05 # 假阳性率(5%)

# 计算检测为阳性的总概率
p_positive <- (p_true_positive * p_disease) + (p_false_positive * (1 - p_disease))

# 计算后验概率(贝叶斯定理)
p_disease_given_positive <- (p_true_positive * p_disease) / p_positive

# 打印结果
cat(sprintf("检测为阳性时实际患病的概率: %.2f%%", p_disease_given_positive * 100))

# 可视化先验、似然和后验
results <- data.frame(
  Scenario = c("Prior", "Likelihood (Disease)", "Likelihood (No Disease)", "Posterior"),
  Probability = c(p_disease, p_true_positive, p_false_positive, p_disease_given_positive)
)

ggplot(results, aes(x = Scenario, y = Probability, fill = Scenario)) +
  geom_bar(stat = "identity") +
  geom_text(aes(label = scales::percent(Probability)), vjust = -0.5) +
  scale_y_continuous(labels = scales::percent) +
  labs(title = "贝叶斯分析: 疾病检测案例",
       subtitle = "先验、似然和后验概率比较",
       y = "概率") +
  theme_minimal()

# 更详细的模拟分析
set.seed(123)
population_size <- 10000

# 生成模拟人群
sim_data <- data.frame(
  id = 1:population_size,
  disease = rbinom(population_size, 1, p_disease),
  test_result = NA
)

# 模拟检测结果
sim_data <- sim_data %>%
  mutate(
    test_result = case_when(
      disease == 1 ~ rbinom(n(), 1, p_true_positive),
      disease == 0 ~ rbinom(n(), 1, p_false_positive)
    )
  )

# 计算实际阳性中患病的比例
empirical_posterior <- sim_data %>%
  filter(test_result == 1) %>%
  summarise(
    p_disease_given_positive = mean(disease)
  ) %>%
  pull()

cat(sprintf("\n模拟数据中阳性结果实际患病的比例: %.2f%%", empirical_posterior * 100))

# 创建列联表展示结果
table(sim_data$disease, sim_data$test_result) %>%
  addmargins() %>%
  knitr::kable(caption = "疾病与检测结果的列联表")

输出:

检测为阳性时实际患病的概率: 16.67%> 
模拟数据中阳性结果实际患病的比例: 17.09%> 
Table: 疾病与检测结果的列联表

|    |    0|   1|   Sum|
|:---|----:|---:|-----:|
|0   | 9414| 485|  9899|
|1   |    1| 100|   101|
|Sum | 9415| 585| 10000|

从输出中我们可以看到,模拟数据中的实际患病率与我们用贝叶斯统计计算得出的差不多,而输出的列表里也能观察到阳性预测的正确率只有100/585≈17.09%,这说明在实际上的应用场景中,我们在看到类似于“99%的准确率”时不要被误导,要验证过基础概率后再下结论,相对的,假如在信用评分模型中,用户的违约率很低,那么如果出现很多高风险决策,这也可能是误报,当然最好是要查看背后的解释。

结构方程模型

事实上,并不是所有的数据都会跟我们平常做练习和书上的一样是每个都很好理解的,很多时候我们想研究的变量,是分散在不同的数据上的,这个时候探讨出变量间的关系,就需要我们潜在的变量了,是这些变量在无形之中担当了变量之间的桥梁。

例如我们知道不同的老师对学生成绩有不同的影响,也许我们会去分析不同的教学方法是怎么影响成绩的,但如果直接对这两者做回归分析,很多时候效果都会很差,究其原因是因为教学方法很可能是间接影响了B,再由B去影响学生成绩,所以这里的重点工作就变成了如何找出B。当然,这种潜在变量并不是唯一的,需要我们提出假设去验证。

但往往潜在变量都是抽象的(毕竟明显的变量在拿到数据的时候,就会看出来),这个时候就需要我们自己人为设计一些测量的工具来代替潜在变量的影响,假设我们认为潜在变量是教学质量,那么我们就可以从教师评分,课程设计、作业反馈上得出来(前提是控制同一批能力差不多的学生),通过设计问题和学生对知识掌握过程中的进步,我们可以得出该教师方法的质量处于什么水平;亦或者认为潜在变量是学生的学习动机,那么可能就需要用学生平时的学习时间,在课堂上的互动程度,对自我的评价,通过学生对学习的感兴趣程度和是否积极的态度来判断。

依然是举一个例子来说明:

library(lavaan)    # 结构方程建模
library(tidySEM)   # 可视化
library(ggplot2)   # 基础绘图
library(dplyr)


set.seed(123)
n_students <- 300

# 生成潜变量(不可直接观测的真实值)
true_motivation <- rnorm(n_students, mean = 0, sd = 1)
true_teaching <- rnorm(n_students, mean = 0, sd = 1)

# 生成观测指标(带测量误差)
student_data <- data.frame(
  # 学习动机的观测指标
  study_time = round(2 + 0.5*true_motivation + rnorm(n_students, 0, 0.3), 1),
  participation = round(3 + 0.7*true_motivation + rnorm(n_students, 0, 1)),
  self_eval = round(5 + 0.6*true_motivation + rnorm(n_students, 0, 1)),
  
  # 教学质量的观测指标
  teach_rating = round(6 + 0.8*true_teaching + rnorm(n_students, 0, 1)),
  course_diff = round(4 - 0.5*true_teaching + rnorm(n_students, 0, 1)), # 难度与质量负相关
  hw_feedback = round(5 + 0.7*true_teaching + rnorm(n_students, 0, 1)),
  
  # 成绩变量(受潜变量影响)
  math_score = round(75 + 5*true_motivation + 3*true_teaching + rnorm(n_students, 0, 5)),
  chinese_score = round(78 + 6*true_motivation + rnorm(n_students, 0, 4))
)


model <- '
  # 测量模型 (Latent Variable Definitions)
  motivation =~ 1*study_time + participation + self_eval  # 固定第一个载荷为1
  teaching_quality =~ 1*teach_rating + course_diff + hw_feedback
  
  # 结构模型 (Regressions)
  math_score ~ motivation + teaching_quality
  chinese_score ~ motivation
  
  # 潜变量间关系
  motivation ~ teaching_quality
  
  # 残差相关(可选)
  # study_time ~~ participation
'

fit <- sem(model, 
           data = student_data,
           estimator = "MLR")  # 使用稳健最大似然估计


summary(fit, standardized = TRUE, fit.measures = TRUE)


# 路径图
graph_sem(fit, 
          layout = get_layout(
            "teaching_quality", "motivation",
            "math_score", "chinese_score",
            rows = 2),
          edge.label.cex = 0.9,
          sizeMan = 6, sizeLat = 8,
          curvature = 60)

# 因子载荷可视化(ggplot2版)
parameterEstimates(fit, standardized = TRUE) %>%
  filter(op == "=~") %>%
  ggplot(aes(x = rhs, y = std.all, fill = lhs)) +
  geom_col(position = "dodge") +
  geom_hline(yintercept = 0.5, linetype = "dashed", color = "red") +
  labs(title = "标准化因子载荷",
       y = "标准化系数", x = "观测指标") +
  coord_flip() +
  theme_minimal()

输出:

lavaan 0.6-19 ended normally after 77 iterations

  Estimator                                         ML
  Optimization method                           NLMINB
  Number of model parameters                        19

  Number of observations                           300

Model Test User Model:
                                              Standard      Scaled
  Test Statistic                                20.142      20.875
  Degrees of freedom                                17          17
  P-value (Chi-square)                           0.267       0.232
  Scaling correction factor                                  0.965
    Yuan-Bentler correction (Mplus variant)                       

Model Test Baseline Model:

  Test statistic                               544.259     551.226
  Degrees of freedom                                28          28
  P-value                                        0.000       0.000
  Scaling correction factor                                  0.987

User Model versus Baseline Model:

  Comparative Fit Index (CFI)                    0.994       0.993
  Tucker-Lewis Index (TLI)                       0.990       0.988
                                                                  
  Robust Comparative Fit Index (CFI)                         0.993
  Robust Tucker-Lewis Index (TLI)                            0.988

Loglikelihood and Information Criteria:

  Loglikelihood user model (H0)              -4424.677   -4424.677
  Scaling correction factor                                  1.023
      for the MLR correction                                      
  Loglikelihood unrestricted model (H1)      -4414.606   -4414.606
  Scaling correction factor                                  0.995
      for the MLR correction                                      
                                                                  
  Akaike (AIC)                                8887.354    8887.354
  Bayesian (BIC)                              8957.726    8957.726
  Sample-size adjusted Bayesian (SABIC)       8897.470    8897.470

Root Mean Square Error of Approximation:

  RMSEA                                          0.025       0.028
  90 Percent confidence interval - lower         0.000       0.000
  90 Percent confidence interval - upper         0.060       0.063
  P-value H_0: RMSEA <= 0.050                    0.857       0.828
  P-value H_0: RMSEA >= 0.080                    0.002       0.004
                                                                  
  Robust RMSEA                                               0.027
  90 Percent confidence interval - lower                     0.000
  90 Percent confidence interval - upper                     0.061
  P-value H_0: Robust RMSEA <= 0.050                         0.845
  P-value H_0: Robust RMSEA >= 0.080                         0.003

Standardized Root Mean Square Residual:

  SRMR                                           0.038       0.038

Parameter Estimates:

  Standard errors                             Sandwich
  Information bread                           Observed
  Observed information based on                Hessian

Latent Variables:
                      Estimate  Std.Err  z-value  P(>|z|)   Std.lv  Std.all
  motivation =~                                                            
    study_time           1.000                               0.473    0.859
    participation        1.435    0.159    9.022    0.000    0.679    0.547
    self_eval            1.013    0.159    6.355    0.000    0.480    0.435
  teaching_quality =~                                                      
    teach_rating         1.000                               0.812    0.618
    course_diff         -0.680    0.135   -5.046    0.000   -0.552   -0.473
    hw_feedback          0.840    0.174    4.838    0.000    0.683    0.580

Regressions:
                   Estimate  Std.Err  z-value  P(>|z|)   Std.lv  Std.all
  math_score ~                                                          
    motivation        9.753    1.191    8.190    0.000    4.616    0.623
    teaching_qulty    3.792    0.726    5.224    0.000    3.080    0.415
  chinese_score ~                                                       
    motivation       11.243    1.080   10.411    0.000    5.322    0.750
  motivation ~                                                          
    teaching_qulty   -0.073    0.055   -1.333    0.183   -0.125   -0.125

Covariances:
                   Estimate  Std.Err  z-value  P(>|z|)   Std.lv  Std.all
 .math_score ~~                                                         
   .chinese_score     3.446    2.361    1.460    0.144    3.446    0.139

Variances:
                   Estimate  Std.Err  z-value  P(>|z|)   Std.lv  Std.all
   .study_time        0.079    0.017    4.801    0.000    0.079    0.262
   .participation     1.079    0.099   10.941    0.000    1.079    0.701
   .self_eval         0.986    0.085   11.611    0.000    0.986    0.811
   .teach_rating      1.070    0.146    7.344    0.000    1.070    0.619
   .course_diff       1.062    0.098   10.879    0.000    1.062    0.777
   .hw_feedback       0.918    0.111    8.297    0.000    0.918    0.663
   .math_score       27.738    3.708    7.480    0.000   27.738    0.504
   .chinese_score    22.022    2.557    8.613    0.000   22.022    0.437
   .motivation        0.221    0.028    7.806    0.000    0.984    0.984
    teaching_qulty    0.660    0.168    3.935    0.000    1.000    1.000

从输出中我们可以看到,CFI,RMSEA等指标都在标准范围内,说明模型与数据的匹配性较高;观察学习动机和教学质量中的指标,可以发现学习时间最能反映学习动机,而自我评价对模型的影响则较小;观察教学质量中的指标,可以发现课程难度与教学质量是负相关的,这也在我们的预期之内,毕竟越难的课,学生听懂的难度也会加大,而不是好的老师就能一下子把学生教会。综合来看,学习动机对语文和数学的成绩都较大,而教学质量会直接影响数学成绩,而对学习动机的影响反倒不大((β=-0.125, p=0.183))。

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

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

相关文章

[Python] 如何使用 Python 调用 Dify 工作流服务实现自动化翻译

在实际项目中,自动化工作流服务可以大大简化复杂任务的处理流程。本文将介绍如何通过 Python 脚本调用 Dify 提供的工作流 API,实现文本翻译的自动化操作。该流程包括设置 API 接口、构造请求体并处理返回结果。 一、背景介绍:什么是 Dify 工作流服务? Dify 是一款支持多种…

PTA-根据已有类Worker,使用LinkedList编写一个WorkerList类,实现计算所有工人总工资的功能。

目录 1.问题描述 2.函数接口定义&#xff1a; 3.裁判测试程序样例&#xff1a; 4.输入和输出样例 输入样例&#xff1a; 输出样例&#xff1a; 5.实现代码 1.问题描述 Main类&#xff1a;在main方法中&#xff0c;调用constructWorkerList方法构建一个Worker对象链表…

微软markitdown PDF/WORD/HTML文档转Markdown格式软件整合包下载

本次和大家分享另一个微软发布的非常热门的文件文档转Markdown格式文档的软件markitdown&#xff0c;软件可以将PDF&#xff0c;word&#xff0c;ppt&#xff0c;Excel等十几种格式文档转换为markdown格式文档&#xff0c;我基于当前最新0.1.2版本制作了免安装一键启动整合包。…

BayesFlow:基于神经网络的摊销贝叶斯推断框架

贝叶斯推断为不确定性条件下的推理、复杂系统建模以及基于观测数据的预测提供了严谨且功能强大的理论框架。尽管贝叶斯建模在理论上具有优雅性&#xff0c;但在实际应用中经常面临显著的计算挑战&#xff1a;后验分布通常缺乏解析解&#xff0c;模型验证和比较需要进行重复的推…

基于FPGA的DES加解密系统verilog实现,包含testbench和开发板硬件测试

目录 1.课题概述 2.系统测试效果 3.核心程序与模型 4.系统原理简介 5.完整工程文件 1.课题概述 基于FPGA的DES加解密系统verilog实现,包含testbench和开发板硬件测试。输入待加密数据&#xff0c;密钥&#xff0c;输出加密数据&#xff0c;然后通过解密模块输出解密后的原…

Python----目标检测(《用于精确目标检测和语义分割的丰富特征层次结构》和R-CNN)

一、《用于精确目标检测和语义分割的丰富特征层次结构》 1.1、基本信息 原文标题&#xff1a;Rich feature hierarchies for accurate object detection and semantic segmentation 中文译名&#xff1a;用于精确目标检测与语义分割的丰富特征层次结构 版本&#xff1a;第5版技…

极简以太彩光网络解决方案4.0正式发布,“彩光”重构园区网络极简之道

5月28日下午,锐捷网络在京举办以“光,本该如此‘简单’”为主题的发布会,正式发布极简以太彩光网络解决方案4.0。作为“彩光”方案的全新进化版本,极简以太彩光4.0从用户需求出发,聚焦场景洞察,开启了一场从底层基因出发的极简革命,通过架构、部署、运维等多维度的创新升级,以强…

国芯思辰| 霍尔电流传感器AH811为蓄电池负载检测系统安全护航

在电动车、储能电站、不间断电源&#xff08;UPS&#xff09;等设备中&#xff0c;蓄电池作为关键的储能单元&#xff0c;其运行状态直接关系到设备的稳定性和使用寿命。而准确监测蓄电池的负载情况&#xff0c;是保障其安全、高效运行的关键。霍尔电流传感器 AH811凭借独特的技…

TortoiseSVN账号切换

SVN登录配置及账号切换 本文主要为了解答svn客户端如何进行账号登录及切换不同权限账号的方式。 一、环境准备与客户端安装 安装TortoiseSVN客户端 ​​下载地址​​&#xff1a;TortoiseSVN官网 ​​安装步骤​​&#xff1a; 双击安装包&#xff0c;按向导完成安装后&#x…

2025年05月28日Github流行趋势

项目名称&#xff1a;agenticSeek 项目地址url&#xff1a;https://github.com/Fosowl/agenticSeek项目语言&#xff1a;Python历史star数&#xff1a;10352今日star数&#xff1a;2444项目维护者&#xff1a;Fosowl, steveh8758, klimentij, ganeshnikhil, apps/copilot-pull-…

篇章五 数据结构——链表(一)

目录 1.ArrayList的缺陷 2. 链表 2.1 链表的概念及结构 2.2 链表结构 1. 单向或者双向 2.带头或者不带头 3.循环或者非循环 2.3 链表的实现 1.完整代码 2.图解 3.显示方法 4.链表大小 5. 链表是否存在 key 值 6.头插法 7.尾插法 8.中间插入 9.删除key值节点 10.…

一文清晰理解目标检测指标计算

一、核心概念 1.交并比IoU 预测边界框与真实边界框区域的重叠比&#xff0c;取值范围为[0,1] 设预测边界框为&#xff0c;真实边界框为 公式&#xff1a; IoU计算为两个边界框交集面积与并集面积之比&#xff0c;图示如下 IoU值越高&#xff0c;表示预测边界框与真实边界框的对…

Artificial Analysis2025年Q1人工智能发展六大趋势总结

2025年第一季度人工智能发展六大趋势总结 ——基于《Artificial Analysis 2025年Q1人工智能报告》 趋势一&#xff1a;AI持续进步&#xff0c;竞争格局白热化 前沿模型竞争加剧&#xff1a;OpenAI凭借“o4-mini&#xff08;高智能版&#xff09;”保持领先&#xff0c;但谷歌&…

高效管理 Python 项目的 UV 工具指南

&#x1f49d;&#x1f49d;&#x1f49d;欢迎莅临我的博客&#xff0c;很高兴能够在这里和您见面&#xff01;希望您在这里可以感受到一份轻松愉快的氛围&#xff0c;不仅可以获得有趣的内容和知识&#xff0c;也可以畅所欲言、分享您的想法和见解。 持续学习&#xff0c;不断…

初识vue3(vue简介,环境配置,setup语法糖)

一&#xff0c;前言 今天学习vue3 二&#xff0c;vue简介及如何创建vue工程 Vue 3 简介 Vue.js&#xff08;读音 /vjuː/&#xff0c;类似 “view”&#xff09;是一款流行的渐进式 JavaScript 框架&#xff0c;用于构建用户界面。Vue 3 是其第三代主要版本&#xff0c;于 …

LeetCode-链表操作题目

虚拟头指针&#xff0c;在当前head的前面建立一个虚拟头指针&#xff0c;然后哪怕当前的head的val等于提供的val也能进行统一操作 203移除链表元素简单题 /*** Definition for singly-linked list.* public class ListNode {* int val;* ListNode next;* ListNode(…

【ARM】MDK浏览信息的生成对于构建时间的影响

1、 文档目标 用于了解MDK的代码浏览信息的生成对于工程的构建是否会产生影响。 2、 问题场景 客户在MDK中使用Compiler 5对于工程进行构建过程中发现&#xff0c;对于是否产生浏览信息会对于构建时间产生一定的影响。在Options中Output栏中勾选了Browse Information后&#…

py爬虫的话,selenium是不是能完全取代requests?

selenium适合动态网页抓取&#xff0c;因为它可以控制浏览器去点击、加载网页&#xff0c;requests则比较适合静态网页采集&#xff0c;它非常轻量化速度快&#xff0c;没有浏览器开销&#xff0c;占用资源少。当然如果不考虑资源占用和速度&#xff0c;selenium是可以替代requ…

docker B站学习

镜像是一个只读的模板&#xff0c;用来创建容器 容器是docker的运行实例&#xff0c;提供了独立可移植的环境 https://www.bilibili.com/video/BV11L411g7U1?spm_id_from333.788.videopod.episodes&vd_sourcee60c804914459274157197c4388a4d2f&p3 目录挂载 尚硅谷doc…

SpringBoot高校宿舍信息管理系统小程序

概述 基于SpringBoot的高校宿舍信息管理系统小程序项目&#xff0c;这是一款非常适合高校使用的信息化管理工具。该系统包含了完整的宿舍管理功能模块&#xff0c;采用主流技术栈开发&#xff0c;代码结构清晰&#xff0c;非常适合学习和二次开发。 主要内容 这个宿舍管理系…