从箱线图到统计指标表

生信技能树

发布于 2023-09-04 15:37:59

3370

发布于 2023-09-04 15:37:59

文章被收录于专栏：生信技能树

最近有粉丝提问到如何从表达量差异分析后的某个基因或者蛋白质或者其它元素在两个分组的差异情况的箱线图到其相关的一系列统计指标表，出处是2023年4月的一个文章：《Saliva biopsy: Detecting the difference of EBV DNA methylation in the diagnosis of nasopharyngeal carcinoma》，如下所示的箱线图是为了说明肿瘤病人队列里面的这5个指标都是统计学显著的高于正常人队列的：

箱线图

关于箱线图，大家都不陌生了，可以问一下chatGPT它的定义，概念，以及应用。

在差异表达基因分析后，我们通常会选择一些显著差异表达的基因进行进一步的可视化分析，例如箱线图。箱线图是一种用于显示一组数据分散情况资料的统计图，包括最大值、最小值、中位数、上四分位数（Q3，75th percentile）和下四分位数（Q1，25th percentile）。

在基因表达分析中，箱线图可以用来显示不同组（例如，疾病组和对照组）中基因的表达水平。箱线图的中位线表示基因在每个组中的中位表达水平，箱子的上下边界表示上四分位数和下四分位数，即表达水平的分布范围。箱线图的“须”（即线的部分）则表示数据的整体分布范围，通常定义为1.5倍的四分位距（IQR，即Q3-Q1），超过这个范围的点通常被视为异常值。在比较两组数据时，我们通常关注的是中位数（即箱线图的中线）是否有显著差异，以及数据的分布是否不同。如果两个组的箱线图有显著的差异，这可能表明基因在这两个组中的表达有显著差异。

为了确定差异是否显著，我们通常会进行统计测试，例如t检验或Mann-Whitney U检验。这些测试可以提供一个p值，用于量化观察到的差异是否可能仅仅是由随机变异引起的。如果p值小于某个阈值（例如0.05），我们通常会认为差异是显著的。然而，需要注意的是，这些统计测试假设数据是独立的，且在t检验的情况下，还假设数据是正态分布的。在实际应用中，这些假设可能不成立，因此在解释结果时需要谨慎。

其实从箱线图到如上所示的各种统计指标表，比较陌生的是最后两列统计学指标（Sensitivity和Specificity）而已。前面的两个分组的内部的样品数量，每个分组的平均值中位值最大值最小值都很容易计算，最后两列统计学指标其实涉及到了一个机器学习的概念了，同样的，可以先问一下chatGPT它的定义，概念，以及应用。

敏感性（Sensitivity）和特异性（Specificity）是用于评估医学测试或生物信息学预测模型性能的两个重要参数。

敏感性：也被称为真阳性率，是指测试对疾病阳性个体的识别能力。换句话说，它是测试正确识别出的阳性结果（真阳性）占所有实际阳性样本（真阳性+假阴性）的比例。如果一个测试的敏感性很高，那么它错过真正的阳性结果（即产生假阴性结果）的可能性就很小。
特异性：也被称为真阴性率，是指测试对疾病阴性个体的识别能力。换句话说，它是测试正确识别出的阴性结果（真阴性）占所有实际阴性样本（真阴性+假阳性）的比例。如果一个测试的特异性很高，那么它错误地将阴性样本识别为阳性（即产生假阳性结果）的可能性就很小。

在理想情况下，我们希望一个测试的敏感性和特异性都能达到100%，但在实际情况中，这两者往往需要进行权衡。例如，如果我们希望尽可能少地错过阳性结果（即提高敏感性），那么我们可能需要接受更多的假阳性结果（即降低特异性）。反之亦然。因此，选择哪种测试或预测模型，以及如何设定阈值，往往需要根据具体的应用场景和接受的错误类型进行决定。

这两个统计学指标（Sensitivity和Specificity）其实绝大部分情况下会通过ROC曲线（Receiver Operating Characteristic curve）和AUC（Area Under the Curve），这两个用于评估分类模型性能的重要工具体现：

ROC曲线：ROC曲线是一种图形化表示分类模型性能的工具。在ROC曲线中，横轴表示假阳性率（False Positive Rate，FPR），纵轴表示真阳性率（True Positive Rate，TPR，也就是敏感性）。ROC曲线通过改变分类阈值来描绘出所有可能的TPR和FPR组合，从而展示出模型在所有阈值下的性能。ROC曲线的左上角（即FPR=0，TPR=1）表示理想的分类模型，该模型能够在没有产生任何假阳性的情况下，完全正确地识别出所有的阳性样本。
AUC：AUC是ROC曲线下的面积，用于量化分类模型的整体性能。AUC的值介于0和1之间，值越接近1，表示模型的性能越好。具体来说，AUC可以被解释为：在随机选取一个阳性样本和一个阴性样本的情况下，分类模型将阳性样本的得分排在阴性样本之前的概率。因此，AUC为0.5表示模型的性能等同于随机猜测，而AUC为1表示模型在所有情况下都能完美地区分阳性样本和阴性样本。

需要注意的是，虽然AUC是一个有用的性能指标，但它并不能完全反映出模型的性能，特别是在类别不平衡的情况下。因此，在评估模型性能时，通常需要结合其他指标，如精确率、召回率、F1分数等，来进行全面的评估。

起码从R的角度来说，箱线图直接到ROC曲线，顺便计算得到AUC值是很容易的。同样的，我也是让chatGPT做了一下：使用R代码举例一个差异分析，并且绘制ROC曲线和表达量差异箱线图

以下是一个使用R进行差异分析、绘制ROC曲线和箱线图的示例。这个示例使用了pROC包进行ROC分析和绘图，使用ggplot2包进行箱线图的绘制。注意，这只是一个示例，实际的分析可能需要根据你的数据和问题进行调整。

# 安装必要的包
install.packages(c("pROC", "ggplot2"))

# 加载必要的包
library(pROC)
library(ggplot2)

# 假设我们有一个数据框，其中包含两列：表达量（expression）和分组（group）
data <- data.frame(
  expression = rnorm(200),  # 这里我们只是生成了一些随机数据
  group = rep(c("A", "B"), each = 100)
)

# 进行t检验，得到p值
t_test_result <- t.test(expression ~ group, data = data)
p_value <- t_test_result$p.value

# 打印p值
print(paste("p-value:", p_value))

# 绘制箱线图
ggplot(data, aes(x = group, y = expression)) +
  geom_boxplot() +
  ggtitle("Boxplot of expression by group")

# 假设我们有一些真实的类别标签和一些预测的概率
labels <- sample(c(0, 1), 200, replace = TRUE)  # 这里我们只是生成了一些随机数据
predictions <- runif(200)  # 这里我们只是生成了一些随机数据

# 计算ROC曲线
roc_result <- roc(labels, predictions)

# 打印AUC
print(paste("AUC:", auc(roc_result)))

# 绘制ROC曲线
plot(roc_result, main="ROC curve")

这个脚本首先进行了一个t检验来检查两个组的表达量是否有显著差异，然后绘制了一个箱线图来可视化这些差异。然后，它计算了一个ROC曲线，并打印了AUC值，最后绘制了ROC曲线。这只是一个基本的示例，实际的分析可能需要更复杂的统计测试和更复杂的图形。

让我非常吃惊，chatGPT居然给出来的案例是可以运行的，但是它的问题在于，没有体现出来差异，因为是同样的分布函数得到的数据分成两组，本来就不应该是有差异，后面的一切计算都看不出来差异。而且前面的表达量差异分析和后面的ROC曲线没有关系。。。。

其实R语言提供了一系列的函数来处理各种统计分布，包括正态分布、二项分布、泊松分布等。这些函数通常有四种形式，分别用于生成密度函数（d）、累积分布函数（p）、生成随机变量（r）和分位数函数（q）。以下是一些常见的统计分布函数及其使用方法：

正态分布（Normal Distribution）：
- dnorm(x, mean = 0, sd = 1)：正态分布的密度函数。
- pnorm(q, mean = 0, sd = 1)：正态分布的累积分布函数。
- qnorm(p, mean = 0, sd = 1)：正态分布的分位数函数。
- rnorm(n, mean = 0, sd = 1)：生成正态分布的随机变量。
二项分布（Binomial Distribution）：
- dbinom(x, size, prob)：二项分布的密度函数。
- pbinom(q, size, prob)：二项分布的累积分布函数。
- qbinom(p, size, prob)：二项分布的分位数函数。
- rbinom(n, size, prob)：生成二项分布的随机变量。
泊松分布（Poisson Distribution）：
- dpois(x, lambda)：泊松分布的密度函数。
- ppois(q, lambda)：泊松分布的累积分布函数。
- qpois(p, lambda)：泊松分布的分位数函数。
- rpois(n, lambda)：生成泊松分布的随机变量。

在上述函数中，x 和 q 是向量，n 是要生成的随机变量的数量，p 是概率，mean 和 sd 分别是正态分布的均值和标准差，size 和 prob 分别是二项分布的试验次数和成功概率，lambda 是泊松分布的参数。这些函数可以用于各种统计分析和模拟实验，是R语言中非常重要的工具。

一个比较好的案例

为了体现两个分组的数值差异，我这里人为的创作了两个不同统计学分布的各自100个分组，然后去比较差异，所以很明显他们应该是统计学显著的。

set.seed(123) 
library(pROC)
library(ggplot2) 
data <- data.frame(
  expression =c( rnorm(100,-1), rpois(100,1)),  # 这里我们只是生成了一些随机数据
  group = rep(c("A", "B"), each = 100)
) 
t_test_result <- t.test(expression ~ group, data = data)
p_value <- t_test_result$p.value  
p1 = ggplot(data, aes(x = group, y = expression)) +
  geom_boxplot() + 
  ggtitle(paste("p-value:", p_value)) +theme_bw()
 
# 根据表达量计算ROC曲线
roc_result <- pROC::roc(data$group, data$expression)
# 打印AUC
print(paste("AUC:", auc(roc_result)))
# 绘制ROC曲线
p2=pROC::ggroc(roc_result)+ggtitle(auc(roc_result))+theme_bw()
library(patchwork)
p1+p2

如下所示：

统计学显著的

当然了，前面的文章：《Saliva biopsy: Detecting the difference of EBV DNA methylation in the diagnosis of nasopharyngeal carcinoma》，给出来了两列统计学指标（Sensitivity和Specificity），所以还需要继续计算；

coords <- coords(roc_result, "best", ret = c("threshold", "sensitivity", "specificity"))
cut_off <- coords[1]
coords
###    threshold     sensitivity specificity
###     -0.003248072     1        0.83

可以看到这两个分组的数值，以0为分界线，很容易就“近乎完美”的区分开来了，而且这个分类模型效果非常棒。。。