CI.RMU: Calculate the confidence interval of the results of fitRMU()

In phenology: Tools to Manage a Parametric Function that Describes Phenology and More

Calculate the confidence interval of the results of fitRMU()

Description

The data must be a data.frame with the first column being years
and two columns for each beach: the average and the se for the estimate.
The correspondence between mean, se and density for each rookery are given in the RMU.names data.frame.
This data.frame must have a column named mean, another named se and a third named density. If no sd column exists, no sd will be considered for the series and is no density column exists, it will be considered as being "dnorm".
In the result list, the mean proportions for each rookeries are in $proportions.
The names of beach columns must not begin by T_, SD_, a0_, a1_ or a2_ and cannot be r.
A RMU is the acronyme for Regional Managment Unit. See:
Wallace, B.P., DiMatteo, A.D., Hurley, B.J., Finkbeiner, E.M., Bolten, A.B., Chaloupka, M.Y., Hutchinson, B.J., Abreu-Grobois, F.A., Amorocho, D., Bjorndal, K.A., Bourjea, J., Bowen, B.W., Dueñas, R.B., Casale, P., Choudhury, B.C., Costa, A., Dutton, P.H., Fallabrino, A., Girard, A., Girondot, M., Godfrey, M.H., Hamann, M., López-Mendilaharsu, M., Marcovaldi, M.A., Mortimer, J.A., Musick, J.A., Nel, R., Seminoff, J.A., Troëng, S., Witherington, B., Mast, R.B., 2010. Regional management units for marine turtles: a novel framework for prioritizing conservation and research across multiple scales. PLoS One 5, e15465.
Variance for each value is additive based on both the observed SE (in the RMU.data object) and a constant value dependent on the rookery when model.SD is equal to "Rookery-constant". The value is a global constant when model.SD is "global-constant". The value is proportional to the observed number of nests when model.SD is "global-proportional" with aSD_*observed+SD_ with aSD_ and SD_ being fitted values. This value is fixed to zero when model.SD is "Zero".
If replicate.CI is 0, no CI is estimated, and only point estimation is returned.

Usage

CI.RMU(
  result = stop("A result obtained from fitRMU is necessary"),
  resultMCMC = NULL,
  chain = 1,
  replicate.CI = 10000,
  regularThin = TRUE,
  probs = c(0.025, 0.5, 0.975),
  silent = FALSE
)

Arguments

result	A result of fitRMu()
resultMCMC	A resuts of fitRMU_MHmcmc()
chain	Number of MCMC chain to be used
replicate.CI	Number of replicates
regularThin	If TRUE, use regular thin for MCMC
probs	The probabilities to return for quantiles
silent	If TRUE does not display anything

Details

CI.RMU calculates the confidence interval of the results of fitRMU()

Value

Return a list with Total, Proportions, and Numbers

Author(s)

Marc Girondot

See Also

Other Fill gaps in RMU: fitRMU_MHmcmc_p(), fitRMU_MHmcmc(), fitRMU(), logLik.fitRMU(), plot.fitRMU()

Examples

## Not run: 
library("phenology")
RMU.names.AtlanticW <- data.frame(mean=c("Yalimapo.French.Guiana", 
                                         "Galibi.Suriname", 
                                         "Irakumpapy.French.Guiana"), 
                                 se=c("se_Yalimapo.French.Guiana", 
                                      "se_Galibi.Suriname", 
                                      "se_Irakumpapy.French.Guiana"), 
                                 density=c("density_Yalimapo.French.Guiana", 
                                           "density_Galibi.Suriname", 
                                           "density_Irakumpapy.French.Guiana"), 
                                           stringsAsFactors = FALSE)
data.AtlanticW <- data.frame(Year=c(1990:2000), 
      Yalimapo.French.Guiana=c(2076, 2765, 2890, 2678, NA, 
                               6542, 5678, 1243, NA, 1566, 1566),
      se_Yalimapo.French.Guiana=c(123.2, 27.7, 62.5, 126, NA, 
                                 230, 129, 167, NA, 145, 20),
      density_Yalimapo.French.Guiana=rep("dnorm", 11), 
      Galibi.Suriname=c(276, 275, 290, NA, 267, 
                       542, 678, NA, 243, 156, 123),
      se_Galibi.Suriname=c(22.3, 34.2, 23.2, NA, 23.2, 
                           4.3, 2.3, NA, 10.3, 10.1, 8.9),
      density_Galibi.Suriname=rep("dnorm", 11), 
      Irakumpapy.French.Guiana=c(1076, 1765, 1390, 1678, NA, 
                               3542, 2678, 243, NA, 566, 566),
      se_Irakumpapy.French.Guiana=c(23.2, 29.7, 22.5, 226, NA, 
                                 130, 29, 67, NA, 15, 20), 
      density_Irakumpapy.French.Guiana=rep("dnorm", 11), stringsAsFactors = FALSE
      )
                           
cst <- fitRMU(RMU.data=data.AtlanticW, RMU.names=RMU.names.AtlanticW, 
               colname.year="Year", model.trend="Constant", 
               model.SD="Zero")
               
out.CI.Cst <- CI.RMU(result=cst)



cst <- fitRMU(RMU.data=data.AtlanticW, RMU.names=RMU.names.AtlanticW, 
               colname.year="Year", model.trend="Constant", 
               model.SD="Zero", 
               control=list(trace=1, REPORT=100, maxit=500, parscale = c(3000, -0.2, 0.6)))
               
cst <- fitRMU(RMU.data=data.AtlanticW, RMU.names=RMU.names.AtlanticW, 
               colname.year="Year", model.trend="Constant", 
               model.SD="Zero", method=c("Nelder-Mead","BFGS"), 
               control = list(trace = 0, REPORT = 100, maxit = 500, 
               parscale = c(3000, -0.2, 0.6)))
expo <- fitRMU(RMU.data=data.AtlanticW, RMU.names=RMU.names.AtlanticW, 
               colname.year="Year", model.trend="Exponential", 
               model.SD="Zero", method=c("Nelder-Mead","BFGS"), 
               control = list(trace = 0, REPORT = 100, maxit = 500, 
               parscale = c(6000, -0.05, -0.25, 0.6)))
YS <- fitRMU(RMU.data=data.AtlanticW, RMU.names=RMU.names.AtlanticW, 
             colname.year="Year", model.trend="Year-specific", method=c("Nelder-Mead","BFGS"), 
             model.SD="Zero")
YS1 <- fitRMU(RMU.data=data.AtlanticW, RMU.names=RMU.names.AtlanticW, 
             colname.year="Year", model.trend="Year-specific", method=c("Nelder-Mead","BFGS"), 
             model.SD="Zero", model.rookeries="First-order")
YS1_cst <- fitRMU(RMU.data=data.AtlanticW, RMU.names=RMU.names.AtlanticW, 
             colname.year="Year", model.trend="Year-specific", 
             model.SD="Constant", model.rookeries="First-order", 
             parameters=YS1$par, method=c("Nelder-Mead","BFGS"))
YS2 <- fitRMU(RMU.data=data.AtlanticW, RMU.names=RMU.names.AtlanticW, 
             colname.year="Year", model.trend="Year-specific",
             model.SD="Zero", model.rookeries="Second-order", 
             parameters=YS1$par, method=c("Nelder-Mead","BFGS"))
YS2_cst <- fitRMU(RMU.data=data.AtlanticW, RMU.names=RMU.names.AtlanticW, 
             colname.year="Year", model.trend="Year-specific",
             model.SD="Constant", model.rookeries="Second-order", 
             parameters=YS1_cst$par, method=c("Nelder-Mead","BFGS"))
               
compare_AIC(Constant=cst, 
            Exponential=expo, 
            YearSpecific=YS)

compare_AIC(YearSpecific_ProportionsFirstOrder_Zero=YS1,
YearSpecific_ProportionsFirstOrder_Constant=YS1_cst)

compare_AIC(YearSpecific_ProportionsConstant=YS,
           YearSpecific_ProportionsFirstOrder=YS1,
           YearSpecific_ProportionsSecondOrder=YS2)
           
compare_AIC(YearSpecific_ProportionsFirstOrder=YS1_cst,
           YearSpecific_ProportionsSecondOrder=YS2_cst)

plot(cst, main="Use of different beaches along the time", what="total")
plot(expo, main="Use of different beaches along the time", what="total")
plot(YS2_cst, main="Use of different beaches along the time", what="total")

plot(YS1, main="Use of different beaches along the time")
plot(YS1_cst, main="Use of different beaches along the time")
plot(YS1_cst, main="Use of different beaches along the time", what="numbers")

# Gamma distribution should be used for MCMC outputs

RMU.names.AtlanticW <- data.frame(mean=c("Yalimapo.French.Guiana", 
                                         "Galibi.Suriname", 
                                         "Irakumpapy.French.Guiana"), 
                                 se=c("se_Yalimapo.French.Guiana", 
                                      "se_Galibi.Suriname", 
                                      "se_Irakumpapy.French.Guiana"), 
                                 density=c("density_Yalimapo.French.Guiana", 
                                           "density_Galibi.Suriname", 
                                           "density_Irakumpapy.French.Guiana"))
                                           
data.AtlanticW <- data.frame(Year=c(1990:2000), 
      Yalimapo.French.Guiana=c(2076, 2765, 2890, 2678, NA, 
                               6542, 5678, 1243, NA, 1566, 1566),
      se_Yalimapo.French.Guiana=c(123.2, 27.7, 62.5, 126, NA, 
                                 230, 129, 167, NA, 145, 20),
      density_Yalimapo.French.Guiana=rep("dgamma", 11), 
      Galibi.Suriname=c(276, 275, 290, NA, 267, 
                       542, 678, NA, 243, 156, 123),
      se_Galibi.Suriname=c(22.3, 34.2, 23.2, NA, 23.2, 
                           4.3, 2.3, NA, 10.3, 10.1, 8.9),
      density_Galibi.Suriname=rep("dgamma", 11), 
      Irakumpapy.French.Guiana=c(1076, 1765, 1390, 1678, NA, 
                               3542, 2678, 243, NA, 566, 566),
      se_Irakumpapy.French.Guiana=c(23.2, 29.7, 22.5, 226, NA, 
                                 130, 29, 67, NA, 15, 20), 
      density_Irakumpapy.French.Guiana=rep("dgamma", 11)
      )
cst <- fitRMU(RMU.data=data.AtlanticW, RMU.names=RMU.names.AtlanticW, 
               colname.year="Year", model.trend="Constant", 
               model.SD="Zero")

## End(Not run)

phenology documentation built on Oct. 16, 2023, 9:06 a.m.