1. 背景 之前寫(xiě)過(guò)一篇文章,介紹多重比較的可視化:方差分析多重比較可視化 。
高廣雄同學(xué)問(wèn)我如何增加趨勢(shì)線,以及調(diào)整顯示的順序。
我不了解如何加趨勢(shì)線,調(diào)整順序的話可以設(shè)置factor的level,然后建議他研究一下,如果成功的話,可以把經(jīng)驗(yàn)寫(xiě)出來(lái),我轉(zhuǎn)發(fā)一下。
今天他給我了一個(gè)鏈接,上面是他實(shí)現(xiàn)的過(guò)程,效果如下:
本來(lái)直接就想用他給的代碼,水一篇。但是?。?!
? 分享的是有道云筆記的小程序,沒(méi)有辦法在電腦端打開(kāi),只能在手機(jī)端查看,我太難了。^[所以,發(fā)給別人文章時(shí),極客用md文件,正規(guī)的是word或者pdf,不靠譜的就是有道小程序了。。。]
? 容易的路不好走,那我就走一條難走的路 ^[然后,又做起了編程的大招:粘貼復(fù)制,修改代碼]。另外,我覺(jué)得上面的示例不夠美觀,好的例子至關(guān)重要,壞的例子南轅北轍。
我想到了數(shù)據(jù)分析師的格言:沒(méi)有數(shù)據(jù),創(chuàng)造數(shù)據(jù)也要上!沒(méi)有理想結(jié)果,多試幾次就能達(dá)到,所以,關(guān)于圖片誤用,一切皆有可能,哈哈。
如果一個(gè)程序員,表情平靜,鍵盤(pán)飛快
? 他一定在寫(xiě)注釋
? 如果一個(gè)程序員,鍵盤(pán)很少動(dòng),鼠標(biāo)飛快
? 他一定在粘貼復(fù)制別人的代碼
? 如果一個(gè)程序員,抓耳撓腮,時(shí)而百度,時(shí)而有道云筆記,好久才動(dòng)一次鍵盤(pán)
? 他應(yīng)該在寫(xiě)程序,磕磕絆絆,連抄帶借鑒
? 2. 模擬數(shù)據(jù) set.seed(123) library(tidyverse) dat = data.frame("Mon" = rnorm(10)+1,"Tue" = rnorm(10)+10,"Wed" = rnorm(10) + 5, "Thu" = rnorm(10)+15,"Fri" = rnorm(10) +4) dat
數(shù)據(jù)如下:
> dat Mon Tue Wed Thu Fri 1 0.4395244 11.224082 3.932176 15.42646 3.305293 2 0.7698225 10.359814 4.782025 14.70493 3.792083 3 2.5587083 10.400771 3.973996 15.89513 2.734604 4 1.0705084 10.110683 4.271109 15.87813 6.168956 5 1.1292877 9.444159 4.374961 15.82158 5.207962 6 2.7150650 11.786913 3.313307 15.68864 2.876891 7 1.4609162 10.497850 5.837787 15.55392 3.597115 8 -0.2650612 8.033383 5.153373 14.93809 3.533345 9 0.3131471 10.701356 3.861863 14.69404 4.779965 10 0.5543380 9.527209 6.253815 14.61953 3.916631
這個(gè)數(shù)據(jù),需要變?yōu)殚L(zhǎng)數(shù)據(jù) ^[就是ggplot2作圖用的數(shù)據(jù)]
dd = pivot_longer(dat,1:5,values_to = "y",names_to = "Day") dd
長(zhǎng)數(shù)據(jù)來(lái)到:
> dd # A tibble: 50 x 2 Day y <fct> <dbl> 1 Mon 1.38 2 Tue 9.51 3 Wed 5.01 4 Thu 16.0 5 Fri 3.29 6 Mon 0.498 7 Tue 7.69 8 Wed 5.39 9 Thu 15.5 10 Fri 4.26 # ... with 40 more rows
3. 方差分析 因?yàn)槟M的數(shù)據(jù),我們直接用aov
進(jìn)行方差分析,直接用agricolae
進(jìn)行多重比較吧,就不整那些正態(tài)性檢驗(yàn),齊性檢驗(yàn)了。
注意,這里要講Day
轉(zhuǎn)化為因子!
str(dd) dd$Day = as.factor(dd$Day) mod = aov(y ~ Day, data=dd) summary(mod)
方差分析結(jié)果:
> summary(mod) Df Sum Sq Mean Sq F value Pr(>F) Day 4 1270.5 317.6 417.1 <2e-16 *** Residuals 45 34.3 0.8 --- Signif. codes: 0 '***’ 0.001 '**’ 0.01 '*’ 0.05 '.’ 0.1 ' ’ 1
可以看到,Day達(dá)到極顯著,可以進(jìn)行多重比較。
4. 多重比較 這里,對(duì)Day
進(jìn)行多重比較,方法是LSD
,水平是0.05
library(agricolae) re = LSD.test(mod,"Day",alpha = 0.05) re1 = re$groups re1
多重比較結(jié)果:
> re1 y groups Thu 15.437094 a Tue 9.637082 b Wed 5.313095 c Fri 3.604152 d Mon 1.123084 e
可以看到,周四最高,齊次是周二,它們之間均達(dá)到顯著性水平!
5. 多重比較作圖 5.1 數(shù)據(jù)準(zhǔn)備為何還要準(zhǔn)備數(shù)據(jù)?因?yàn)樾枰獪?zhǔn)備數(shù)據(jù)作圖。
? 作為數(shù)據(jù)分析師,我們80%的時(shí)間,都是在準(zhǔn)備數(shù)據(jù)。。。
? 整體思路:
xx = aggregate(y ~ Day, dd,sd) names(xx) = c("Day","sd") xx re2 = re1 %>% mutate(Day = rownames(re1)) %>% inner_join(.,xx,by="Day") re2
結(jié)果:
> re2 y groups Day sd 1 15.437094 a Thu 1.0688850 2 9.637082 b Tue 0.9955075 3 5.313095 c Wed 0.5475116 4 3.604152 d Fri 0.7055487 5 1.123084 e Mon 0.9363923
5.2 做柱形圖## 做直方圖 p1 = re2 %>% ggplot(aes(Day,y)) + geom_col(aes(fill = Day), width=.4) p1
5.3 在每個(gè)柱子上加標(biāo)準(zhǔn)誤## 加上標(biāo)準(zhǔn)誤 p2 = p1 + geom_errorbar(aes(ymax = y + sd, ymin = y - sd),width = .1,size=.5) p2
5.4 添加顯著性## 加上顯著性 p3 = p2 + geom_text(aes(label = groups,y = y + sd +0.5)) p3
5.5 背景設(shè)置為空白為何要設(shè)置空白背景?ggplot默認(rèn)的風(fēng)格不香嗎?
? 空白的背景顯得比較專業(yè)。。。
? ## 背景設(shè)置為空白 p4 = p3 + theme(panel.grid = element_blank(), panel.background = element_rect(color = "black",fill = "transparent")) p4
5.6 設(shè)置標(biāo)簽? 廣告無(wú)處不在,有一種廣告叫做鑲嵌到圖中!
? ## 增加標(biāo)簽 p5 = p4 + labs(x = "星期幾?",y = "表現(xiàn)如何?",title = "公眾號(hào):育種數(shù)據(jù)分析之放飛自我") + theme(plot.title = element_text(hjust = 0.5)) p5
6. 調(diào)整順序 本來(lái),應(yīng)該是盡善盡美的,但是如果我想將順序按照:Mon,Tue,Wed,Thu,F(xiàn)ri進(jìn)行排布,應(yīng)該如何處理?
用factor函數(shù),重新定義一下水平:
str(re2) re2$Day = factor(re2$Day,levels = c("Mon","Tue","Wed","Thu","Fri")) re2
「看一下作圖效果:」
p1 = re2 %>% ggplot(aes(Day,y)) + geom_col(aes(fill = Day), width=.4) p2 = p1 + geom_errorbar(aes(ymax = y + sd, ymin = y - sd),width = .1,size=.5) p3 = p2 + geom_text(aes(label = groups,y = y + sd +0.5)) p4 = p3 + theme(panel.grid = element_blank(), panel.background = element_rect(color = "black",fill = "transparent")) p5 = p4 + labs(x = "星期幾?",y = "表現(xiàn)如何?",title = "公眾號(hào):育種數(shù)據(jù)分析之放飛自我") + theme(plot.title = element_text(hjust = 0.5)) p5
搞定!
7. 增加趨勢(shì)線 為何要增加趨勢(shì)線呢?
? 因?yàn)榧由馅厔?shì)線,比較學(xué)術(shù)!
? p6 = p5 + geom_line(aes(group=""),color = "red") p6