PlantWise adaptation logistic model

Logistic Clean

#lOGIT

setwd("c:/Users/machariam/Desktop/Logit")
list.files()
Dat=read.csv("Dat.csv")
names(Dat)
str(Dat)
hist(LandPercapita)
hist(MaizeLand)
hist(PPI)
hist(Labour)
hist(InputRecieved)
hist(Yield)
hist(InputCosts)
hist(TotalShocks)
hist(OwnExperience)
hist(RadioTv)
hist(Extension)
hist(Plantclinics)
hist(Dat$CreditSource)
hist(Dat$RespAge)

summary(LandPercapita)
summary(MaizeLand)
summary(Yield) 
summary(InputCosts)

PYield <- quantile(Yield, probs=c(0.10, 0.9))
PYield
PInputCost<-quantile(InputCosts, prob=c(0.10,0.9))
PInputCost

class(Dat)

Dat.sub=subset(Dat,(Yield<= 4000) & (Yield >= 330) & (InputCosts<=207) & (InputCosts>=7))

View(Dat.sub)


write.csv(Dat.sub, file = ("C:/users/machariam/Desktop/PlotData.csv"))

length(Dat.sub$HHID)

hist(Dat.sub$Yield)
hist(Dat.sub$InputCosts)
hist(Dat.sub$Labour)
hist(Dat.sub$EducCode)
hist(Dat.sub$GenderCode)
hist(Dat.sub$CreditSource)
str(Dat.sub)

cor(Dat.sub[sapply(Dat.sub, is.numeric)])

cor(subset(Dat.sub, select = c(LandPercapita,MaizeLand,Yield, InputCosts,Labour)))

cor(Dat.sub$LandPercapita,Dat.sub$Labour)

fCulturalPA=as.factor(Dat.sub$CulturalPA)

fGender=as.factor(Dat.sub$GenderCode)

class(fCulturalPA)
class(fGender)
class(EducCode)

#Cultural prevention/ avoidance

reg0 = glm(fCulturalPA~LandPercapita+MaizeLand+Yield+HHMembers+CreditSource+InputCosts+fGender+RespAge+EducCode+PPI+Labour+InputRecieved+TotalShocks+OwnExperience+RadioTv+Extension+Plantclinics,family=binomial,data=Dat.sub)
summary(reg0)

anova(reg0, test="Chisq")

exp(coef(reg0))

exp(confint(reg1))

#Cultural planting / cropping patterns

names(Dat.sub)

fCulturalPCP=as.factor(Dat.sub$CulturalPCP)
class(fCulturalPCP)

reg2 = glm(fCulturalPCP~LandPercapita+MaizeLand+Yield+HHMembers+RespAge+CreditSource+InputCosts+Gender+EducCode+PPI+Labour+InputRecieved+TotalShocks+OwnExperience+RadioTv+Extension+Plantclinics,family=binomial,data=Dat.sub)
summary(reg2)

exp(coef(reg2))

predict(reg2)
x=predict(reg2, type="response")# Predicted probabilities
summary(x)


#Cultural resistance / Tolerant varieties

names(Dat.sub)
fCulturalRTV=as.factor(Dat.sub$CulturalRTV)

reg4 = glm(fCulturalRTV~LandPercapita+MaizeLand+HHMembers+RespAge+Yield+CreditSource+InputCosts+Gender+EducCode+PPI+Labour+InputRecieved+TotalShocks+OwnExperience+RadioTv+Extension+Plantclinics,family=binomial,data=Dat.sub)
summary(reg4)

exp(coef(reg4))

#Mechanical late stage intervention

fMechanicalLS=as.factor(Dat.sub$MechanicalLS)

names(Dat.sub)

reg6 = glm(fMechanicalLS~LandPercapita+MaizeLand+HHMembers+RespAge+Yield+CreditSource+InputCosts+Gender+EducCode+PPI+Labour+InputRecieved+TotalShocks+OwnExperience+RadioTv+Extension+Plantclinics,family=binomial,data=Dat.sub)
summary(reg6)

exp(coef(reg6))


#Mechanical early stage intervention

names(Dat.sub)

fMechanicalES=as.factor(Dat.sub$MechanicalES)

reg8 = glm(fMechanicalES~LandPercapita+MaizeLand+HHMembers+RespAge+Yield+CreditSource+InputCosts+Gender+EducCode+PPI+Labour+InputRecieved+TotalShocks+OwnExperience+RadioTv+Extension+Plantclinics,family=binomial,data=Dat.sub)
summary(reg8)

exp(coef(reg8))

#Biological -NA

#Pesticides

names(Dat.sub)

fPesticides=as.factor(Dat.sub$Pesticides)

reg10 = glm(fPesticides~LandPercapita+MaizeLand+HHMembers+RespAge+Yield+CreditSource+InputCosts+Gender+EducCode+PPI+Labour+InputRecieved+TotalShocks+OwnExperience+RadioTv+Extension+Plantclinics,family=binomial,data=Dat.sub)

summary(reg10)

exp(coef(reg10))


#Monitoring

names(Dat.sub)

fMonitoring=as.factor(Dat.sub$Monitoring)

reg12 = glm(fMonitoring~LandPercapita+MaizeLand+Yield+HHMembers+RespAge+CreditSource+InputCosts+Gender+EducCode+PPI+Labour+InputRecieved+TotalShocks+OwnExperience+RadioTv+Extension+Plantclinics,family=binomial,data=Dat.sub)

summary(reg12)

exp(coef(reg12))



####Option 2

#Application of pesticides

names(Dat.sub)

fPesticides2=as.factor(Dat.sub$Pesticides2)

reg13 = glm(fPesticides2~LandPercapita+MaizeLand+Yield+InputCosts+Gender+EducCode+PPI+Labour+InputRecieved+TotalShocks+OwnExperience+RadioTv+Extension+Plantclinics,family=binomial,data=Dat.sub)

summary(reg13)


reg14 = glm(fPesticides2~MaizeLand+Yield+InputCosts+Gender+EducCode+PPI+InputRecieved+OwnExperience+RadioTv+Extension+Plantclinics,family=binomial,data=Dat.sub)

summary(reg14)

#Remove crop residues

names(Dat.sub)

fCropResidues2=as.factor(Dat.sub$CropResidues2)

reg15 = glm(fCropResidues2~LandPercapita+MaizeLand+Yield+InputCosts+Gender+EducCode+PPI+Labour+InputRecieved+TotalShocks+OwnExperience+RadioTv+Extension+Plantclinics,family=binomial,data=Dat.sub)

summary(reg15)


reg16 = glm(fCropResidues2~MaizeLand+Gender+EducCode+InputRecieved+OwnExperience+Extension+Plantclinics,family=binomial,data=Dat.sub)

summary(reg16)

#Uproot and burn

names(Dat.sub)

fUproot2=as.factor(Dat.sub$Uproot2)

reg17 = glm(fCropResidues2~LandPercapita+MaizeLand+Yield+InputCosts+Gender+EducCode+PPI+Labour+InputRecieved+TotalShocks+OwnExperience+RadioTv+Extension+Plantclinics,family=binomial,data=Dat.sub)

summary(reg17)

reg18 = glm(fCropResidues2~LandPercapita+MaizeLand+Gender+EducCode+InputRecieved+TotalShocks+OwnExperience+RadioTv+Extension+Plantclinics,family=binomial,data=Dat.sub)

summary(reg18)

#Hand picking

names(Dat.sub)

fHandPick2=as.factor(Dat.sub$HandPick2)

reg19 = glm(fHandPick2~LandPercapita+MaizeLand+Yield+InputCosts+Gender+EducCode+PPI+Labour+InputRecieved+TotalShocks+OwnExperience+RadioTv+Extension+Plantclinics,family=binomial,data=Dat.sub)

summary(reg19)

reg20 = glm(fHandPick2~Yield+InputCosts+PPI+TotalShocks+OwnExperience+RadioTv+Extension+Plantclinics,family=binomial,data=Dat.sub)

summary(reg20)

#Monitoring

names(Dat.sub)

fMonitor2=as.factor(Dat.sub$Monitor2)

reg21 = glm(fMonitor2~LandPercapita+MaizeLand+Yield+InputCosts+Gender+EducCode+PPI+Labour+InputRecieved+TotalShocks+OwnExperience+RadioTv+Extension+Plantclinics,family=binomial,data=Dat.sub)

summary(reg21)

reg22 = glm(fMonitor2~MaizeLand+Yield+InputCosts+Gender+EducCode+InputRecieved+OwnExperience+RadioTv+Extension+Plantclinics,family=binomial,data=Dat.sub)

summary(reg22)


#Regression

names(Dat.sub)
str(Dat.sub)

summary(Dat.sub$Yield)


lm0 = lm(Yield~LandPercapita+MaizeLand+InputCosts+Gender+HHMembers+RespAge+CreditSource+EducCode+PPI+Labour+InputRecieved+TotalShocks+OwnExperience+RadioTv+Extension+Plantclinics, data=Dat.sub)

summary(lm0)

install.packages("car")

library(car)

outlierTest(lm0)

qqPlot(lm0, main="QQ Plot")

#Plot the siignificant variables

plot(Dat.sub$InputCosts, Dat.sub$Yield)

plot(Dat.sub$InputCosts~Dat.sub$Yield,ylab="Yield (kg/ha)", xlab="Input costs (USD)",xlim=c(0,200), ylim=c(300,4000),las=1,cex.lab=1,cex.axis =0.9)

A=1325.213
B=8.077

lines(x<-c(0:200), A+B*(x),col='black ',lwd=2, lty=1)

legend(50, 3600,legend=c("Y=1325.21+8.077x"), bty="n", cex=1.2)