東方水蠊感染強度分析
===
## 感染強度
---
###

rm(list=ls())
setwd("C:/Users/max96/Desktop/實驗")

dat <- read.csv("food.csv", head = T)

dat <- na.omit(dat) #刪除樣本缺失值(NA)#

dat$totalinfection <- apply(dat[,4:6], 1, sum)#算出蟯蟲總感染隻數#
dat$rate <- as.numeric(dat$totalinfection>0)#算出感染率#
dat$infectionstrength <- dat$totalinfection/dat$rate#算出每隻的感染強度#
dat$F.infstr <- dat$F
dat$M.infstr <- dat$M
dat$L.infstr <- dat$L


nymphfiber <- subset(dat , stage=="nymph" & food=="fiber")#建立獨立資料#
nymphprotein <- subset(dat , stage=="nymph" & food=="protein")
malefiber <- subset(dat , stage=="male" & food=="fiber")
maleprotein <- subset(dat , stage=="male" & food=="protein")
femalefiber <- subset(dat , stage=="female" & food=="fiber")
femaleprotein <- subset(dat , stage=="female" & food=="protein")


nf.ratio=sum(nymphfiber$totalinfection)/sum(nymphfiber$rate)#感染強度#
np.ratio=sum(nymphprotein$totalinfection)/sum(nymphprotein$rate)
mf.ratio=sum(malefiber$totalinfection)/sum(malefiber$rate)
mp.ratio=sum(maleprotein$totalinfection)/sum(maleprotein$rate)
ff.ratio=sum(femalefiber$totalinfection)/sum(femalefiber$rate)
fp.ratio=sum(femaleprotein$totalinfection)/sum(femaleprotein$rate)



F.model <- lm (F.infstr ~ stage * food, 
                    data = dat)
anova(F.model)

M.model <- lm (M.infstr ~ stage * food, 
               data = dat)
anova(F.model)

L.model <- lm (L.infstr ~ stage * food,
               data = dat)
anova(L.model)

infstr.model <- lm (infectionstrength ~ stage * food, 
                    data = dat)
anova(infstr.model)

t.test(nymphfiber$infectionstrength , malefiber$infectionstrength, paired=F)
t.test(nymphfiber$infectionstrength , femalefiber$infectionstrength, paired=F)
t.test(malefiber$infectionstrength , femalefiber$infectionstrength, paired=F)