library(readr)
GMCCoffstation <- read_csv("GMCCoffstation.csv")
## Parsed with column specification:
## cols(
## season = col_double(),
## ward = col_double(),
## farmer = col_character(),
## crop = col_character(),
## treat = col_character(),
## rep = col_double(),
## biomass = col_double(),
## `biomass/1000` = col_double(),
## grain = col_double(),
## `grain/1000` = col_double()
## )
View(GMCCoffstation)
attach(GMCCoffstation)
library(lme4)
## Loading required package: Matrix
library(lmerTest)
##
## Attaching package: 'lmerTest'
## The following object is masked from 'package:lme4':
##
## lmer
## The following object is masked from 'package:stats':
##
## step
library(ggplot2)
library(psych)
##
## Attaching package: 'psych'
## The following objects are masked from 'package:ggplot2':
##
## %+%, alpha
library(agricolae)
require(MASS)
## Loading required package: MASS
library(Hmisc)
## Loading required package: lattice
## Loading required package: survival
## Loading required package: Formula
##
## Attaching package: 'Hmisc'
## The following object is masked from 'package:psych':
##
## describe
## The following objects are masked from 'package:base':
##
## format.pval, units
library(lubridate) # for working with dates
##
## Attaching package: 'lubridate'
## The following object is masked from 'package:base':
##
## date
library(scales) # to access breaks/formatting functions
##
## Attaching package: 'scales'
## The following objects are masked from 'package:psych':
##
## alpha, rescale
## The following object is masked from 'package:readr':
##
## col_factor
library(gridExtra) # for arranging plots
library(ggthemes)
library(ggrepel)
library(tidyr)
##
## Attaching package: 'tidyr'
## The following object is masked from 'package:Matrix':
##
## expand
library(dplyr)
##
## Attaching package: 'dplyr'
## The following object is masked from 'package:gridExtra':
##
## combine
## The following objects are masked from 'package:lubridate':
##
## intersect, setdiff, union
## The following objects are masked from 'package:Hmisc':
##
## src, summarize
## The following object is masked from 'package:MASS':
##
## select
## The following objects are masked from 'package:stats':
##
## filter, lag
## The following objects are masked from 'package:base':
##
## intersect, setdiff, setequal, union
library(dplyr) # for data manipulation
library(tidyr) # for data tidying
library(ggplot2) # for generating the visualizations
library(treemap)
library(reshape2)
##
## Attaching package: 'reshape2'
## The following object is masked from 'package:tidyr':
##
## smiths
library(ggalluvial)
season=as.numeric(season)
ward=as.factor(ward)
farmer=as.factor(farmer)
treat=as.factor(treat)
rep=as.factor(rep)
######################ANOVA###############
mod<-aov(grain/1000~rep+farmer+ward+season:rep+treat+season+farmer:ward+farmer:treat:ward+season:ward,data=GMCCoffstation)
anova(mod)
## Analysis of Variance Table
##
## Response: grain/1000
## Df Sum Sq Mean Sq F value Pr(>F)
## rep 1 4.287 4.2872 4.5834 0.033816 *
## farmer 20 162.584 8.1292 8.6908 < 2.2e-16 ***
## treat 5 26.772 5.3545 5.7244 6.976e-05 ***
## season 1 27.098 27.0982 28.9704 2.612e-07 ***
## rep:season 1 5.067 5.0671 5.4172 0.021208 *
## ward:season 1 9.128 9.1284 9.7591 0.002123 **
## farmer:ward:treat 100 81.343 0.8134 0.8696 0.774441
## Residuals 158 147.789 0.9354
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
summary(mod)
## Df Sum Sq Mean Sq F value Pr(>F)
## rep 1 4.29 4.287 4.583 0.03382 *
## farmer 20 162.58 8.129 8.691 < 2e-16 ***
## treat 5 26.77 5.354 5.724 6.98e-05 ***
## season 1 27.10 27.098 28.970 2.61e-07 ***
## rep:season 1 5.07 5.067 5.417 0.02121 *
## ward:season 1 9.13 9.128 9.759 0.00212 **
## farmer:ward:treat 100 81.34 0.813 0.870 0.77444
## Residuals 158 147.79 0.935
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
y<-HSD.test(mod,"treat",group = TRUE)
y
## $statistics
## MSerror Df Mean CV MSD
## 0.9353752 158 2.549194 37.93937 0.5695539
##
## $parameters
## test name.t ntr StudentizedRange alpha
## Tukey treat 6 4.080022 0.05
##
## $means
## grain/1000 std r Min Max Q25 Q50
## COWP_INT 2.505654 1.221901 48 0.3836 6.0609 1.767800 2.54715
## JCKBN_INT 2.418983 1.008070 48 0.2537 4.9269 1.800750 2.36910
## LABLAB_INT 2.418181 1.292090 48 0.4825 6.5697 1.510825 2.39455
## PGPEA_COWP_INT 2.172323 1.112348 48 0.5073 5.1913 1.366975 1.98030
## PGPEA_INT 2.618533 1.249121 48 0.6387 6.0833 1.640575 2.35075
## SOLE_MZ 3.161487 1.525724 48 0.5355 6.6848 2.090525 2.71220
## Q75
## COWP_INT 3.162225
## JCKBN_INT 2.917725
## LABLAB_INT 3.011300
## PGPEA_COWP_INT 2.917450
## PGPEA_INT 3.269675
## SOLE_MZ 4.337925
##
## $comparison
## NULL
##
## $groups
## grain/1000 groups
## SOLE_MZ 3.161487 a
## PGPEA_INT 2.618533 ab
## COWP_INT 2.505654 b
## JCKBN_INT 2.418983 b
## LABLAB_INT 2.418181 b
## PGPEA_COWP_INT 2.172323 b
##
## attr(,"class")
## [1] "group"
x<-HSD.test(mod,"season",group = TRUE)
x
## $statistics
## MSerror Df Mean CV MSD
## 0.9353752 158 2.549194 37.93937 0.330285
##
## $parameters
## test name.t ntr StudentizedRange alpha
## Tukey season 3 3.346043 0.05
##
## $means
## grain/1000 std r Min Max Q25 Q50 Q75
## 2017 2.843977 1.366597 96 0.5073 6.6848 1.778000 2.69950 3.754375
## 2018 2.738111 1.132469 96 0.9399 6.1279 1.937725 2.52270 3.135900
## 2019 2.065493 1.172418 96 0.2537 5.1913 1.082975 1.98715 2.822875
##
## $comparison
## NULL
##
## $groups
## grain/1000 groups
## 2017 2.843977 a
## 2018 2.738111 a
## 2019 2.065493 b
##
## attr(,"class")
## [1] "group"
v<-HSD.test(mod,"farmer",group = TRUE)
v
## $statistics
## MSerror Df Mean CV
## 0.9353752 158 2.549194 37.93937
##
## $parameters
## test name.t ntr StudentizedRange alpha
## Tukey farmer 21 5.135691 0.05
##
## $means
## grain/1000 std r Min Max Q25 Q50
## E. Maresa 2.478561 1.3101922 18 0.7193 4.8348 1.350725 2.46155
## E. Matanda 2.656083 0.5783633 12 1.9428 3.8277 2.065325 2.70585
## E. Pangwa 1.964733 0.6196082 12 0.7998 2.8411 1.651825 2.03390
## E.Matanda 1.806967 0.3775486 6 1.4516 2.3965 1.540900 1.66680
## E.Pangwa 1.079467 0.6002911 6 0.5384 2.2410 0.782175 0.93520
## F. Masango 2.564500 0.4976239 12 1.8057 3.4873 2.283400 2.52735
## F.Masango 3.014800 1.1298755 6 1.7877 5.0675 2.398425 2.86675
## J.Ackchomwe 2.232250 1.0063666 6 0.7087 3.1801 1.579375 2.66785
## J.Mukandi 1.817194 0.6088046 18 0.6928 3.1896 1.472750 1.82325
## Jakuchomwe 3.018175 0.6930173 12 2.0529 4.0870 2.447650 3.13120
## L. Kambarami 2.994161 1.9116775 18 0.4278 6.5697 1.534750 2.83300
## L. Masango 2.057233 0.8364770 12 0.9399 3.9311 1.637875 1.90390
## L.Chikaka 3.431733 1.5742558 18 0.7895 6.0833 2.346600 3.30985
## L.Masango 2.594133 1.0972136 6 0.5355 3.8251 2.626550 2.81755
## M. Masango 1.864950 0.9205786 18 0.5073 4.3559 1.244450 1.59845
## M.Chiweshe 3.637050 0.7211141 18 2.0781 4.9142 3.205700 3.69995
## S. Chikauro 4.383767 1.1424939 18 1.7859 6.6848 3.948925 4.36210
## S. Phiri 2.003256 1.0100625 18 0.2537 3.8888 1.190375 2.30385
## S.Manhembe 2.161783 0.8010658 18 0.7849 3.7792 1.684750 1.83680
## S.Vambe 2.492867 1.5323083 18 0.3836 5.3343 1.031150 2.45080
## T. Matanda 1.772089 0.3872743 18 1.0686 2.5193 1.547150 1.80515
## Q75
## E. Maresa 3.237825
## E. Matanda 3.026175
## E. Pangwa 2.446125
## E.Matanda 2.037575
## E.Pangwa 1.055675
## F. Masango 2.775075
## F.Masango 3.160025
## J.Ackchomwe 2.881300
## J.Mukandi 2.190225
## Jakuchomwe 3.481000
## L. Kambarami 4.379225
## L. Masango 2.413125
## L.Chikaka 4.440850
## L.Masango 2.959050
## M. Masango 2.297050
## M.Chiweshe 4.161150
## S. Chikauro 4.985475
## S. Phiri 2.602075
## S.Manhembe 2.398525
## S.Vambe 3.030775
## T. Matanda 1.950200
##
## $comparison
## NULL
##
## $groups
## grain/1000 groups
## S. Chikauro 4.383767 a
## M.Chiweshe 3.637050 ab
## L.Chikaka 3.431733 ab
## Jakuchomwe 3.018175 bc
## F.Masango 3.014800 bcd
## L. Kambarami 2.994161 bcd
## E. Matanda 2.656083 bcd
## L.Masango 2.594133 bcd
## F. Masango 2.564500 bcd
## S.Vambe 2.492867 bcd
## E. Maresa 2.478561 bcd
## J.Ackchomwe 2.232250 bcd
## S.Manhembe 2.161783 cd
## L. Masango 2.057233 cd
## S. Phiri 2.003256 cd
## E. Pangwa 1.964733 cd
## M. Masango 1.864950 cd
## J.Mukandi 1.817194 cd
## E.Matanda 1.806967 cd
## T. Matanda 1.772089 cd
## E.Pangwa 1.079467 d
##
## attr(,"class")
## [1] "group"
#################density plots#################
ggplot(data = GMCCoffstation, mapping = aes(x = grain, fill = ward)) + geom_density(alpha = 0.5)

#################density plots#################
ggplot(data = GMCCoffstation, mapping = aes(x = grain, fill = farmer)) + geom_density(alpha = 0.5)

###############boxplots##############
theme_set(theme_gray(base_size =12))
ggplot(GMCCoffstation, aes(x = treat, y =grain/1000, color = treat)) +
geom_boxplot(size=1.2,varwidth = TRUE) +
geom_point(data = GMCCoffstation, aes(y = mean(grain/1000))) +
geom_line(data = GMCCoffstation, aes(y =mean(grain/1000)))+ylab("Grain yield [t/ha] ") + xlab("treatment")+
theme( axis.text.x = element_text( angle = 90, hjust = 1))+ scale_y_continuous(breaks = c(0,500,1000,1500,2000,2500,3000,3500,4000,4500,5000,5500,6000,6500,7000))

###############facets#################
theme_set(theme_gray(base_size =12))
ggplot(GMCCoffstation, aes(x = treat, y =grain/1000, color = treat)) +
geom_boxplot(size=1.2,varwidth = TRUE) +
geom_point(data = GMCCoffstation, aes(y = mean(grain/1000))) +
geom_line(data = GMCCoffstation, aes(y =mean(grain/1000)))+ylab("Grain yield [t/ha] ") + xlab("treatment") + facet_wrap(.~season)+
theme( axis.text.x = element_text( angle = 90, hjust = 1))+ theme(legend.position = "none")+ scale_y_continuous(breaks = c(0,1,1.5,2,2.5,3,3.5,4,4.5,5,5.5,6,6.5,7))

################################################
theme_set(theme_gray(base_size =12))
ggplot(GMCCoffstation, aes(x = treat, y =grain/1000, color = treat)) +
geom_boxplot(size=1.2,varwidth = TRUE) +
geom_point(data = GMCCoffstation, aes(y = mean(grain/1000))) +
geom_line(data = GMCCoffstation, aes(y =mean(grain/1000)))+ylab("Grain yield [t/ha] ") + xlab("treatment") + facet_wrap(.~farmer)+
theme( axis.text.x = element_text( angle = 90, hjust = 1))+ theme(legend.position = "none")+ scale_y_continuous(breaks = c(0,2.5,4.5,6.5,7.5))
## geom_path: Each group consists of only one observation. Do you need to
## adjust the group aesthetic?
## geom_path: Each group consists of only one observation. Do you need to
## adjust the group aesthetic?
## geom_path: Each group consists of only one observation. Do you need to
## adjust the group aesthetic?
## geom_path: Each group consists of only one observation. Do you need to
## adjust the group aesthetic?
## geom_path: Each group consists of only one observation. Do you need to
## adjust the group aesthetic?

######################################################
theme_set(theme_gray(base_size =12))
ggplot(GMCCoffstation, aes(x = treat, y =grain/1000, color = treat)) +
geom_boxplot(size=1.2,varwidth = TRUE) +
geom_point(data = GMCCoffstation, aes(y = mean(grain/1000))) +
geom_line(data = GMCCoffstation, aes(y =mean(grain/1000)))+ylab("Grain yield [t/ha] ") + xlab("treatment") + facet_wrap(.~season+ward)+
theme( axis.text.x = element_text( angle = 90, hjust = 1))+ theme(legend.position = "none")+ scale_y_continuous(breaks = c(0,1.5,2.5,3.5,4.5,5.5,6.5,7.5))

############biomass/1000###########
theme_set(theme_gray(base_size =15))
m <- ggplot(data=GMCCoffstation,aes(x=treat, y=biomass))
m + geom_violin(size=1.3,shape=8) + geom_boxplot(width=.2, outlier.size=0,fill=c("red","yellow","grey", "blue", "gold", "green"))+ylab("biomass [t/ha]") + xlab("treatment")+ theme(axis.title = element_text(size=20, face="bold"), axis.text.x = element_text(size=12, face="bold", angle = 90, hjust = 1), axis.text.y = element_text(size=16, face="bold"))+ scale_y_continuous(breaks = c(0,500,1000,1500,2000,2500,3000,3500,4000,4500,5000,5500,6000,6500,7000))
## Warning: Ignoring unknown parameters: shape

theme_set(theme_gray(base_size =15))
m <- ggplot(data=GMCCoffstation,aes(x=treat, y=grain))
m + geom_violin(size=1.3,shape=8) + geom_boxplot(width=.2, outlier.size=0,fill=c("red","yellow","grey", "blue", "gold", "green"))+ylab("grain [t/ha]") +stat_summary(fun.y=median, geom="point", size=2, color="red")+ xlab("treatment")+ geom_jitter(shape=16, position=position_jitter(0.2))+ theme(axis.title = element_text(size=20, face="bold"), axis.text.x = element_text(size=12, face="bold", angle = 90, hjust = 1), axis.text.y = element_text(size=16, face="bold"))+ scale_y_continuous(breaks = c(0,500,1000,1500,2000,2500,3000,3500,4000,4500,5000,5500,6000,6500,7000))
## Warning: Ignoring unknown parameters: shape

#####################facet###################
theme_set(theme_gray(base_size =12))
ggplot(GMCCoffstation, aes(x = treat, y =biomass/1000, color = treat)) +
geom_boxplot(size=1.2,varwidth = TRUE) +
geom_point(data = GMCCoffstation, aes(y = mean(biomass/1000))) +
geom_line(data = GMCCoffstation, aes(y =mean(biomass/1000)))+ylab("Biomass [t/ha] ") + xlab("treatment") + facet_wrap(.~season)+
theme( axis.text.x = element_text( angle = 90, hjust = 1))+ theme(legend.position = "none")+ scale_y_continuous(breaks = c(0,1,1.5,2,2.5,3,3.5,4,4.5,5,5.5,6,6.5,7))

ggplot(GMCCoffstation, aes(treat, grain/1000)) + geom_boxplot(fill = "springgreen3") + theme_bw() +
xlab("") + ylab("Yield (t/ha)") +
theme(axis.title = element_text(size=24, face="bold"),
axis.text.x = element_text(size=12, face="bold", angle = 90, hjust = 1),
axis.text.y = element_text(size=20, face="bold"))+ scale_y_continuous(breaks = c(0,500,1000,1500,2000,2500,3000,3500,4000,4500,5000,5500,6000,6500,7000))

ggplot(GMCCoffstation, aes(x = treat, y = grain, fill = treat)) +
geom_violin() + scale_y_log10() + geom_boxplot(width = 0.2) + theme(axis.title = element_text(size=20, face="bold"), axis.text.x = element_text(size=12, face="bold", angle = 90, hjust = 1), axis.text.y = element_text(size=16, face="bold")) + theme(legend.position = "none")+ scale_y_continuous(breaks = c(0,500,1000,1500,2000,2500,3000,3500,4000,4500,5000,5500,6000,6500,7000))
## Scale for 'y' is already present. Adding another scale for 'y', which
## will replace the existing scale.

#############basic ggplot###################
theme_set(theme_gray(base_size =15))
p<-ggplot(GMCCoffstation, aes(x = treat, y = grain, fill = treat)) + geom_boxplot(aes(fill = treat) + scale_y_log10())+ theme( axis.text.x = element_text( angle = 90, hjust = 1))+ theme(legend.position = "none")+ scale_y_continuous(breaks = c(0,500,1000,1500,2000,2500,3000,3500,4000,4500,5000,5500,6000,6500,7000))
p

#################using facet grid##################
###########split in vertical direction##################
theme_set(theme_gray(base_size =15))
p + facet_grid(rows = vars(season)) + theme(legend.position = "none")+ scale_y_continuous(breaks = c(500,1500,2500,3500,4500,5500,6500))
## Scale for 'y' is already present. Adding another scale for 'y', which
## will replace the existing scale.

theme_set(theme_gray(base_size =15))
p + facet_grid(rows = vars(ward)) + theme(legend.position = "none")+ scale_y_continuous(breaks = c(500,1500,2500,3500,4500,5500,6500))
## Scale for 'y' is already present. Adding another scale for 'y', which
## will replace the existing scale.

###############split in horizontal direction#############
p + facet_grid(cols = vars(season), scales= "free", space= "free")+ theme(legend.position = "none")+ scale_y_continuous(breaks = c(0,500,1000,1500,2000,2500,3000,3500,4000,4500,5000,5500,6000,6500,7000))
## Scale for 'y' is already present. Adding another scale for 'y', which
## will replace the existing scale.

p + facet_grid(cols = vars(ward), scales= "free", space= "free")+ theme(legend.position = "none")+ scale_y_continuous(breaks = c(0,500,1000,1500,2000,2500,3000,3500,4000,4500,5000,5500,6000,6500,7000))
## Scale for 'y' is already present. Adding another scale for 'y', which
## will replace the existing scale.

p + facet_grid(cols = vars(season), scales= "free", space= "free")+ theme(legend.position = "none")+ scale_y_continuous(breaks = c(0,500,1000,1500,2000,2500,3000,3500,4000,4500,5000,5500,6000,6500,7000))
## Scale for 'y' is already present. Adding another scale for 'y', which
## will replace the existing scale.

#############basic biomass ggplot###################
theme_set(theme_gray(base_size =15))
p<-ggplot(GMCCoffstation, aes(x = treat, y = biomass, fill = treat)) + geom_boxplot(aes(fill = treat) + scale_y_log10())+ theme( axis.text.x = element_text( angle = 90, hjust = 1))+ theme(legend.position = "none")+ scale_y_continuous(breaks = c(0,500,1000,1500,2000,2500,3000,3500,4000,4500,5000,5500,6000,6500,7000))
p

#################using facet grid##################
###########split in vertical direction##################
theme_set(theme_gray(base_size =15))
p + facet_grid(rows = vars(season)) + theme(legend.position = "none")+ scale_y_continuous(breaks = c(500,1500,2500,3500,4500,5500,6500))
## Scale for 'y' is already present. Adding another scale for 'y', which
## will replace the existing scale.

theme_set(theme_gray(base_size =15))
p + facet_grid(rows = vars(ward)) + theme(legend.position = "none")+ scale_y_continuous(breaks = c(500,1500,2500,3500,4500,5500,6500))
## Scale for 'y' is already present. Adding another scale for 'y', which
## will replace the existing scale.

###############split in horizontal direction#############
p + facet_grid(cols = vars(season), scales= "free", space= "free")+ theme(legend.position = "none")+ scale_y_continuous(breaks = c(0,500,1000,1500,2000,2500,3000,3500,4000,4500,5000,5500,6000,6500,7000))
## Scale for 'y' is already present. Adding another scale for 'y', which
## will replace the existing scale.

p + facet_grid(cols = vars(ward), scales= "free", space= "free")+ theme(legend.position = "none")+ scale_y_continuous(breaks = c(0,500,1000,1500,2000,2500,3000,3500,4000,4500,5000,5500,6000,6500,7000))
## Scale for 'y' is already present. Adding another scale for 'y', which
## will replace the existing scale.

p + facet_grid(vars(season), vars(ward))+ scale_y_continuous(breaks = c(500,1500,2500,3500,4500,5500,6500))
## Scale for 'y' is already present. Adding another scale for 'y', which
## will replace the existing scale.

p + facet_grid(cols = vars(season), scales= "free", space= "free")+ theme(legend.position = "none")+ scale_y_continuous(breaks = c(0,500,1000,1500,2000,2500,3000,3500,4000,4500,5000,5500,6000,6500,7000))
## Scale for 'y' is already present. Adding another scale for 'y', which
## will replace the existing scale.

ggplot() + geom_col(data = GMCCoffstation, aes(x = treat, y = grain))+ theme_classic(base_size = 20)+ theme(axis.title = element_text(size=20, face="bold"), axis.text.x = element_text(size=12, face="bold", angle = 75, hjust = 1), axis.text.y = element_text(size=16, face="bold"))

ggplot(GMCCoffstation, aes(x=treat, y=grain, fill=treat))+
geom_bar(stat='identity', fill="forest green")+
ylab("grain")+ theme_classic(base_size = 20)+ theme(axis.title = element_text(size=20, face="bold"), axis.text.x = element_text(size=12, face="bold", angle = 75, hjust = 1), axis.text.y = element_text(size=16, face="bold"))+
facet_wrap(~season)

#########################################3
ggplot(GMCCoffstation, aes(x=treat, y=grain, fill=treat))+
geom_bar(stat='identity', fill="forest green")+
ylab("grain")+ theme_classic(base_size = 20)+ theme(axis.title = element_text(size=20, face="bold"), axis.text.x = element_text(size=12, face="bold", angle = 75, hjust = 1), axis.text.y = element_text(size=16, face="bold"))+
facet_wrap(~season, ncol=1)

###############################3
ggplot(GMCCoffstation, aes(x=treat, y=grain, fill=treat))+
geom_bar(stat='identity', fill="forest green")+
ylab("grain")+ theme_classic(base_size = 20)+ theme(axis.title = element_text(size=20, face="bold"), axis.text.x = element_text(size=12, face="bold", angle = 75, hjust = 1), axis.text.y = element_text(size=16, face="bold"))+
facet_wrap(~season, ncol=1, strip.position = "left")

###################GEOM PLOTS #########################
ggplot(data = GMCCoffstation, aes(x=season,y=grain/1000,color=treat))+
geom_point()+ geom_smooth()
## `geom_smooth()` using method = 'loess' and formula 'y ~ x'
## Warning in simpleLoess(y, x, w, span, degree = degree, parametric =
## parametric, : pseudoinverse used at 2017
## Warning in simpleLoess(y, x, w, span, degree = degree, parametric =
## parametric, : neighborhood radius 2.01
## Warning in simpleLoess(y, x, w, span, degree = degree, parametric =
## parametric, : reciprocal condition number 6.7241e-017
## Warning in simpleLoess(y, x, w, span, degree = degree, parametric =
## parametric, : There are other near singularities as well. 4.0401
## Warning in predLoess(object$y, object$x, newx = if
## (is.null(newdata)) object$x else if (is.data.frame(newdata))
## as.matrix(model.frame(delete.response(terms(object)), : pseudoinverse used
## at 2017
## Warning in predLoess(object$y, object$x, newx = if
## (is.null(newdata)) object$x else if (is.data.frame(newdata))
## as.matrix(model.frame(delete.response(terms(object)), : neighborhood radius
## 2.01
## Warning in predLoess(object$y, object$x, newx = if
## (is.null(newdata)) object$x else if (is.data.frame(newdata))
## as.matrix(model.frame(delete.response(terms(object)), : reciprocal
## condition number 6.7241e-017
## Warning in predLoess(object$y, object$x, newx = if
## (is.null(newdata)) object$x else if (is.data.frame(newdata))
## as.matrix(model.frame(delete.response(terms(object)), : There are other
## near singularities as well. 4.0401
## Warning in simpleLoess(y, x, w, span, degree = degree, parametric =
## parametric, : pseudoinverse used at 2017
## Warning in simpleLoess(y, x, w, span, degree = degree, parametric =
## parametric, : neighborhood radius 2.01
## Warning in simpleLoess(y, x, w, span, degree = degree, parametric =
## parametric, : reciprocal condition number 6.7241e-017
## Warning in simpleLoess(y, x, w, span, degree = degree, parametric =
## parametric, : There are other near singularities as well. 4.0401
## Warning in predLoess(object$y, object$x, newx = if
## (is.null(newdata)) object$x else if (is.data.frame(newdata))
## as.matrix(model.frame(delete.response(terms(object)), : pseudoinverse used
## at 2017
## Warning in predLoess(object$y, object$x, newx = if
## (is.null(newdata)) object$x else if (is.data.frame(newdata))
## as.matrix(model.frame(delete.response(terms(object)), : neighborhood radius
## 2.01
## Warning in predLoess(object$y, object$x, newx = if
## (is.null(newdata)) object$x else if (is.data.frame(newdata))
## as.matrix(model.frame(delete.response(terms(object)), : reciprocal
## condition number 6.7241e-017
## Warning in predLoess(object$y, object$x, newx = if
## (is.null(newdata)) object$x else if (is.data.frame(newdata))
## as.matrix(model.frame(delete.response(terms(object)), : There are other
## near singularities as well. 4.0401
## Warning in simpleLoess(y, x, w, span, degree = degree, parametric =
## parametric, : pseudoinverse used at 2017
## Warning in simpleLoess(y, x, w, span, degree = degree, parametric =
## parametric, : neighborhood radius 2.01
## Warning in simpleLoess(y, x, w, span, degree = degree, parametric =
## parametric, : reciprocal condition number 6.7241e-017
## Warning in simpleLoess(y, x, w, span, degree = degree, parametric =
## parametric, : There are other near singularities as well. 4.0401
## Warning in predLoess(object$y, object$x, newx = if
## (is.null(newdata)) object$x else if (is.data.frame(newdata))
## as.matrix(model.frame(delete.response(terms(object)), : pseudoinverse used
## at 2017
## Warning in predLoess(object$y, object$x, newx = if
## (is.null(newdata)) object$x else if (is.data.frame(newdata))
## as.matrix(model.frame(delete.response(terms(object)), : neighborhood radius
## 2.01
## Warning in predLoess(object$y, object$x, newx = if
## (is.null(newdata)) object$x else if (is.data.frame(newdata))
## as.matrix(model.frame(delete.response(terms(object)), : reciprocal
## condition number 6.7241e-017
## Warning in predLoess(object$y, object$x, newx = if
## (is.null(newdata)) object$x else if (is.data.frame(newdata))
## as.matrix(model.frame(delete.response(terms(object)), : There are other
## near singularities as well. 4.0401
## Warning in simpleLoess(y, x, w, span, degree = degree, parametric =
## parametric, : pseudoinverse used at 2017
## Warning in simpleLoess(y, x, w, span, degree = degree, parametric =
## parametric, : neighborhood radius 2.01
## Warning in simpleLoess(y, x, w, span, degree = degree, parametric =
## parametric, : reciprocal condition number 6.7241e-017
## Warning in simpleLoess(y, x, w, span, degree = degree, parametric =
## parametric, : There are other near singularities as well. 4.0401
## Warning in predLoess(object$y, object$x, newx = if
## (is.null(newdata)) object$x else if (is.data.frame(newdata))
## as.matrix(model.frame(delete.response(terms(object)), : pseudoinverse used
## at 2017
## Warning in predLoess(object$y, object$x, newx = if
## (is.null(newdata)) object$x else if (is.data.frame(newdata))
## as.matrix(model.frame(delete.response(terms(object)), : neighborhood radius
## 2.01
## Warning in predLoess(object$y, object$x, newx = if
## (is.null(newdata)) object$x else if (is.data.frame(newdata))
## as.matrix(model.frame(delete.response(terms(object)), : reciprocal
## condition number 6.7241e-017
## Warning in predLoess(object$y, object$x, newx = if
## (is.null(newdata)) object$x else if (is.data.frame(newdata))
## as.matrix(model.frame(delete.response(terms(object)), : There are other
## near singularities as well. 4.0401
## Warning in simpleLoess(y, x, w, span, degree = degree, parametric =
## parametric, : pseudoinverse used at 2017
## Warning in simpleLoess(y, x, w, span, degree = degree, parametric =
## parametric, : neighborhood radius 2.01
## Warning in simpleLoess(y, x, w, span, degree = degree, parametric =
## parametric, : reciprocal condition number 6.7241e-017
## Warning in simpleLoess(y, x, w, span, degree = degree, parametric =
## parametric, : There are other near singularities as well. 4.0401
## Warning in predLoess(object$y, object$x, newx = if
## (is.null(newdata)) object$x else if (is.data.frame(newdata))
## as.matrix(model.frame(delete.response(terms(object)), : pseudoinverse used
## at 2017
## Warning in predLoess(object$y, object$x, newx = if
## (is.null(newdata)) object$x else if (is.data.frame(newdata))
## as.matrix(model.frame(delete.response(terms(object)), : neighborhood radius
## 2.01
## Warning in predLoess(object$y, object$x, newx = if
## (is.null(newdata)) object$x else if (is.data.frame(newdata))
## as.matrix(model.frame(delete.response(terms(object)), : reciprocal
## condition number 6.7241e-017
## Warning in predLoess(object$y, object$x, newx = if
## (is.null(newdata)) object$x else if (is.data.frame(newdata))
## as.matrix(model.frame(delete.response(terms(object)), : There are other
## near singularities as well. 4.0401
## Warning in simpleLoess(y, x, w, span, degree = degree, parametric =
## parametric, : pseudoinverse used at 2017
## Warning in simpleLoess(y, x, w, span, degree = degree, parametric =
## parametric, : neighborhood radius 2.01
## Warning in simpleLoess(y, x, w, span, degree = degree, parametric =
## parametric, : reciprocal condition number 6.7241e-017
## Warning in simpleLoess(y, x, w, span, degree = degree, parametric =
## parametric, : There are other near singularities as well. 4.0401
## Warning in predLoess(object$y, object$x, newx = if
## (is.null(newdata)) object$x else if (is.data.frame(newdata))
## as.matrix(model.frame(delete.response(terms(object)), : pseudoinverse used
## at 2017
## Warning in predLoess(object$y, object$x, newx = if
## (is.null(newdata)) object$x else if (is.data.frame(newdata))
## as.matrix(model.frame(delete.response(terms(object)), : neighborhood radius
## 2.01
## Warning in predLoess(object$y, object$x, newx = if
## (is.null(newdata)) object$x else if (is.data.frame(newdata))
## as.matrix(model.frame(delete.response(terms(object)), : reciprocal
## condition number 6.7241e-017
## Warning in predLoess(object$y, object$x, newx = if
## (is.null(newdata)) object$x else if (is.data.frame(newdata))
## as.matrix(model.frame(delete.response(terms(object)), : There are other
## near singularities as well. 4.0401

#Simple alluvial diagram
ggplot(GMCCoffstation, aes(kg = grain/1000, axis1 = treat, axis2 = season)) +
geom_stratum() + geom_alluvium(aes(fill = treat), color = "black")

#####################
ggplot(GMCCoffstation, aes(kg = grain/1000, axis1 = treat, axis2 = season)) +
geom_alluvium(aes(fill = treat), color = "black") + geom_stratum()

#Add text
ggplot(GMCCoffstation, aes(kg = grain/1000, axis1 = treat, axis2 = season)) +
geom_alluvium(aes(fill = treat), color = "black") + geom_stratum() +
geom_text(stat = "stratum", label.strata = TRUE) +
theme(legend.position = "true",
axis.text = element_text(size=14))

###################GEOM PLOTS #########################
ggplot(data = GMCCoffstation, aes(x=season,y=grain/1000,color=treat))+
geom_point()+ geom_smooth()
## `geom_smooth()` using method = 'loess' and formula 'y ~ x'
## Warning in simpleLoess(y, x, w, span, degree = degree, parametric =
## parametric, : pseudoinverse used at 2017
## Warning in simpleLoess(y, x, w, span, degree = degree, parametric =
## parametric, : neighborhood radius 2.01
## Warning in simpleLoess(y, x, w, span, degree = degree, parametric =
## parametric, : reciprocal condition number 6.7241e-017
## Warning in simpleLoess(y, x, w, span, degree = degree, parametric =
## parametric, : There are other near singularities as well. 4.0401
## Warning in predLoess(object$y, object$x, newx = if
## (is.null(newdata)) object$x else if (is.data.frame(newdata))
## as.matrix(model.frame(delete.response(terms(object)), : pseudoinverse used
## at 2017
## Warning in predLoess(object$y, object$x, newx = if
## (is.null(newdata)) object$x else if (is.data.frame(newdata))
## as.matrix(model.frame(delete.response(terms(object)), : neighborhood radius
## 2.01
## Warning in predLoess(object$y, object$x, newx = if
## (is.null(newdata)) object$x else if (is.data.frame(newdata))
## as.matrix(model.frame(delete.response(terms(object)), : reciprocal
## condition number 6.7241e-017
## Warning in predLoess(object$y, object$x, newx = if
## (is.null(newdata)) object$x else if (is.data.frame(newdata))
## as.matrix(model.frame(delete.response(terms(object)), : There are other
## near singularities as well. 4.0401
## Warning in simpleLoess(y, x, w, span, degree = degree, parametric =
## parametric, : pseudoinverse used at 2017
## Warning in simpleLoess(y, x, w, span, degree = degree, parametric =
## parametric, : neighborhood radius 2.01
## Warning in simpleLoess(y, x, w, span, degree = degree, parametric =
## parametric, : reciprocal condition number 6.7241e-017
## Warning in simpleLoess(y, x, w, span, degree = degree, parametric =
## parametric, : There are other near singularities as well. 4.0401
## Warning in predLoess(object$y, object$x, newx = if
## (is.null(newdata)) object$x else if (is.data.frame(newdata))
## as.matrix(model.frame(delete.response(terms(object)), : pseudoinverse used
## at 2017
## Warning in predLoess(object$y, object$x, newx = if
## (is.null(newdata)) object$x else if (is.data.frame(newdata))
## as.matrix(model.frame(delete.response(terms(object)), : neighborhood radius
## 2.01
## Warning in predLoess(object$y, object$x, newx = if
## (is.null(newdata)) object$x else if (is.data.frame(newdata))
## as.matrix(model.frame(delete.response(terms(object)), : reciprocal
## condition number 6.7241e-017
## Warning in predLoess(object$y, object$x, newx = if
## (is.null(newdata)) object$x else if (is.data.frame(newdata))
## as.matrix(model.frame(delete.response(terms(object)), : There are other
## near singularities as well. 4.0401
## Warning in simpleLoess(y, x, w, span, degree = degree, parametric =
## parametric, : pseudoinverse used at 2017
## Warning in simpleLoess(y, x, w, span, degree = degree, parametric =
## parametric, : neighborhood radius 2.01
## Warning in simpleLoess(y, x, w, span, degree = degree, parametric =
## parametric, : reciprocal condition number 6.7241e-017
## Warning in simpleLoess(y, x, w, span, degree = degree, parametric =
## parametric, : There are other near singularities as well. 4.0401
## Warning in predLoess(object$y, object$x, newx = if
## (is.null(newdata)) object$x else if (is.data.frame(newdata))
## as.matrix(model.frame(delete.response(terms(object)), : pseudoinverse used
## at 2017
## Warning in predLoess(object$y, object$x, newx = if
## (is.null(newdata)) object$x else if (is.data.frame(newdata))
## as.matrix(model.frame(delete.response(terms(object)), : neighborhood radius
## 2.01
## Warning in predLoess(object$y, object$x, newx = if
## (is.null(newdata)) object$x else if (is.data.frame(newdata))
## as.matrix(model.frame(delete.response(terms(object)), : reciprocal
## condition number 6.7241e-017
## Warning in predLoess(object$y, object$x, newx = if
## (is.null(newdata)) object$x else if (is.data.frame(newdata))
## as.matrix(model.frame(delete.response(terms(object)), : There are other
## near singularities as well. 4.0401
## Warning in simpleLoess(y, x, w, span, degree = degree, parametric =
## parametric, : pseudoinverse used at 2017
## Warning in simpleLoess(y, x, w, span, degree = degree, parametric =
## parametric, : neighborhood radius 2.01
## Warning in simpleLoess(y, x, w, span, degree = degree, parametric =
## parametric, : reciprocal condition number 6.7241e-017
## Warning in simpleLoess(y, x, w, span, degree = degree, parametric =
## parametric, : There are other near singularities as well. 4.0401
## Warning in predLoess(object$y, object$x, newx = if
## (is.null(newdata)) object$x else if (is.data.frame(newdata))
## as.matrix(model.frame(delete.response(terms(object)), : pseudoinverse used
## at 2017
## Warning in predLoess(object$y, object$x, newx = if
## (is.null(newdata)) object$x else if (is.data.frame(newdata))
## as.matrix(model.frame(delete.response(terms(object)), : neighborhood radius
## 2.01
## Warning in predLoess(object$y, object$x, newx = if
## (is.null(newdata)) object$x else if (is.data.frame(newdata))
## as.matrix(model.frame(delete.response(terms(object)), : reciprocal
## condition number 6.7241e-017
## Warning in predLoess(object$y, object$x, newx = if
## (is.null(newdata)) object$x else if (is.data.frame(newdata))
## as.matrix(model.frame(delete.response(terms(object)), : There are other
## near singularities as well. 4.0401
## Warning in simpleLoess(y, x, w, span, degree = degree, parametric =
## parametric, : pseudoinverse used at 2017
## Warning in simpleLoess(y, x, w, span, degree = degree, parametric =
## parametric, : neighborhood radius 2.01
## Warning in simpleLoess(y, x, w, span, degree = degree, parametric =
## parametric, : reciprocal condition number 6.7241e-017
## Warning in simpleLoess(y, x, w, span, degree = degree, parametric =
## parametric, : There are other near singularities as well. 4.0401
## Warning in predLoess(object$y, object$x, newx = if
## (is.null(newdata)) object$x else if (is.data.frame(newdata))
## as.matrix(model.frame(delete.response(terms(object)), : pseudoinverse used
## at 2017
## Warning in predLoess(object$y, object$x, newx = if
## (is.null(newdata)) object$x else if (is.data.frame(newdata))
## as.matrix(model.frame(delete.response(terms(object)), : neighborhood radius
## 2.01
## Warning in predLoess(object$y, object$x, newx = if
## (is.null(newdata)) object$x else if (is.data.frame(newdata))
## as.matrix(model.frame(delete.response(terms(object)), : reciprocal
## condition number 6.7241e-017
## Warning in predLoess(object$y, object$x, newx = if
## (is.null(newdata)) object$x else if (is.data.frame(newdata))
## as.matrix(model.frame(delete.response(terms(object)), : There are other
## near singularities as well. 4.0401
## Warning in simpleLoess(y, x, w, span, degree = degree, parametric =
## parametric, : pseudoinverse used at 2017
## Warning in simpleLoess(y, x, w, span, degree = degree, parametric =
## parametric, : neighborhood radius 2.01
## Warning in simpleLoess(y, x, w, span, degree = degree, parametric =
## parametric, : reciprocal condition number 6.7241e-017
## Warning in simpleLoess(y, x, w, span, degree = degree, parametric =
## parametric, : There are other near singularities as well. 4.0401
## Warning in predLoess(object$y, object$x, newx = if
## (is.null(newdata)) object$x else if (is.data.frame(newdata))
## as.matrix(model.frame(delete.response(terms(object)), : pseudoinverse used
## at 2017
## Warning in predLoess(object$y, object$x, newx = if
## (is.null(newdata)) object$x else if (is.data.frame(newdata))
## as.matrix(model.frame(delete.response(terms(object)), : neighborhood radius
## 2.01
## Warning in predLoess(object$y, object$x, newx = if
## (is.null(newdata)) object$x else if (is.data.frame(newdata))
## as.matrix(model.frame(delete.response(terms(object)), : reciprocal
## condition number 6.7241e-017
## Warning in predLoess(object$y, object$x, newx = if
## (is.null(newdata)) object$x else if (is.data.frame(newdata))
## as.matrix(model.frame(delete.response(terms(object)), : There are other
## near singularities as well. 4.0401
