embryogrowth (version 6.4)

searchR: Fit the parameters that best represent nest incubation data.

Description

Fit the parameters that best represent data. test can be a data.frame with two columns Mean and SD and rownames with the nest name. Function to fit thermal reaction norm can be also expressed as a Weibull function with k (shape), lambda (scale) and theta parameters.

Usage

searchR(parameters = stop("Initial set of parameters must be provided"), fixed.parameters = NULL, temperatures = stop("Formated temperature must be provided !"), derivate = dydt.Gompertz, test = c(Mean = 39.33, SD = 1.92), M0 = 1.7, saveAtMaxiter = FALSE, fileName = "intermediate", weight = NULL, SE = FALSE, control = list(trace = 1, REPORT = 100, maxit = 500))

Arguments

parameters

A set of parameters used as initial point for searching

fixed.parameters

A set of parameters that will not be changed

temperatures

Timeseries of temperatures after formated using FormatNests()

derivate

Function used to fit embryo growth: dydt.Gompertz, dydt.exponential or dydt.linear

test

A vector with Mean and SD of size of hatchlings, ex. test=c(Mean=39, SD=3). Can be a data.frame also. See description

M0

Measure of hatchling size or mass proxi at laying date

saveAtMaxiter

If True, each time number of interation reach maxiter, current data are saved in file with filename name

fileName

The intermediate results are saved in file with fileName.Rdata name

weight

A named vector of the weight for each nest for likelihood estimation

SE

If TRUE, the SE of parameters are estimated.

control

List for control parameters for optimx

Value

A result object

Details

searchR fits the parameters that best represent nest incubation data.

Examples

## Not run: 
# library(embryogrowth)
# data(nest)
# formated <- FormatNests(nest)
# # The initial parameters value can be:
# # "T12H", "DHA",  "DHH", "Rho25"
# # Or
# # "T12L", "DT", "DHA",  "DHH", "DHL", "Rho25"
# # K for Gompertz must be set as fixed parameter or being a constant K  
# # or relative to the hatchling size rK
# x <- structure(c(106.59891311201, 614.181133951497, 306.267053513175, 
# 120.327257089974), .Names = c("DHA", "DHH", "T12H", "Rho25"))
# # pfixed <- c(K=82.33) or rK=82.33/39.33
# pfixed <- c(rK=2.093313)
# resultNest_4p <- searchR(parameters=x, fixed.parameters=pfixed, 
# 	temperatures=formated, derivate=dydt.Gompertz, M0=1.7, 
# 	test=c(Mean=39.33, SD=1.92))
# data(resultNest_4p)
# plot(resultNest_4p, xlim=c(0,70), ylimT=c(22, 32), ylimS=c(0,45), series=1)
# x <- structure(c(106.567809092008, 527.359011254683, 614.208632495199, 
# 2720.94506457237, 306.268259715624, 120.336791245212), .Names = c("DHA", 
# "DHH", "DHL", "DT", "T12L", "Rho25"))
# pfixed <- c(rK=2.093313)
# 
# # exemple of data.frame for test
# ttest <- data.frame(Mean=rep(25.5, formated$IndiceT["NbTS"]), 
#                      SD=rep(0.75, formated$IndiceT["NbTS"]), 
#                      row.names=names(formated)[1:formated$IndiceT["NbTS"]])
# resultNest_6p <- searchR(parameters=x, fixed.parameters=pfixed, 
# 	                        temperatures=formated, derivate=dydt.Gompertz, M0=1.7, 
# 	                        test=ttest)
# 	                        
# data(resultNest_6p)
# pMCMC <- TRN_MHmcmc_p(resultNest_6p, accept=TRUE)
# # Take care, it can be very long, sometimes several days
# result_mcmc_6p <- GRTRN_MHmcmc(result=resultNest_6p,  
# 	parametersMCMC=pMCMC, n.iter=10000, n.chains = 1, n.adapt = 0,  
# 	thin=1, trace=TRUE)
# data(result_mcmc_6p)
# # compare_AIC() is a function from the package "HelpersMG"
# compare_AIC(test1=resultNest_4p, test2=resultNest_6p)
# ############ with new parametrization
# data(resultNest_4p)
# x0 <- resultNest_4p$par
# t <- hist(resultNest_4p, plot=FALSE)
# x <- c(3.4, 3.6, 5.4, 5.6, 7.6, 7.5, 3.2)
# names(x) <- seq(from=range(t$temperatures)[1], to=range(t$temperatures)[2], 
#      length.out=7)
# newx <- ChangeSSM(temperatures = (200:350)/10, parameters = x0, 
#        initial.parameters = x, 
#        control=list(maxit=5000))
# pfixed <- c(rK=2.093313)
# resultNest_newp <- searchR(parameters=newx, fixed.parameters=pfixed,
#  temperatures=formated, derivate=dydt.Gompertz, M0=1.7,
#  test=c(Mean=39.33, SD=1.92))
# plotR_hist(resultNest_newp, ylim=c(0,0.3), xlimR=c(23, 34), ylimH=c(0, 0.3))
# compare_AIC(test4p=resultNest_4p, 
#             test6p=resultNest_6p, 
#             testAnchor=resultNest_newp)
# ############################################
# # example with thermal reaction norm fitted from Weibull function
# ############################################
# x <- structure(c(298890.11996796, 1229465.06811278, -1229160.54956529, 
# 27.2843254770591), .Names = c("k", "lambda", "theta", "scale"))
# pfixed <- c(rK=2.093313)
# x <- ChangeSSM(temperatures = (200:350)/10,
#                parameters = resultNest_4p$par,
#                initial.parameters = x, 
#                control=list(maxit=1000))
# resultNest_4p_Weibull <- searchR(parameters=x, fixed.parameters=pfixed, 
#                          temperatures=formated, derivate=dydt.Gompertz, M0=1.7, 
#                          test=c(Mean=39.33, SD=1.92))
# plotR(list(resultNest_4p, resultNest_4p_Weibull), ylim=c(0,0.3), col=c("Black", "red"))
# compare_AIC(SSM=resultNest_4p, Weibull=resultNest_4p_Weibull)
# ## End(Not run)