-
2021-03-07 22:16:46
科研工作,需要扎实的统计学知识做基础。众所周知,GraphPad Prism不只是一个强大的绘图工具,还是一个强大的统计学工具!
只有掌握了统计分析方法的原理,才会清楚“什么类型的数据该选用什么分析方法”,绘制出的图表才会讲出明确的故事。
现在,GraphPad官方中文用户指南已经上线!今天,我们就先从一个基本的‘P值’概念带着大家学习统计!
什么是P值?
假设您从两个用不同药物治疗的动物样本中收集了数据。您已测量每种动物血浆中的一种酶,且测量方法不同。您想知道这种差异是否是由于药物的作用 - 这两种群体是否有不同的方法。
观察不同样本方法并不足以说服您得出群体有不同的方法的结论。有可能群体具有相同的平均值(即药物对您正在测量的酶无影响),且您观察到的样本平均值之间的差异只会偶然发生。您永远无法确定您观察到的差异是否反映了真正的差异,或只是在随机抽样过程中出现。您所能做的就是计算概率。
第一步是注明零假设。事实上,治疗并不影响您测量的结果(因此,所有差异均因随机抽样所致)。
P值是一个概率值,范围为0 - 1,可回答该问题( 可能您从未想过要问):
在这种规模的试验中,如果群体真的具有相同的平均值,则观察到样本平均值之间的差异至少与实际观察到的差异一样大的概率是多少?
对P值的最常见错误解读
许多人误解了P值的含义。让我们假设您比较了两个平均值,得到等于0.03的P值。
该P值的正确定义:
观察到的差异与您所观察到的差异一样大的几率为3%,即使两个群体平均值相同(零假设为真)。
或者
从相同群体中随机抽样会导致比您在97%的实验中观察到的差异更小,比您在3%的实验中观察到的差异更大。
错误:
您所观察到的差异反映出群体之间的真实差异的几率为97%,差异偶然发生的几率为3%。
对P值的更多错误解读
line(1)列出了普遍认为的P值谬误,我们在这里总结如下:谬论:p值是由于抽样误差而造成结果的概率
假设零假设为真,则计算P值。换言之,P值的计算是基于由于抽样误差而造成差异的假设。因此,P值无法表示由于抽样误差而造成结果的概率。谬论:P值是零假设为真的概率
没有。假设零假设为真,则计算P值,因此P值不可能是真的值。谬论:1 - P是替代假设为真的概率
如果P值为0.03,则很容易想到:如果差异只有3%的概率由随机因素造成,则差异肯定有97%的概率由真实因素造成。但这是一种错误的想法!
您可以说,如果零假设为真,则97%的实验会导致比您观察到的差异更小的差异,3%的实验会导致与您观察到的差异一样大的差异或比您观察到的差异更大的差异。
P值的计算基于零假设正确的假设。P值无法表明该假设是否正确。P值表明,如果零假设为真,则很少能观察到与您观察到的差异一样大的差异或比您观察到的差异更大的差异。
科学家必须回答的问题是,结果是否不太可能导致放弃零假设。谬论:1 - P是重复实验后结果保持不变的概率
如果P值为0.03,则很容易认为这意味着有97%的几率能在重复实验中得到“相似”的结果。并非如此。谬论:高P值证明零假设为真。
否。高P值意味着,如果零假设为真,则在本实验中观察到的治疗效果便不足为奇了。但这并不能证明零假设为真。谬论:P值是拒绝零假设的概率
在一个特定实验的P值小于显著性水平α(您(应该)将α作为实验设计的一部分)时,您会拒绝零假设(并认为结果具有统计学显著性)。所以如果零假设为真,则α是拒绝零假设的概率。
P值和α不一样。P值是从每次比较中计算得出,并且是对证据强度的一种度量。将显著性水平α设置为实验设计的一部分。
GraphPad官方中文用户指南现已上线,大家可以登录中国官网找到入口:https://www.graphpad-prism.cn
GraphPad中文用户指南是一个线上知识库平台,支持关键词搜索,不论是统计学小白想从基础开始学习,还是资深用户遇到了个别棘手的问题,都可以在这个平台找到自己的方向。希望GraphPad中文用户指南的上线能给大家的科研学习和工作带来帮助!
大家如果有任何问题,欢迎给我们评论或者发私信哦!
更多相关内容 -
python 相关性检验怎么计算p值_数据分析---用Python进行相关性分析(兼谈假设检验)...
2020-11-20 23:00:53一、相关关系和相关系数世上除了因果关系,还有相关关系。有一个叫“相关系数”的指标来量化两个事物之间的相关程度。其中相关系数用“r”表示,取值范围介于-1和1之间。当(X,Y)正相关的时候,r=1;当(X,Y)负相关的...一、相关关系和相关系数
世上除了因果关系,还有相关关系。
有一个叫“相关系数”的指标来量化两个事物之间的相关程度。
其中相关系数用“r”表示,取值范围介于-1和1之间。
当(X,Y)正相关的时候,r=1;当(X,Y)负相关的时候,r=-1;当(X,Y)不相关的时候,r=0。
当然一般的线性相关有更严格的划分:
- r|<0.3 不存在线性关系
- 0.3<|r|<0.5 低度线性关系
- 0.5<|r|<0.8 显著线性关系
- |r|>0.8 高度线性关系
二、相关性和假设检验
有指标来衡量两者之间的相关程度,不代表能够去衡量相关程度。
因引入两个概念:
现在,针对我们分析的两组数据(X,Y)(两组数据被称为抽样),我们的疑问来了:
抽样的(X,Y)是否可以正确反应总体的情况呢?
这里涉及:假设检验。
具体操作如下:
零假设H0:总体的数据不呈相关性(相关系数为0),并先认为H0正确 备选假设H1:总体的数据呈现相关性(相关系数不为0) 引入一个指标:显著性水平p,一般将其设定为0.05或者0.01 当p<0.05,拒绝原假设,备选假设正确; 当p>0.05,原假设正确。 所以,在进行相关性分析实验的之前,我们需要分两步走: 1.进行假设检验,获得p值<0.05,得到结论:总体的数据呈现相关性 2.进行相关性分析,得到r值 如果p值>0.05(或者0.01),则实验失败,抽样数据无法反应整体情况。不管r值表现如何都是偶然事件。 只有在p值<0.05(或者0.01)的前提下,才可以参考r值,进而判断相关程度。
三、兼谈假设检
假设检验的3种类型
假设检验的套路:
四、利用Python进行相关性分析
判定两者相关的方式有两种:
- 图形观测法:通过绘制散点图判断两者是否存在一定相关关系
- 科学计算法:通过计算相关性系数r
我们用第二种
import numpy as np import pandas as pd import matplotlib.pyplot as plt import seaborn as sns from scipy.stats import kstest from scipy import stats #读入数据 data=pd.read_csv('http://jse.amstat.org/datasets/normtemp.dat.txt',header=None,sep='s+',names=['Temperature','sex','heart']) print(data.describe()) Temperature_data = data['Temperature'] u = data['Temperature'].mean() std = data['Temperature'].std() r,p = stats.pearsonr(data.Temperature, data.heart) print('相关系数r为 = %6.3f,p值为 = %6.3f'%(r,p)) 相关系数r为 = 0.254,p值为 = 0.004
得到:相关系数r=0.021,p值为=0.004
结论:总体的数据呈相关性,且相关系数为:0.021,但不是线性相关。
五、拓展
对于多维数据,需要计算两两之间的相关性。
比如是思维数据,列名分别为:A、B、C、D
就需要计算:
A:B、C、D
B:A、C、D
C:A、B、D
D:A、B、C
代码如下:
import numpy as np import pandas as pd import matplotlib.pyplot as plt import scipy.stats as stats # 导入数据 data = pd.DataFrame(数据地址) #或者 data=pd.read_csv(数据地址) # 相关性计算 print(data.corr()) # 绘图 fig = pd.plotting.scatter_matrix(data,figsize=(6,6),c ='blue',marker = 'o',diagonal='',alpha = 0.8,range_padding=0.2) # diagonal只能为'hist'/'kde' plt.show()
假设检验
Python真香:用Python作假设检验zhuanlan.zhihu.com -
数据可视化——R语言多种方式绘制分组数据(箱形图、小提琴图、散点图、误差条图等及相互组合并自动添加P值...
2020-05-27 14:17:11R语言多种方式呈现分组数据(箱形图、小提琴图、散点图、误差条图等)准备分组的离散变量分组的连续变量箱形图(Box plots)小提琴图(Violin plots)点图(Dot plots)一维散点图(Stripcharts)Sinaplot带有误差条的均值和...数据可视化——R语言多种方式绘制分组数据(箱形图、小提琴图、散点图、误差条图等及相互组合并自动添加P值)
本文内容为一篇英文博客的翻译内容,博客原文链接:Plot Grouped Data: Box plot, Bar Plot and More
概述:本文重点介绍各种用于分组数据(基于类别变量分组)的绘图方式。包括箱形图(Box plots)、小提琴图(Violin plots)、点图(Dot plots)、一维散点图(Stripcharts)、Sinaplot、条形图及折线图等多种形式及各种图形的相互组合。类型多样,灵活易用。算作是对之前博文很好的总结。另外,本文还实现了对分组变量组间比较自动添加P值(T检验或wilcoxon检验的P值)和显著性水平。
使用工具:R语言中的ggplot2包和ggpubr包
准备
加载所需的包,并使用函数theme_set()设定主题样式为theme_pubclean()
library(dplyr) library(ggplot2) library(ggpubr) theme_set(theme_pubclean())
分组的离散变量
- 绘图类型:条形图,用于表示每个组别数据的频数。关键函数:geom_bar()
- 示例数据集:ggplot2包中的diamonds,是一个表示钻石性质的数据集
- cut:分类变量,表示钻石切割的质量,由坏到好共5个等级:Fair, Good, Very Good, Premium, Ideal
- color:分类变量,表示钻石的颜色,“J”表示较差,“D”表示较好
仅采用diamonds中的一个子集数据用于图形绘制,子集数据选择流程如下:
- 仅保留颜色为“J”或“D”的数据
- 按cut和color两个属性对数据分组
- 统计每组的数量
- 采用geom_bar()绘制条形图
- 筛选数据并统计每个分组的数量
df <- diamonds %>% filter(color %in% c("J", "D")) %>% group_by(cut, color) %>% summarise(counts = n()) head(df, 4)
数据形式如下:
## # A tibble: 4 x 3 ## # Groups: cut [2] ## cut color counts ## ## 1 Fair D 163 ## 2 Fair J 119 ## 3 Good D 662 ## 4 Good J 307
- 绘制各分组的条形图
- 关键函数:geom_bar();关键参数:stat = “identity”,按原始方式绘制图形
- scale_color_manual()和scale_fill_manual()可以设置条形图的边框颜色和填充颜色
# Stacked bar plots of y = counts by x = cut, # colored by the variable color ggplot(df, aes(x = cut, y = counts)) + geom_bar( aes(color = color, fill = color), stat = "identity", position = position_stack() ) + scale_color_manual(values = c("#0073C2FF", "#EFC000FF"))+ scale_fill_manual(values = c("#0073C2FF", "#EFC000FF")) # Use position = position_dodge() p <- ggplot(df, aes(x = cut, y = counts)) + geom_bar( aes(color = color, fill = color), stat = "identity", position = position_dodge(0.8), width = 0.7 ) + scale_color_manual(values = c("#0073C2FF", "#EFC000FF"))+ scale_fill_manual(values = c("#0073C2FF", "#EFC000FF")) p
注意:设置position = position_stack(reverse = TRUE)后可改变条形图的堆叠顺序,即可改变color为”D“和”J“的两个组堆叠顺序,reverse = FALSE时”D“在上,reverse = TRUE时”J“在上。另外,也可以采用点图替代条形图。
ggplot(df, aes(cut, counts)) + geom_linerange( aes(x = cut, ymin = 0, ymax = counts, group = color), color = "lightgray", size = 1.5, position = position_dodge(0.3) )+ geom_point( aes(color = color), position = position_dodge(0.3), size = 3 )+ scale_color_manual(values = c("#0073C2FF", "#EFC000FF"))+ theme_pubclean()
3. 为分离条形图添加数值p + geom_text( aes(label = counts, group = color), position = position_dodge(0.8), vjust = -0.3, size = 3.5 )
4. 为堆叠条形图添加数值包括三步:
- 将数据按cut和color排序,由于position_stack()会颠倒顺序,color按降序排列。
- 计算每个cut分类的累计数量和,用于作为数值显示时的Y轴位置。将数值显示于条形图的中间位置,即cumsum(counts) - 0.5 * counts处。
- 创建条形图并添加数值
# Arrange/sort and compute cumulative summs df <- df %>% arrange(cut, desc(color)) %>% mutate(lab_ypos = cumsum(counts) - 0.5 * counts) head(df, 4)
数据形式如下:
## # A tibble: 4 x 4 ## # Groups: cut [2] ## cut color counts lab_ypos ## ## 1 Fair J 119 59.5 ## 2 Fair D 163 200.5 ## 3 Good J 307 153.5 ## 4 Good D 662 638.0
# Create stacked bar graphs with labels ggplot(df, aes(x = cut, y = counts)) + geom_bar(aes(color = color, fill = color), stat = "identity") + geom_text( aes(y = lab_ypos, label = counts, group = color), color = "white" ) + scale_color_manual(values = c("#0073C2FF", "#EFC000FF"))+ scale_fill_manual(values = c("#0073C2FF", "#EFC000FF"))
- 或者,可以使用ggpubr包中的ggbarplot()函数快速创建上图
ggbarplot(df, x = "cut", y = "counts", color = "color", fill = "color", palette = c("#0073C2FF", "#EFC000FF"), label = TRUE, lab.pos = "in", lab.col = "white", ggtheme = theme_pubclean() )
- 条形图的替代方案
用点的数量表示各组中的记录数。X轴和Y轴均表示分类变量,每一个点隶属于一个分组。对于给定的组,点数与该组中的记录数相对应。
- 关键函数:geom_jitter()
- 关键参数:alpha, color, fill, shape 和size
下列示例仅绘制了一小部分数据(整个diamonds 中的1/5)。
diamonds.frac <- dplyr::sample_frac(diamonds, 1/5) ggplot(diamonds.frac, aes(cut, color)) + geom_jitter(aes(color = cut), size = 0.3)+ ggpubr::color_palette("jco")+ ggpubr::theme_pubclean()
分组的连续变量
采用箱形图、小提琴图、点图等方式呈现分组连续变量的数据分布;
同时,还描述了如何自动添加组间比较的P值。设置主题样式为theme_bw()
theme_set(theme_bw())
- 示例数据集:ToothGrowth
- len:连续变量,表示牙齿的长度,用作Y轴
- dose:分类变量,表示维生素C的剂量水平,取值分别为0.5, 1, 和 2 毫克/天,用作X轴
首先将变量dose从数值型转换为离散因子型。
data("ToothGrowth") ToothGrowth$dose <- as.factor(ToothGrowth$dose) head(ToothGrowth)
数据形式如下:
## len supp dose ## 1 4.2 VC 0.5 ## 2 11.5 VC 0.5 ## 3 7.3 VC 0.5 ## 4 5.8 VC 0.5 ## 5 6.4 VC 0.5 ## 6 10.0 VC 0.5
箱形图(Box plots)
- 关键函数:geom_boxplot()
- 可自定义的关键参数
- width:指定箱形图的宽度。
- notch:逻辑变量,如果为TRUE,创建具有凹槽的箱形图,凹槽范围表示围绕中位数的置信区间,median ± 1.58*IQR/sqrt(n)。凹槽可用于组间比较:如果两个箱形图的凹槽不重叠,说明这两组间的中位数存在差异。
- color, size, linetype: 指定箱形图边框的颜色,粗细及线型。
- fill:指定箱形图的填充颜色。
- outlier.colour, outlier.shape, outlier.size:指定异常点的颜色、形状和大小。
1.创建最基本的箱形图
- 标准的带凹槽的箱形图
# Default plot e <- ggplot(ToothGrowth, aes(x = dose, y = len)) e + geom_boxplot() # Notched box plot with mean points e + geom_boxplot(notch = TRUE, fill = "lightgray")+ stat_summary(fun.y = mean, geom = "point", shape = 18, size = 2.5, color = "#FC4E07")
- 依据组别改变箱形图的颜色
# Color by group (dose) e + geom_boxplot(aes(color = dose))+ scale_color_manual(values = c("#00AFBB", "#E7B800", "#FC4E07")) # Change fill color by group (dose) e + geom_boxplot(aes(fill = dose)) + scale_fill_manual(values = c("#00AFBB", "#E7B800", "#FC4E07"))
注意:scale_x_discrete()的灵活使用- 如设置scale_x_discrete(limits = c(“0.5”, “2”))仅显示dose为"0.5"和"2"的两个分组
- 改变scale_x_discrete(limits = c(“2”, “0.5”, “1”))可改变X轴分类变量的显示顺序
例如:
# Choose which items to display: group "0.5" and "2" e + geom_boxplot() + scale_x_discrete(limits=c("0.5", "2")) # Change the default order of items e + geom_boxplot() + scale_x_discrete(limits=c("2", "0.5", "1"))
2. 创建具有多个分组的箱形图- 将两个分类变量用于数据分组,dose和supp
- 使用position_dodge()可以调整多个分组的箱形图之间的间隔
e2 <- e + geom_boxplot( aes(fill = supp), position = position_dodge(0.9) ) + scale_fill_manual(values = c("#999999", "#E69F00")) e2
- 还可以依据facet_wrap()函数将不同分组显示于不同面板
e2 + facet_wrap(~supp)
小提琴图(Violin plots)
小提琴图类似于箱形图,不同之处在于小提琴图还显示了数据的核概率密度估计。通常,小提琴图还包含数据分布的中位数及四分位数范围的框,与标准箱形图类似。
关键函数:
- geom_violin(),用于创建小提琴图。关键参数:
- color, size, linetype:指定边框的颜色,粗细及线型。
- fill:指定填充颜色。
- trim:逻辑变量,如果为TRUE(默认),将小提琴的尾部修整到与数据范围一致,如果为FALSE,将不修整尾部。
- stat_summary(),用于在绘制小提琴图时添加汇总统计信息,如添加均值,中值等。
- 使用汇总统计创建基本的小提琴图
# Add mean points +/- SD # Use geom = "pointrange" or geom = "crossbar" e + geom_violin(trim = FALSE) + stat_summary( fun.data = "mean_sdl", fun.args = list(mult = 1), geom = "pointrange", color = "black" ) # Combine with box plot to add median and quartiles # Change color by groups e + geom_violin(aes(fill = dose), trim = FALSE) + geom_boxplot(width = 0.2)+ scale_fill_manual(values = c("#00AFBB", "#E7B800", "#FC4E07"))+ theme(legend.position = "none")
注意:函数mean_sdl()用于添加平均值和标准差,以误差条形式显示均值 ± 标准差。参数 mult需被设置为常数,在上面的R代码中, mult设为1,表示均值±(1倍)标准差;默认情况下,mult = 2,表示均值±(2倍)标准差- 创建具有多个分组的小提琴图
e + geom_violin( aes(color = supp), trim = FALSE, position = position_dodge(0.9) ) + geom_boxplot( aes(color = supp), width = 0.15, position = position_dodge(0.9) ) + scale_color_manual(values = c("#00AFBB", "#E7B800"))
点图(Dot plots)
- 关键函数:geom_dotplot(),创建堆叠的点,每个点代表一个观测值
- 关键参数:
- stackdir:表示向哪个方向堆放点,默认值为“up”,向上;“down”,向下;“center”, 居中;“centerwhole”,居中,但点与点间也需对齐。
- stackratio:点之间堆放的间距。默认为1,表示点与点之间只是接触,设置更小的值将使得点之间靠得更近,点间将存在重叠。
- color, fill:设置点的边缘颜色或填充颜色。
- dotsize:点相对于binwidth的直径,默认为1。
也可以将汇总统计与点图叠加。
- 创建基本的点图
# Violin plots with mean points +/- SD e + geom_dotplot( binaxis = "y", stackdir = "center", fill = "lightgray" ) + stat_summary( fun.data = "mean_sdl", fun.args = list(mult=1), geom = "pointrange", color = "red" ) # Combine with box plots e + geom_boxplot(width = 0.5) + geom_dotplot( binaxis = "y", stackdir = "center", fill = "white" ) # Dot plot + violin plot + stat summary e + geom_violin(trim = FALSE) + geom_dotplot( binaxis='y', stackdir='center', color = "black", fill = "#999999" ) + stat_summary( fun.data="mean_sdl", fun.args = list(mult=1), geom = "pointrange", color = "#FC4E07", size = 0.4 )
2. 创建具有多个分组的点图# Color dots by groups e + geom_boxplot(width = 0.5, size = 0.4) + geom_dotplot( aes(fill = supp), trim = FALSE, binaxis='y', stackdir='center' )+ scale_fill_manual(values = c("#00AFBB", "#E7B800")) # Change the position : interval between dot plot of the same group e + geom_boxplot( aes(color = supp), width = 0.5, size = 0.4, position = position_dodge(0.8) ) + geom_dotplot( aes(fill = supp, color = supp), trim = FALSE, binaxis='y', stackdir='center', dotsize = 0.8, position = position_dodge(0.8) )+ scale_fill_manual(values = c("#00AFBB", "#E7B800"))+ scale_color_manual(values = c("#00AFBB", "#E7B800"))
一维散点图(Stripcharts)
当样本量较小时,可绘制一维散点图,将比箱形图更合适
- 关键函数:geom_jitter()
- 关键参数:color, fill, size, shape:指定点的边缘颜色、填充颜色、大小及形状。
- 创建基本的一维散点图
- 按组别来设定点的形状和颜色
- 调整抖动程度:position_jitter(0.2)
- 添加汇总统计
e + geom_jitter( aes(shape = dose, color = dose), position = position_jitter(0.2), size = 1.2 ) + stat_summary( aes(color = dose), fun.data="mean_sdl", fun.args = list(mult=1), geom = "pointrange", size = 0.4 )+ scale_color_manual(values = c("#00AFBB", "#E7B800", "#FC4E07"))
2. 创建具有多个分组的一维散点图代码类似于点图的绘制,但需要注意的是,使用函数position_jitterdodge()来调整点之间的抖动程度,而不是position_dodge()。
e + geom_jitter( aes(shape = supp, color = supp), position = position_jitterdodge(jitter.width = 0.2, dodge.width = 0.8), size = 1.2 ) + stat_summary( aes(color = supp), fun.data="mean_sdl", fun.args = list(mult=1), geom = "pointrange", size = 0.4, position = position_dodge(0.8) )+ scale_color_manual(values = c("#00AFBB", "#E7B800"))
Sinaplot
- Sinaplot能结合一维散点图和小提琴图。点沿X轴抖动,点的抖动密度与数据分布的密度对应,该图显示的轮廓与小提琴图相似,但又类似于少量数据点的一维散点图 (Sidiropoulos et al. 2015)。
- sinaplot以精简和易于理解的方式呈现数据分布信息,能表示点数量、密度分布、离群值和散布信息等。
- 关键函数:geom_sina()
library(ggforce) # Create some data d1 <- data.frame( y = c(rnorm(200, 4, 1), rnorm(200, 5, 2), rnorm(400, 6, 1.5)), group = rep(c("Grp1", "Grp2", "Grp3"), c(200, 200, 400)) ) # Sinaplot ggplot(d1, aes(group, y)) + geom_sina(aes(color = group), size = 0.7)+ scale_color_manual(values = c("#00AFBB", "#E7B800", "#FC4E07"))
带有误差条的均值和中值图
将展示如何绘制具有一个组别或多个组别的连续变量的汇总统计信息。
ggpubr软件包提供了一种简单的方法,只需较少的输入即可创建均值/中位数图。请参考以下相关文章:ggpubr-Plot Means/Medians and Error Bars设置主题样式为:theme_pubr()
theme_set(ggpubr::theme_pubr())
- 基本的均值/中位数图。一个连续变量和一个分组变量的情况:
- 示例数据集:ToothGrowth。
df <- ToothGrowth df$dose <- as.factor(df$dose) head(df, 3)
数据形式如下:
## len supp dose ## 1 4.2 VC 0.5 ## 2 11.5 VC 0.5 ## 3 7.3 VC 0.5
- 计算按dose分组的len的汇总统计量
library(dplyr) df.summary <- df %>% group_by(dose) %>% summarise( sd = sd(len, na.rm = TRUE), len = mean(len) ) df.summary
数据形式如下:
## # A tibble: 3 x 3 ## dose sd len ## ## 1 0.5 4.50 10.6 ## 2 1 4.42 19.7 ## 3 2 3.77 26.1
- 使用汇总统计数据创建误差条。关键函数:
- geom_crossbar():用于绘制中间带有水平线的hollow bar,
- geom_errorbar():用于绘制误差条
- geom_errorbarh():用于绘制水平误差条
- geom_linerange():用于绘制用垂直线表示的线段(区间)
- geom_pointrange():用于绘制用垂直线表示的线段(区间),但中间有一个点
首先使用汇总统计数据初始化ggplot,指定x、y、ymin和ymax,通常指定ymin = len-sd,ymax = len+sd 来添加向下和向上的误差条(均值 ± 标准差)。如果仅需要向上的误差条而不需要向下的误差条,可设置 ymin = len ,ymax = len+sd
# Initialize ggplot with data f <- ggplot( df.summary, aes(x = dose, y = len, ymin = len-sd, ymax = len+sd) ) f + geom_crossbar() +ggtitle("f + geom_crossbar()") f + geom_errorbar(width = 0.4) +ggtitle("f + geom_errorbar()") f + geom_linerange() +ggtitle("f + geom_linerange()") f + geom_pointrange() +ggtitle("f + geom_pointrange()")
各种误差条类型:
创建简单的误差条:
# Vertical line with point in the middle f + geom_pointrange() # Standard error bars f + geom_errorbar(width = 0.2) + geom_point(size = 1.5)
创建水平的误差条,Y轴表示dose ,X轴表示len,还需指定xmin和xmax。
# Horizontal error bars with mean points # Change the color by groups ggplot( df.summary, aes(x = len, y = dose, xmin = len-sd, xmax = len+sd) ) + geom_point(aes(color = dose)) + geom_errorbarh(aes(color = dose), height=.2)+ theme_light()
- 添加一维散点图和小提琴图。为此,需要使用df初始化ggplot中的数据data(绘制一维散点图时使用);将df.summary汇总统计信息用作函数geom_pointrange()的输入(绘制小提琴图时使用)
# Combine with jitter points ggplot(df, aes(dose, len)) + geom_jitter( position = position_jitter(0.2), color = "darkgray" ) + geom_pointrange( aes(ymin = len-sd, ymax = len+sd), data = df.summary ) # Combine with violin plots ggplot(df, aes(dose, len)) + geom_violin(color = "darkgray", trim = FALSE) + geom_pointrange( aes(ymin = len-sd, ymax = len+sd), data = df.summary )
- 创建均值 ± 误差的基本折线图/条形图,仅需使用df.summary 数据
- 为折线图添加向下或向上的误差条:设置ymin = len-sd 和ymax = len+sd
- 只为条形图添加向上的误差条:设置ymin = len和ymax = len+sd
值得注意的是:绘制折线图时,如果仅有一个分组时,应始终在aes()中设置group = 1
# (1) Line plot ggplot(df.summary, aes(dose, len)) + geom_line(aes(group = 1)) + geom_errorbar( aes(ymin = len-sd, ymax = len+sd),width = 0.2) + geom_point(size = 2) # (2) Bar plot ggplot(df.summary, aes(dose, len)) + geom_bar(stat = "identity", fill = "lightgray", color = "black") + geom_errorbar(aes(ymin = len, ymax = len+sd), width = 0.2)
对于折线图,可以设置X轴为数值型变量,而非因子变量df.sum2 <- df.summary df.sum2$dose <- as.numeric(df.sum2$dose) ggplot(df.sum2, aes(dose, len)) + geom_line() + geom_errorbar( aes(ymin = len-sd, ymax = len+sd),width = 0.2) + geom_point(size = 2)
- 折线图 + 一维散点图:将df作为一维散点图的输入,将df.summary 作为其他geom的输入
- 对于折线图,首先添加一维散点,再添加折线,再添加误差条,再添加均值点
- 对于条形图,首先添加条形图,再添加一维散点,再添加误差条
# (1) Create a line plot of means + # individual jitter points + error bars ggplot(df, aes(dose, len)) + geom_jitter( position = position_jitter(0.2), color = "darkgray") + geom_line(aes(group = 1), data = df.summary) + geom_errorbar( aes(ymin = len-sd, ymax = len+sd), data = df.summary, width = 0.2) + geom_point(data = df.summary, size = 2) # (2) Bar plots of means + individual jitter points + errors ggplot(df, aes(dose, len)) + geom_bar(stat = "identity", data = df.summary, fill = NA, color = "black") + geom_jitter( position = position_jitter(0.2), color = "black") + geom_errorbar( aes(ymin = len-sd, ymax = len+sd), data = df.summary, width = 0.2)
2. 具有多个分组的均值/中位数图。一个连续变量(len),两个分组变量(dose, supp)的情况- 计算len按dose 和supp分组的汇总统计:
library(dplyr) df.summary2 <- df %>% group_by(dose, supp) %>% summarise( sd = sd(len), len = mean(len) ) df.summary2
示例数据如下:
## # A tibble: 6 x 4 ## # Groups: dose [?] ## dose supp sd len ## ## 1 0.5 OJ 4.46 13.23 ## 2 0.5 VC 2.75 7.98 ## 3 1 OJ 3.91 22.70 ## 4 1 VC 2.52 16.77 ## 5 2 OJ 2.66 26.06 ## 6 2 VC 4.80 26.14
- 创建具有多个分组的误差条
- geom_pointrange()按组别supp设置不同的颜色
- 误差条和均值点按组别supp设置不同的颜色
# (1) Pointrange: Vertical line with point in the middle ggplot(df.summary2, aes(dose, len)) + geom_pointrange( aes(ymin = len-sd, ymax = len+sd, color = supp), position = position_dodge(0.3) )+ scale_color_manual(values = c("#00AFBB", "#E7B800")) # (2) Standard error bars ggplot(df.summary2, aes(dose, len)) + geom_errorbar( aes(ymin = len-sd, ymax = len+sd, color = supp), position = position_dodge(0.3), width = 0.2 )+ geom_point(aes(color = supp), position = position_dodge(0.3)) + scale_color_manual(values = c("#00AFBB", "#E7B800"))
-
创建具有多个分组的折线图/条形图
- 折线图:按组别supp设置不同的线型
- 条形图:按组别supp设置不同的填充颜色
# (1) Line plot + error bars ggplot(df.summary2, aes(dose, len)) + geom_line(aes(linetype = supp, group = supp))+ geom_point()+ geom_errorbar( aes(ymin = len-sd, ymax = len+sd, group = supp), width = 0.2 ) # (2) Bar plots + upper error bars. ggplot(df.summary2, aes(dose, len)) + geom_bar(aes(fill = supp), stat = "identity", position = position_dodge(0.8), width = 0.7)+ geom_errorbar( aes(ymin = len, ymax = len+sd, group = supp), width = 0.2, position = position_dodge(0.8) )+ scale_fill_manual(values = c("grey80", "grey30"))
- 创建具有多个分组的均值 ± 标准差图。使用ggpubr包,将自动计算汇总统计信息并创建图形。
library(ggpubr) # Create line plots of means ggline(ToothGrowth, x = "dose", y = "len", add = c("mean_sd", "jitter"), color = "supp", palette = c("#00AFBB", "#E7B800")) # Create bar plots of means ggbarplot(ToothGrowth, x = "dose", y = "len", add = c("mean_se", "jitter"), color = "supp", palette = c("#00AFBB", "#E7B800"), position = position_dodge(0.8))
- 使用ggplot2包绘制上图
# Create line plots ggplot(df, aes(dose, len)) + geom_jitter( aes(color = supp), position = position_jitter(0.2) ) + geom_line( aes(group = supp, color = supp), data = df.summary2 ) + geom_errorbar( aes(ymin = len-sd, ymax = len+sd, color = supp), data = df.summary2, width = 0.2 )+ scale_color_manual(values = c("#00AFBB", "#E7B800"))
添加P值和显著性水平
介绍如何轻松地 i)比较两个或多个组的均值; ii)并将p值和显著性水平自动添加到ggplot图中。
关键函数:
- compare_means():计算单次或多次均值比较的结果
- stat_compare_means:可将P值和显著性水平自动添加到ggplot图中
比较均值的最常用方法包括:
方法 R实现函数 描述 T-test t.test() 比较两组(参数检验) Wilcoxon test wilcox.test() 比较两组(非参数检验) ANOVA aov() or anova() 比较多组(参数检验) Kruskal-Wallis kruskal.test() 比较多组(非参数检验) - 比较两个独立组
- T检验
library(ggpubr) compare_means(len ~ supp, data = ToothGrowth, method = "t.test")
示例结果如下:
## # A tibble: 1 x 8 ## .y. group1 group2 p p.adj p.format p.signif method ## ## 1 len OJ VC 0.0606 0.0606 0.061 ns T-test
- 创建箱形图并添加P值。操作时设置method = "t.test或method = “wilcox.test”(默认值)。
# Create a simple box plot and add p-values p <- ggplot(ToothGrowth, aes(supp, len)) + geom_boxplot(aes(color = supp)) + scale_color_manual(values = c("#00AFBB", "#E7B800")) p + stat_compare_means(method = "t.test") # Display the significance level instead of the p-value # Adjust label position p + stat_compare_means( aes(label = ..p.signif..), label.x = 1.5, label.y = 40 )
2. 比较两个配对样本ggpaired(ToothGrowth, x = "supp", y = "len", color = "supp", line.color = "gray", line.size = 0.4, palette = "jco")+ stat_compare_means(paired = TRUE)
3. 比较多组(超过两组)如果分类变量包含两个以上组别时,则将自动执行成对测试(pairwise tests)。 默认方法是“ wilcox.test”。 也可以将其更改为“t.test”。
# Perorm pairwise comparisons compare_means(len ~ dose, data = ToothGrowth)
示例结果如下:
## # A tibble: 3 x 8 ## .y. group1 group2 p p.adj p.format p.signif method ## ## 1 len 0.5 1 7.02e-06 1.40e-05 7.0e-06 **** Wilcoxon ## 2 len 0.5 2 8.41e-08 2.52e-07 8.4e-08 **** Wilcoxon ## 3 len 1 2 1.77e-04 1.77e-04 0.00018 *** Wilcoxon
# Visualize: Specify the comparisons you want my_comparisons <- list( c("0.5", "1"), c("1", "2"), c("0.5", "2") ) ggboxplot(ToothGrowth, x = "dose", y = "len", color = "dose", palette = "jco")+ stat_compare_means(comparisons = my_comparisons)+ stat_compare_means(label.y = 50)
- 多个分组变量
- (1) 创建一个按组划分的多面板框图(此处为dose)
# Use only p.format as label. Remove method name. ggplot(ToothGrowth, aes(supp, len)) + geom_boxplot(aes(color = supp))+ facet_wrap(~dose) + scale_color_manual(values = c("#00AFBB", "#E7B800")) + stat_compare_means(label = "p.format")
- (2) 创建一个包含所有框图的单一面板。X表示dose,Y表示len,颜色表示supp,stat_compare_means()中的group可设置分组。
ggplot(ToothGrowth, aes(dose, len)) + geom_boxplot(aes(color = supp))+ scale_color_manual(values = c("#00AFBB", "#E7B800")) + stat_compare_means(aes(group = supp), label = "p.signif")
注意出现如下信息时,说明显著性绘制失败,需要更新ggpubr包
Warning message: Computation failed in `stat_compare_means()`: Column `p` must be length 1 (the group size), not 3
- 多组配对比较
# Box plot facetted by "dose" p <- ggpaired(ToothGrowth, x = "supp", y = "len", color = "supp", palette = "jco", line.color = "gray", line.size = 0.4, facet.by = "dose", short.panel.labs = FALSE) # Use only p.format as label. Remove method name. p + stat_compare_means(label = "p.format", paired = TRUE)
阅读更多: Add P-values and Significance Levels to ggplots总结
- 可视化分组连续变量的数据分布:X轴表示分组变量,Y轴表示连续变量。
ggplot2中包含的函数:
- geom_boxplot() 用于绘制箱形图
- geom_violin() 用于绘制小提琴图
- geom_dotplot() 用于绘制点图
- geom_jitter()用于绘制一维散点图
- geom_line() 用于绘制折线图
- geom_bar() 用于绘制条形图
R示例代码:首先创建一个名为 e 的图,然后添加一个图层:
ToothGrowth$dose <- as.factor(ToothGrowth$dose) e <- ggplot(ToothGrowth, aes(x = dose, y = len)) e + geom_boxplot() + ggtitle("e + geom_boxplot()") e + geom_violin(trim = FALSE) + ggtitle("e + geom_violin()") e + geom_dotplot(binaxis = "y", stackdir = "center", fill = "lightgray") + ggtitle("e + geom_dotplot()") e + geom_jitter(position = position_jitter(0.2)) + ggtitle("e + geom_jitter()") library(dplyr) df <- ToothGrowth df$dose <- as.factor(df$dose) df.summary <- df %>% group_by(dose) %>% summarise( sd = sd(len, na.rm = TRUE), len = mean(len) ) df.summary e = ggplot(df.summary, aes(dose, len)) e + geom_line(aes(group = 1)) + geom_point(size = 2)+ ggtitle("e + geom_line() ") e + geom_bar(stat = "identity", color = "black") + ggtitle("e + geom_bar() ")
2. 创建均值和中值的误差条图X轴表示分组变量;Y轴表示连续变量的汇总统计(均值/中值)。
- 计算汇总统计数据并使用汇总数据初始化ggplot
# Summary statistics library(dplyr) df.summary <- ToothGrowth %>% group_by(dose) %>% summarise( sd = sd(len, na.rm = TRUE), len = mean(len) ) # Initialize ggplot with data f <- ggplot( df.summary, aes(x = dose, y = len, ymin = len-sd, ymax = len+sd) ) f + geom_crossbar() +ggtitle("f + geom_crossbar()") f + geom_linerange() +ggtitle("f + geom_linerange()") f + geom_errorbar(width = 0.4) +ggtitle("f + geom_errorbar()") f + geom_pointrange() +ggtitle("f + geom_pointrange()")
- 各种误差条类型:
- 将误差条图与小提琴图、点图、折线图及条形图结合。
# Combine with violin plots ggplot(ToothGrowth, aes(dose, len))+ geom_violin(trim = FALSE) + geom_pointrange(aes(ymin = len-sd, ymax = len + sd), data = df.summary) # Combine with dot plots ggplot(ToothGrowth, aes(dose, len))+ geom_dotplot(stackdir = "center", binaxis = "y", fill = "lightgray", dotsize = 1) + geom_pointrange(aes(ymin = len-sd, ymax = len + sd), data = df.summary) # Combine with line plot ggplot(df.summary, aes(dose, len))+ geom_line(aes(group = 1)) + geom_pointrange(aes(ymin = len-sd, ymax = len + sd)) # Combine with bar plots ggplot(df.summary, aes(dose, len))+ geom_bar(stat = "identity", fill = "lightgray") + geom_pointrange(aes(ymin = len-sd, ymax = len + sd))
扩展阅读
- ggpubr: Publication Ready Plots. https://goo.gl/7uySha
- Facilitating Exploratory Data Visualization: Application to TCGA
Genomic Data. https://goo.gl/9LNsRi - Add P-values and Significance Levels to ggplots.
https://goo.gl/VH7Yq7 - Plot Means/Medians and Error Bars. https://goo.gl/zRwAeV
References
原文链接:Plot Grouped Data: Box plot, Bar Plot and More
Sidiropoulos, Nikos, Sina Hadi Sohi, Nicolas Rapin, and Frederik Otzen Bagger. 2015. “SinaPlot: An Enhanced Chart for Simple and Truthful Representation of Single Observations over Multiple Classes.” bioRxiv. Cold Spring Harbor Laboratory. doi:10.1101/028191.
-
数据可视化——R语言为ggplot图形添加P值和显著性水平
2020-05-28 15:54:24数据可视化——R语言为ggplot图形添加P值和显著性水平准备安装和加载R包示例数据均值比较的方法用于添加P值的R函数compare_means()stat_compare_means()独立双样本组间比较配对双样本组间比较多组样本的组间比较多个...数据可视化——R语言为ggplot图形添加P值和显著性水平
本文对一篇英文博客进行翻译,博客原文链接:Add P-values and Significance Levels to ggplots
概述:本文介绍如何轻松地为ggplot图形添加P值和显著性水平:
- 比较两组或多组的均值
- 自动地将P值和显著性水平添加到ggplot图形中,如箱形图,点图,条形图和折线图等
使用工具: R语言中的ggplot2包和ggpubr包
准备
安装和加载R包
本文使用ggpubr包,要求版本高于0.1.3。ggpubr是一个基于ggplot2的计算工具包。
- 直接输入以下命令从CARN中下载安装
install.packages("ggpubr")
- 也可以从GitHub中下载安装最新版本
if(!require(devtools)) install.packages("devtools") devtools::install_github("kassambara/ggpubr")
- 加载ggpubr包
library(ggpubr)
ggpubr的官方文档可在以下位置获得:http://www.sthda.com/english/rpkgs/ggpubr
示例数据
示例数据集:ToothGrowth
data("ToothGrowth") head(ToothGrowth)
示例数据如下:
## len supp dose ## 1 4.2 VC 0.5 ## 2 11.5 VC 0.5 ## 3 7.3 VC 0.5 ## 4 5.8 VC 0.5 ## 5 6.4 VC 0.5 ## 6 10.0 VC 0.5
均值比较的方法
R中用于两组或多组间均值比较的标准统计方法在之前的文章也有描述:comparing means in R
均值比较的常见方法:
方法 R实现函数 描述 T-test t.test() 比较两组(参数检验) Wilcoxon test wilcox.test() 比较两组(非参数检验) ANOVA aov()或anova() 比较多组(参数检验) Kruskal-Wallis kruskal.test() 比较多组(非参数检验) 以下链接提供了各种方法的详细介绍:
-
单样本均值与已知均值的比较
-
独立双样本均值比较
-
配对双样本均值比较
-
两组以上的均值比较
- 方差分析(ANOVA, 参数检验)
- Kruskal-Wallis检验(替代单因素方差分析的非参数检验)
用于添加P值的R函数
介绍两个ggpubr包中的函数
- compare_means():用于执行均值比较
- stat_compare_means():用于在ggplot图形中自动添加P值和显著性水平
compare_means()
该函数用于执行均值比较。该函数与标准的R函数相比,灵活性更强。
简化形式如下:
compare_means(formula, data, method = "wilcox.test", paired = FALSE, group.by = NULL, ref.group = NULL, ...)
- formula:指定一个公式,公式形式为 x ~ group,其中,x 表示一个数值型变量,group 表示一个因子型变量,包含一个或多个水平。例如,一个示例公式为 formula = TP53 ~ cancer_group,表示在 cancer_group 对应的各水平间比较TP53的表达水平;也可以同时指定多个响应变量,如 formula = c(TP53, PTEN) ~ cancer_group。
- data:指定一个数据框(data.frame),数据框需包含formula中的变量。
- method:指定统计检验的方法。默认为“wilcox.test”,即Wilcoxon检验(非参数检验);也可指定其他统计方法:
- “t.test”,即T检验(参数检验)。“t.test”和“wilcox.test”用于两组样本间的比较。当超过两组时,将会执行两两比较(pairwise comparison)。
- “anova”(参数检验)或 “kruskal.test”(非参数检验),用于执行多组间的单因素方差分析。
- paired:指定一个逻辑变量,表示是否需要执行配对检验,仅适用于t.test 和wilcox.test。
- group.by:指定一个分组变量的字符名,用于在统计检验之前对数据进行分组。当存在group.by指定的变量时,均值比较将在不同水平的各个子集数据中执行。
- ref.group:指定一个组别的字符名,作为对照组(reference group)。如果指定,各个分组水平将与对照组水平进行比较。也可指定ref.group为“.all.”,表示每个分组水平将于所有分组水平(如base-mean)进行比较。
stat_compare_means()
该函数是对ggplot2的扩展,可将均值比较后的P值添加到ggplot图形中,如箱形图、点图、条形图和折线图等。
简化形式如下:
stat_compare_means(mapping = NULL, comparisons = NULL hide.ns = FALSE, label = NULL, label.x = NULL, label.y = NULL, ...)
-
mapping:通过 aes() 设置绘图时的aesthetic
-
comparisons:指定一个列表(list),每个列表元素需为长度等于2的向量。向量的内容可以为X轴的两个组别名(字符型),也可以是两个感兴趣组的组别索引(整数值),表示采用指定的两个组别进行比较。
-
hide.ns:逻辑变量,如果设为TRUE,显示显著性水平时将隐藏 ns 字样,即组间差异不显著时不显示 ns 字样。
-
label:指定一个字符串,表示标签类型。可为:“p.signif”(显示显著性水平),“p.format”(显示格式化的P值)。
-
label.x, label.y:指定一个数值,表示显示标签的绝对坐标位置。
-
…:传递给函数compare_means()的参数,如method、paired、ref.group。
独立双样本组间比较
执行统计检验
compare_means(len ~ supp, data = ToothGrowth)
示例结果如下:
## # A tibble: 1 x 8 ## .y. group1 group2 p p.adj p.format p.signif method ## ## 1 len OJ VC 0.0645 0.0645 0.064 ns Wilcoxon
method默认为“wilcox.test”(非参数检验),可指定method = “t.test”,表示T检验(参数检验)
返回值为具有以下列的数据框:
-
.y.:用于统计检验的数值变量
-
p:P值
-
p.adj:调整后的P值,调整P值的默认方法为p.adjust.method = “holm”
-
p.format: 格式化的P值
-
p.signif:显著性水平,即用不同数量的 * 表示显著性水平
-
method:用于组间比较的统计方法
创建添加P值的箱形图
p <- ggboxplot(ToothGrowth, x = "supp", y = "len", color = "supp", palette = "jco", add = "jitter") # Add p-value p + stat_compare_means() # Change method p + stat_compare_means(method = "t.test")
注意:显示P值的标签位置可以通过如下参数来调整:label.x, label.y, hjust 和vjust。显示P值的标签默认为 compare_means() 返回值中的 method 和 p 的组合。也可以通过 aes() 函数指定为其他显示形式。例如:
aes(label = ..p.format..) 或 aes(label = paste0(“p =”, ..p.format..)) 表示只显示格式化的P值,而不显示method aes(label = ..p.signif..) 表示仅显示显著性水平 aes(label = paste0(..method.., “\n”, “p =”, ..p.format..)) 表示在method名和P值之间添加换行符(“\n”)
示例如下:
p + stat_compare_means( aes(label = ..p.signif..), label.x = 1.5, label.y = 40)
另外,也可以将参数label指定为字符向量:p + stat_compare_means( label = "p.signif", label.x = 1.5, label.y = 40)
配对双样本组间比较
执行统计检验
compare_means(len ~ supp, data = ToothGrowth, paired = TRUE)
示例结果如下:
## # A tibble: 1 x 8 ## .y. group1 group2 p p.adj p.format p.signif method ## ## 1 len OJ VC 0.00431 0.00431 0.0043 ** Wilcoxon
使用函数 ggpaired() 可视化配对数据
ggpaired(ToothGrowth, x = "supp", y = "len", color = "supp", line.color = "gray", line.size = 0.4, palette = "jco")+ stat_compare_means(paired = TRUE)
多组样本的组间比较
- 全局检验(所有组的均值比较)
# Global test compare_means(len ~ dose, data = ToothGrowth, method = "anova")
示例结果如下:
## # A tibble: 1 x 6 ## .y. p p.adj p.format p.signif method ## ## 1 len 9.53e-16 9.53e-16 9.5e-16 **** Anova
添加全局检验的P值(所有组比较总的P值)
# Default method = "kruskal.test" for multiple groups ggboxplot(ToothGrowth, x = "dose", y = "len", color = "dose", palette = "jco")+ stat_compare_means() # Change method to anova ggboxplot(ToothGrowth, x = "dose", y = "len", color = "dose", palette = "jco")+ stat_compare_means(method = "anova")
- 两两比较(Pairwise comparisons)
如果分组变量包含两个以上的水平,两两比较的检验(pairwise test)将自动执行。默认方法为“wilcox.test”,也可设置为“t.test”。
# Perorm pairwise comparisons compare_means(len ~ dose, data = ToothGrowth)
示例结果如下:
## # A tibble: 3 x 8 ## .y. group1 group2 p p.adj p.format p.signif method ## ## 1 len 0.5 1 7.02e-06 1.40e-05 7.0e-06 **** Wilcoxon ## 2 len 0.5 2 8.41e-08 2.52e-07 8.4e-08 **** Wilcoxon ## 3 len 1 2 1.77e-04 1.77e-04 0.00018 *** Wilcoxon
# Visualize: Specify the comparisons you want my_comparisons <- list( c("0.5", "1"), c("1", "2"), c("0.5", "2") ) ggboxplot(ToothGrowth, x = "dose", y = "len", color = "dose", palette = "jco")+ stat_compare_means(comparisons = my_comparisons)+ # Add pairwise comparisons p-value stat_compare_means(label.y = 50) # Add global p-value
如果需要指定标签显示的Y轴位置,可使用参数label.y
ggboxplot(ToothGrowth, x = "dose", y = "len", color = "dose", palette = "jco")+ stat_compare_means(comparisons = my_comparisons, label.y = c(29, 35, 40))+ stat_compare_means(label.y = 45)
注:ggsignif 包也可以很方便的为条形图添加组间比较的P值- 相对于对照组的多重两两比较的检验(Multiple pairwise tests)
# Pairwise comparison against reference compare_means(len ~ dose, data = ToothGrowth, ref.group = "0.5", method = "t.test")
示例结果如下:
## # A tibble: 2 x 8 ## .y. group1 group2 p p.adj p.format p.signif method ## ## 1 len 0.5 1 6.70e-09 6.70e-09 6.7e-09 **** T-test ## 2 len 0.5 2 1.47e-16 2.94e-16 < 2e-16 **** T-test
# Visualize ggboxplot(ToothGrowth, x = "dose", y = "len", color = "dose", palette = "jco")+ stat_compare_means(method = "anova", label.y = 40)+ # Add global p-value stat_compare_means(label = "p.signif", method = "t.test", ref.group = "0.5") # Pairwise comparison against reference
- 相对于所有组(base-mean)的多重两两比较的检验
# Comparison of each group against base-mean compare_means(len ~ dose, data = ToothGrowth, ref.group = ".all.", method = "t.test")
示例结果如下:
## # A tibble: 3 x 8 ## .y. group1 group2 p p.adj p.format p.signif method ## ## 1 len .all. 0.5 1.24e-06 3.73e-06 1.2e-06 **** T-test ## 2 len .all. 1 5.67e-01 5.67e-01 0.57 ns T-test ## 3 len .all. 2 1.37e-05 2.74e-05 1.4e-05 **** T-test
# Visualize ggboxplot(ToothGrowth, x = "dose", y = "len", color = "dose", palette = "jco")+ stat_compare_means(method = "anova", label.y = 40)+ # Add global p-value stat_compare_means(label = "p.signif", method = "t.test", ref.group = ".all.") # Pairwise comparison against all
下面使用Github中可用的骨髓瘤数据集展示一些典型情况,其中,与“.all.”的比较将很有用。
将依据患者molecular进行分组,绘制各个组别DEPDC1基因的表达水平。目的是比较各个组别之间是否存在差别,如果有差别,差别又在哪里。
要回答以上问题,可以在所有7个组之间进行两两比较(pairwise comparison)。 由于组别较多,将导致很多种组别组合,这将很难解释。
一个简单的解决办法是将7组中的每一组与“.all.”(如base-mean)比较。当检验结果显著时,可以得出结论:与所有组相比,xxx组的DEPDC1的表达水平显著降低或显著升高。
# Load myeloma data from GitHub myeloma <- read.delim("https://raw.githubusercontent.com/kassambara/data/master/myeloma.txt") # Perform the test compare_means(DEPDC1 ~ molecular_group, data = myeloma, ref.group = ".all.", method = "t.test")
示例结果如下:
## # A tibble: 7 x 8 ## .y. group1 group2 p p.adj p.format p.signif method ## ## 1 DEPDC1 .all. Cyclin D-1 0.149690 0.44907 0.14969 ns T-test ## 2 DEPDC1 .all. Cyclin D-2 0.523143 1.00000 0.52314 ns T-test ## 3 DEPDC1 .all. Hyperdiploid 0.000282 0.00169 0.00028 *** T-test ## 4 DEPDC1 .all. Low bone disease 0.005084 0.02542 0.00508 ** T-test ## 5 DEPDC1 .all. MAF 0.086107 0.34443 0.08611 ns T-test ## 6 DEPDC1 .all. MMSET 0.576291 1.00000 0.57629 ns T-test ## # ... with 1 more rows
注意:上面的R代码中myeloma数据的下载地址存在错误,可从以下链接找到myeloma数据:https://github.com/kassambara/data/blob/master/myeloma.txt
然后自己将全部数据复制并保存为myeloma.txt,从本地读取myeloma.txt。# Visualize the expression profile ggboxplot(myeloma, x = "molecular_group", y = "DEPDC1", color = "molecular_group", add = "jitter", legend = "none") + rotate_x_text(angle = 45)+ geom_hline(yintercept = mean(myeloma$DEPDC1), linetype = 2)+ # Add horizontal line at base mean stat_compare_means(method = "anova", label.y = 1600)+ # Add global annova p-value stat_compare_means(label = "p.signif", method = "t.test", ref.group = ".all.") # Pairwise comparison against all
根据上图的结果,可以得出proliferation组的DEPDC1表达水平显著升高;Hyperdiploid组和Low bone disease组中DEPDC1表达水平显著降低。
注意:想要隐藏 ns 标志,可以设置参数 hide.ns = TRUE
# Visualize the expression profile ggboxplot(myeloma, x = "molecular_group", y = "DEPDC1", color = "molecular_group", add = "jitter", legend = "none") + rotate_x_text(angle = 45)+ geom_hline(yintercept = mean(myeloma$DEPDC1), linetype = 2)+ # Add horizontal line at base mean stat_compare_means(method = "anova", label.y = 1600)+ # Add global annova p-value stat_compare_means(label = "p.signif", method = "t.test", ref.group = ".all.", hide.ns = TRUE) # Pairwise comparison against all
多个分组变量
- 使用另一个变量进行分组后再执行独立双样本比较
执行统计检验:
compare_means(len ~ supp, data = ToothGrowth, group.by = "dose")
示例结果如下:
## # A tibble: 3 x 9 ## dose .y. group1 group2 p p.adj p.format p.signif method ## ## 1 0.5 len OJ VC 0.02319 0.0464 0.023 * Wilcoxon ## 2 1.0 len OJ VC 0.00403 0.0121 0.004 ** Wilcoxon ## 3 2.0 len OJ VC 1.00000 1.0000 1.000 ns Wilcoxon
在上面的示例中,对于分类变量dose的每一个水平,分类变量supp又将数据分为两个水平:OJ和VC,然后在这两个水平上对数值变量len进行均值比较。
可视化:创建一个按组划分的多面板框图(此处为“dose”)
# Box plot facetted by "dose" p <- ggboxplot(ToothGrowth, x = "supp", y = "len", color = "supp", palette = "jco", add = "jitter", facet.by = "dose", short.panel.labs = FALSE) # Use only p.format as label. Remove method name. p + stat_compare_means(label = "p.format")
# Or use significance symbol as label p + stat_compare_means(label = "p.signif", label.x = 1.5)
注意:想要隐藏 ns 标志,可以设置参数hide.ns = TRUE
可视化:创建一个包含所有箱形图的单一面板。X表示dose,Y表示len,颜色表示supp
p <- ggboxplot(ToothGrowth, x = "dose", y = "len", color = "supp", palette = "jco", add = "jitter") p + stat_compare_means(aes(group = supp))
# Show only p-value p + stat_compare_means(aes(group = supp), label = "p.format")
# Use significance symbol as label p + stat_compare_means(aes(group = supp), label = "p.signif")
- 使用另一个变量进行分组后再执行配对双样本比较
执行统计检验:
compare_means(len ~ supp, data = ToothGrowth, group.by = "dose", paired = TRUE)
实例结果如下:
## # A tibble: 3 x 9 ## dose .y. group1 group2 p p.adj p.format p.signif method ## ## 1 0.5 len OJ VC 0.0330 0.0659 0.033 * Wilcoxon ## 2 1.0 len OJ VC 0.0191 0.0572 0.019 * Wilcoxon ## 3 2.0 len OJ VC 1.0000 1.0000 1.000 ns Wilcoxon
可视化:创建一个按组划分的多面板框图(此处为“dose”)
# Box plot facetted by "dose" p <- ggpaired(ToothGrowth, x = "supp", y = "len", color = "supp", palette = "jco", line.color = "gray", line.size = 0.4, facet.by = "dose", short.panel.labs = FALSE) # Use only p.format as label. Remove method name. p + stat_compare_means(label = "p.format", paired = TRUE)
其他绘图方式
- 条形图和折线图(一个分组变量)
# Bar plot of mean +/-se ggbarplot(ToothGrowth, x = "dose", y = "len", add = "mean_se")+ stat_compare_means() + # Global p-value stat_compare_means(ref.group = "0.5", label = "p.signif", label.y = c(22, 29)) # compare to ref.group # Line plot of mean +/-se ggline(ToothGrowth, x = "dose", y = "len", add = "mean_se")+ stat_compare_means() + # Global p-value stat_compare_means(ref.group = "0.5", label = "p.signif", label.y = c(22, 29))
- 条形图和折线图(两个分组变量)
ggbarplot(ToothGrowth, x = "dose", y = "len", add = "mean_se", color = "supp", palette = "jco", position = position_dodge(0.8))+ stat_compare_means(aes(group = supp), label = "p.signif", label.y = 29) ggline(ToothGrowth, x = "dose", y = "len", add = "mean_se", color = "supp", palette = "jco")+ stat_compare_means(aes(group = supp), label = "p.signif", label.y = c(16, 25, 29))
注意:经过实际测试,笔者发现R语言中的统计方法计算结果的P值与SPSS中的P值存在差异。如,常规的方差分析(ANOVA) + 事后两两组间比较(如Bonferroni校正)使用上述R函数就很难得出与SPSS中一致的结果。如果需要使用SPSS的统计P值,建议对生成的图形进行后期修改。
References
-
excle如何计算t的值,TTEST 在EXCEL计算出的结果是t还是p值,用哪一个公式在excel中计算出t值和p值是多少
2021-01-17 12:56:47EXCEL的TTEST求出t值,怎么算P?TTEST(array1,array2,tails,type) 返回与 t 检验相概率。Array1 为第一个数据集。Array2 为第二个数据集。Tails 指示分布曲线的尾数。如果 tails = 1, TTEST 使用单尾分布果 tails =... -
终于有人把 p 值讲明白了!
2021-09-11 00:30:32导读:p值(P value)就是当原假设为真时,比所得到的样本观察结果更极端的结果出现的概率,是用来判定假设检验结果的一个参数。p值是根据实际统计量计算出的显著性水平。本文带你了解p值和对... -
t检验的p值对照表_第十讲 R-两独立样本t检验
2020-11-21 02:29:09两独立样本t检验用于比较的平均两个独立的组是否存在差异。例如,假设我们测量了100个人的体重:50名女性(A组...注意,仅当数据呈正态分布时,才可以使用两独立样本t检验。可以使用Shapiro-Wilk test进行检查。(请... -
显著性检验:P值和置信度
2021-05-26 16:21:17显著性差异(ρ,Statistical significance) 是统计学上对数据差异性的评价。...在统计学中,显著性检验是“统计假设检验”(Statistical hypothesis tesing)的一种,显著性检验是检测科学实验中的实验组与对照. -
统计学中的P值应该怎么计算
2020-12-20 22:02:58谁浮夸了年华2019-12-10 18:363370P 值即概率,反映某一事件发生的可能性大小。统计学根据显著性检验方法所得到的P 值,一般以P < 0.05 为显著, P <0.01 为非常显著,其含义是样本间的差异由抽样误差所致的... -
c语言利用指针求一组数的最大值,最小值。平均值
2018-01-06 18:04:25最近一直在学算法,c语言的指针在算法中用的不多,所以就没怎么学,直到后来帮我同学做课程设计的时候,才学了一段时间,为了防止过段时间忘了指针,所以就把课程设计保存下来,以后忘了还可以看看这个回忆一下指针... -
python 如何判断一组数据是否符合正态分布
2020-08-06 17:51:26正态分布: 若随机变量x服从有个数学期望为μ,方差为σ2的正态分布,记为N(μ,σ) 其中期望值决定密度函数的位置,标准差决定分布的幅度,当υ=0,σ=0 时的...#构造一组随机数据 s = pd.DataFrame(np.random.randn. -
P值 卡方值
2019-10-17 13:01:48P值即概率,反映某一事件发生的可能性大小。 不同的P数值所表达的含义也是不一样的。 统计学根据显著性检验方法所得到的P 值,一般以P < 0.05为有统计学差异, P<0.01 为有显著统计学差异,P<0.001为有... -
通过给定一组数据点并反求控制点的NURBS曲线插值生成Matlab编程实例
2020-08-02 23:10:12非均匀有理B样条(NURBS)是几何造型的一个十分有用的工具,以其良好的特性越来越受到工程界的重视。它可以用统一的数学形式表示规则曲面与自由曲面;有操纵控制顶点及权因子,为各种形状设计提供了充分的灵活性,有强... -
【Python计算检验值】一元线性回归拟合,t值与p值,显著性检验
2018-04-24 12:27:14本文主要讨论Python实现一元回归的线性拟合、最小二乘法估计回归参数和显著性检验(t检验和p值). 一元线性回归模型是描述两个变量之间相关关系的最简单的回归模型. 通常人们对所要研究的问题首先要收集与它有关的n... -
如何判断一组数据是否符合正态分布呢?
2021-03-04 14:21:28在很多模型及假设检验中都需要满足一个假设条件:数据需服从正态分布。这篇文章主要讲讲如何判断数据是否符合正态分布。主要分为两种方法:描述统计方法和统计检验方法。 描述统计方法 描述统计就是用描述的数字或... -
用MATLAB判断一组数据是否符合正态分布
2021-04-18 03:49:47问题描述有一组已知的数据,例如为:1006.1,1014,1001.6,996.4,997.8,981.6,996.4,991.9,993.3,1000.6,987.3,1015.6,981.6,996.2,999.2,994.5,1005.9,1001.9,986.4,1007.6,1001.4,1014.6,1010.2,993.9,1001.4这组... -
6-6 求一组数中的最大值、最小值和平均值 (10分)
2019-12-25 15:38:30编写函数,求一组数中的最大值、最小值和平均值。 函数接口定义: float fun(int a[],int n,int *max,int *min); 其中 a、n、max 和 min 都是用户传入的参数。函数求a数组中n个元素的最大值、最小值和平均值... -
S7-1200PLC求数组里数据最大值最小值FB块
2019-12-13 08:27:38最近做一个检测设备,需要知道任意大小的数组里的数据最大值和最小值,如果数组里元素个数少可以用梯形图一个个比较,但是如果元素个数多了这种方式就比较麻烦了,所以我用SCL语言开发了一个标准的FB块如下,需要... -
p值小于0.05拒绝还是接受_干货:关乎你的实验成败,0.05这个值不容小觑!
2020-11-21 03:54:09大部分科研搬砖者们都会进行假设检验,求算出P值,如果P值小于0.05,我们就说两者之间有显著性差异。那么你真的了解P值君,到底是啥吗?下面和小编一起走进统计学的世界吧,让你的数据分析地更有理有据,文章看起来... -
【特征工程】判断一组数据的分布形态
2019-11-12 17:42:31思考:输入到NN模型中的特征要做归一化处理...言归正传,airbnb根据不同特征做不一样的归一化,因为他们对数据进行了观察,发现了部分长尾数据,因此做了log的归一化处理[1],这点很惊喜。在我刚工作的时候,也有人... -
qRT-PCR差异分析及P值计算
2020-12-30 08:07:39qRT-PCR是一种相对表达定量的方法,他的计算方法有很多,常用的相对定量数据分析方法是KJ Livak(Applied Biosystems)等人在2001年提出的“比较Ct法相对定量”,即:利用ΔCt值差异来推算基因表达差异(Ct目的基因 – ... -
两组数据相关性分析,找到其中一组的最优值
2018-08-24 02:36:13现有两组数据长度L和弹性E,分析其相关性,想找到长度L为多少时,弹性值E最大,在graphpad中得到结果是线性相关,却无法给出想要的L最优值。用什么统计方法能解决此问题呢? Number of XY Pairs 5332 Pearson r -0.... -
如何使用 一行代码 搞定一组数据的(极值、平均值、中位数、四分位数、数量统计和标准差)
2020-02-28 21:33:34或许当你看到一行代码的...)就像进行数据处理的时候,有时会遇到求极值(最大值、最小值)、平均值、中位数和四分位数(25%、 75%)的情况。 这一篇博客就是你的福音,让你绝对0基础使用python 进行数据分析。 ... -
p值 统计学意义_什么是统计意义? P值定义以及如何计算
2020-08-16 17:43:45p值 统计学意义P values are one of the most widely used concepts in statistical analysis. They are used by researchers, analysts and statisticians to draw insights from data and make informed decisions... -
p值还是 FDR
2018-03-24 22:15:44p值:在零假设下,观测到某一特定实验结果的概率称为p值。假阳性:得到了阳性结果,但这个阳性结果是假的。假阴性:得到了阴性结果,但这个阴性结果是假的单次检验:针对单个基因(蛋白),采用统计检验,假设采用的... -
手把手教你用Graphpad做单因素方差分析
2021-03-07 22:16:450 2首先打开Graphpad软件,输入数据,更改组名,更改Data1为Total distance0 3进行数据正态性检验:选中数据 – Analyze - Column analyses - Normality and Log normality Tests -选中组别-OK0 4进行正态性检验软件... -
Python数据分析基础: 数据缺失值处理
2020-10-31 21:56:01作者:东哥起飞 公众号:Python数据科学 圣人曾说过:数据和特征决定了机器学习的上限,而模型和算法只是逼近这个上限而已。 再好的模型,如果没有好...本篇我们来说说面对数据的缺失值,我们该如何处理。文末有. -
【数理统计】神奇的P值
2019-06-29 13:46:21工作中经常会通过AB Test帮助做产品决策,简单说就是为产品制作两个(A/B)或多个(A/B/C/...)版本,在同一时间维度,分别让不同组的用户群随机的访问这些版本,收集各群组的用户的数据,最后分析评估出最好版本... -
python数据预处理之异常值、缺失值处理方法
2020-05-03 20:23:27数据预处理是明确分析目标与思路之后进行数据分析的第一步,也是整个项目中最基础、花费时间较长的工作。除了互联网埋点的数据或企业内部的业务数据之外,往往我们拿到的,比如说网上采集的数据并不是那样规整,这类...