# NOT RUN {
library(embryogrowth)
totalIncubation_Cc <- subset(DatabaseTSD,
Species=="Caretta caretta" &
Note != "Sinusoidal pattern" &
!is.na(Total) & Total != 0)
par <- c(S.low=0.5, S.high=0.3,
P.low=25, deltaP=10, MaxHS=0.8)
g.logistic <- HatchingSuccess.fit(par=par, data=totalIncubation_Cc)
HatchingSuccess.lnL(par=g.logistic$par, data=totalIncubation_Cc)
plot(g.logistic)
par <- c(S.low=0.5, S.high=0.3,
P.low=25, deltaP=10,
K1.low=1, K2.low=-1, K1.high=1, K2.high=-1,
MaxHS=0.8)
g.flexit <- HatchingSuccess.fit(par=par, data=totalIncubation_Cc)
HatchingSuccess.lnL(par=g.flexit$par, data=totalIncubation_Cc)
compare_AICc(logistic=g.logistic, flexit=g.flexit)
plot(x=g.logistic, what = c("observations", "ML", "CI"), replicates=10000)
pMCMC <- HatchingSuccess.MHmcmc_p(result = g.logistic, accept = TRUE)
MCMC <- HatchingSuccess.MHmcmc(result = g.logistic, parametersMCMC = pMCMC,
n.iter = 100000,
adaptive = TRUE)
plot(MCMC, parameters = "S.low")
plot(MCMC, parameters = "S.high")
plot(MCMC, parameters = "P.low")
plot(MCMC, parameters = "deltaP")
plot(MCMC, parameters = "MaxHS")
plot(x=g.logistic, what = c("observations", "ML", "CI"),
replicates=10000, resultmcmc=MCMC)
# }
Run the code above in your browser using DataLab