#STATISTICS SECTION # compare the mean CPUE of legal-sized crab between baited traps and tested traps using Non-parametric Wilcoxon Rank-Sum Test > Catch=read.csv("C:/KHANH/PhD/snow crab/2017 with Bob/field data/CPUE_MannwhitneyUtest.csv", header=T) > attach(Catch) > names(Catch) > wilcox.test(Control,HU) > wilcox.test(Control,HUD) > wilcox.test(Control,LU) > wilcox.test(Control,LUD) > wilcox.test(HU,HUD) > wilcox.test(HU,LU) > wilcox.test(HU,LUD) > wilcox.test(HUD,LU) > wilcox.test(HUD,LUD) > wilcox.test(LU,LUD) #COMPARE ON sublegal > subCatch=read.csv("C:/KHANH/PhD/snow crab/2017 with Bob/field data/CPUE_MannwhitneyUtest_sublegal.csv", header=T) > attach(subCatch) > wilcox.test(Control,HU) > wilcox.test(Control,HUD) > wilcox.test(Control,LU) > wilcox.test(Control,LUD) > wilcox.test(HU,HUD) > wilcox.test(HU,LU) > wilcox.test(HU,LUD) > wilcox.test(HUD,LU) > wilcox.test(HUD,LUD) > wilcox.test(LU,LUD) # Test mean CW on diff treatment using MannwhitneyUtest Lengthdis=read.csv("C:/KHANH/PhD/snow crab/2017 with Bob/field data/length.csv", header=T) attach(Lengthdis) > wilcox.test(Control,HU) > wilcox.test(Control,HUD) > wilcox.test(Control,LU) > wilcox.test(Control,LUD) > wilcox.test(HU,HUD) > wilcox.test(HU,LU) > wilcox.test(HU,LUD) > wilcox.test(HUD,LU) > wilcox.test(HUD,LUD) > wilcox.test(LU,LUD) # Compare diff soak time > Difsoak=read.csv("C:/KHANH/PhD/snow crab/2017 with Bob/field data/CPUE_MannwhitneyUtest_difday.csv", header=T) > attach(Difsoak) > names(Difsoak) > wilcox.test(Control_lagal4day,Control_lagal15day) > wilcox.test(light_lagal4day,light_lagal15day) > wilcox.test(Control_sublagal4day,Control_sublagal15day) > wilcox.test(light_sublagal4day,light_sublagal15day) # Compare CPUE between light and no light trap > LvsnoL=read.csv("C:/KHANH/PhD/snow crab/2017 with Bob/field data/CPUE_MannwhitneuUtest_light and nolight.csv", header=T) > attach(LvsnoL) > names(LvsnoL) #for legal sized > wilcox.test(CPUE_legal_light,CPUE_legal_nolight) #for sublegal sized > wilcox.test(CPUE_sublegal_light,CPUE_sublegal_nolight) # Compare CPUE between positions > pos=read.csv("C:/KHANH/PhD/snow crab/2017 with Bob/field data/CPUE_MannwhitneuUtest_positon.csv", header=T) > attach(pos) > names(pos) #for legal sised > wilcox.test(CPUE_legal_high,CPUE_legal_low) #for sublegal sized > wilcox.test(CPUE_sublegal_high,CPUE_sublegal_low) # Compare CPUE between orientations > orientation=read.csv("C:/KHANH/PhD/snow crab/2017 with Bob/field data/CPUE_MannwhitneuUtest_orientation.csv", header=T) > attach(orientation) > names(orientation) #for legal sized > wilcox.test(CPUE_legal_upright,CPUE_legal_upsidedown) # for sublegal sized > wilcox.test(CPUE_sublegal_upright,CPUE_sublegal_upsidedown) # The crab size frequency distributions between the experimental treatments were compared using Kolmogorov-Smirnov two-sample Z test > Lengthdis=read.csv("C:/KHANH/PhD/snow crab/2017 with Bob/field data/length.csv", header=T) > attach(Lengthdis) > ks.test(Control,HU) > ks.test(Control,HUD) > ks.test(Control,LU) > ks.test(Control,LUD) > ks.test(HU,HUD) > ks.test(HU,LU) > ks.test(HU,LUD) > ks.test(HUD,LU) > ks.test(HUD,LUD) > ks.test(LU,LUD) # test overdispersion for the legal sized crab > Data=read.csv("C:/KHANH/PhD/snow crab/2017 with Bob/field data/CPUE.csv", header=T) > attach(Data) > names(Data) #check overdispersion > model_total = glm(CPUE_total ~ Treatment*Soak_Poisson, quasipoisson) > summary(model_total) # Dispersion parameter for quasipoisson family taken to be 3.77. thus we can conclude that overdispersion. We have to use negative binomila distribution > model = glm(CPUE_legal ~ Treatment+Soak_Poisson, quasipoisson) > summary(model) # for sublegal-sized > model2 = glm(CPUE_sublegal ~ Treatment+Soak_Poisson, quasipoisson) > summary(model2) # Dispersion parameter for quasipoisson family taken to be 2.399. thus we can conclude that overdispersion. We have to use negative binomila distribution > require(MASS) #legal > m1 = glm.nb(CPUE_legal ~ Treatment*Soak_Poisson, data = Data) > summary(m1) # No interaction term > m1 = glm.nb(CPUE_legal ~ Treatment+Soak_Poisson, data = Data) > summary(m1) #sublegal > m2 = glm.nb(CPUE_sublegal ~ Treatment*Soak_Poisson, data = Data) > summary(m2) # PLOT SECTION > Data=read.csv("C:/KHANH/PhD/snow crab/2017 with Bob/field data/CPUE.csv", header=T) > attach(Data) > names(Data) > library(easyGgplot2) # plot CPUE versus diff Treatments > ggplot2.boxplot(data=Data, xName="Treatment",yName="CPUE_legal", groupName="Treatment",groupColors=c("white","white","white", "white","white"),xTickLabelFont=c(12, "plain", "black"),ytitle="CPUE",ytitleFont=c(12,"plain", "black"),yTickLabelFont=c(12, "plain", "black"),axisLine=c(1, "solid", "black"),backgroundColor="white",showLegend=FALSE,xtitle="Treatment", xtitleFont=c(12,"plain", "black")) #sublegal sized > ggplot2.boxplot(data=Data, xName="Treatment",yName="CPUE_sublegal", groupName="Treatment",groupColors=c("white","white","white", "white","white"),xTickLabelFont=c(12, "plain", "black"),ytitle="CPUE",ytitleFont=c(12,"plain", "black"),yTickLabelFont=c(12, "plain", "black"),axisLine=c(1, "solid", "black"),backgroundColor="white",showLegend=FALSE,xtitle="") # Light vs nolight > LightvsNo=read.csv("C:/KHANH/PhD/snow crab/2017 with Bob/field data/CPUE_lightvsnolight.csv", header=T) > attach(LightvsNo) > ggplot2.boxplot(data=LightvsNo, xName="light_vs_nolight",yName="CPUE_light_vs_nolight", groupName="size_light_vs_nolight",groupColors=c("white","grey"), showLegend=T, xtitle="", xTickLabelFont=c(12, "plain", "black"),ytitle="CPUE",ytitleFont=c(12,"plain", "black"),yTickLabelFont=c(12, "plain", "black"), axisLine=c(0.8, "solid", "black"),backgroundColor="white",legendTitle="") # plot CPUE for diff possitions > position=read.csv("C:/KHANH/PhD/snow crab/2017 with Bob/field data/CPUE_position.csv", header=T) > attach(position) > ggplot2.boxplot(data=position, xName="position",yName="CPUE_position", groupName="size_position",groupColors=c("white","grey"), showLegend=T, xtitle="", xTickLabelFont=c(12, "plain", "black"),ytitle="CPUE",ytitleFont=c(12,"plain", "black"),yTickLabelFont=c(12, "plain", "black"), axisLine=c(0.8, "solid", "black"),backgroundColor="white",legendTitle="") # plot CPUE for diff orientaions > orientation=read.csv("C:/KHANH/PhD/snow crab/2017 with Bob/field data/CPUE_orientation.csv", header=T) > attach(orientation) > ggplot2.boxplot(data=orientation, xName="orientation",yName="CPUE_orientation", groupName="size_orientation",groupColors=c("white","grey"), showLegend=T, xtitle="", xTickLabelFont=c(12, "plain", "black"),ytitle="CPUE",ytitleFont=c(12,"plain", "black"),yTickLabelFont=c(12, "plain", "black"), axisLine=c(0.8, "solid", "black"),backgroundColor="white",legendTitle="") # Plot CPUE for diff soak time > soak=read.csv("C:/KHANH/PhD/snow crab/2017 with Bob/field data/CPUE_difday.csv", header=T) > attach(soak) > names(soak) # for legal size > p=ggplot2.boxplot(data=soak, xName="Soak",yName="CPUE_legal", groupName="Treatment",groupColors=c("white","white"), showLegend=F, xtitle="Soak time (day)",xtitleFont=c(12,"plain", "black"), xTickLabelFont=c(12, "plain", "black"),ytitle="CPUE for legal-sized",ytitleFont=c(12,"plain", "black"),yTickLabelFont=c(12, "plain", "black"), axisLine=c(0.8, "solid", "black"),backgroundColor="white",legendTitle="",faceting=T, facetingVarNames="Treatment",facetingDirection="horizontal") > p+scale_x_discrete(limits=c("4", "6", "15")) # for sublegal size crab > p=ggplot2.boxplot(data=soak, xName="Soak",yName="CPUE_sublegal", groupName="Treatment",groupColors=c("white","white"), showLegend=F, xtitle="Soak time (day)",xtitleFont=c(12,"plain", "black"), xTickLabelFont=c(12, "plain", "black"),ytitle="CPUE for sublegal-sized",ytitleFont=c(12,"plain", "black"),yTickLabelFont=c(12, "plain", "black"), axisLine=c(0.8, "solid", "black"),backgroundColor="white",legendTitle="",faceting=T, facetingVarNames="Treatment",facetingDirection="horizontal") > p+scale_x_discrete(limits=c("4", "6", "15")) # plot CW versus diff Treatment > lengthmean=read.csv("C:/KHANH/PhD/snow crab/2017 with Bob/field data/CW.csv", header=T) > attach(lengthmean) > ggplot(lengthmean, aes(x=CW, group=Treatment))+geom_density(aes(colour=Treatment,linetype=Treatment),size=1)+scale_linetype_manual(values=c("solid","longdash","dotted","dashed","twodash"))+scale_color_manual(values=c("blue","purple","black","red","green"))+xlab("Carapace width (mm)")+ylab("Frequency in proportion")+theme_bw()+theme(axis.text.x = element_text(face="plain", color="black",size=12),axis.text.y = element_text(face="plain", color="black",size=12))+theme(axis.line = element_line(colour = "black", size = 0.5, linetype = "solid"))+theme(legend.title = element_text(colour="black", size=12, face="plain"))+xlim(70,130)+ theme(axis.text=element_text(size=12), axis.title=element_text(size=12,face="plain"))+ theme(legend.title=element_blank(),legend.key = element_rect(fill = "white", size = 0.5, linetype="solid"), legend.key.size = unit(1.5, "lines")) # economic feasibility > eco=read.csv("C:/KHANH/PhD/snow crab/2017 with Bob/field data/economic_feasibility.csv", header=T) > attach(eco) > names(eco) > library(easyGgplot2) > ggplot2.lineplot(data=eco, xName="Year", yName="Cost_new", groupName="Fishery", groupColors=c("black","black"), xtitle="Year", xtitleFont=c(12,"plain", "black"),xTickLabelFont=c(12, "plain", "black"),ytitle="Cumulative variable costs (CND 1000)",ytitleFont=c(12,"plain", "black"),yTickLabelFont=c(12, "plain", "black"), axisLine=c(0.5, "solid", "black"), backgroundColor="white") +theme(legend.title=element_blank()) # Income of LED light fishing method after 14 years vs diff crab price > procrabprice=read.csv("C:/KHANH/PhD/snow crab/2017 with Bob/field data/income after 14 years of diff crab price.csv", header=T) > attach(procrabprice) > ggplot(procrabprice, aes(x=Crab_price, y=Profit_light)) + geom_line()+ theme_bw()+theme(axis.text.x = element_text(face="plain", color="black",size=12),axis.text.y = element_text(face="plain", color="black",size=12)) + theme( axis.line = element_line(colour = "black", size = 0.5, linetype = "solid"))+xlab("Crab price (CND/kg)")+ ylab("Income (CND 1000)")+ theme(axis.text=element_text(size=12), axis.title=element_text(size=12,face="plain"))+ theme(legend.title = element_text(colour="black", size=12, face="plain")) # Income of traditional fishery method after 14 years vs diff crab price > procrabprice=read.csv("C:/KHANH/PhD/snow crab/2017 with Bob/field data/income after 14 years of diff crab price.csv", header=T) > attach(procrabprice) > ggplot(procrabprice, aes(x=Crab_price, y=Income, group=Fishery)) + geom_line(aes(linetype=Fishery))+ scale_linetype_manual(values=c("solid","dashed"))+ theme_bw()+theme(axis.text.x = element_text(face="plain", color="black",size=12),axis.text.y = element_text(face="plain", color="black",size=12)) + theme( axis.line = element_line(colour = "black", size = 0.5, linetype = "solid"))+xlab("Crab price (CND/kg)")+ ylab("Income (CND 1000)")+ theme(axis.text=element_text(size=12), axis.title=element_text(size=12,face="plain"))+ theme(legend.title = element_text(colour="black", size=12, face="plain"))+theme(legend.title=element_blank()) # D at different quota allocation > difquota=read.csv("C:/KHANH/PhD/snow crab/2017 with Bob/field data/Different quota.csv", header=T) > attach(difquota) > ggplot2.scatterplot(data=difquota, xName="Quota", yName="Diff", size=0, color="white", addRegLine=TRUE, regLineColor="black", regLineSize=1,xtitle="Quota allocation (Ton/year)", xtitleFont=c(12,"plain", "black"),xTickLabelFont=c(12, "plain", "black"),ytitle="D (CND 1000)",ytitleFont=c(12,"plain", "black"),yTickLabelFont=c(12, "plain", "black"), axisLine=c(0.5, "solid", "black"), backgroundColor="white")