此文內(nèi)容首發(fā)于微信公眾號：R語言搬運工，關(guān)注公眾號瀏覽更多精彩內(nèi)容

原文鏈接

image

今天在Github上看到了一個有意思的包墓猎，作者IndrajeetPatil集合開發(fā)的ggstatsplot包楚午，在繪圖方面將不同的圖和統(tǒng)計檢驗結(jié)果整合在了一起刚盈，使得視覺效果更加直觀，實現(xiàn)了一圖勝多言的效果。以下內(nèi)容主要基于作者IndrajeetPatil發(fā)布的關(guān)于包ggstatsplot的詳細介紹說明丈屹，在此表示感謝。
Patil, I. (2021). Visualizations with statistical details: The
'ggstatsplot' approach. Journal of Open Source Software, 6(61), 3167,
doi:10.21105/joss.03167

由于包里集合多種統(tǒng)計方法伶棒，內(nèi)容非常豐富旺垒，現(xiàn)簡單介紹其中第一部分：ggbetweenstats函數(shù)（一個整合了箱型圖、小提琴圖和統(tǒng)計檢驗方法結(jié)果的函數(shù)肤无。

PS：想要弄清楚這個包返回的結(jié)果先蒋，必須要對目前主要的統(tǒng)計檢驗分析方法有一個相對較熟悉的了解。

使用的數(shù)據(jù)為包內(nèi)所含電影分類時長的數(shù)據(jù)：

（https://github.com/IndrajeetPatil/ggstatsplot）

image.png

數(shù)據(jù)大體是說不同類型電影的時長

通過繪圖來看Action和Comedy類型電影中宛渐，mpaa類型電影的時長分布模式竞漾，和mpaa分組之間的統(tǒng)計檢驗差異性。

先上代碼和圖窥翩，隨后解釋~~~~~~~

# for reproducibility
set.seed(123)

# plot
grouped_ggbetweenstats(
  data = dplyr::filter(movies_long, genre %in% c("Action", "Comedy")),
  x = mpaa,
  y = length,
  grouping.var = genre, # grouping variable
  outlier.tagging = TRUE, # whether outliers need to be tagged
  outlier.label = title, # variable to be used for tagging outliers
  outlier.coef = 2,
  ggsignif.args = list(textsize = 4, tip_length = 0.01),
  p.adjust.method = "hommel", # method for adjusting p-values for multiple comparisons
  # adding new components to `ggstatsplot` default
  ggplot.component = list(ggplot2::scale_y_continuous(sec.axis = ggplot2::dup_axis())),
  caption = substitute(paste(italic("Source"), ": IMDb (Internet Movie Database)")),
  palette = "default_jama",
  package = "ggsci",
  plotgrid.args = list(nrow = 1),
  annotation.args = list(title = "Differences in movie length by mpaa ratings for different genres")
)

image.png

圖形中业岁，一眼便可以看到是使用了箱型圖和小提琴圖繪制了數(shù)據(jù)的分布模式，這個在ggplot中可以很容易實現(xiàn)寇蚊，但是統(tǒng)計檢驗方法結(jié)果在圖中給出了詳細的展示叨襟。（看左邊半圖，紅線框出）

1. 首先進行了參數(shù)檢驗幔荒，看一下三個組之間均值的差異性糊闽，并給出了檢驗參 數(shù)、置信區(qū)間等爹梁。

2. 然后對三個組進行兩兩比較右犹，比較方法為Games-Howell，矯正p值使用Games-Howell方法姚垃。通過添加倒扣的中括號進行詮釋（只展示顯著性的組）念链。

3. 標記了異常值等數(shù)據(jù)。

4.給出了貝葉斯假設(shè)檢驗結(jié)果积糯。

是不是非常的全面掂墓，可能對于新手還不是很能理解】闯桑看一下具體的展示說明

image.png

圖形中以一個公式的形式展示統(tǒng)計檢驗結(jié)果君编，具體的參數(shù)意義包的開發(fā)者在網(wǎng)站上也已給出（如上圖）。

image.png

image.png

分組的數(shù)量不同，在檢驗差異性時使用的方法是不一樣的兑燥，具體的可以查看上圖亮瓷。

其次看一下數(shù)據(jù)符不符合統(tǒng)計假設(shè)條件：

image.png

在數(shù)理統(tǒng)計中，為檢驗不同分組間數(shù)據(jù)的差異性通常根據(jù)數(shù)據(jù)的分布形式使用不同的方法降瞳，當數(shù)據(jù)符合正態(tài)分布模式且方差齊性條件滿足的時候嘱支，通常選擇ANOVA分析方法進行檢驗。如果條件不滿足就要根據(jù)具體情況再選擇方法了挣饥，比如使用非參數(shù)檢驗除师。但是，也不一定亮靴！

當數(shù)據(jù)不符合正態(tài)分布且齊性條件滿足的時候其實也可以使用方差檢驗馍盟，因為方差分析對數(shù)據(jù)的正態(tài)性并不是十分的敏感，同時可以考慮使用非參數(shù)檢驗茧吊。但是當數(shù)據(jù)也不符合方差齊性的時候這時候可以使用另一個檢驗方法welck檢驗贞岭，當你的實驗設(shè)計為混合設(shè)計時，數(shù)據(jù)不符合正態(tài)分布和方差齊性搓侄，你可以考慮使用robust分析.

在ggstatsplot包中瞄桨，作者也都集成了不同的參數(shù)用來設(shè)置：

ggbetweenstats(
data,
x,
y,
plot.type = "boxviolin",
type = "parametric",
pairwise.comparisons = TRUE,
pairwise.display = "significant",
p.adjust.method = "holm",
effsize.type = "unbiased",
bf.prior = 0.707,
bf.message = TRUE,
results.subtitle = TRUE,
xlab = NULL,
ylab = NULL,
caption = NULL,
title = NULL,
subtitle = NULL,
k = 2L,
var.equal = FALSE,
conf.level = 0.95,
nboot = 100L,
tr = 0.2,
centrality.plotting = TRUE,
centrality.type = type,
centrality.point.args = list(size = 5, color = "darkred"),
centrality.label.args = list(size = 3, nudge_x = 0.4, segment.linetype = 4,
min.segment.length = 0),
outlier.tagging = FALSE,
outlier.label = NULL,
outlier.coef = 1.5,
outlier.shape = 19,
outlier.color = "black",
outlier.label.args = list(size = 3),
point.args = list(position = ggplot2::position_jitterdodge(dodge.width = 0.6), alpha
= 0.4, size = 3, stroke = 0),
violin.args = list(width = 0.5, alpha = 0.2),
ggsignif.args = list(textsize = 3, tip_length = 0.01),
ggtheme = ggstatsplot::theme_ggstatsplot(),
package = "RColorBrewer",
palette = "Dark2",
ggplot.component = NULL,
output = "plot",
...
)

通過查看幫助文件，可以根據(jù)自己的要求選擇不同的方法讶踪，比如type參數(shù)：

image.png

image.png

當測試兩個分組之間差異性的時候同樣可以使用該包進行繪制芯侥，這時候涉及的統(tǒng)計方法將更加豐富：

使用的數(shù)據(jù)為包gapminder提供的數(shù)據(jù)集，該數(shù)據(jù)集提供了從 1952 年到 2007 年 142 個國家中每一個國家的預(yù)期壽命乳讥、人均 GDP 和人口的值柱查，每 5 年一次 （Gapminder Foundation）。

image.png

假設(shè)我們要檢查的第一件事是 2007 年一個大陸國家的預(yù)期壽命分布云石。我們還想知道大陸之間預(yù)期壽命的平均差異是否具有統(tǒng)計顯著性唉工。

# since the confidence intervals for the effect sizes are computed using
# bootstrapping, important to set a seed for reproducibility
set.seed(123)

# function call
ggstatsplot::ggbetweenstats(
  data = dplyr::filter(gapminder::gapminder, year == 2007, continent != "Oceania"),
  x = continent,
  y = lifeExp
)

image.png

從圖中可以看出，該函數(shù)默認返回貝葉斯因子進行測試汹忠。如果零假設(shè)顯著性檢驗 (NHST) 方法無法拒絕零假設(shè)淋硝，則貝葉斯方法可以幫助索引支持零假設(shè)的證據(jù)（即BF01）

到目前為止，我們只使用了經(jīng)典的參數(shù)測試和箱形小提琴圖宽菜，但我們也可以使用其他可用選項：

type參數(shù)包含了下列縮寫："p"（為參數(shù)檢驗）谣膳， "np"（為非參數(shù)檢驗）， "r"（robust檢驗）铅乡， "bf"（貝葉斯因子）继谚。

顯示的圖型（"box"，"violin"隆判，或"boxviolin"）犬庇。

可以修改調(diào)色板僧界。
首先侨嘀，使用參數(shù)檢驗方法檢驗兩組之間的差異性臭挽，使用的是welch統(tǒng)計方法

# for reproducibility
set.seed(123)
library(ggstatsplot)
library(gapminder)

# selecting subset of the data
df_year <- dplyr::filter(gapminder::gapminder, year == 2007 | year == 1957)

# parametric t-test and box plot
p1 <- ggstatsplot::ggbetweenstats(
  data = df_year,
  x = year,
  y = lifeExp,
  xlab = "Year",
  ylab = "Life expectancy",
  plot.type = "box",
  type = "p",
  conf.level = 0.99,
  title = "Parametric test",
  package = "ggsci",
  palette = "nrc_npg"
)

image.png

其次，使用非參數(shù)檢驗方法咬腕，使用的是曼-惠特尼U檢驗方法

p2 <- ggstatsplot::ggbetweenstats(
  data = df_year,
  x = year,
  y = lifeExp,
  xlab = "Year",
  ylab = "Life expectancy",
  plot.type = "violin",
  type = "np",
  conf.level = 0.99,
  title = "Non-parametric Test (violin plot)",
  package = "ggsci",
  palette = "uniform_startrek"
)

image.png

再次欢峰，使用穩(wěn)健性分析（robust），檢驗方法為Yuen 的修整均值檢驗

p3 <- ggstatsplot::ggbetweenstats(
  data = df_year,
  x = year,
  y = lifeExp,
  xlab = "Year",
  ylab = "Life expectancy",
  plot.type = "boxviolin",
  type = "r",
  conf.level = 0.99,
  title = "Robust Test (box & violin plot)",
  tr = 0.005,
  package = "wesanderson",
  palette = "Royal2",
  k = 3
)

image.png

最后使用貝葉斯方法涨共，

p4 <- ggstatsplot::ggbetweenstats(
  data = df_year,
  x = year,
  y = lifeExp,
  xlab = "Year",
  ylab = "Life expectancy",
  type = "bayes",
  plot.type = "box",
  title = "Bayesian Test (box plot)",
  package = "ggsci",
  palette = "nrc_npg"
)

image.png

image.png

精彩推薦：

R語言繪制散點圖geom_point
R語言添加擬合曲線geom_smooth
R語言箱線圖boxplot
R語言線圖geom_line