用R语言创建样本空间，计算事件发生概率

起初想把标题定为“大一概率习题，难倒了10年R语言大佬”，虽然主题明确，但再三考虑读者的阅读情景，还是用这个“用R语言创建样本空间，计算事件发生概率”标题吧：
1> 写代码着急找函数，标题写一些关键词，易于读者检索；
2> 本篇重点在于讲解题目的R编程实现，而不是讨论中西教育差异；
3> 不要看起来博主技术很菜。

表弟今年去英国读书，前几天发给我一份概率论习题，“姐，这个题教授让拿RMarkdown实现，求指导。”
看到题目的那一刻，I am soooooo fang（慌）！！！
这能是大一学习2个月的练习题？？？完全可以拿来考R专业证书好嘛！！！

【题目】
随机试验：从n+2张card中随机抽取，1,2,3,...,n是绿牌，n+1,n+2是红牌，取出后不放回，每张牌拿到的概率相等。
1> 抽到第一张红牌之前，所有card都放在X桌上，包括第一张红牌；
2> 接下来抽到的card放在Y桌上，直到抽到第二张红牌，第二张红牌也放在Y桌上；
3> 剩余抽到的card放在Z桌上。
下图为当n=4时，XYZ桌上card组合的一种结果示例。

【要求】
1> 使用R语言编程计算概率事件，且使用规范的数学符号；
2> 假设n≥2，给出这些问题的答案：
    a) 定义随机实验的样本空间 Ω
    b) 计算样本空间的大小 |Ω|
    c) 定义事件A（3张桌子上都至少有1张绿牌），定义事件B（X桌上至少有1张绿牌）
    d) 计算事件A、事件B发生的概率（验证：当n=4时，P(A)=1/5）

全网搜了2小时资料后，得出一个结论：这个题，“拿来主义”，它不好用了！
为啥呢？？？
1> 与这个题型相似的资料非常多，但是没有R语言解题的连贯逻辑；
2> 许多技术博客要付费才能查看全文，且提及的函数不适用。
想想自己写一个吧。

接下来分3个部分来说明题目的实现过程：
1. 解题逻辑&函数
2. Rcode示例
3. RMarkdown交作业

一、解题逻辑&函数

1. 解题逻辑

假设：n=4，有n+2张牌，即6张牌，1,2,3,4是绿牌，5,6是红牌。
1> 创建样本空间大矩阵Ω
矩阵Ω的行数为随机试验所有可能结果的总数，列数为card数量；
当取出第一张card时有6种可能，第二张有5种可能，...，第六张有1种可能，也就是总共有 6*5*4*3*2*1=6! 种card排列组合结果，即行数为6! ；
第1列记录第一张card值，第2列记录第二张card值，以此类推，总共六张card，即列数为6；
2> 查找每个事件中两张红牌的位置；
3> 根据红牌位置将随机事件结果赋值给X、Y、Z；
4> 计算每个事件中X、Y、Z分别包含几张绿牌；
5> 定义目标事件&计算概率
图中所示为步骤1-2-3的矩阵示例。

2. R语言相关函数

对于上述解题步骤，涉及到这些R语言操作以及函数使用：
1> 计算N的阶乘：prod(1:N)
2> 生成1:N所有排列组合，调用gtools包中的函数：permutations()
（注意：这里不可以使用sample()函数，sample()只能用来模拟随机事件，不能生成样本空间。）
3> 字符串拼接：paste()
4> 在字符串中查找指定字符的位置，调用stringr包中的函数：str_locate()
5> 字符串截取：str_sub()
6> 计算字符串中指定字符的个数：str_count()
7> 对数据框筛选指定条件的行：subset()

二、Rcode示例

【1. 创建样本空间Ω】

n<-4# 定义n
N<-n+2 # card数量
cards<-1:N # 牌中包含的数字
library(gtools)
Omega<-permutations(n=N, r=N, v=cards,repeats.allowed=FALSE)# 全排列组合矩阵，不允许重复
Omega<-as.data.frame(Omega)# 矩阵转为数据框格式

【2. 查找2张红牌位置】
分成3段code来实现

# 1. 对每一行生成一个拼接字符串
for (i in 1:prod(1:N)) {
      Omega[i,N+1] <- do.call(paste,c(Omega[i,1:N],sep=''))
 }

# 2. 计算指定字符：n+1、n+2所在位置
library(stringi)
library(stringr)
n1<-as.character(n+1)
n2<-as.character(n+2)
Omega$label1<-str_locate(Omega$V7,n1)# n+1的位置
Omega$label2<-str_locate(Omega$V7,n2)# n+2的位置

# 3. 计算第一次出现红牌、第二次出现红牌的位置
for (i in 1:prod(1:N)) {
 Omega$label_1st_red[i]<-min(Omega$label1[i],Omega$label2[i])# 第一张红牌的位置
 Omega$label_sec_red[i]<-max(Omega$label1[i],Omega$label2[i])# 第二张红牌的位置
 }

【3. 将随机事件结果赋值给X、Y、Z】

# 根据红牌位置截取字符串，生成X、Y、Z的数据集
Omega$table_X<-str_sub(Omega$V7,1,Omega$label_1st_red)
Omega$table_Y<-str_sub(Omega$V7,Omega$label_1st_red+1,Omega$label_sec_red)
Omega$table_Z<-str_sub(Omega$V7,Omega$label_sec_red+1,N)

【4. 计算每个事件中X、Y、Z分别包含几张绿牌】

for (i in 1:prod(1:N)) {
 Omega$table_X_green_count[i]<-sum(str_count(Omega$table_X[i],as.character(c(1:n))))
 Omega$table_Y_green_count[i]<-sum(str_count(Omega$table_Y[i],as.character(c(1:n))))
 Omega$table_Z_green_count[i]<-sum(str_count(Omega$table_Z[i],as.character(c(1:n))))
 }

【5. 定义目标事件&计算概率】

# 事件A：每张桌子都至少有1张绿牌
# 事件B：X桌子上至少有1张绿牌
event_A<-(Omega$table_X_green_count>=1 & Omega$table_Y_green_count>=1 & Omega$table_Z_green_count>=1)
event_B<-(Omega$table_X_green_count>=1)

# 计算事件概率
prob_event_A<-nrow(subset(Omega,event_A))/prod(1:N)
prob_event_B<-nrow(subset(Omega,event_B))/prod(1:N)

# 输出结果
head(Omega)
print(paste0('|Ω|：',prod(1:N)))
print(paste0('P(A)：',prob_event_A))
print(paste0('P(B)：',prob_event_B))

三、RMarkdown交作业

将Rcode嵌入RMarkdown原文件的代码框中，可以将程序结果输出为HTML、docx、PDF等文件，这部分对rmd文件格式做个简单介绍。

rmd代码主要包括3个模块：
1. 文件信息设置，如：标题、作者、输出文件的格式；
2. 全局设置，如：是否在文件中显示代码部分、警告信息是否写入结果文件；
3. Rcode部分，即代码主体。

1. 文件信息设置

---
title: "用R语言创建样本空间，计算事件发生概率"
author: "千行"
date: "2023/11/11"
output:
    html_document: default
---

2. 全局设置

```{r setup, include=FALSE}
knitr::opts_chunk$set(echo = FALSE)# 文件中不显示代码行
```

3. Rcode部分

```{r results='hold',warning=FALSE}
n<-4# 定义n
N<-n+2 # card数量
cards<-1:N # 牌中包含的数字

library(gtools)
Omega<-permutations(n=N, r=N, v=cards,repeats.allowed=FALSE)# 全排列组合矩阵，不允许重复
Omega<-as.data.frame(Omega)# 矩阵转为数据框格式
head(Omega,10)
```

把R代码主体替换到RMarkdown的Rcode中，点击Knit，就可以交作业啦~

文末：
1> 这个题最难的部分在于样本空间的创建，博主真的是不知道如何快速生成全排列组合，就先用sample()模拟了一个大样本矩阵；
2> 从逻辑到code实现总共用了6小时，后来问了ChatGPT，它1分钟内给出了结果，且代码逻辑精炼、计算过程节省内存，最惊喜的是它检索出了全排列组合函数permutations()，帮助解决了第一步的问题，真棒！
3> 当文心一言慢悠悠走来，让我觉得AI还很远的时候，ChatGPT出现就是一道闪电。