## This shows how the get_workflow_design() may be used,
## but in most cases one does not need to manipulate it.
if (FALSE) {
blurred <- function(mu,s2,sample.size) {
s <- rnorm(n=sample.size,mean=mu,sd=sqrt(s2))
s <- exp(s/4)
return(c(mean=mean(s),var=var(s)))
}
## Simulated data:
set.seed(123)
dSobs <- blurred(mu=4,s2=1,sample.size=40)
workflow_design <- get_workflow_design(npar=2L)
## Construct initial reference table:
# Sample its parameters:
if (IMPLICIT <- TRUE) { # use implicit control by get_workflow_design()
parsp_j <- init_reftable(lower=c(mu=2.5,s2=0.25,sample.size=40),
upper=c(mu=5.2,s2=2.4,sample.size=40))
# => get_workflow_design() has been called internally with default values
# to provide the dimension of the initial reference table.
# The following syntax provides a more explicit control:
} else {
parsp_j <- init_reftable(lower=c(mu=2.5,s2=0.25,sample.size=40),
upper=c(mu=5.2,s2=2.4,sample.size=4,
nUnique=workflow_design$init_reft_size))
}
# and yet another way to the same result could be
#
# parsp_j <- data.frame(mu=runif(init_reft_size,min=2.5,max=5.2),
# s2=runif(init_reft_size,min=0.25,max=2.4),
# sample.size=40)
# Generate the initial simulations:
dsimuls <- add_reftable(, Simulate="blurred", parsTable=parsp_j, verbose=FALSE)
## Construct projections
mufit <- project("mu",stats=c("mean","var"),data=dsimuls,verbose=FALSE)
s2fit <- project("s2",stats=c("mean","var"),data=dsimuls,verbose=FALSE)
dprojectors <- list(MEAN=mufit,VAR=s2fit)
## Apply projections on simulated statistics and 'data':
dprojSimuls <- project(dsimuls,projectors=dprojectors,verbose=FALSE)
dprojSobs <- project(dSobs,projectors=dprojectors)
## Summary-likelihood inference:
# Initial Inference of log-likelihood surface
slik_j <- infer_SLik_joint(dprojSimuls, stat.obs=dprojSobs, verbose=TRUE)
# Find maximum, confidence intervals...
slik_j <- MSL(slik_j, eval_RMSEs=FALSE, CIs=FALSE)
## Refinements over iterations
# Here, with only two estimated parameters, workflow_design$final_reft_size
# suggests a final reference table of 5000 simulations, attained through
# 6 refine() calls with intermediate sizes given by
(workflow_design$reftable_sizes)
# here 500, 1000, 2000, 3000, 4000, 5000.
#
if (IMPLICIT) { # Again using implicit control by get_workflow_design()
# Essentially, it suffices to call
for (it in seq(6)) slik_j <-
refine(slik_j, eval_RMSEs= it==6L, CIs= it==6L)
# to run the default workflow. Again, the following syntax,
# showing how successive table sizes are controlled internally,
# provides a more explicit control:
} else {
reftable_sizes <- workflow_design$reftable_sizes
init_reft_size <- workflow_design$init_reft_size
refine_sizes <- diff(c(init_reft_size, reftable_sizes))
maxit <- workflow_design$subblock_nbr
for(it in seq_len(length(refine_sizes)-1L)) {
add_it <- refine_sizes[it]
slik_j <- refine(slik_j, ntot=add_it, maxit=maxit,
eval_RMSEs=FALSE, CIs=FALSE)
}
add_it <- tail(refine_sizes,1L)
slik_j <- refine(slik_j, ntot=add_it, maxit=maxit,
CIs=add_it, eval_RMSEs=add_it)
}
}
Run the code above in your browser using DataLab