最近考种数据回来了,今年测得单株,本想用箱线图表示,后来一并把小提琴图也做了,顺便做个笔记吧。
利用SUMIF对相同rowname数据合并
首先遇到了一个问题,由于之前在收花时分了好几批收的,结果在最后考种统计时存在一个编号对应了好几个数据,所以第一步要把这些数据求和。
一定要相信推动科技进步的永远是那些“懒汉”们,如果你比较勤快,可能会一个一个加起来在总结,但是懒汉们不会,于是就有了SUMIF
函数。 首先看下SUMIF
函数是个什么玩意。
SUMIF函数是Excel常用函数。使用 SUMIF 函数可以对报表范围中符合指定条件的值求和。Excel中sumif函数的用法是根据指定条件对若干单元格、区域或引用求和。
SUMIF(range,criteria,sum_range)
1)range 为用于条件判断的单元格区域。
2)criteria 为确定哪些单元格将被相加求和的条件,其形式可以为数字、文本、表达式或单元格内容。例如,条件可以表示为 32、"32"、">32" 、"apples"或A1。条件还可以使用通配符:问号 (?) 和星号 (* ),如需要求和的条件为第二个数字为2的,可表示为"?2*",从而简化公式设置。
3)sum_range 是需要求和的实际单元格。
- 第一步就是先复制编号到另外一列,(不要和已有的数据挨着,中间空几列),然后去重,就剩下唯一的编号。
- 然后就按照上面
SUMIF
的语法写自己的公式了,我的公式:=SUMIF($A$2:$A$802,$F2,B$2:B$504)
。dolor符号$
在excel中表示绝对引用,因为后面要直接往下拉,同时要找出其他条目还要横向拉,excel在拉的时候对应的单元格也会自动调整,但是我要找的范围其实是固定的从A2单元格到A802单元格,这个范围是我们全部的编号区域,
- 所以为了避免在下拉或者向右拉的时候这个区域改变我们限定不管你函数跑到哪个单元格
range
的范围总是在A2到A802,所以range
是$A$2:$A$802
; - 然后是
criteria
也就是我们要找的唯一值,也就是刚才去重的编号,在这里我们需要纵向拉的时候单元格随着变化,但是横向拉不变,所以用$F
限定criteria数据永远在F列,当你向下拉的时候随着你公式位置的变化,criteria随着公式的位置变化成F2,F3,F...。注意这里不能用$限定2! 原因..没必要解释了吧; - 最后是
sum_range
,意思就是你要求和的数据,因为我这里统计了铃数,籽棉重,皮棉重三个数据,所以在做完铃数后向右一拉籽棉重和皮棉重的数据也出来了,但是前提是范围仍旧是2-504这个大范围,向右拉的时候要求求和数据列从B变为C再到D,所以我用$
限定了2-504的范围,而没有限定列,让列名随着函数位置的改变随着改变,所以sum_range
部分为B$2:B$504
。
自己胡编了个数据做了个动图:
用ggplot2画箱线图和小提琴🎻图
-
ggplot2要求的数据输入方式
额,数据没发表所以有码图如下,
首先看下ggplot的图形语法。
尽管qplot作为ggplot2的快速作图(quick plot)函数, 能够极大的简化作图步骤, 容易入门和上手, 但是qplot却不是泛型函数, 而ggplot()作为泛型函数, 能对任意类型的R对象进行可视化操作, 是ggplot2的精髓所在, 因而在本文中主要的绘图都是通过ggplot()来完成的。有关于qplot的介绍可以细看Hadley的官方介绍。
在Hadley的ggplot2官方文档中, Hadely这样对Wilkinson的图形语法进行了描述:“一张统计图形就是从数据到几何对象(geometric object, 缩写为geom, 包括点、线、条形等)的图形属性(aesthetic attributes, 缩写为aes, 包括颜色、形状、大小等)的一个映射。此外, 图形中还可能包含数据的统计变换(statistical transformation, 缩写为stats), 最后绘制在某个特定的坐标系(coordinate system, 缩写为coord)中, 而分面(facet, 指将绘图窗口划分为若干个子窗口)则可以用来生成数据中不同子集的图形。”因此在ggplot2中, 图形语法中至少包括了如下几个图形部件:
- 数据(data)
- 映射(mapping)
- 几何对象(geom)
- 统计变换(stats)
- 标度(scale)
- 坐标系(coord)
- 分面(facet)
这些组件之间是通过“+”, 以图层(layer)的方式来粘合构图的, 所以图层是ggplot2中一个重要的概念。当然, 在掌握基本的图形部件基础上, 要完成一幅高质量的统计绘图, 仍然需要其他图形部件来进一步扩展, 这包括了:- 主题(theme)
- 存储和输出
对于数据,ggplot规定了必须是dataframe的形式,
#数据输入,将全部数据先导进来
phenotype = read.table("2018phenotype.txt", header = T, sep = "\t")
phenotype[1:4,1:4]
rownames(phenotype)=phenotype[,1]
phenotype=phenotype[,-1]
head(phenotype)
#将铃数的数据新建一个数据框
MBN = data.frame(BN = phenotype$Ball_Num.,
cultivar = phenotype$cultivar)
这里将所有数据导入后想用铃数画箱线图,那么要把铃数的数据先拎出来,此外还得“标注”清楚那个这些铃重数据来自哪个品种。那么就新建一个数据框,一列为ball number,BN
,另一列就是cultivar
。
看下ggplotmanual中的例子
base <- ggplot(mpg, aes(displ, hwy)) + geom_point()
base + geom_smooth()
# To override the data, you must use %+%
base %+% subset(mpg, fl == "p")
# Alternatively, you can add multiple components with a list.
# This can be useful to return from a function.
base + list(subset(mpg, fl == "p"), geom_smooth())
- Description
Aesthetic mappings describe how variables in the data are mapped to visual properties (aesthetics) of geoms. Aesthetic mappings can be set in ggplot2() and in individual layers.- Usage
aes(x, y, ...)- Arguments
x, y, ...- Details
List of name value pairs giving aesthetics to map to variables. The names for x and y aesthetics are typically omitted because they are so common; all other aesthetics must be named.
This function also standardises aesthetic names by converting color to colour (also in substrings, e.g. point_color to point_colour) and translating old style R names to ggplot names (eg. pch to shape, cex to size).- Value
A list with class uneval. Components of the list are either quosures or constants. Quasiquotation
aes() is a quoting function. This means that its inputs are quoted to be evaluated in the context of the data. This makes it easy to work with variables from the data frame because you can name those directly. The flip side is that you have to use quasiquotation to program with aes(). See a tidy evaluation tutorial such as the dplyr programming vignette to learn more about these techniques.- See Also
vars() for another quoting function designed for faceting specifications.
其实主体就是:
base <- ggplot(mpg, aes(displ, hwy))
其中mpg
可以理解为我们画图数据新建的data.frame,displ
为x轴,hwy
为y轴数据。所以我们图的主体就是:
pMBN <- ggplot(MBN, aes(x=cultivar, y=BN, fill=cultivar))
# fill = cultivar 表示颜色按照cultivar分类填充。
然后通过+
对图形进行调整。
library(RColorBrewer)
colourCount = length(unique(mtcars$hp))
getPalette = colorRampPalette(brewer.pal(9, "Set3"))
pMBN <- ggplot(MBN, aes(x=cultivar, y=BN, fill=cultivar))+
geom_boxplot()+
scale_fill_manual(values = getPalette(colourCount))+
labs(title = "Ball number per plant of each cultivar")+
theme(plot.title = element_text(hjust = 0.5),
axis.title = element_text(size = 10, color = "blue", face = "bold"),
axis.text.x = element_text(size = 8, angle = 45, hjust = 1, vjust = 1, face = "bold"))
其中发生了一点小意外,就是在用scale_fill_brewer()
函数调用调色板给box上色的时候发现超过8个后都成了灰色,布吉岛为啥,百度后发现:
有3类调色板(palettes)- sequential, diverging, and qualitative - 每一类调色板包含8到12种颜色(可以利用
brewer.pal.info
或者?RColorBrewer
看到更多的细节)。
好奇的读者可能注意到如果柱状图包含13或者更多的柱子,就可能出现显示不出更多颜色,
的确,length(unique(mtcars$hp))
有 22 个唯一值,然而调色板 Set2 只有 8 中颜色。调色板中的颜色缺乏导致了ggplot如下的警告:
- Warning message:
In RColorBrewer::brewer.pal(n, pal) : n too large, allowed maximum for palette Set2 is 8
Returning the palette you asked for with that many colors- RColorBrewer为我们提供了一种通过使用构造函数colorRampPalette插入现有调色板来生成更大调色板的方法。它生成实际工作的函数:它们通过插入现有的调色板来构建具有任意数量颜色的调色板。
总结就是用Set2
或者Set3
中的颜色超过8个可能就GG了,所以用RColorBrewer
包扩充Set2
或者Set3
的颜色,所以相应的用scale_fill_manual()
代替了scale_fill_brewer()
。前面的geom_boxplot()
直接简单理解成画了个箱线图把,而后面的labs()
可以看出来是加标题的,同样用xlab()
加x轴标题,ylab()
加y轴标题,theme()
则是对图的美化,比如字体,颜色啥的,偶然发现一个思路可以事先新建好一个自己喜欢的theme,然后每次作图的时候直接+mytheme
于是就模仿着画了个。
##自己的主题1
mytheme <- theme(plot.title = element_text(face = "bold.italic",
size = "14", color = "brown"),
axis.title = element_text(face = "bold.italic",
size = "10",color = "blue"),
axis.text.x = element_text(face = "bold",
size = 9, angle = 45, hjust = 1, vjust = 1),
panel.background = element_rect(fill = "white",
color = "black"),
panel.grid.major.y = element_line(color = "black",
linetype = 2),
panel.grid.minor.y = element_line(color = "black",
linetype = 2),
panel.grid.minor.x = element_blank())
#自己的主题2
mytheme2 <- theme(plot.title = element_text(face = "bold.italic",
size = "14", color = "brown"),
axis.title = element_text(face = "bold.italic",
size = "10",color = "blue"),
axis.text.x = element_text(face = "bold",
size = 8, angle = 45, hjust = 1, vjust = 1),
axis.text.y = element_text(face = "bold",size = 8),
panel.background = element_rect(fill = "white",
color = "black"),
panel.grid.major.y = element_blank(),
panel.grid.minor.y = element_blank(),
legend.text = element_text(size = 8),
legend.title = element_text(size = 10,
face = "bold"),
panel.grid.minor.x = element_blank())
mythemepolor <- theme(plot.title = element_text(face = "bold.italic",
size = "14", color = "brown", hjust = 0.5),
axis.title = element_text(face = "bold.italic",
size = "10",color = "blue"),
axis.text.x = element_text(face = "bold", color = "purple",
size = 9),
panel.background = element_rect(fill = "white",
color = "black"),
panel.grid.major.y = element_line(color = "black",
linetype = 2),
panel.grid.minor.y = element_line(color = "black",
linetype = 2),
panel.grid.minor.x = element_blank())
一共搞了3个theme,第一个属于有y轴标线的,第二个为空白背景的,第三个是画极坐标的,有了这些,就开始画图
## 直接用新主题加前面的图,覆盖掉前面的主题
pMBNfinal = pMBN_ylab + mytheme2
## 保存成pdf
ggsave(file="ball number per plant of each cultivar boxplot.pdf",plot=pMBNfinal, width = 12, height = 6)
看看加主题和没加主题的效果吧:
- 小提琴图的画法
根据上面画boxplot的思路,小提琴的画法如下
vMBN <- ggplot(MBN, aes(x=cultivar, y=BN)) +
geom_violin(fill="lightblue") +
geom_boxplot(fill="lightgreen", width= .2) +
labs(title = "Ball number per plant of each cultivar") +
ylab("Number")
vMBNfinal = vMBN + mytheme2
当然这里不是一个单纯的小提琴图,因为geom_boxplot(fill="lightgreen", width= .2)
加入了箱线图,如果不想要箱线图,直接去掉这一行,并且补充fill
参数因为单纯的箱线图可以提供的信息是最大值,最小值,上四分位数,下四分位数,但是无法反映数据分布和数据密度,此次单株性状统计按道理来说应该符合正态分布,当然可以用夏皮罗检验做正态分布检验,但是也可以从小提琴图看数据是否符合正态分布。
结果显而易见,有部分的数据缺失不太符合正态分布,究其原因可能是因为边际效用等作用造成的,同一行中在两边的单株通常长势较好,但是中间的有些差,所以铃数这个性状可能不会优先考虑。认为比较有说服力的还是铃重。
- 极坐标下的小提琴图画法
纯属娱乐..不过效果不错。
#polor violin
polor2 = ggplot(lwpb, aes(x = cultivar,
y = lw,
fill = cultivar))+
geom_violin(alpha = 0.95, width = 1) +
theme_bw()+
coord_polar()+
scale_fill_discrete(c=100, l=100)+
labs(title = "lint weight per ball")+
xlab("cultivar")+
ylab("weight(g)")
polor2finnal = polor2 +
mythemepolor
ggsave(file = "lint weight per ball_polor.pdf", plot = polor2finnal)
效果:
这个是皮棉重,相比铃重的数据,这个数据的分布基本上都会趋向于正态分布了。
就这么多吧。Reference太多了,后面会整理补上。
最后用李老师曾经说过的一句名言结束吧:生命在于寻找,想搞科研,先学会百度google吧,在准备问别人之前先问问自己,自己努力寻找过了吗?