library(ggplot2)
library(extrafontdb)
library(MASS)
library(pscl)
## Classes and Methods for R developed in the
## Political Science Computational Laboratory
## Department of Political Science
## Stanford University
## Simon Jackman
## hurdle and zeroinfl functions by Achim Zeileis
library(psych)
##
## Attaching package: 'psych'
## The following objects are masked from 'package:ggplot2':
##
## %+%, alpha
library(gridExtra)
library(repr) ### adjusting the length and width of your plot
library(beanplot)
library("devtools")
## Loading required package: usethis
library("yarrr")
## Loading required package: jpeg
## Loading required package: BayesFactor
## Loading required package: coda
## Loading required package: Matrix
## ************
## Welcome to BayesFactor 0.9.12-4.2. If you have questions, please contact Richard Morey (richarddmorey@gmail.com).
##
## Type BFManual() to open the manual.
## ************
## Loading required package: circlize
## ========================================
## circlize version 0.4.8
## CRAN page: https://cran.r-project.org/package=circlize
## Github page: https://github.com/jokergoo/circlize
## Documentation: http://jokergoo.github.io/circlize_book/book/
##
## If you use it in published research, please cite:
## Gu, Z. circlize implements and enhances circular visualization
## in R. Bioinformatics 2014.
## ========================================
## yarrr v0.1.5. Citation info at citation('yarrr'). Package guide at yarrr.guide()
## Email me at Nathaniel.D.Phillips.is@gmail.com
##
## Attaching package: 'yarrr'
## The following object is masked from 'package:ggplot2':
##
## diamonds
library(agricolae)
library(easynls)
library(MVN)
## Registered S3 method overwritten by 'GGally':
## method from
## +.gg ggplot2
## sROC 0.1-2 loaded
library(lme4)
## Registered S3 methods overwritten by 'lme4':
## method from
## cooks.distance.influence.merMod car
## influence.merMod car
## dfbeta.influence.merMod car
## dfbetas.influence.merMod car
library(ggsignif)
library(ggpubr)
## Loading required package: magrittr
library(multcomp)
## Loading required package: mvtnorm
## Loading required package: survival
## Loading required package: TH.data
##
## Attaching package: 'TH.data'
## The following object is masked from 'package:MASS':
##
## geyser
library(multcompView)
library(Rmisc)
## Loading required package: lattice
## Loading required package: plyr
##
## Attaching package: 'plyr'
## The following object is masked from 'package:ggpubr':
##
## mutate
library(readr)
yld <- read_csv("~/Dr Nyagumbo/Chemical analysis/final yield with socio-analysis/yld.csv")
## Warning: Duplicated column names deduplicated: 'biomass' => 'biomass_1' [33]
## Parsed with column specification:
## cols(
## .default = col_character(),
## age = col_double(),
## hhsize = col_double(),
## experience_CA = col_double(),
## weeding = col_double(),
## row_spacing = col_double(),
## plant_spacing = col_double(),
## number_plants = col_double(),
## plant_population = col_double(),
## days_seeding = col_double(),
## biomass = col_double(),
## grain = col_double(),
## biomass_1 = col_double()
## )
## See spec(...) for full column specifications.
View(yld)
attach(yld)
names(yld)
## [1] "farmer_name" "tech" "season" "district"
## [5] "village" "site" "agroecology" "cropping_system"
## [9] "class" "age" "sex" "hhsize"
## [13] "experience_CA" "local_soil_name" "soil_colour" "texture"
## [17] "field_treatment" "weeding" "crop_name" "variety"
## [21] "row_spacing" "plant_spacing" "number_plants" "tillage"
## [25] "residue_applied" "residue_type" "plant_population" "planting_rains"
## [29] "seeding" "days_seeding" "biomass" "grain"
## [33] "biomass_1"
tech=as.factor(tech)
season=as.factor(season)
district=as.factor(district)
agroecology=as.factor(agroecology)
site=as.factor(site)
#age
sex=as.factor(sex)
## hhsze
local_soil_name=as.factor(local_soil_name)
texture=as.factor(texture)
soil_colour=as.factor(soil_colour)
variet=factor(variety,levels = c("DKC 8033","MH33","SC627","Kanyani","Kanthochi","SC503","Mapasa","Syngenta","ZM523","MH18","Mphangala","Peacock","DKC 9089","Local","Masika","MH44","Nkango","Pannar","QPM","SC403"),ordered = TRUE)
###############################################
cropping_systems=factor(class,levels = c("CA","CONV","Intercrop maize","Monocrop maize","Mbeya fert","Normal fert"))
residue_applied=as.factor(residue_applied)
residue_type=as.factor(residue_type)
class=as.factor(class)
##############################################
biod<-lm(grain~tech+cropping_system+agroecology+residue_applied+plant_population+weeding+number_plants,data = yld)
summary(biod)
##
## Call:
## lm(formula = grain ~ tech + cropping_system + agroecology + residue_applied +
## plant_population + weeding + number_plants, data = yld)
##
## Residuals:
## Min 1Q Median 3Q Max
## -3213.5 -1415.4 -344.4 873.5 7142.2
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 2.840e+03 1.278e+03 2.222 0.028302 *
## techNon-CSA -6.962e+02 4.614e+02 -1.509 0.134135
## cropping_systemIntercrop -1.332e+02 4.779e+02 -0.279 0.780949
## cropping_systemMbeya 1.265e+03 6.444e+02 1.963 0.052106 .
## agroecologyLow rainfall -8.937e+01 4.194e+02 -0.213 0.831636
## residue_appliedYes -5.809e+01 5.280e+02 -0.110 0.912592
## plant_population 6.963e-02 2.167e-02 3.213 0.001709 **
## weeding -7.879e+01 3.798e+02 -0.207 0.836030
## number_plants -1.171e+03 3.358e+02 -3.485 0.000701 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 2128 on 113 degrees of freedom
## (2 observations deleted due to missingness)
## Multiple R-squared: 0.2582, Adjusted R-squared: 0.2057
## F-statistic: 4.917 on 8 and 113 DF, p-value: 3.131e-05
anova(biod)
## Analysis of Variance Table
##
## Response: grain
## Df Sum Sq Mean Sq F value Pr(>F)
## tech 1 25687973 25687973 5.6748 0.0188808 *
## cropping_system 2 20505345 10252672 2.2649 0.1085321
## agroecology 1 5852766 5852766 1.2929 0.2579121
## residue_applied 1 6322492 6322492 1.3967 0.2397565
## plant_population 1 64571763 64571763 14.2646 0.0002550 ***
## weeding 1 146504 146504 0.0324 0.8575527
## number_plants 1 54990576 54990576 12.1480 0.0007005 ***
## Residuals 113 511517993 4526708
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
model_1<-aov(grain~tech+agroecology+site+texture+soil_colour+variety+cropping_system+tech*agroecology+tech*site+tech*texture+tech*soil_colour+tech*cropping_system+tech*residue_applied+tech*cropping_system*agroecology+tech*cropping_system*site+tech*cropping_system*soil_colour+cropping_system*soil_colour*residue_applied+cropping_system*soil_colour*variety+residue_applied+residue_type+plant_population+weeding+number_plants+age+sex+hhsize+experience_CA,data = yld)
#############
out<-LSD.test(model_1,"tech",p.adj="hommel",console=TRUE)
##
## Study: model_1 ~ "tech"
##
## LSD t Test for grain
## P value adjustment method: hommel
##
## Mean Square Error: 1623275
##
## tech, means and individual ( 95 %) CI
##
## grain std r LCL UCL Min Max
## CSA 3903.330 2662.174 51 3534.267 4272.392 539.844 11091.00
## Non-CSA 2869.872 2127.760 51 2500.809 3238.934 0.000 8502.05
##
## Alpha: 0.05 ; DF Error: 23
## Critical Value of t: 2.068658
##
## Minimum Significant Difference: 521.9328
##
## Treatments with the same letter are not significantly different.
##
## grain groups
## CSA 3903.330 a
## Non-CSA 2869.872 b
out
## $statistics
## MSerror Df Mean CV t.value MSD
## 1623275 23 3386.601 37.62115 2.068658 521.9328
##
## $parameters
## test p.ajusted name.t ntr alpha
## Fisher-LSD hommel tech 2 0.05
##
## $means
## grain std r LCL UCL Min Max Q25
## CSA 3903.330 2662.174 51 3534.267 4272.392 539.844 11091.00 1817.968
## Non-CSA 2869.872 2127.760 51 2500.809 3238.934 0.000 8502.05 1421.559
## Q50 Q75
## CSA 3103.922 5272.583
## Non-CSA 2028.926 4131.753
##
## $comparison
## NULL
##
## $groups
## grain groups
## CSA 3903.330 a
## Non-CSA 2869.872 b
##
## attr(,"class")
## [1] "group"
###################################################
out<-LSD.test(model_1,"cropping_system",p.adj="hommel",console=TRUE)
##
## Study: model_1 ~ "cropping_system"
##
## LSD t Test for grain
## P value adjustment method: hommel
##
## Mean Square Error: 1623275
##
## cropping_system, means and individual ( 95 %) CI
##
## grain std r LCL UCL Min Max
## CA 3236.625 2224.705 56 2884.424 3588.826 0.000 10082.23
## Intercrop 3179.784 2519.803 28 2681.696 3677.871 539.844 10186.92
## Mbeya 4174.907 2977.093 18 3553.682 4796.131 813.145 11091.00
##
## Alpha: 0.05 ; DF Error: 23
## Critical Value of t: 2.068658
##
## Groups according to probability of means differences and alpha level( 0.05 )
##
## Treatments with the same letter are not significantly different.
##
## grain groups
## Mbeya 4174.907 a
## CA 3236.625 b
## Intercrop 3179.784 b
out
## $statistics
## MSerror Df Mean CV
## 1623275 23 3386.601 37.62115
##
## $parameters
## test p.ajusted name.t ntr alpha
## Fisher-LSD hommel cropping_system 3 0.05
##
## $means
## grain std r LCL UCL Min Max Q25
## CA 3236.625 2224.705 56 2884.424 3588.826 0.000 10082.23 1566.081
## Intercrop 3179.784 2519.803 28 2681.696 3677.871 539.844 10186.92 1681.100
## Mbeya 4174.907 2977.093 18 3553.682 4796.131 813.145 11091.00 1483.844
## Q50 Q75
## CA 2791.701 4219.068
## Intercrop 2216.329 4155.100
## Mbeya 4294.150 6186.985
##
## $comparison
## NULL
##
## $groups
## grain groups
## Mbeya 4174.907 a
## CA 3236.625 b
## Intercrop 3179.784 b
##
## attr(,"class")
## [1] "group"
##################################################
md_social<-lm(grain~age+sex+hhsize+experience_CA,data = yld)
summary(md_social)
##
## Call:
## lm(formula = grain ~ age + sex + hhsize + experience_CA, data = yld)
##
## Residuals:
## Min 1Q Median 3Q Max
## -3341.4 -1680.6 -675.3 1074.8 6872.6
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 179.21 1495.05 0.120 0.9048
## age 70.89 30.51 2.324 0.0222 *
## sexMale 414.13 568.48 0.728 0.4680
## hhsize 13.13 123.84 0.106 0.9158
## experience_CA -78.41 31.49 -2.490 0.0144 *
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 2349 on 99 degrees of freedom
## (20 observations deleted due to missingness)
## Multiple R-squared: 0.1312, Adjusted R-squared: 0.09614
## F-statistic: 3.739 on 4 and 99 DF, p-value: 0.007086
anova(md_social)
## Analysis of Variance Table
##
## Response: grain
## Df Sum Sq Mean Sq F value Pr(>F)
## age 1 32457026 32457026 5.8811 0.01712 *
## sex 1 15804384 15804384 2.8637 0.09374 .
## hhsize 1 50810 50810 0.0092 0.92375
## experience_CA 1 34224125 34224125 6.2013 0.01443 *
## Residuals 99 546363913 5518827
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
### mean and sample size function
stat_box_data <- function(grain, upper_limit = max(grain)*1.5) {
return(
data.frame(
y =11000,
label = paste('N =', length(grain), '\n',
'Mean =', round(mean(grain), 0), '\n')
)
)
}
#######################################################################
##tech in general
theme_set(theme_gray(base_size = 10))
ggplot(yld, aes(x =tech, y = grain))+ geom_boxplot(size=0.7,varwidth = FALSE,outlier.colour = "red",outlier.shape = 2, shape=2,fill=c("grey","grey")) + geom_smooth(method=lm)+ ylab("Grain Yield [kg/ha]") + xlab("")+stat_summary(fun.y=mean, geom="point", shape=10, size=2, color="blue", fill="red")+stat_summary(
fun.data = stat_box_data ,
geom = "text",
hjust = 0.5,
vjust = 0.9,
size=3,
color="blue"
)+ theme_classic(base_size = 12)+ theme(axis.text.x = element_text(face="bold", color="black",
size=12),
axis.text.y = element_text(face="bold", color="black",
size=11))+
scale_y_continuous(breaks = c(0, 1000, 2000,3000,4000,5000,6000,7000,8000,9000))+
geom_text(data = yld, x = 1, y = 4500, label = "a",color="blue",size=4)+
geom_text(data = yld, x = 2, y = 3800, label = "b",color="blue",size=4)
agric<-subset(yld,agroecology=="Low rainfall")
ag<-aov(grain~tech,data=agric)
h<-HSD.test(model_1,"tech",group =TRUE)
h
## $statistics
## MSerror Df Mean CV MSD
## 1623275 23 3386.601 37.62115 521.9327
##
## $parameters
## test name.t ntr StudentizedRange alpha
## Tukey tech 2 2.925523 0.05
##
## $means
## grain std r Min Max Q25 Q50 Q75
## CSA 3903.330 2662.174 51 539.844 11091.00 1817.968 3103.922 5272.583
## Non-CSA 2869.872 2127.760 51 0.000 8502.05 1421.559 2028.926 4131.753
##
## $comparison
## NULL
##
## $groups
## grain groups
## CSA 3903.330 a
## Non-CSA 2869.872 b
##
## attr(,"class")
## [1] "group"
out<-LSD.test(ag,"tech",p.adj="hommel",console=TRUE)
##
## Study: ag ~ "tech"
##
## LSD t Test for grain
## P value adjustment method: hommel
##
## Mean Square Error: 6475313
##
## tech, means and individual ( 95 %) CI
##
## grain std r LCL UCL Min Max
## CSA 3594.786 2950.114 22 2499.928 4689.644 539.844 11091.004
## Non-CSA 2674.928 2060.935 22 1580.071 3769.786 675.177 7786.251
##
## Alpha: 0.05 ; DF Error: 42
## Critical Value of t: 2.018082
##
## Minimum Significant Difference: 1548.363
##
## Treatments with the same letter are not significantly different.
##
## grain groups
## CSA 3594.786 a
## Non-CSA 2674.928 a
out
## $statistics
## MSerror Df Mean CV t.value MSD
## 6475313 42 3134.857 81.17319 2.018082 1548.363
##
## $parameters
## test p.ajusted name.t ntr alpha
## Fisher-LSD hommel tech 2 0.05
##
## $means
## grain std r LCL UCL Min Max Q25
## CSA 3594.786 2950.114 22 2499.928 4689.644 539.844 11091.004 1542.783
## Non-CSA 2674.928 2060.935 22 1580.071 3769.786 675.177 7786.251 1213.626
## Q50 Q75
## CSA 2317.118 4529.520
## Non-CSA 1732.519 4191.446
##
## $comparison
## NULL
##
## $groups
## grain groups
## CSA 3594.786 a
## Non-CSA 2674.928 a
##
## attr(,"class")
## [1] "group"
theme_set(theme_gray(base_size = 10))
ggplot(yld, aes(x =tech, y = grain))+ geom_boxplot(size=0.7,varwidth = FALSE,outlier.colour = "red",outlier.shape = 2, shape=2,fill=c("grey","grey","grey","grey")) + geom_smooth(method=lm)+ ylab("Grain Yield [kg/ha]") + xlab("")+stat_summary(fun.y=mean, geom="point", shape=10, size=2, color="blue", fill="red")+stat_summary(
fun.data = stat_box_data ,
geom = "text",
hjust = 0.5,
vjust = 0.9,
size=3,
color="blue"
)+ theme_classic(base_size = 12)+ theme(axis.text.x = element_text(face="bold", color="black",
size=12),
axis.text.y = element_text(face="bold", color="black",
size=11))+facet_wrap(.~agroecology)+
scale_y_continuous(breaks = c(0, 1000, 2000,3000,4000,5000,6000,7000,8000,9000))+
geom_text(data = yld, x = 1, y = 4400, label = "a",color="blue",size=4)+
geom_text(data = yld, x = 2, y = 3800, label = "b",color="blue",size=4)
theme_set(theme_gray(base_size = 10))
ggplot(yld, aes(x =cropping_system, y = grain))+ geom_boxplot(size=0.7,varwidth = FALSE,outlier.colour = "red",outlier.shape = 2, shape=2,fill=c("grey","grey","grey")) + geom_smooth(method=lm)+ ylab("Grain Yield [kg/ha]") + xlab("")+stat_summary(fun.y=mean, geom="point", shape=10, size=2, color="blue", fill="red")+stat_summary(
fun.data = stat_box_data ,
geom = "text",
hjust = 0.5,
vjust = 0.9,
size=3,
color="blue"
)+ theme_classic(base_size = 12)+ theme(axis.text.x = element_text(face="bold", color="black",
size=12),
axis.text.y = element_text(face="bold", color="black",
size=11))+
scale_y_continuous(breaks = c(0, 1000, 2000,3000,4000,5000,6000,7000,8000,9000))+
geom_text(data = yld, x = 1, y = 3800, label = "b",color="blue",size=4)+
geom_text(data = yld, x = 2, y = 3500, label = "b",color="blue",size=4)+
geom_text(data = yld, x = 3, y = 4800, label = "a",color="blue",size=4)
ca<-subset(yld,cropping_system=="CA")
a1<-aov(grain~class,data=ca)
out<-LSD.test(a1,"class",p.adj="hommel",console=TRUE)
##
## Study: a1 ~ "class"
##
## LSD t Test for grain
## P value adjustment method: hommel
##
## Mean Square Error: 5037591
##
## class, means and individual ( 95 %) CI
##
## grain std r LCL UCL Min Max
## CA 4106.249 2617.562 32 3313.121 4899.376 624.155 10082.226
## CONV 2713.485 1795.425 32 1920.358 3506.613 0.000 6777.729
##
## Alpha: 0.05 ; DF Error: 62
## Critical Value of t: 1.998972
##
## Minimum Significant Difference: 1121.652
##
## Treatments with the same letter are not significantly different.
##
## grain groups
## CA 4106.249 a
## CONV 2713.485 b
out
## $statistics
## MSerror Df Mean CV t.value MSD
## 5037591 62 3409.867 65.82245 1.998972 1121.652
##
## $parameters
## test p.ajusted name.t ntr alpha
## Fisher-LSD hommel class 2 0.05
##
## $means
## grain std r LCL UCL Min Max Q25 Q50
## CA 4106.249 2617.562 32 3313.121 4899.376 624.155 10082.226 2073.808 3392.720
## CONV 2713.485 1795.425 32 1920.358 3506.613 0.000 6777.729 1439.203 2196.515
## Q75
## CA 5472.482
## CONV 3842.888
##
## $comparison
## NULL
##
## $groups
## grain groups
## CA 4106.249 a
## CONV 2713.485 b
##
## attr(,"class")
## [1] "group"
####
biod<-lm(grain~tech+class+agroecology+residue_applied+plant_population+weeding+number_plants,data = ca)
summary(biod)
##
## Call:
## lm(formula = grain ~ tech + class + agroecology + residue_applied +
## plant_population + weeding + number_plants, data = ca)
##
## Residuals:
## Min 1Q Median 3Q Max
## -3461.4 -1134.1 -185.9 804.2 6000.7
##
## Coefficients: (1 not defined because of singularities)
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) -995.74131 1980.45835 -0.503 0.61712
## techNon-CSA -658.82627 1188.67884 -0.554 0.58165
## classCONV NA NA NA NA
## agroecologyLow rainfall -55.37537 571.52509 -0.097 0.92317
## residue_appliedYes 547.00610 1201.05082 0.455 0.65059
## plant_population 0.08740 0.02927 2.986 0.00421 **
## weeding 917.61604 415.54828 2.208 0.03142 *
## number_plants -293.43188 549.25733 -0.534 0.59533
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 1960 on 55 degrees of freedom
## (2 observations deleted due to missingness)
## Multiple R-squared: 0.3434, Adjusted R-squared: 0.2718
## F-statistic: 4.794 on 6 and 55 DF, p-value: 0.0005361
anova(biod)
## Analysis of Variance Table
##
## Response: grain
## Df Sum Sq Mean Sq F value Pr(>F)
## tech 1 30738005 30738005 8.0044 0.0065003 **
## agroecology 1 51500 51500 0.0134 0.9082284
## residue_applied 1 3647145 3647145 0.9497 0.3340504
## plant_population 1 53565167 53565167 13.9488 0.0004478 ***
## weeding 1 21356511 21356511 5.5614 0.0219428 *
## number_plants 1 1095992 1095992 0.2854 0.5953327
## Residuals 55 211206536 3840119
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##################################
mbeya<-subset(yld,cropping_system=="Mbeya")
m1<-aov(grain~class,data=mbeya)
out<-LSD.test(m1,"class",p.adj="hommel",console=TRUE)
##
## Study: m1 ~ "class"
##
## LSD t Test for grain
## P value adjustment method: hommel
##
## Mean Square Error: 8266295
##
## class, means and individual ( 95 %) CI
##
## grain std r LCL UCL Min Max
## Mbeya fert 4535.499 3188.240 10 2625.358 6445.640 1339.578 11091.004
## Normal fert 3792.518 2523.434 10 1882.377 5702.659 813.145 7786.251
##
## Alpha: 0.05 ; DF Error: 18
## Critical Value of t: 2.100922
##
## Minimum Significant Difference: 2701.347
##
## Treatments with the same letter are not significantly different.
##
## grain groups
## Mbeya fert 4535.499 a
## Normal fert 3792.518 a
out
## $statistics
## MSerror Df Mean CV t.value MSD
## 8266295 18 4164.008 69.04685 2.100922 2701.347
##
## $parameters
## test p.ajusted name.t ntr alpha
## Fisher-LSD hommel class 2 0.05
##
## $means
## grain std r LCL UCL Min Max Q25
## Mbeya fert 4535.499 3188.240 10 2625.358 6445.640 1339.578 11091.004 1569.714
## Normal fert 3792.518 2523.434 10 1882.377 5702.659 813.145 7786.251 1502.000
## Q50 Q75
## Mbeya fert 4761.506 5797.551
## Normal fert 3800.357 5789.862
##
## $comparison
## NULL
##
## $groups
## grain groups
## Mbeya fert 4535.499 a
## Normal fert 3792.518 a
##
## attr(,"class")
## [1] "group"
######
biod<-lm(grain~tech+class+agroecology+residue_applied+plant_population+weeding+number_plants,data = mbeya)
summary(biod)
##
## Call:
## lm(formula = grain ~ tech + class + agroecology + residue_applied +
## plant_population + weeding + number_plants, data = mbeya)
##
## Residuals:
## Min 1Q Median 3Q Max
## -3393.1 -1293.6 -234.2 827.5 4896.5
##
## Coefficients: (1 not defined because of singularities)
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 1.100e+04 4.803e+03 2.289 0.0394 *
## techNon-CSA -7.127e+02 1.181e+03 -0.604 0.5565
## classNormal fert NA NA NA NA
## agroecologyLow rainfall 7.874e+01 1.585e+03 0.050 0.9611
## residue_appliedYes -1.822e+03 1.429e+03 -1.275 0.2247
## plant_population 7.451e-02 6.007e-02 1.240 0.2368
## weeding -2.004e+03 1.372e+03 -1.461 0.1677
## number_plants -2.698e+03 1.168e+03 -2.310 0.0379 *
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 2582 on 13 degrees of freedom
## Multiple R-squared: 0.428, Adjusted R-squared: 0.164
## F-statistic: 1.621 on 6 and 13 DF, p-value: 0.2185
anova(biod)
## Analysis of Variance Table
##
## Response: grain
## Df Sum Sq Mean Sq F value Pr(>F)
## tech 1 2760105 2760105 0.4139 0.53119
## agroecology 1 1143890 1143890 0.1715 0.68551
## residue_applied 1 72085 72085 0.0108 0.91878
## plant_population 1 5446893 5446893 0.8168 0.38256
## weeding 1 19846756 19846756 2.9761 0.10818
## number_plants 1 35591268 35591268 5.3371 0.03793 *
## Residuals 13 86692425 6668648
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
######################################
intercrop<-subset(yld,cropping_system=="Intercrop")
int<-aov(grain~class,data=intercrop)
out<-LSD.test(int,"class",p.adj="hommel",console=TRUE)
##
## Study: int ~ "class"
##
## LSD t Test for grain
## P value adjustment method: hommel
##
## Mean Square Error: 5138851
##
## class, means and individual ( 95 %) CI
##
## grain std r LCL UCL Min Max
## Intercrop maize 3047.227 2329.012 20 2021.071 4073.382 539.844 10186.92
## Monocrop maize 2802.389 2203.044 20 1776.234 3828.544 884.808 8502.05
##
## Alpha: 0.05 ; DF Error: 38
## Critical Value of t: 2.024394
##
## Minimum Significant Difference: 1451.203
##
## Treatments with the same letter are not significantly different.
##
## grain groups
## Intercrop maize 3047.227 a
## Monocrop maize 2802.389 a
out
## $statistics
## MSerror Df Mean CV t.value MSD
## 5138851 38 2924.808 77.50606 2.024394 1451.203
##
## $parameters
## test p.ajusted name.t ntr alpha
## Fisher-LSD hommel class 2 0.05
##
## $means
## grain std r LCL UCL Min Max
## Intercrop maize 3047.227 2329.012 20 2021.071 4073.382 539.844 10186.92
## Monocrop maize 2802.389 2203.044 20 1776.234 3828.544 884.808 8502.05
## Q25 Q50 Q75
## Intercrop maize 1726.993 2322.294 3911.369
## Monocrop maize 1629.169 1931.782 3184.842
##
## $comparison
## NULL
##
## $groups
## grain groups
## Intercrop maize 3047.227 a
## Monocrop maize 2802.389 a
##
## attr(,"class")
## [1] "group"
######
biod<-lm(grain~tech+class+agroecology+residue_applied+plant_population+weeding+number_plants,data = intercrop)
summary(biod)
##
## Call:
## lm(formula = grain ~ tech + class + agroecology + residue_applied +
## plant_population + weeding + number_plants, data = intercrop)
##
## Residuals:
## Min 1Q Median 3Q Max
## -3176.3 -1247.0 -43.6 783.3 4743.7
##
## Coefficients: (1 not defined because of singularities)
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 1.044e+04 3.028e+03 3.448 0.00156 **
## techNon-CSA -3.705e+02 5.485e+02 -0.675 0.50407
## classMonocrop maize NA NA NA NA
## agroecologyLow rainfall -1.017e+03 5.998e+02 -1.695 0.09942 .
## residue_appliedYes -8.617e+01 8.326e+02 -0.103 0.91820
## plant_population 5.548e-02 3.974e-02 1.396 0.17208
## weeding -3.806e+03 1.087e+03 -3.502 0.00135 **
## number_plants -1.099e+03 4.003e+02 -2.746 0.00970 **
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 1703 on 33 degrees of freedom
## Multiple R-squared: 0.5116, Adjusted R-squared: 0.4228
## F-statistic: 5.762 on 6 and 33 DF, p-value: 0.0003457
anova(biod)
## Analysis of Variance Table
##
## Response: grain
## Df Sum Sq Mean Sq F value Pr(>F)
## tech 1 599456 599456 0.2068 0.652270
## agroecology 1 21621126 21621126 7.4586 0.010055 *
## residue_applied 1 1604797 1604797 0.5536 0.462115
## plant_population 1 23890176 23890176 8.2413 0.007097 **
## weeding 1 30644062 30644062 10.5712 0.002648 **
## number_plants 1 21854922 21854922 7.5392 0.009696 **
## Residuals 33 95661247 2898826
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
################################# facet wrap by cropping system
theme_set(theme_gray(base_size = 10))
ggplot(ca, aes(x =class, y = grain))+ geom_boxplot(size=0.7,varwidth = FALSE,outlier.colour = "red",outlier.shape = 2, shape=2,fill=c("grey","grey")) + geom_smooth(method=lm)+ ylab("Grain Yield [kg/ha]") + xlab("")+stat_summary(fun.y=mean, geom="point", shape=10, size=2, color="blue", fill="red")+stat_summary(
fun.data = stat_box_data ,
geom = "text",
hjust = 0.5,
vjust = 0.9,
size=4,
color="blue"
)+ theme_bw(base_size = 12)+ theme(axis.text.x = element_text(face="bold", color="black",
size=12),
axis.text.y = element_text(face="bold", color="black",
size=12))+
scale_y_continuous(breaks = c(0, 1000, 2000,3000,4000,5000,6000,7000,8000,9000)) +
geom_text(data = ca, x = 1, y = 5000, label = "a",color="blue",size=5)+
geom_text(data = ca, x = 2, y = 3500, label = "b",color="blue",size=5)
#+facet_wrap(.~agroecology)
####################################################################################
theme_set(theme_gray(base_size = 10))
ggplot(mbeya, aes(x =class, y = grain))+ geom_boxplot(size=0.7,varwidth = FALSE,outlier.colour = "red",outlier.shape = 2, shape=2,fill=c("grey","grey")) + geom_smooth(method=lm)+ ylab("Grain Yield [kg/ha]") + xlab("")+stat_summary(fun.y=mean, geom="point", shape=10, size=2, color="blue", fill="red")+stat_summary(
fun.data = stat_box_data ,
geom = "text",
hjust = 0.5,
vjust = 0.9,
size=4,
color="blue"
)+ theme_bw(base_size = 12)+ theme(axis.text.x = element_text(face="bold", color="black",
size=12),
axis.text.y = element_text(face="bold", color="black",
size=12))+
scale_y_continuous(breaks = c(0, 1000, 2000,3000,4000,5000,6000,7000,8000,9000)) +
geom_text(data = mbeya, x = 1, y = 3500, label = "a",color="blue",size=5)+
geom_text(data = mbeya, x = 2, y = 5000, label = "a",color="blue",size=5)
#+facet_wrap(.~agroecology)
####################################################################################
theme_set(theme_gray(base_size = 10))
ggplot(intercrop, aes(x =class, y = grain))+ geom_boxplot(size=0.7,varwidth = FALSE,outlier.colour = "red",outlier.shape = 2, shape=2,fill=c("grey","grey")) + geom_smooth(method=lm)+ ylab("Grain Yield [kg/ha]") + xlab("")+stat_summary(fun.y=mean, geom="point", shape=10, size=2, color="blue", fill="red")+stat_summary(
fun.data = stat_box_data ,
geom = "text",
hjust = 0.5,
vjust = 0.9,
size=4,
color="blue"
)+ theme_bw(base_size = 12)+ theme(axis.text.x = element_text(face="bold", color="black",
size=12),
axis.text.y = element_text(face="bold", color="black",
size=12))+
scale_y_continuous(breaks = c(0, 1000, 2000,3000,4000,5000,6000,7000,8000,9000)) +
geom_text(data = intercrop, x = 1, y = 7000, label = "a",color="blue",size=6)+
geom_text(data = intercrop, x = 2, y = 7000, label = "a",color="blue",size=6)
#+facet_wrap(.~agroecology)
########################## suggested error bar plots
ggplot(intercrop, aes(x=class, y=grain)) +
stat_summary(fun.y=mean,
geom="bar",position=position_dodge(),colour="black",width=.5,size=.5) +
stat_summary(fun.ymin=min,fun.ymax=max,geom="errorbar",
color="blue",position=position_dodge(.7), width=.1) +
scale_fill_manual("Legend")+
xlab("")+
ylab("Grain Yield [kg/ha] ") +
theme(panel.background = element_rect(fill = 'white', colour = 'black'))+theme_bw()+
scale_y_continuous(breaks = c(0, 1000, 2000,3000,4000,5000,6000,7000,8000,9000))+ theme(axis.text.x = element_text(face="bold", color="black",
size=12),
axis.text.y = element_text(face="bold", color="black",
size=10))
#################################################
ggplot(ca, aes(x=class, y=grain)) +
stat_summary(fun.y=mean,
geom="bar",position=position_dodge(),colour="black",width=.5,size=.5) +
stat_summary(fun.ymin=min,fun.ymax=max,geom="errorbar",
color="blue",position=position_dodge(.7), width=.1) +
scale_fill_manual("Legend")+
xlab("")+
ylab("Grain Yield [kg/ha] ") +
theme(panel.background = element_rect(fill = 'white', colour = 'black'))+theme_bw()+
scale_y_continuous(breaks = c(0, 1000, 2000,3000,4000,5000,6000,7000,8000,9000))+ theme(axis.text.x = element_text(face="bold", color="black",
size=12),
axis.text.y = element_text(face="bold", color="black",
size=10))
######################################
ggplot(mbeya, aes(x=class, y=grain)) +
stat_summary(fun.y=mean,
geom="bar",position=position_dodge(),colour="black",width=.5,size=.5) +
stat_summary(fun.ymin=min,fun.ymax=max,geom="errorbar",
color="blue",position=position_dodge(.7), width=.1) +
scale_fill_manual("Legend")+
xlab("")+
ylab("Grain Yield [kg/ha] ") +
theme(panel.background = element_rect(fill = 'white', colour = 'black'))+theme_bw()+
scale_y_continuous(breaks = c(0, 1000, 2000,3000,4000,5000,6000,7000,8000,9000))+ theme(axis.text.x = element_text(face="bold", color="black",
size=12),
axis.text.y = element_text(face="bold", color="black",
size=10))
######################### alll
##################
tgc <- summarySE(mbeya, measurevar="grain", groupvars=c("class","cropping_system"))
tgc
## class cropping_system N grain sd se ci
## 1 Mbeya fert Mbeya 10 4535.499 3188.240 1008.2099 2280.729
## 2 Normal fert Mbeya 10 3792.518 2523.434 797.9799 1805.156
# Error bars represent standard error of the mean
ggplot(tgc, aes(x=class, y=grain)) +
geom_bar(position=position_dodge(), stat="identity",width = 0.3) +
geom_errorbar(aes(ymin=grain-se, ymax=grain+se),
width=.1,color="blue" , # Width of the error bars
position=position_dodge(.9))+
ylab("Grain Yield [kg/ha] ")+theme_bw()+
scale_y_continuous(breaks = c(0, 1000, 2000,3000,4000,5000,6000,7000,8000,9000))+ theme(axis.text.x = element_text(face="bold", color="black",
size=13),
axis.text.y = element_text(face="bold",color="black",
size=10))+
geom_text(data =ca, x = 1, y = 3300, label = "a",color="blue",size=6)+
geom_text(data = ca, x = 2, y = 2500, label = "a",color="blue",size=6)
stat_summary(geom = ‘text’, fun.y = max, position = position_dodge(.7), label = c(“a”,“a”))+
out<-LSD.test(model_1,"soil_colour",p.adj="hommel",console=TRUE)
##
## Study: model_1 ~ "soil_colour"
##
## LSD t Test for grain
## P value adjustment method: hommel
##
## Mean Square Error: 1623275
##
## soil_colour, means and individual ( 95 %) CI
##
## grain std r LCL UCL Min Max
## Black 3285.136 2469.0849 56 2932.935 3637.337 622.803 11091.004
## Brown 3888.476 2593.3389 34 3436.469 4340.483 0.000 10186.923
## Grey 2831.054 1544.7883 10 1997.594 3664.514 539.844 4591.708
## Red 473.463 213.1107 2 -1390.210 2337.136 322.771 624.155
##
## Alpha: 0.05 ; DF Error: 23
## Critical Value of t: 2.068658
##
## Groups according to probability of means differences and alpha level( 0.05 )
##
## Treatments with the same letter are not significantly different.
##
## grain groups
## Brown 3888.476 a
## Black 3285.136 a
## Grey 2831.054 ab
## Red 473.463 b
out
## $statistics
## MSerror Df Mean CV
## 1623275 23 3386.601 37.62115
##
## $parameters
## test p.ajusted name.t ntr alpha
## Fisher-LSD hommel soil_colour 4 0.05
##
## $means
## grain std r LCL UCL Min Max Q25
## Black 3285.136 2469.0849 56 2932.935 3637.337 622.803 11091.004 1504.778
## Brown 3888.476 2593.3389 34 3436.469 4340.483 0.000 10186.923 1915.004
## Grey 2831.054 1544.7883 10 1997.594 3664.514 539.844 4591.708 1390.146
## Red 473.463 213.1107 2 -1390.210 2337.136 322.771 624.155 398.117
## Q50 Q75
## Black 2302.666 4794.970
## Brown 3397.566 4700.572
## Grey 3417.270 4084.907
## Red 473.463 548.809
##
## $comparison
## NULL
##
## $groups
## grain groups
## Brown 3888.476 a
## Black 3285.136 a
## Grey 2831.054 ab
## Red 473.463 b
##
## attr(,"class")
## [1] "group"
### soil color
theme_set(theme_gray(base_size = 10))
ggplot(yld, aes(x =soil_colour, y = grain))+ geom_boxplot(size=0.7,varwidth = FALSE,outlier.colour = "red",outlier.shape = 2, shape=2,fill=c("grey","grey","grey","grey")) + geom_smooth(method=lm)+ ylab("Grain Yield [kg/ha]") + xlab("")+stat_summary(fun.y=mean, geom="point", shape=10, size=2, color="blue", fill="red")+stat_summary(
fun.data = stat_box_data ,
geom = "text",
hjust = 0.5,
vjust = 0.9,
size=3,
color="blue"
)+ theme_bw(base_size = 12)+ theme(axis.text.x = element_text(face="bold", color="black",
size=10),
axis.text.y = element_text(face="bold", color="black",
size=8))+
scale_y_continuous(breaks = c(0, 1000, 2000,3000,4000,5000,6000,7000,8000,9000))+
geom_text(data = yld, x = 1, y = 4000, label = "b",color="blue",size=4)+
geom_text(data = yld, x = 2, y = 4000, label = "a",color="blue",size=4)+
geom_text(data = yld, x = 3, y = 3800, label = "b",color="blue",size=4)+
geom_text(data = yld, x = 4, y = 3000, label = "c",color="blue",size=4)
############## facet wrapping soil color by agroecology
theme_set(theme_gray(base_size = 10))
ggplot(yld, aes(x =soil_colour, y = grain))+ geom_boxplot(size=0.7,varwidth = FALSE,outlier.colour = "red",outlier.shape = 2, shape=2) + geom_smooth(method=lm)+ ylab("Grain Yield [kg/ha]") + xlab("")+stat_summary(fun.y=mean, geom="point", shape=10, size=2, color="blue", fill="red")+stat_summary(
fun.data = stat_box_data ,
geom = "text",
hjust = 0.5,
vjust = 0.9,
size=3,
color="blue"
)+ theme_bw(base_size = 12)+ theme(axis.text.x = element_text(face="bold", color="black",
size=10),
axis.text.y = element_text(face="bold", color="black",
size=8))+
scale_y_continuous(breaks = c(0, 1000, 2000,3000,4000,5000,6000,7000,8000,9000))+facet_wrap(.~agroecology)
################## wrapping by soil color
theme_set(theme_gray(base_size = 10))
ggplot(yld, aes(x =tech, y = grain))+ geom_boxplot(size=0.7,width = 0.5,outlier.colour = "red",outlier.shape = 2, shape=2,fill=c("grey","grey","grey","grey","grey","grey","grey","grey")) + geom_smooth(method=lm)+ ylab("Grain Yield [kg/ha]") + xlab("")+stat_summary(fun.y=mean, geom="point", shape=10, size=2, color="blue", fill="red")+stat_summary(
fun.data = stat_box_data ,
geom = "text",
hjust = 0.5,
vjust = 0.9,
size=3,
color="blue"
)+ theme_bw(base_size = 12)+ theme(axis.text.x = element_text(face="bold", color="black",
size=10),
axis.text.y = element_text(face="bold", color="black",
size=8))+
scale_y_continuous(breaks = c(0, 1000,3000,5000,7000,9000))+facet_wrap(.~soil_colour)+
geom_text(data = yld, x = 1, y = 8000, label = "a",color="blue",size=4)+
geom_text(data = yld, x = 2, y = 8000, label = "a",color="blue",size=4)
#### grain yield and texture
theme_set(theme_gray(base_size = 10))
ggplot(yld, aes(x =tech, y = grain))+ geom_boxplot(size=0.7,varwidth = FALSE,outlier.colour = "red",outlier.shape = 2, shape=2,fill=c("grey","grey","grey","grey","grey","grey")) + geom_smooth(method=lm)+ ylab("Grain Yield [kg/ha]") + xlab("")+stat_summary(fun.y=mean, geom="point", shape=10, size=2, color="blue", fill="red")+stat_summary(
fun.data = stat_box_data ,
geom = "text",
hjust = 0.5,
vjust = 0.9,
size=3,
color="blue"
)+ theme_bw(base_size = 12)+ theme(axis.text.x = element_text(face="bold", color="black",
size=7),
axis.text.y = element_text(face="bold", color="black",
size=8))+
scale_y_continuous(breaks = c(0, 1000, 2000,3000,4000,5000,6000,7000,8000,9000))+
facet_wrap(.~texture)+
geom_text(data = yld, x = 1, y = 9000, label = "a",color="blue",size=6)+
geom_text(data = yld, x = 2, y = 9000, label = "b",color="blue",size=6)
p <- ggplot(yld, aes(x=variet))
p + geom_bar(aes(y = (..count..)/sum(..count..)))+ theme_bw(base_size = 12)+ theme(axis.text.x = element_text(face="bold", color="black",
size=8,angle = 90),
axis.text.y = element_text(face="bold", color="black",
size=8))+
scale_y_continuous() +
geom_text(aes( label = scales::percent(..count../sum(..count..)),
y= ..count../sum(..count..)), stat= "count", vjust = -.5,size=2.5)+ ylab("Proportion of Varieties Grown")+xlab("Varieties Grown")
#### anova model
resiapp<-aov(grain~residue_applied+agroecology,data = yld)
out<-LSD.test(resiapp,"agroecology",p.adj="hommel",console=TRUE)
##
## Study: resiapp ~ "agroecology"
##
## LSD t Test for grain
## P value adjustment method: hommel
##
## Mean Square Error: 5588037
##
## agroecology, means and individual ( 95 %) CI
##
## grain std r LCL UCL Min Max
## High rainfall 3507.128 2322.774 80 2983.892 4030.364 0.000 10186.92
## Low rainfall 3134.857 2557.573 44 2429.325 3840.389 539.844 11091.00
##
## Alpha: 0.05 ; DF Error: 121
## Critical Value of t: 1.979764
##
## Groups according to probability of means differences and alpha level( 0.05 )
##
## Treatments with the same letter are not significantly different.
##
## grain groups
## High rainfall 3507.128 a
## Low rainfall 3134.857 a
out
## $statistics
## MSerror Df Mean CV
## 5588037 121 3375.032 70.0409
##
## $parameters
## test p.ajusted name.t ntr alpha
## Fisher-LSD hommel agroecology 2 0.05
##
## $means
## grain std r LCL UCL Min Max Q25
## High rainfall 3507.128 2322.774 80 2983.892 4030.364 0.000 10186.92 1732.171
## Low rainfall 3134.857 2557.573 44 2429.325 3840.389 539.844 11091.00 1332.043
## Q50 Q75
## High rainfall 2986.055 4741.526
## Low rainfall 2067.860 4300.314
##
## $comparison
## NULL
##
## $groups
## grain groups
## High rainfall 3507.128 a
## Low rainfall 3134.857 a
##
## attr(,"class")
## [1] "group"
### effect of residue application
theme_set(theme_gray(base_size = 10))
ggplot(yld, aes(x =residue_applied, y = grain))+ geom_boxplot(size=0.7,varwidth = FALSE,outlier.colour = "red",outlier.shape = 2, shape=2,fill=c("grey","grey")) + geom_smooth(method=lm)+ ylab("Grain Yield [kg/ha]") + xlab("Residue Application")+stat_summary(fun.y=mean, geom="point", shape=10, size=2, color="blue", fill="red")+stat_summary(
fun.data = stat_box_data ,
geom = "text",
hjust = 0.5,
vjust = 0.9,
size=3.3,
color="blue"
)+ theme_bw(base_size = 12)+ theme(axis.text.x = element_text(face="bold", color="black",
size=10),
axis.text.y = element_text(face="bold", color="black",
size=8))+
scale_y_continuous(breaks = c(0, 1000, 2000,3000,4000,5000,6000,7000,8000,9000))+
geom_text(data = yld, x = 1, y = 3500, label = "b",color="blue",size=5)+
geom_text(data = yld, x = 2, y = 4500, label = "a",color="blue",size=5)
################## gender by technology
theme_set(theme_gray(base_size = 10))
ggplot(yld, aes(x =sex, y = grain))+ geom_boxplot(size=0.7,varwidth = FALSE,outlier.colour = "red",outlier.shape = 2, shape=2,fill=c("grey","grey","grey","grey")) + geom_smooth(method=lm)+ ylab("Grain Yield [kg/ha]") + xlab("")+stat_summary(fun.y=mean, geom="point", shape=10, size=2, color="blue", fill="red")+stat_summary(
fun.data = stat_box_data ,
geom = "text",
hjust = 0.5,
vjust = 0.9,
size=3,
color="blue"
)+ theme_bw(base_size = 12)+ theme(axis.text.x = element_text(face="bold", color="black",
size=10),
axis.text.y = element_text(face="bold", color="black",
size=8))+
scale_y_continuous(breaks = c(0, 1000, 2000,3000,4000,5000,6000,7000,8000,9000))+facet_wrap(.~tech)+
geom_text(data = yld, x = 1, y = 9000, label = "b",color="blue",size=5)+
geom_text(data = yld, x = 2, y = 9000, label = "a",color="blue",size=5)
socio
################################## experience by sex
ggplot(yld, aes(x =experience_CA ,y =grain)) +
geom_point(size=3) + geom_smooth(method=lm, position = "jitter", level = 0.95, formula = y ~ x + I(x^2))+ylab("Grain Yield [kg/ha]") +
xlab("Experience in CA [yrs]")+theme_bw()+ theme(axis.text.x = element_text(face="bold", color="black", size=10),
axis.text.y = element_text(face="bold", color="black",
size=10))+ theme( axis.line = element_line(colour = "black",
size = 0.8, linetype = "solid"))+xlim(0,10)+
stat_cor(label.x = 6) +
scale_y_continuous(breaks = c(0, 1000, 2000,3000,4000,5000,6000,7000,8000,9000))
## Warning: Removed 20 rows containing non-finite values (stat_smooth).
## Warning: Removed 20 rows containing non-finite values (stat_cor).
## Warning: Removed 20 rows containing missing values (geom_point).
####################################age in general
ggplot(yld, aes(x =age ,y =grain)) +
geom_point(size=3,shape=2,color="black") + geom_smooth(method=lm, position = "jitter", level = 0.95, formula = y ~ x + I(x^2)+I(x^3))+ylab("Grain Yield [kg/ha]") +
xlab("Age[yrs]")+theme_bw()+ theme(axis.text.x = element_text(color="black", size=10),
axis.text.y = element_text(color="black",
size=10))+ theme( axis.line = element_line(colour = "black",
size = 0.8, linetype = "solid"))+ stat_cor(method = "pearson", label.x = 55, label.y = 9000)+xlim(30,70)+
scale_y_continuous(breaks = c(0, 1000, 2000,3000,4000,5000,6000,7000,8000,9000))
## Warning: Removed 8 rows containing non-finite values (stat_smooth).
## Warning: Removed 8 rows containing non-finite values (stat_cor).
## Warning: Removed 8 rows containing missing values (geom_point).
#######################################
################################################# age by sex
ggplot(yld, aes(x =age ,y =grain,color=sex,shape=sex)) +
geom_point(size=3) + geom_smooth(method=lm, position = "jitter", level = 0.95, formula = y ~ x + I(x^2))+ylab("Grain Yield [kg/ha]") +
xlab("Age[yrs]")+theme_bw()+ theme(axis.text.x = element_text(face="bold", color="black", size=10),
axis.text.y = element_text(face="bold", color="black",
size=10))+ theme( axis.line = element_line(colour = "black",
size = 0.8, linetype = "solid"))+
stat_cor(aes(color = sex), label.x = 45) +
scale_y_continuous(breaks = c(0, 1000, 2000,3000,4000,5000,6000,7000,8000,9000))
## Warning: Removed 8 rows containing non-finite values (stat_smooth).
## Warning: Removed 8 rows containing non-finite values (stat_cor).
## Warning: Removed 8 rows containing missing values (geom_point).
plant population
## effect of plant population
ggplot(yld, aes(x =plant_population,y =grain)) +
geom_point(size=3) + geom_smooth(method=lm, position = "jitter", level = 0.95)+ylab("Grain Yield [kg/ha]") + xlab("Plant population")+theme_bw()+ stat_cor(method = "pearson", label.x = 5000, label.y = 9000)+ theme(axis.text.x = element_text(face="bold", color="black", size=10),
axis.text.y = element_text(face="bold", color="black",
size=10))+
scale_y_continuous(breaks = c(0, 1000, 2000,3000,4000,5000,6000,7000,8000,9000))+
scale_x_continuous(breaks = c(0, 10000, 20000,30000,40000,50000,60000,70000))
################################################
#### facet wrap by technology
ggplot(yld, aes(x =plant_population,y =grain,color=tech,shape=tech)) +
geom_point(size=3) + geom_smooth(method=lm, position = "jitter", level = 0.95, formula = y ~ x + I(x^2))+ylab("Grain Yield [kg/ha]") + xlab("Plant population")+theme_bw()+ theme(axis.text.x = element_text(face="bold", color="black", size=10),
axis.text.y = element_text(face="bold", color="black",
size=10))+
scale_y_continuous(breaks = c(0, 1000, 2000,3000,4000,5000,6000,7000,8000,9000))+
scale_x_continuous(breaks = c(0, 10000, 20000,30000,40000,50000,60000,70000))+
stat_cor(aes(color = tech), label.x = 5000)
##################################################
### using the facet wrap approach
ggplot(yld, aes(x =plant_population,y =grain,color=tech,shape=tech)) +
geom_point(size=3) + geom_smooth(method=lm, position = "jitter", level = 0.95)+ylab("Grain Yield [kg/ha]") + xlab("Plant population")+theme_bw()+ theme(axis.text.x = element_text(color="black", size=7),
axis.text.y = element_text(face="bold", color="black",
size=10))+
scale_y_continuous(breaks = c(0, 1000, 2000,3000,4000,5000,6000,7000,8000,9000))+
scale_x_continuous(breaks = c(0, 10000,30000,50000,70000))+facet_wrap(.~cropping_system)+
stat_cor(aes(color = tech), label.x = 5000)
Effects of weeding
## effects of weeding
ggplot(ca, aes(x =weeding,y =grain)) +
geom_point(size=3,shape=2) + geom_smooth(method=lm, position = "jitter", level = 0.95)+ylab("Grain Yield [kg/ha]") + xlab("Weeding")+theme_bw()+ stat_cor(method = "pearson", label.x = 2, label.y = 9000)+ theme(axis.text.x = element_text(face="bold", color="black", size=10),
axis.text.y = element_text(face="bold", color="black",
size=10))+
scale_y_continuous(breaks = c(0, 1000, 2000,3000,4000,5000,6000,7000,8000,9000))
## Warning: Removed 2 rows containing non-finite values (stat_smooth).
## Warning: Removed 2 rows containing non-finite values (stat_cor).
## Warning: Removed 2 rows containing missing values (geom_point).
#######################################################
#######################################
ggplot(ca, aes(x =weeding,y =grain,color=tech,shape=tech)) +
geom_point(size=3) + geom_smooth(method=lm, position = "jitter", level = 0.95, formula = y ~ x + I(x^2))+ylab("Grain Yield [kg/ha]") + xlab("Weeding")+theme_bw()+ stat_cor(method = "pearson", label.x = 1, label.y = 9000)+ theme(axis.text.x = element_text(face="bold", color="black", size=10),
axis.text.y = element_text(face="bold", color="black",
size=10))+
scale_y_continuous(breaks = c(0, 1000, 2000,3000,4000,5000,6000,7000,8000,9000))+facet_wrap(.~tech) + theme(legend.position="none")
## Warning: Removed 2 rows containing non-finite values (stat_smooth).
## Warning: Removed 2 rows containing non-finite values (stat_cor).
## Warning: Removed 2 rows containing missing values (geom_point).
Number of plants
###### number of plants
## number_plants
ggplot(yld, aes(x =number_plants,y =grain)) +
geom_bar(stat="identity")+ylab("Grain Yield [kg/ha]") + xlab("Number of plants")+theme_bw()+ theme(axis.text.x = element_text(face="bold", color="black", size=10),
axis.text.y = element_text(face="bold", color="black",
size=10))
####################################
ggplot(yld, aes(x =number_plants,y =grain))+
geom_bar(stat="identity")+ylab("Grain Yield [kg/ha]") + xlab("Number of plants")+theme_bw()+ theme(axis.text.x = element_text(face="bold", color="black", size=10),
axis.text.y = element_text(face="bold", color="black",
size=10))+facet_wrap(.~tech)+ scale_y_continuous(breaks = c(0, 1000, 2000,3000,4000,5000,6000,7000,8000,9000))
library(readr)
rainfall <- read_csv("~/Dr Nyagumbo/Chemical analysis/final yield with socio-analysis/rainfall.csv")
## Parsed with column specification:
## cols(
## month = col_character(),
## label = col_character(),
## rfall = col_double(),
## cumm = col_double()
## )
attach(rainfall)
month=as.factor(month)
label=as.factor(label)
mo=factor(month,levels =c("Sep","Oct","Nov","Dec","Jan","Feb","Mar","Apr","May","Jun","Jul"),ordered = TRUE)
ggplot(rainfall, aes(x =mo,y =rfall,color=label,shape=label)) + geom_smooth(method=lm, position = "jitter", level = 0.95, formula = y ~ x + I(x^2))+ylab("Average Rainfall [mm]") + xlab("Month")+theme_bw()+ theme(axis.text.x = element_text(face="bold", color="black", size=10),
axis.text.y = element_text(face="bold", color="black",
size=10))+ theme(legend.position = c(0.8, 0.75))
###########################################################
ggplot(data=rainfall, aes(x=mo, y=rfall, group=label,color=label,shape=label)) +
geom_line(aes(linetype=label),size=1)+
geom_point(aes(shape=label))+theme_bw()+ theme(axis.text.x = element_text(face="bold", color="black", size=10),
axis.text.y = element_text(face="bold", color="black",
size=10))+ theme(legend.position = c(0.8, 0.75))+ylab(" Average Rainfall [mm]") + xlab("Month")
#######################################################################
ggplot(data=rainfall, aes(x=mo, y=cumm, group=label,color=label,shape=label)) +
geom_line(aes(linetype=label),size=1)+
geom_point(aes(shape=label),size=2)+theme_bw()+ theme(axis.text.x = element_text(face="bold", color="black", size=10),
axis.text.y = element_text(face="bold", color="black",
size=10))+ theme(legend.position = c(0.8, 0.45))+
scale_y_continuous(breaks = c(0,200,400,600,800,1000,1200,1400,1600))+ylab("Commulative Average Rainfall [mm]") + xlab("Month")+
geom_text(data = rainfall, x = 7, y = 1200, label = "Cyclone IDAI (March 14 and 15, 2019) ",color="blue",size=3)
##############################################################
ggplot(data=rainfall, aes(x=mo, y=cumm, group=label,color=label,shape=label))+
geom_smooth(method = "loess",level = 0.95)+theme_bw()+ theme(axis.text.x = element_text(face="bold", color="black", size=10),
axis.text.y = element_text(face="bold", color="black",
size=10))+ theme(legend.position = c(0.8, 0.45))+
scale_y_continuous(breaks = c(0,200,400,600,800,1000,1200,1400,1600))+ylab("Commulative Average Rainfall [mm]") + xlab("Month")
geom_text(data = rainfall, x = 7, y = 1200, label = “Cyclone IDAI (March 14 and 15, 2019)”,color=“blue”,size=3)
geom_segment(aes(x=7,xend=9,y=1400,yend=1500)
, arrow=arrow(length=unit(0.6,“cm”) ),color=“blue” )