Nothing
R/MonteCarloFunction.R
In TempStable: A Collection of Methods to Estimate Parameters of Different Tempered Stable Distributions
Defines functions
getTime_
updateOutputFile
get_filename
writeCheckPoint
deleteCheckPoint
updateCheckPointValues
get_filename_checkPoint_Temp
readCheckPoint
Estim_Des_Temp
initOutputFile
getTempEstimation
ComputeMCSimForTempered
getSeedVector
parallelizeMCsimulation
TemperedEstim_Simulation
Documented in
parallelizeMCsimulation
TemperedEstim_Simulation
#' Monte Carlo Simulation
#'
#' @description
#' Runs Monte Carlo simulation for a selected estimation method. The function
#' can save results in a file.
#'
#' @details
#' \strong{TemperedTyp} With the parameter 'TemperedTyp' you can choose the
#' tempered stable distribution you want to use. Here is a list of distribution
#' you can choose from:
#' \describe{
#' \item{TSS}{Tempered stabel subordinator: See [charTSS()] for details.}
#' \item{CTS}{Classical tempered stable distribution: See [charCTS()] for
#' details.}
#' \item{GTS}{Generalized classical tempered stable distribution: See
#' [charGTS()] for details.}
#' \item{NTS}{Normal tempered stable distribution: See [charNTS()] for
#' details.}
#' \item{MTS}{Modified tempered stable distribution: See [charMTS()] for
#' details.}
#' \item{RDTS}{Rapid decreasing tempered stable distribution: See [charRDTS()]
#' for details.}
#' \item{KRTS}{Kim-Rachev tempered stable distribution: See [charKRTS()] for
#' details.}
#' }
#'
#' \strong{Error Handling} It is advisable to set it to TRUE when user is
#' planning to launch long simulations as it will prevent the procedure to stop
#' if an error occurs for one sample data. The estimation function will produce
#' a vector of NA as estimated parameters related to this (error generating)
#' sample data and move on to the next Monte Carlo step.
#'
#' \strong{Output file} Setting \code{saveOutput} to \code{TRUE} will have the
#' side effect of saving a csv file in the working directory. This file will
#' have \code{MCparam*length(SampleSizes)} lines and its columns will be:
#' \describe{
#' \item{alphaT, ...:}{the true value of the parameters.}
#' \item{data size:}{the sample size used to generate the simulated data.}
#' \item{seed:}{the seed value used to generate the simulated data.}
#' \item{alphaE, ...:}{the estimate of the parameters.}
#' \item{failure:}{binary: 0 for success, 1 for failure.}
#' \item{time:}{estimation running time in seconds.}
#' }
#' The file name is informative to let the user identify the value of the true
#' parameters, the MC parameters as well as the options selected for the
#' estimation method. The csv file is updated after each MC estimation which is
#' useful when the simulation stops before it finishes.
#'
#' \strong{SeedOptions} If users does not want to control the seed generation,
#' they could ignore this argument (default value NULL). This argument can be
#' more useful when they wants to cut the simulation (even for one parameter
#' value) into pieces. In that case, they can control which part of the seed
#' vector they want to use.
#' \describe{
#' \item{MCtot:}{total values of MC simulations in the entire process.}
#' \item{seedStart:}{starting index in the seed vector. The vector extracted
#' will be of size MCparam.}
#' }
#'
#' \strong{Estimfct} Additional parameters are needed for different estimation
#' functions. These are listed below for each function. The list of additional
#' parameters starts after the parameter \code{eps} in the parameter list.
#' \describe{
#' \item{For ML:}{ See usage of Maximum likelihood estimation in Kim et al.
#' (2008).No additional parameters are needed.}
#' \item{For GMM:}{Generalized Method of Moments by Feuerverger (1981).
#' The parameters \code{algo, alphaReg, regularization, WeightingMatrix, and
#' t_scheme} must be specified.
#'
#' Parameter \code{t_scheme}: One of the most important features of this
#' method is that it allows the user to choose how to place the points where
#' the moment conditions are evaluated. One can choose among 6 different
#' options. Depending on the option, further parameters have to be passed.
#' \describe{
#' \item{"equally":}{equally placed points in \code{min_t,max_t}. When
#' provided, user's \code{min_t} and \code{max_t} will be used (when
#' \code{Coinstrained == FALSE}).
#' }
#' \item{"NonOptAr":}{non optimal arithmetic placement.
#' }
#' \item{"uniformOpt":}{uniform optimal placement.
#' }
#' \item{"ArithOpt":}{arithmetic optimal placement.
#' }
#' \item{"Var Opt":}{optimal variance placement as explained above.
#' }
#' \item{"free":}{user needs to pass own set of points in \code{t_free}.
#' }
#' }
#'
#' Parameter \code{WeightingMatrix}: One can choose among 3 different options:
#' \describe{
#' \item{"OptAsym":}{the optimal asymptotic choice.
#' }
#' \item{"DataVar":}{the covariance matrix of the data provided.
#' }
#' \item{"Id":}{the identity matrix.
#' }
#' }
#' }
#' \item{For Cgmm:}{Continuum Generalized Methods of Moments by Carrasco &
#' Kotchoni (2017). The parameters \code{algo, alphaReg, subdivisions,
#' IntegrationMethod, randomIntegrationLaw, s_min, and s_max} must be
#' specified.
#' }
#' \item{For GMC:}{Generalized Method of Cumulants (GMC) by Massing, T.
#' (2022). The parameters \code{algo, alphaReg, regularization,
#' WeightingMatrix, and ncond} must be specified.
#' }
#' }
#'
#' \strong{IterationControl} If \code{algo = "IT..."} or \code{algo =
#' "Cue..."} the user can control each iteration by setting up the list
#' IterationControl which contains the following elements:
#' \describe{
#' \item{NbIter}{maximum number of iteration. The loop stops when NBIter is
#' reached; default = 10.}
#' \item{PrintIterlogical}{if set to TRUE, the value of the current parameter
#' estimation is printed to the screen at each iteration; default = TRUE.}
#' \item{RelativeErrMax}{the loop stops if the relative error between two
#' consecutive estimation steps is smaller than RelativeErrMax;
#' default = 1e-3.}
#' }
#'
#' \strong{methodR} Random numbers must be generated for each MC study. For
#' each distribution, different methods are available for this (partly also
#' depending on alpha). For more information, the documentation of the
#' respective \code{r...()} distribution can be called up. By default, the fastest
#' method is selected. Since the deviation error can amplify to the edges of
#' alpha depending on the method, it is recommended to check the generated
#' random numbers once for each distribution using the density function before
#' starting the simulation.
#'
#' \strong{Parallelization} Parallelization of the function is possible with
#' using [parallelizeMCsimulation()]. If someone wants to parallelize the
#' function manually, the parameter \code{MCparam} must be set to \code{1} and
#' the parameter \code{SeedOption} must be changed for each iteration.
#'
#' Since this package is structurally based on the \strong{"StableEstim"
#' package by Tarak Kharrat and Georgi N. Boshnakov}, more detailed
#' documentation can be found in their documentation.
#'
#' @seealso
#' \url{https://github.com/GeoBosh/StableEstim/blob/master/R/Simulation.R}
#'
#' @references
#' Massing, T. (2023), 'Parametric Estimation of Tempered Stable Laws'
#'
#' Kim, Y. s.; Rachev, S. T.; Bianchi, M. L. & Fabozzi, F. J. (2008), 'Financial
#' market models with lévy processes and time-varying volatility'
#' \doi{10.1016/j.jbankfin.2007.11.004}
#'
#' Hansen, L. P. (1982), 'Large sample properties of generalized method of
#' moments estimators' \doi{10.2307/1912775}
#'
#' Hansen, L. P.; Heaton, J. & Yaron, A. (1996), 'Finite-Sample Properties of
#' Some Alternative GMM Estimators' \doi{10.1080/07350015.1996.10524656}
#'
#' Feuerverger, A. & McDunnough, P. (1981), 'On the efficiency of empirical
#' characteristic function procedures'
#' \doi{10.1111/j.2517-6161.1981.tb01143.x}
#'
#' Carrasco, M. & Kotchoni, R. (2017), 'Efficient estimation using the
#' characteristic function' \doi{10.1017/S0266466616000025};
#'
#' Kuechler, U. & Tappe, S. (2013), 'Tempered stable distribution and processes'
#' \doi{10.1016/j.spa.2013.06.012}
#'
#' @param ParameterMatrix The matrix is to be composed of vectors, row by row.
#' Each vector must fit the pattern of theta of the \code{TemperedType}.
#' @param SampleSizes Sample sizes to be used to simulate the data. By default,
#' we use 200 (small sample size) and 1600 (large sample size);
#' vector of integer.
#' @param MCparam Number of Monte Carlo simulation for each couple of parameter,
#' default=100; integer
#' @param TemperedType A String. Either "CTS", "TSS", "NTS", "MTS", "GTS",
#' "KRTS", "RDTS".
#' @param Estimfct The estimation function to be used. A String.
#' Either "ML", "GMM", "Cgmm", or "GMC".
#' @param HandleError Logical flag: if set to TRUE, the simulation doesn't stop
#' when an error in the estimation function is encountered. A vector of
#' (size 4) NA is saved and the the simulation carries on. See details.
#' @param saveOutput Logical flag: if set to TRUE, a csv file (for each couple
#' of parameter) with the the estimation
#' information is saved in the current directory. See details.
#' @param SeedOptions List to control the seed generation. See details.
#' @param eps Numerical error tolerance. \code{1e-06} by default.
#' @param algo algorithm: For GMM: \code{"2SGMM"} is the two step GMM proposed
#' by Hansen (1982). \code{"CueGMM"} and \code{"ITGMM"} are respectively the
#' continuous updated and the iterative GMM proposed by Hansen, Eaton et Yaron
#' (1996) and adapted to the continuum case. For GMC: \code{"2SGMC", "CueGMC"}.
#' For Cgmm: \code{"2SCgmm", "CueCgmm", ...}.
#' @param regularization regularization scheme to be used for moment methods,
#' one of \code{"Tikhonov"} (Tikhonov), \code{"LF"} (Landweber-Fridmann) and
#' \code{"cut-off"} (spectral cut-off).
#' @param WeightingMatrix type of weighting matrix used to compute the
#' objective function for the GMM and GMC methods, one of \code{"OptAsym"} (the
#' optimal asymptotic), \code{"DataVar"} (the data driven, only for GMM) and
#' \code{"Id"} (the identity matrix).
#' @param t_scheme scheme used to select the points for the GMM method where the
#' moment conditions are evaluated, one of \code{"equally"} (equally placed),
#' \code{"NonOptAr"} (non optimal arithmetic placement), \code{"uniformOpt"}
#' (uniform optimal placement), \code{"ArithOpt"} (arithmetic optimal placement)
#' , \code{"Var Opt"} (optimal variance placement) and \code{"free"} (users need
#' to pass their own set of points in ...).
#' @param alphaReg value of the regularisation parameter; numeric. Example Value
#' could be ==0.01.
#' @param t_free sequence, if \code{t_scheme=="free"}.
#' @param subdivisions Number of subdivisions used to compute the different
#' integrals involved in the computation of the objective function for the Cgmm
#' method (to minimise); numeric.
#' @param IntegrationMethod Numerical integration method to be used to
#' approximate the (vectorial) integrals for the Cgmm method. Users can choose
#' between "Uniform" discretization or the "Simpson"'s rule (the 3-point
#' Newton-Cotes quadrature rule).
#' @param randomIntegrationLaw Probability measure associated to the Hilbert
#' space spanned by the moment conditions for the Cgmm method.
#' @param s_min,s_max Lower and Upper bounds of the interval where the moment
#' conditions are considered for the Cgmm method; numeric.
#' @param ncond Integer. Number of moment conditions (until order \code{ncond})
#' for the GMC method. Must not be less than 3 for TSS, 6 for CTS, 5 for NTS.
#' @param IterationControl only used if algo = "IT..." or algo = "Cue..."
#' to control the iterations. See Details.
#' @param nb_t integer, if you set \code{t_scheme <- "equally"}. nb_t could be
#' == 20 for example.
#' @param methodR A string. Method generates random variates of TS distribution.
#' "TM" by default. Switches automatically if the method is not applicable in
#' this way.
#' @param ... Other arguments to be passed to the estimation function.
#'
#' @return If \code{saveOutput == FALSE}, the return object is a list of 2.
#' Results of the simulation are listed in \code{$outputMat}. If \code{
#' saveOutput == TRUE}, only a csv file is saved and nothing is returned.
#'
#' @examples
#' \donttest{
#' TemperedEstim_Simulation(ParameterMatrix = rbind(c(1.5,1,1,1,1,0),
#' c(0.5,1,1,1,1,0)),
#' SampleSizes = c(4), MCparam = 4,
#' TemperedType = "CTS", Estimfct = "ML",
#' saveOutput = FALSE)
#'
#' TemperedEstim_Simulation(ParameterMatrix = rbind(c(1.5,1,1,1,1,0)),
#' SampleSizes = c(4), MCparam = 4,
#' TemperedType = "CTS", Estimfct = "GMM",
#' saveOutput = FALSE, algo = "2SGMM",
#' regularization = "cut-off",
#' WeightingMatrix = "OptAsym", t_scheme = "free",
#' alphaReg = 0.01,
#' t_free = seq(0.1,2,length.out=12))
#'
#' TemperedEstim_Simulation(ParameterMatrix = rbind(c(1.45,0.55,1,1,1,0)),
#' SampleSizes = c(4), MCparam = 4,
#' TemperedType = "CTS", Estimfct = "Cgmm",
#' saveOutput = FALSE, algo = "2SCgmm",
#' alphaReg = 0.01, subdivisions = 50,
#' IntegrationMethod = "Uniform",
#' randomIntegrationLaw = "unif",
#' s_min = 0, s_max= 1)
#'
#' TemperedEstim_Simulation(ParameterMatrix = rbind(c(1.45,0.55,1,1,1,0)),
#' SampleSizes = c(4), MCparam = 4,
#' TemperedType = "CTS", Estimfct = "GMC",
#' saveOutput = FALSE, algo = "2SGMC",
#' alphaReg = 0.01, WeightingMatrix = "OptAsym",
#' regularization = "cut-off", ncond = 8)
#' }
#'
#' @export
TemperedEstim_Simulation <- function(ParameterMatrix,
SampleSizes = c(200, 1600),
MCparam = 100,
TemperedType = c("CTS", "TSS",
"NTS", "MTS", "GTS",
"KRTS", "RDTS"),
Estimfct = c("ML", "GMM", "Cgmm", "GMC"),
HandleError = TRUE, saveOutput = FALSE,
SeedOptions = NULL, eps = 1e-06,
algo = NULL, regularization = NULL,
WeightingMatrix = NULL, t_scheme = NULL,
alphaReg = NULL, t_free = NULL,
nb_t = NULL, subdivisions = NULL,
IntegrationMethod = NULL,
randomIntegrationLaw = NULL, s_min = NULL,
s_max = NULL, ncond = NULL,
IterationControl = NULL,
methodR = "TM", ...) {
#seeAlso: https://github.com/GeoBosh/StableEstim/blob/master/R/Simulation.R
SeedVector <- getSeedVector(MCparam, SeedOptions)
Estimfct <- match.arg(arg = Estimfct,
choices = c("ML", "GMM", "Cgmm", "GMC"))
TemperedType <- match.arg(arg = TemperedType, choices =
c("CTS", "TSS","NTS", "MTS", "GTS",
"KRTS", "RDTS"))
nab <- nrow(ParameterMatrix)
npar <- ncol(ParameterMatrix)
lS <- length(SampleSizes)
nRowOutput <- nab * lS
OutputCollection <- empty_list <- vector(mode = "list", length = nab)
returnList <- empty_list <- vector(mode = "list")
#returnList <- matrix(data = NA, ncol = npar, nrow = nab*MCparam)
indexStatOutput <- 1
#Values that needs to be set in the beginning as checkpoint files are not
# allowed
ab <- 1
sample <- 1
mc <- 0
if(MCparam != 1){
CheckPointValues <- readCheckPoint(ParameterMatrix, TemperedType,
Estimfct, nab, npar, lS, MCparam,
eps = eps,
algo = algo,
regularization = regularization,
WeightingMatrix =
WeightingMatrix,
t_scheme = t_scheme,
alphaReg = alphaReg,
t_free = t_free,
nb_t = nb_t,
subdivisions = subdivisions,
IntegrationMethod =
IntegrationMethod,
randomIntegrationLaw =
randomIntegrationLaw,
s_min = s_min,
s_max = s_max,
ncond = ncond,
IterationControl = IterationControl,
ab = ab, sample = sample, mc = mc,
...)
}
else{
CheckPointValues <- list(ab = ab, nab = nab, npar = npar, sample = sample,
nSS = lS, mc = mc, MCparam = 1)
}
updatedCheckPointValues <- updateCheckPointValues(CheckPointValues, MCparam,
lS, nab)
# Ist im Zielordner bereits eine csv-Datei mit dem gleichen Namen (thetaT
# und MCparam sind gleich), wird die Datei aktuell um die weiteren
# Ergebnisse aktualisiert. Mention in Details
# Wird dieser Code aktiviert, würde das Überschreiben einer Datei verboten
# werden.
#
#if (updatedCheckPointValues$mc_start != 1){
# print("'Can't Compute Stat summary when the process doesn't
# start from the beginning!!")
#}
for (ab in updatedCheckPointValues$ab_start:nab) {
thetaT <- ParameterMatrix[ab, ]
outputString <- switch(TemperedType,
CTS = paste("Alpha=", thetaT[1] ,
" *** DeltaP=", thetaT[2],
" *** DeltaM=", thetaT[3],
" *** LambdaP=", thetaT[4],
" *** LambdaM=", thetaT[5],
" *** mu=", thetaT[6], sep = ""),
TSS = paste("Alpha=", thetaT[1] ,
" *** Delta=", thetaT[2],
" *** Lambda=", thetaT[3], sep = ""),
NTS = paste("Alpha=", thetaT[1] ,
" *** Beta=", thetaT[2],
" *** Delta=", thetaT[3],
" *** Lambda=", thetaT[4],
" *** mu=", thetaT[5], sep = ""),
MTS = paste("Alpha=", thetaT[1] ,
" *** Delta=", thetaT[2],
" *** LambdaP=", thetaT[3],
" *** LambdaM=", thetaT[4],
" *** mu=", thetaT[5], sep = ""),
GTS = paste("AlphaP=", thetaT[1] ,
" *** AlphaM=", thetaT[2],
" *** DeltaP=", thetaT[3],
" *** DeltaM=", thetaT[4],
" *** LambdaP=", thetaT[5],
" *** LambdaM=", thetaT[6],
" *** mu=", thetaT[7], sep = ""),
KRTS = paste("Alpha=", thetaT[1] ,
" *** kP=", thetaT[2],
" *** kM=", thetaT[3],
" *** rP=", thetaT[4],
" *** rM=", thetaT[5],
" *** pP=", thetaT[6],
" *** pM=", thetaT[7],
" *** mu=", thetaT[8], sep = ""),
RDTS = paste("Alpha=", thetaT[1] ,
" *** Delta=", thetaT[2],
" *** LambdaP=", thetaT[3],
" *** LambdaM=", thetaT[4],
" *** mu=", thetaT[5], sep = "")
# ,CGMY = paste("C=", thetaT[1] ,
# " *** G=", thetaT[2],
# " *** M=", thetaT[3],
# " *** Y=", thetaT[4], sep = "")
)
cat("---------------- ", outputString, " --------------- \n", sep = "")
if (saveOutput) initOutputFile(thetaT, MCparam, TemperedType,
Estimfct, ...)
EstimOutput <- ComputeMCSimForTempered(thetaT = thetaT,
MCparam = MCparam,
SampleSizes =
as.vector(SampleSizes),
SeedVector = SeedVector,
TemperedType = TemperedType,
Estimfct = Estimfct,
HandleError = HandleError,
ab_current = ab,
nab = nab,
npar = npar, ParameterMatrix,
CheckPointValues =
updatedCheckPointValues,
saveOutput = saveOutput,
eps = eps,
algo = algo,
regularization = regularization,
WeightingMatrix =
WeightingMatrix,
t_scheme = t_scheme,
alphaReg = alphaReg,
t_free = t_free,
nb_t = nb_t,
subdivisions = subdivisions,
IntegrationMethod =
IntegrationMethod,
randomIntegrationLaw =
randomIntegrationLaw,
s_min = s_min,
s_max = s_max,
ncond = ncond,
IterationControl =
IterationControl,
methodR = methodR,
...)
OutputCollection <- EstimOutput
if (saveOutput == FALSE){
if (length(returnList) == 0) returnList <- EstimOutput
else returnList <- Map(f = rbind, x = returnList, init = EstimOutput)
}
}
# if(MCparam != 1){
# deleteCheckPoint(ParameterMatrix, TemperedType, Estimfct, nab, npar, lS,
# MCparam,
# eps = eps,
# algo = algo,
# regularization = regularization,
# WeightingMatrix =
# WeightingMatrix,
# t_scheme = t_scheme,
# alphaReg = alphaReg,
# t_free = t_free,
# subdivisions = subdivisions,
# IntegrationMethod =
# IntegrationMethod,
# randomIntegrationLaw =
# randomIntegrationLaw,
# s_min = s_min,
# s_max = s_max,
# ncond = ncond,
# IterationControl = IterationControl,
# methodR = methodR,
# ...)
# }
if (saveOutput == FALSE){
return(returnList)
}
}
#' Function to parallelize the Monte Carlo Simulation
#'
#' Since the Monte Carlo Simulation is very computationally intensive, it may
#' be worthwhile to split it across all available processor cores. To do this,
#' simply pass all the parameters from the [TemperedEstim_Simulation()]
#' function to this function in the same way.
#'
#' In this function exactly the arguments must be passed, which are also needed
#' for the function [TemperedEstim_Simulation()]. However, a few functions of
#' [TemperedEstim_Simulation()] are not possible here. The restrictions are
#' described in more detail for the individual arguments.
#'
#' In addition to the arguments of function [TemperedEstim_Simulation()], the
#' argument "cores" can be assigned an integer value. This value determines how
#' many different processes are to be parallelized. If value is \code{NULL}, R
#' tries to read out how many cores the processor has and passes this
#' value to "cores".
#'
#' During the simulation, the progress of the simulation can be viewed in a
#' file in the workspace named "IterationControlForParallelization.txt".
#'
#' @param ParameterMatrix The matrix is to be composed of vectors, row by row.
#' Each vector must fit the pattern of theta of the \code{TemperedType}.
#' Compared to the function [TemperedEstim_Simulation()], the matrix here may
#' contain only one parameter vector.
#' @param MCparam Number of Monte Carlo simulation for each couple of parameter,
#' default=100; integer
#' @param SampleSizes Sample sizes to be used to simulate the data. By default,
#' we use 200 (small sample size). Vector of integer. Compared to the function
#' [TemperedEstim_Simulation()], the vector here may contain only one integer.
#' @param saveOutput Logical flag: In the function [TemperedEstim_Simulation()]
#' the argument can be true. Then an external csv file is created. Here the
#' argument must be false. The output of the values works in this function
#' exclusively via the return of the function.
#' @param SeedOptions is an argument what can be used in
#' [TemperedEstim_Simulation()] but must be NULL here.
#' @param cores size of cluster for parallelization. Positive Integer.
#' @param iterationDisplayToFileSystem creates a text file in your file system
#' that displays the current iteration of the simulation.
#' @param ... The function works only if all necessary arguments from the
#' function [TemperedEstim_Simulation()] are passed. See description and
#' details.
#'
#' @return The return object is a list of 2. Results of the simulation are
#' listed in \code{$outputMat}.
#'
#' @export
#' @importFrom foreach %dopar%
parallelizeMCsimulation <- function(
ParameterMatrix,
MCparam = 10000,
SampleSizes = c(200),
saveOutput = FALSE,
cores = 2,
SeedOptions = NULL,
iterationDisplayToFileSystem = FALSE,
...){
mc <- NULL
if (!is.null(SeedOptions)){
stop("SeedOptions is used by the function and cannot be passed as an
argument")
}
if(is.null(ParameterMatrix) || nrow(ParameterMatrix) != 1){
stop("Compared to the function TemperedEstim_Simulation(), the matrix here
may contain only one parameter vector.")
}
if(is.null(SampleSizes) || length(SampleSizes) != 1){
stop("Compared to the function TemperedEstim_Simulation(), the SampleSizes
vector here may contain only one integer.")
}
if(is.null(saveOutput) || saveOutput){
stop("Compared to the function TemperedEstim_Simulation(), saveOutput
must be FALSE.")
}
if (is.null(cores) || cores < 1){
cores <- parallel::detectCores()
}
if(MCparam < cores){
cores <- MCparam
}
R <- MCparam #MonteCarloRuns
cl <- parallel::makeCluster(cores)
doParallel::registerDoParallel(cl)
resultOfSimulation <- foreach::foreach(mc = 1:R, .combine = rbind,
.export = ls("package:TempStable",
all.names = TRUE)
)%dopar%{
returnValue <- TemperedEstim_Simulation(
ParameterMatrix = ParameterMatrix,
MCparam = 1,
SampleSizes = SampleSizes,
saveOutput = saveOutput,
SeedOptions = list(MCtot = R, seedStart = mc),
... = ...
)
if (iterationDisplayToFileSystem){
utils::write.table(
paste("Last Monte Carlo run: ", mc) ,
file = base::paste("IterationControlForParallelization.txt"),
sep = "\t", row.names = FALSE)
}
returnValue$outputMat
}
parallel::stopCluster(cl)
attr(resultOfSimulation, "rng") <- NULL
attr(resultOfSimulation, "doRNG_version") <- NULL
#Delete txt file
if (iterationDisplayToFileSystem){
base::unlink(x = base::paste("IterationControlForParallelization.txt"),
force = TRUE)
}
return(resultOfSimulation)
}
# No export.
getSeedVector <- function(Outputsize, SeedOptions = NULL) {
set.seed(345)
if (is.null(SeedOptions))
vec <- as.vector(sample.int(n = 3 * Outputsize, size = Outputsize))
else {
MCtot <- SeedOptions$MCtot
seedStart <- SeedOptions$seedStart
seedEnd <- seedStart + Outputsize
vec <-
as.vector(sample.int(n = 3 * MCtot, size = MCtot))[seedStart:seedEnd]
}
vec
}
# No export.
ComputeMCSimForTempered <- function(thetaT, MCparam, SampleSizes, SeedVector,
TemperedType, Estimfct, HandleError,
ab_current, nab, npar, ParameterMatrix,
CheckPointValues = NULL, saveOutput, eps,
algo, regularization, WeightingMatrix,
t_scheme, alphaReg, t_free,nb_t,
subdivisions, IntegrationMethod,
randomIntegrationLaw, s_min, s_max, ncond,
IterationControl, methodR, ...)
{
if (TemperedType == "CTS") {
Ncol <- 16
} else if (TemperedType == "TSS") {
Ncol <- 10
} else if (TemperedType == "NTS" || TemperedType == "MTS" ||
TemperedType == "RDTS") {
Ncol <- 14
} else if (TemperedType == "GTS") {
Ncol <- 18
} else if (TemperedType == "KRTS") {
Ncol <- 20
}
# else {
# Ncol <- 12
# }
nSS <- length(SampleSizes)
Nrow <- nSS * MCparam
Output <- matrix(data = NA, ncol = Ncol, nrow = Nrow)
if (TemperedType == "CTS") {
colnames(Output) <- c("alphaT", "delta+T", "delta-T", "lambda+T",
"lambda-T", "muT", "data size", "seed", "alphaE",
"delta+E", "delta-E", "lambda+E", "lambda-E",
"muE", "failure", "time")
} else if (TemperedType == "TSS") {
colnames(Output) <- c("alphaT", "deltaT", "lambdaT", "data size",
"seed", "alphaE", "deltaE", "lambdaE",
"failure", "time")
} else if (TemperedType == "NTS") {
colnames(Output) <- c("alphaT", "betaT", "deltaT", "lambdaT", "muT",
"data size", "seed", "alphaE", "betaE", "deltaE",
"lambdaE", "muE", "failure", "time")
} else if (TemperedType == "MTS") {
colnames(Output) <- c("alphaT", "delta", "lambda+T", "lambda-T", "muT",
"data size", "seed", "alphaE", "deltaE", "lambda+E",
"lambda-E", "muE", "failure", "time")
} else if (TemperedType == "GTS") {
colnames(Output) <- c("alpha+T", "alpha-T", "delta+T", "delta-T",
"lambda+T", "lambda-T", "muT", "data size", "seed",
"alpha+E", "alpha-E", "delta+E", "delta-E",
"lambda+E", "lambda-E", "muE", "failure", "time")
} else if (TemperedType == "KRTS") {
colnames(Output) <- c("alphaT", "k+T", "k-T", "r+T", "r-T",
"p+T", "p-T", "muT", "data size", "seed",
"alphaE", "k+E", "k-E", "r+E", "r-E",
"p+E", "p-E", "muE", "failure", "time")
} else if (TemperedType == "RDTS") {
colnames(Output) <- c("alphaT", "delta", "lambda+T", "lambda-T", "muT",
"data size", "seed", "alphaE", "deltaE", "lambda+E",
"lambda-E", "muE", "failure", "time")
}
# else {
# colnames(Output) <- c("C.T", "G.T", "M.T", "Y.T", "data size", "seed",
# "C.E", "G.E", "M.E", "Y.E", "failure", "time")
# }
if (ab_current == CheckPointValues$ab_start) {
sample_start = CheckPointValues$sample_start
mc_start = CheckPointValues$mc_start
} else {
sample_start = 1
mc_start = 1
}
for (sample in sample_start:nSS) {
size <- SampleSizes[sample]
if (sample != sample_start) mc_start = 1
for (mc in mc_start:MCparam) {
tIter <- getTime_()
iter <- mc + (sample - 1) * MCparam
set.seed(seed <- SeedVector[mc])
if (TemperedType == "CTS") {
x <- rCTS(n = size, alpha = thetaT[1], deltap = thetaT[2],
deltam = thetaT[3], lambdap = thetaT[4],
lambdam = thetaT[5], mu = thetaT[6], methodR = methodR, ...)
} else if (TemperedType == "TSS") {
x <- rTSS(n = size, alpha = thetaT[1], delta = thetaT[2],
lambda = thetaT[3], methodR = methodR, ...)
} else if (TemperedType == "NTS") {
x <- rNTS(n = size, alpha = thetaT[1], beta = thetaT[2],
delta = thetaT[3], lambda = thetaT[4], mu = thetaT[5],
methodR = methodR, ...)
} else if (TemperedType == "MTS") {
x <- rMTS(n = size, theta = thetaT, methodR = methodR, ...)
} else if (TemperedType == "GTS") {
x <- rGTS(n = size, theta = thetaT, methodR = methodR, ...)
} else if (TemperedType == "KRTS") {
x <- rKRTS(n = size, alpha = thetaT[1], kp = thetaT[2], km = thetaT[3],
rp = thetaT[4], rm = thetaT[5], pp = thetaT[6],
pm = thetaT[7], mu = thetaT[8], methodR = methodR, ...)
} else if (TemperedType == "RDTS") {
x <- rRDTS(n = size, theta = thetaT, methodR = methodR, ...)
}
# else {
# x <- rCGMY(n = size, C = thetaT[1], M = thetaT[2],
# G = thetaT[3], Y = thetaT[4])
# }
Estim <- getTempEstimation(thetaT = thetaT, x = x, seed = seed,
size = size, Ncol = Ncol,
TemperedType = TemperedType,
Estimfct = Estimfct,
HandleError = HandleError, eps = eps,
algo = algo, regularization = regularization,
WeightingMatrix =
WeightingMatrix,
t_scheme = t_scheme,
alphaReg = alphaReg,
t_free = t_free,
nb_t = nb_t,
subdivisions = subdivisions,
IntegrationMethod =
IntegrationMethod,
randomIntegrationLaw =
randomIntegrationLaw,
s_min = s_min,
s_max = s_max,
ncond = ncond,
IterationControl = IterationControl,
...)
Output[iter, ] <- Estim$outputMat
file <- Estim$file
#When checkpoint file should be availabe again
# if (!is.null(CheckPointValues) && MCparam != 1) {
# writeCheckPoint(ParameterMatrix, TemperedType, Estimfct, ab_current,
# nab, npar, sample, nSS, mc,
# MCparam,
# eps,
# algo,
# regularization,
# WeightingMatrix,
# t_scheme,
# alphaReg,
# t_free,
# subdivisions,
# IntegrationMethod,
# randomIntegrationLaw,
# s_min,
# s_max,
# ncond,
# IterationControl,
# ...)
# }
if (saveOutput) updateOutputFile(thetaT, MCparam, TemperedType,
Estim)
StableEstim::PrintEstimatedRemainingTime(iter, tIter, Nrow)
}
}
#End Sample
if (isFALSE(saveOutput)){
return(list(outputMat = Output))
}
else{
return(list(outputMat = Output, file = file))
}
}
# No export.
getTempEstimation <- function(thetaT, x, seed, size, Ncol, TemperedType,
Estimfct, HandleError, eps,
algo,
regularization,
WeightingMatrix,
t_scheme,
alphaReg,
t_free,
nb_t,
subdivisions,
IntegrationMethod,
randomIntegrationLaw,
s_min,
s_max,
ncond,
IterationControl,
...) {
output <- vector(length = Ncol)
if (TemperedType == "CTS") {
output[1:8] <- c(thetaT, size, seed)
} else if (TemperedType == "TSS") {
output[1:5] <- c(thetaT, size, seed)
} else if (TemperedType == "NTS") {
output[1:7] <- c(thetaT, size, seed)
} else if (TemperedType == "MTS") {
output[1:7] <- c(thetaT, size, seed)
} else if (TemperedType == "GTS") {
output[1:9] <- c(thetaT, size, seed)
} else if (TemperedType == "KRTS") {
output[1:10] <- c(thetaT, size, seed)
} else if (TemperedType == "RDTS") {
output[1:7] <- c(thetaT, size, seed)
}
# else { #CGMY
# output[1:6] <- c(thetaT, size, seed)
# }
theta0 <- thetaT - 0.1 #noise
EstimRes <- TemperedEstim(TemperedType = TemperedType,
EstimMethod = Estimfct,
data = x, theta0 = theta0, ComputeCov = FALSE,
HandleError = HandleError, eps = eps,
algo = algo, regularization = regularization,
WeightingMatrix =
WeightingMatrix,
t_scheme = t_scheme,
alphaReg = alphaReg,
t_free = t_free,
nb_t = nb_t,
subdivisions = subdivisions,
IntegrationMethod =
IntegrationMethod,
randomIntegrationLaw =
randomIntegrationLaw,
s_min = s_min,
s_max = s_max,
ncond = ncond,
IterationControl = IterationControl, ...)
if (TemperedType == "CTS") {
output[9:14] <- EstimRes@par
} else if (TemperedType == "TSS") {
output[6:8] <- EstimRes@par
} else if (TemperedType == "NTS") {
output[8:12] <- EstimRes@par
} else if (TemperedType == "MTS") {
output[8:12] <- EstimRes@par
} else if (TemperedType == "GTS") {
output[10:16] <- EstimRes@par
} else if (TemperedType == "KRTS") {
output[11:18] <- EstimRes@par
} else if (TemperedType == "RDTS") {
output[8:12] <- EstimRes@par
}
# else {
# output[7:10] <- EstimRes@par
# }
if (TemperedType == "CTS") {
output[15:16] <- c(EstimRes@failure, EstimRes@duration)
} else if (TemperedType == "TSS") {
output[9:10] <- c(EstimRes@failure, EstimRes@duration)
} else if (TemperedType == "NTS") {
output[13:14] <- c(EstimRes@failure, EstimRes@duration)
} else if (TemperedType == "MTS") {
output[13:14] <- c(EstimRes@failure, EstimRes@duration)
} else if (TemperedType == "GTS") {
output[17:18] <- c(EstimRes@failure, EstimRes@duration)
} else if (TemperedType == "KRTS") {
output[19:20] <- c(EstimRes@failure, EstimRes@duration)
} else if (TemperedType == "RDTS") {
output[13:14] <- c(EstimRes@failure, EstimRes@duration)
}
# else {
# output[11:12] <- c(EstimRes@failure, EstimRes@duration)
# }
list(outputMat = output, file = EstimRes@method)
}
# Function title
# Merge with EstimSimulation
#
# @examples
# ComputeMCSimForTempered_parallel(1,c(1.5, 1, 1, 1, 1, 0),10,"CTS","ML")
# ComputeMCSimForTempered_parallel(1,c(0.5, 1, 1),10,"TSS","Cgmm",
# IntegrationMethod = "Simpson",
# randomIntegrationLaw = "unif")
#
# No export.
# ComputeMCSimForTempered_parallel <- function(MCparam, thetaT, size,
# TemperedType = c("CTS",
# "TSS",
# "NTS","CGMY"),
# Estimfct = c("ML", "GMM", "Cgmm",
# "GMC"),
# HandleError = TRUE, eps, ...) {
# Estimfct <- match.arg(Estimfct)
# TemperedType <- match.arg(TemperedType)
# Ncol <- ifelse(TemperedType == "CTS", 15, 9)
# if (TemperedType == "CTS") {
# Ncol <- 15
# } else if (TemperedType == "TSS") {
# Ncol <- 9
# } else if (TemperedType == "NTS") {
# Ncol <- 13
# } else {
# Ncol <- 11
# }
# Output <- numeric(Ncol)
# if (TemperedType == "CTS") {
# names(Output) <- c("alphaT", "delta+T", "delta-T", "lambda+T",
# "lambda-T", "muT", "data size", "alphaE", "delta+E",
# "delta-E", "lambda+E", "lambda-E", "muE", "failure",
# "time")
# } else if (TemperedType == "TSS") {
# names(Output) <- c("alphaT", "deltaT", "lambdaT", "data size", "alphaE",
# "deltaE", "lambdaE", "failure", "time")
# } else if (TemperedType == "NTS") {
# names(Output) <- c("alphaT", "betaT", "deltaT", "lambdaT", "muT",
# "data size", "alphaE", "betaE", "deltaE", "lambdaE",
# "muE", "failure", "time")
# } else {
# names(Output) <- c("C.T", "G.T", "M.T", "Y.T", "data size", "C.E",
# "G.E", "M.E", "Y.E", "failure", "time")
# }
#
# if (TemperedType == "CTS") {
# x <- rCTS(n = size, alpha = thetaT[1], deltap = thetaT[2],
# deltam = thetaT[3], lambdap = thetaT[4], lambdam = thetaT[5],
# mu = thetaT[6])
# } else if (TemperedType == "TSS") {
# x <- rTSS(n = size, alpha = thetaT[1], delta = thetaT[2],
# lambda = thetaT[3])
# } else if (TemperedType == "NTS") {
# x <- rNTS(n = size, alpha = thetaT[1], beta = thetaT[2],
# delta = thetaT[3], lambda = thetaT[4], mu = thetaT[5])
# } else {
# x <- rCGMY(n = size, C = theta[1], G = theta[2], G = theta[3],
# Y = theta[4])
# }
# Estim <- getTempEstimation_parallel(thetaT = thetaT, x = x, size = size,
# Ncol = Ncol,
# TemperedType = TemperedType,
# Estimfct = Estimfct,
# HandleError = HandleError,
# eps = eps, ...)
# Output <- c(MCparam, Estim)
#
# return(Output)
# }
# No export.
# getTempEstimation_parallel <- function(thetaT, x, size, Ncol, TemperedType,
# Estimfct, HandleError, eps, ...) {
# output <- vector(length = Ncol)
# if (TemperedType == "CTS") {
# output[1:7] <- c(thetaT, size)
# } else if (TemperedType == "TSS") {
# output[1:4] <- c(thetaT, size)
# } else if (TemperedType == "NTS") {
# output[1:6] <- c(thetaT, size)
# } else {
# output[1:5] <- c(thetaT, size)
# }
# theta0 <- thetaT - 0.1 #noise
# EstimRes <- TemperedEstim_v2(TemperedType = TemperedType,
# EstimMethod = Estimfct, data = x,
# theta0 = theta0, ComputeCov = FALSE,
# HandleError = HandleError, eps = eps, ...)
# if (TemperedType == "CTS") {
# output[8:13] <- EstimRes$par
# } else if (TemperedType == "TSS") {
# output[5:7] <- EstimRes$par
# } else if (TemperedType == "NTS") {
# output[7:11] <- EstimRes$par
# } else {
# output[6:9] <- EstimRes$par
# }
# if (TemperedType == "CTS") {
# output[14:15] <- c(EstimRes$failure, EstimRes$duration)
# } else if (TemperedType == "TSS") {
# output[8:9] <- c(EstimRes$failure, EstimRes$duration)
# } else if (TemperedType == "NTS") {
# output[12:13] <- c(EstimRes$failure, EstimRes$duration)
# } else {
# output[10:11] <- c(EstimRes$failure, EstimRes$duration)
# }
# return(output)
# }
##### for statistical summary#####
# Function title
#
# Gap holder for description.
#
# Gap holder for details.
#
# @param EstimOutput A gap holder.
# @param FctsToApply A gap holder.
# @param SampleSizes A gap holder.
# @param CheckMat A gap holder.
# @param tolFailCheck A gap holder.
# @param MCparam A gap holder.
#
# @return Gap holder for return.
#
# @export
# ComputeStatOutput <- function(EstimOutput, FctsToApply, SampleSizes, CheckMat,
# tolFailCheck, MCparam, ...) {
# list(alpha = ComputeStatOutputPar(EstimOutput = EstimOutput,
# FctsToApply = FctsToApply, par = "alpha",
# SampleSizes = SampleSizes,
# CheckMat = CheckMat,
# tolFailCheck = tolFailCheck,
# MCparam = MCparam, ...),
# beta = ComputeStatOutputPar(EstimOutput = EstimOutput,
# FctsToApply = FctsToApply, par = "beta",
# SampleSizes = SampleSizes,
# CheckMat = CheckMat,
# tolFailCheck = tolFailCheck,
# MCparam = MCparam, ...),
# gamma = ComputeStatOutputPar(EstimOutput = EstimOutput,
# FctsToApply = FctsToApply, par = "gamma",
# SampleSizes = SampleSizes,
# CheckMat = CheckMat,
# tolFailCheck = tolFailCheck,
# MCparam = MCparam, ...),
# delta = ComputeStatOutputPar(EstimOutput = EstimOutput,
# FctsToApply = FctsToApply, par = "delta",
# SampleSizes = SampleSizes,
# CheckMat = CheckMat,
# tolFailCheck = tolFailCheck,
# MCparam = MCparam, ...))
# }
##### for Output File#####
# No export.
initOutputFile <- function(thetaT, MCparam, TemperedType, Estimfct, ...) {
method <- Estim_Des_Temp(TemperedType, Estimfct, ...)
fileName <- get_filename(thetaT, MCparam, TemperedType, method)
if (!file.exists(fileName)) {
x <- switch(TemperedType,
CTS = paste("alphaT", "delta+T", "delta-T", "lambda+T",
"lambda-T", "muT", "data size", "seed",
"alphaE", "delta+E", "delta-E", "lambda+E",
"lambda-E", "muE", "failure", "time",
sep = ","),
TSS = paste("alphaT", "deltaT", "lambdaT",
"data size", "seed", "alphaE", "deltaE",
"lambdaE", "failure", "time",
sep = ","),
NTS = paste("alphaT", "betaT", "deltaT", "lambdaT", "muT",
"data size", "seed", "alphaE", "betaE", "deltaE",
"lambdaE", "muE", "failure", "time",
sep = ","),
MTS = paste("alphaT", "delta", "lambda+T", "lambda-T", "muT",
"data size", "seed", "alphaE", "deltaE", "lambda+E",
"lambda-E", "muE", "failure", "time", sep = ""),
GTS = paste("alpha+T", "alpha-T", "delta+T", "delta-T",
"lambda+T", "lambda-T", "muT", "data size", "seed",
"alpha+E", "alpha-E", "delta+E", "delta-E",
"lambda+E", "lambda-E", "muE", "failure", "time",
sep = ""),
KRTS = paste("alphaT", "k+T", "k-T", "r+T", "r-T",
"p+T", "p-T", "muT", "data size", "seed",
"alphaE", "k+E", "k-E", "r+E", "r-E",
"p+E", "p-E", "muE", "failure", "time", sep = ""),
RDTS = paste("alphaT", "delta", "lambda+T", "lambda-T", "muT",
"data size", "seed", "alphaE", "deltaE",
"lambda+E", "lambda-E", "muE", "failure", "time",
sep = "")
# ,CGMY = paste("CT", "GT", "MT", "YT", "data size", "seed", "CE",
# "GE", "ME", "YE","failure", "time",
# sep = ",")
)
write(x, file = fileName, sep = "\n")
}
}
# No export.
Estim_Des_Temp <- function(TemperedType = c("CTS", "TSS", "NTS", "MTS", "GTS",
"KRTS", "RDTS"),
EstimMethod = c("ML", "GMM", "Cgmm", "GMC"),
eps,
algo,
regularization,
WeightingMatrix,
t_scheme,
alphaReg,
t_free,
nb_t,
subdivisions,
IntegrationMethod,
randomIntegrationLaw,
s_min,
s_max,
ncond,
IterationControl,
...) {
TemperedType <- match.arg(TemperedType)
EstimMethod <- match.arg(EstimMethod)
EstimFcts <- getTempEstimFcts(TemperedType, EstimMethod,
eps = eps,
algo = algo,
regularization = regularization,
WeightingMatrix =
WeightingMatrix,
t_scheme = t_scheme,
alphaReg = alphaReg,
t_free = t_free,
nb_t = nb_t,
subdivisions = subdivisions,
IntegrationMethod =
IntegrationMethod,
randomIntegrationLaw =
randomIntegrationLaw,
s_min = s_min,
s_max = s_max,
ncond = ncond,
IterationControl = IterationControl,
...)
EstimFcts$methodDes(eps = eps,
algo = algo,
regularization = regularization,
WeightingMatrix =
WeightingMatrix,
t_scheme = t_scheme,
alphaReg = alphaReg,
t_free = t_free,
nb_t = nb_t,
subdivisions = subdivisions,
IntegrationMethod =
IntegrationMethod,
randomIntegrationLaw =
randomIntegrationLaw,
s_min = s_min,
s_max = s_max,
ncond = ncond,
IterationControl = IterationControl,
...)
}
##### for Checkpoints#####
# No export.
readCheckPoint <- function(ParameterMatrix, TemperedType, Estimfct, nab, npar,
nSS, MCparam,
eps,
algo,
regularization,
WeightingMatrix,
t_scheme,
alphaReg,
t_free,
nb_t,
subdivisions,
IntegrationMethod,
randomIntegrationLaw,
s_min,
s_max,
ncond,
IterationControl,
ab,
sample,
mc,
...) {
method <- Estim_Des_Temp(TemperedType, Estimfct,
eps = eps,
algo = algo,
regularization = regularization,
WeightingMatrix =
WeightingMatrix,
t_scheme = t_scheme,
alphaReg = alphaReg,
t_free = t_free,
nb_t = nb_t,
subdivisions = subdivisions,
IntegrationMethod =
IntegrationMethod,
randomIntegrationLaw =
randomIntegrationLaw,
s_min = s_min,
s_max = s_max,
ncond = ncond,
IterationControl = IterationControl, ...)
#This code can be used to create a checkpoint text file.
# fileName <- get_filename_checkPoint_Temp(ParameterMatrix, nab, npar,
# MCparam, method)
# if (!file.exists(fileName)) {
# write(x = "## ab;nab;npar;sample;nSS;mc;MCparam", file = fileName,
# sep = "\n")
# ab <- 1
# sample <- 1
# mc <- 0
# write(x = paste("--", ab, nab, npar, sample, nSS, mc, MCparam,
# sep = ";"), file = fileName, sep = "\n", append = TRUE)
# } else {
# tab <- as.numeric(utils::read.table(file = fileName,
# header = F, sep = ";"))
# ab <- tab[2]
# sample <- tab[5]
# mc <- tab[7]
# n_ab <- tab[3]
# n_par <- tab[4]
# n_SS <- tab[6]
# mc_Param <- tab[8]
# stopifnot(nab == n_ab, npar == n_par, nSS == n_SS, mc_Param == MCparam)
# }
list(ab = ab, nab = nab, npar = npar, sample = sample, nSS = nSS, mc = mc,
MCparam = MCparam)
}
# No export.
# This function writes every value of the current checkpoint in a string as
# return.
get_filename_checkPoint_Temp <- function(ParameterMatrix, nab, npar, MCparam,
method) {
# This should be adapted for every Tempered Type.
# BUT: The filenames are getting too long and this results in errors.
# Version 0.1.0 will not feature checkpoints during calculation.
#
# case: StableEstim
# if (npar == 3) {
# MC <- paste(paste("alpha0=", ParameterMatrix[1, 1], sep = ""),
# paste("delta0=", ParameterMatrix[1, 2], sep = ""),
# paste("lambda0=", ParameterMatrix[1, 3], sep = ""),
# paste("alphan=", ParameterMatrix[nab, 1], sep = ""),
# paste("deltan=", ParameterMatrix[nab, 2], sep = ""),
# paste("lambdan=", ParameterMatrix[nab, 3], sep = ""),
# paste("MCparam", MCparam, sep = ""), sep = "_")
# } else {
# case: "CTS"
# MC <- paste(paste("alpha0=", ParameterMatrix[1, 1], sep = ""),
# paste("delta+0=", ParameterMatrix[1, 2], sep = ""),
# paste("delta-0=", ParameterMatrix[1, 3], sep = ""),
# paste("lambda+0=", ParameterMatrix[1, 4], sep = ""),
# paste("lambda-0=",ParameterMatrix[1, 5], sep = ""),
# paste("mu0=", ParameterMatrix[1, 6], sep = ""),
# paste("alphan=", ParameterMatrix[nab,1], sep = ""),
# paste("delta+n=", ParameterMatrix[nab, 2], sep = ""),
# paste("delta-n=", ParameterMatrix[nab,3], sep = ""),
# paste("lambda+n=", ParameterMatrix[nab, 4], sep = ""),
# paste("lambda-n=", ParameterMatrix[nab,5], sep = ""),
# paste("mun=", ParameterMatrix[nab, 5], sep = ""),
# paste("MCparam", MCparam, sep = ""), sep = "_")
# }
MC <- paste("Test", sep = "")
methodTrunc <- method
if (base::nchar(methodTrunc) > 30){
methodTrunc <- substring(methodTrunc,1,30)
}
fileName <- paste(MC, methodTrunc, "_CHECKPOINT.txt", sep = "")
fileName
}
# No export.
updateCheckPointValues <- function(CheckPointValues, MCparam, lS, nab) {
ab_start <- CheckPointValues$ab
sample_start <- CheckPointValues$sample
mc_start <- CheckPointValues$mc
if (CheckPointValues$mc == MCparam) {
mc_start = 1
if (CheckPointValues$sample == lS) {
sample_start = 1
if (CheckPointValues$ab == nab)
stop("Simulation finished already! check your output file")
else ab_start = CheckPointValues$ab + 1
} else sample_start = CheckPointValues$sample + 1
} else mc_start = mc_start + 1
list(ab_start = ab_start, sample_start = sample_start, mc_start = mc_start)
}
# No export.
deleteCheckPoint <- function(ParameterMatrix, TemperedType, Estimfct, nab, npar,
nSS, MCparam,
eps,
algo,
regularization,
WeightingMatrix,
t_scheme,
alphaReg,
t_free,
subdivisions,
IntegrationMethod,
randomIntegrationLaw,
s_min,
s_max,
ncond,
IterationControl,
...) {
method <- Estim_Des_Temp(TemperedType, Estimfct,
eps = eps,
algo = algo,
regularization = regularization,
WeightingMatrix =
WeightingMatrix,
t_scheme = t_scheme,
alphaReg = alphaReg,
t_free = t_free,
subdivisions = subdivisions,
IntegrationMethod =
IntegrationMethod,
randomIntegrationLaw =
randomIntegrationLaw,
s_min = s_min,
s_max = s_max,
ncond = ncond,
IterationControl = IterationControl,
...)
##This code can be used to create a checkpoint text file.
# fileName <- get_filename_checkPoint_Temp(ParameterMatrix, nab, npar,
# MCparam, method)
# unlink(x = fileName, force = TRUE)
}
# No export.
writeCheckPoint <- function(ParameterMatrix, TemperedType, Estimfct, ab, nab,
npar, sample, nSS, mc, MCparam, ...) {
method <- Estim_Des_Temp(TemperedType, Estimfct, ...)
fileName <- get_filename_checkPoint_Temp(ParameterMatrix, nab, npar,
MCparam, method)
line = readLines(fileName, -1)
line[2] = paste("--", ab, nab, npar, sample, nSS, mc, MCparam, sep = ";")
writeLines(line, fileName)
}
#Added by Cedric 20220726
# Currently not necessary. Not even adatpted
#NameStatOutput <- function(FctsToApply, StatOutput) {
# Names <- c("alpha", "beta", "n", names(FctsToApply), "failure", "time")
# lapply(StatOutput, function(x) {
# colnames(x) <- Names
# return(x)
# })
#}
#Added by Cedric 20220729
# No export.
get_filename <- function(thetaT, MCparam, TemperedType, method,
extension = ".csv") {
MC <- switch(TemperedType,
CTS = paste(paste("Alpha=", thetaT[1]),
paste("DeltaP=", thetaT[2]),
paste("DeltaM=", thetaT[3]),
paste("LambdaP=", thetaT[4]),
paste("LambdaM=", thetaT[5]),
paste("mu=", thetaT[6]),
"MCparam", MCparam, sep = "_"),
TSS = paste(paste("Alpha=", thetaT[1]),
paste("Delta=", thetaT[2]),
paste("Lambda=", thetaT[3]),
"MCparam", MCparam, sep = "_"),
NTS = paste(paste("Alpha=", thetaT[1]),
paste("Beta=", thetaT[2]),
paste("Delta=", thetaT[3]),
paste("Lambda=", thetaT[4]),
paste("mu=", thetaT[5]),
"MCparam", MCparam, sep = "_"),
MTS = paste(paste("Alpha=", thetaT[1]),
paste("Delta=", thetaT[2]),
paste("LambdaP=", thetaT[3]),
paste("LambdaM=", thetaT[4]),
paste("mu=", thetaT[5]),
"MCparam", MCparam, sep = "_"),
GTS = paste(paste("AlphaP=", thetaT[1]),
paste("AlphaM=", thetaT[2]),
paste("DeltaP=", thetaT[3]),
paste("DeltaM=", thetaT[4]),
paste("LambdaP=", thetaT[5]),
paste("LambdaM=", thetaT[6]),
paste("mu=", thetaT[7]),
"MCparam", MCparam, sep = "_"),
KRTS = paste(paste("Alpha=", thetaT[1]),
paste("kP=", thetaT[2]),
paste("kM=", thetaT[3]),
paste("rP=", thetaT[4]),
paste("rM=", thetaT[5]),
paste("pP=", thetaT[6]),
paste("pM=", thetaT[7]),
paste("mu=", thetaT[8]),
"MCparam", MCparam, sep = "_"),
RDTS = paste(paste("Alpha=", thetaT[1]),
paste("Delta=", thetaT[2]),
paste("LambdaP=", thetaT[3]),
paste("LambdaM=", thetaT[4]),
paste("mu=", thetaT[5]),
"MCparam", MCparam, sep = "_")
# ,CGMY = paste(paste("C=", thetaT[1]),
# paste("G=", thetaT[2]),
# paste("M=", thetaT[3]),
# paste("Y=", thetaT[4]),
# "MCparam", MCparam, sep = "_")
)
fileName <- paste(MC, method, extension, sep = "")
fileName
}
#Added by Cedric 20220805
# No export.
updateOutputFile <- function(thetaT, MCparam, TemperedType, Output){
method <- Output$file
fileName <- get_filename(thetaT, MCparam, TemperedType, method)
if (!file.exists(fileName)) {
x <- switch(TemperedType,
CTS = paste("alphaT", "delta+T", "delta-T", "lambda+T",
"lambda-T", "muT", "data size", "seed",
"alphaE", "delta+E", "delta-E", "lambda+E",
"lambda-E", "muE", "failure", "time",
sep = ","),
TSS = paste("alphaT", "deltaT", "lambdaT",
"data size", "seed", "alphaE", "deltaE",
"lambdaE", "failure", "time",
sep = ","),
NTS = paste("alphaT", "betaT", "deltaT", "lambdaT", "muT",
"data size", "seed", "alphaE", "betaE", "deltaE",
"lambdaE", "muE", "failure", "time",
sep = ","),
MTS = paste("alphaT", "delta", "lambda+T", "lambda-T", "muT",
"data size", "seed", "alphaE", "deltaE", "lambda+E",
"lambda-E", "muE", "failure", "time", sep = ","),
GTS = paste("alpha+T", "alpha-T", "delta+T", "delta-T",
"lambda+T", "lambda-T", "muT", "data size", "seed",
"alpha+E", "alpha-E", "delta+E", "delta-E",
"lambda+E", "lambda-E", "muE", "failure", "time",
sep = ","),
KRTS = paste("alphaT", "k+T", "k-T", "r+T", "r-T",
"p+T", "p-T", "muT", "data size", "seed",
"alphaE", "k+E", "k-E", "r+E", "r-E",
"p+E", "p-E", "muE", "failure", "time", sep = ","),
RDTS = paste("alphaT", "delta", "lambda+T", "lambda-T", "muT",
"data size", "seed", "alphaE", "deltaE",
"lambda+E", "lambda-E", "muE", "failure", "time",
sep = ",")
# ,CGMY = paste("CT", "GT", "MT", "YT", "data size", "seed", "CE",
# "GE", "ME", "YE","failure", "time",
# sep = ",")
)
write(x, file = fileName, sep = "\n")
}
write(x = paste(as.character(Output$outputMat), collapse=","),
file = fileName, sep="\n", append=TRUE)
}
#Added by Cedric 20221011
# No export.
getTime_ <- function() proc.time()[3]
Try the TempStable package in your browser
Any scripts or data that you put into this service are public.
TempStable documentation built on Oct. 24, 2023, 5:06 p.m.
R Package Documentation
Browse R Packages
We want your feedback!
Note that we can't provide technical support on individual packages. You should contact the package authors for that.
Defines functions getTime_ updateOutputFile get_filename writeCheckPoint deleteCheckPoint updateCheckPointValues get_filename_checkPoint_Temp readCheckPoint Estim_Des_Temp initOutputFile getTempEstimation ComputeMCSimForTempered getSeedVector parallelizeMCsimulation TemperedEstim_Simulation
Documented in parallelizeMCsimulation TemperedEstim_Simulation
#' Monte Carlo Simulation
#'
#' @description
#' Runs Monte Carlo simulation for a selected estimation method. The function
#' can save results in a file.
#'
#' @details
#' \strong{TemperedTyp} With the parameter 'TemperedTyp' you can choose the
#' tempered stable distribution you want to use. Here is a list of distribution
#' you can choose from:
#' \describe{
#' \item{TSS}{Tempered stabel subordinator: See [charTSS()] for details.}
#' \item{CTS}{Classical tempered stable distribution: See [charCTS()] for
#' details.}
#' \item{GTS}{Generalized classical tempered stable distribution: See
#' [charGTS()] for details.}
#' \item{NTS}{Normal tempered stable distribution: See [charNTS()] for
#' details.}
#' \item{MTS}{Modified tempered stable distribution: See [charMTS()] for
#' details.}
#' \item{RDTS}{Rapid decreasing tempered stable distribution: See [charRDTS()]
#' for details.}
#' \item{KRTS}{Kim-Rachev tempered stable distribution: See [charKRTS()] for
#' details.}
#' }
#'
#' \strong{Error Handling} It is advisable to set it to TRUE when user is
#' planning to launch long simulations as it will prevent the procedure to stop
#' if an error occurs for one sample data. The estimation function will produce
#' a vector of NA as estimated parameters related to this (error generating)
#' sample data and move on to the next Monte Carlo step.
#'
#' \strong{Output file} Setting \code{saveOutput} to \code{TRUE} will have the
#' side effect of saving a csv file in the working directory. This file will
#' have \code{MCparam*length(SampleSizes)} lines and its columns will be:
#' \describe{
#' \item{alphaT, ...:}{the true value of the parameters.}
#' \item{data size:}{the sample size used to generate the simulated data.}
#' \item{seed:}{the seed value used to generate the simulated data.}
#' \item{alphaE, ...:}{the estimate of the parameters.}
#' \item{failure:}{binary: 0 for success, 1 for failure.}
#' \item{time:}{estimation running time in seconds.}
#' }
#' The file name is informative to let the user identify the value of the true
#' parameters, the MC parameters as well as the options selected for the
#' estimation method. The csv file is updated after each MC estimation which is
#' useful when the simulation stops before it finishes.
#'
#' \strong{SeedOptions} If users does not want to control the seed generation,
#' they could ignore this argument (default value NULL). This argument can be
#' more useful when they wants to cut the simulation (even for one parameter
#' value) into pieces. In that case, they can control which part of the seed
#' vector they want to use.
#' \describe{
#' \item{MCtot:}{total values of MC simulations in the entire process.}
#' \item{seedStart:}{starting index in the seed vector. The vector extracted
#' will be of size MCparam.}
#' }
#'
#' \strong{Estimfct} Additional parameters are needed for different estimation
#' functions. These are listed below for each function. The list of additional
#' parameters starts after the parameter \code{eps} in the parameter list.
#' \describe{
#' \item{For ML:}{ See usage of Maximum likelihood estimation in Kim et al.
#' (2008).No additional parameters are needed.}
#' \item{For GMM:}{Generalized Method of Moments by Feuerverger (1981).
#' The parameters \code{algo, alphaReg, regularization, WeightingMatrix, and
#' t_scheme} must be specified.
#'
#' Parameter \code{t_scheme}: One of the most important features of this
#' method is that it allows the user to choose how to place the points where
#' the moment conditions are evaluated. One can choose among 6 different
#' options. Depending on the option, further parameters have to be passed.
#' \describe{
#' \item{"equally":}{equally placed points in \code{min_t,max_t}. When
#' provided, user's \code{min_t} and \code{max_t} will be used (when
#' \code{Coinstrained == FALSE}).
#' }
#' \item{"NonOptAr":}{non optimal arithmetic placement.
#' }
#' \item{"uniformOpt":}{uniform optimal placement.
#' }
#' \item{"ArithOpt":}{arithmetic optimal placement.
#' }
#' \item{"Var Opt":}{optimal variance placement as explained above.
#' }
#' \item{"free":}{user needs to pass own set of points in \code{t_free}.
#' }
#' }
#'
#' Parameter \code{WeightingMatrix}: One can choose among 3 different options:
#' \describe{
#' \item{"OptAsym":}{the optimal asymptotic choice.
#' }
#' \item{"DataVar":}{the covariance matrix of the data provided.
#' }
#' \item{"Id":}{the identity matrix.
#' }
#' }
#' }
#' \item{For Cgmm:}{Continuum Generalized Methods of Moments by Carrasco &
#' Kotchoni (2017). The parameters \code{algo, alphaReg, subdivisions,
#' IntegrationMethod, randomIntegrationLaw, s_min, and s_max} must be
#' specified.
#' }
#' \item{For GMC:}{Generalized Method of Cumulants (GMC) by Massing, T.
#' (2022). The parameters \code{algo, alphaReg, regularization,
#' WeightingMatrix, and ncond} must be specified.
#' }
#' }
#'
#' \strong{IterationControl} If \code{algo = "IT..."} or \code{algo =
#' "Cue..."} the user can control each iteration by setting up the list
#' IterationControl which contains the following elements:
#' \describe{
#' \item{NbIter}{maximum number of iteration. The loop stops when NBIter is
#' reached; default = 10.}
#' \item{PrintIterlogical}{if set to TRUE, the value of the current parameter
#' estimation is printed to the screen at each iteration; default = TRUE.}
#' \item{RelativeErrMax}{the loop stops if the relative error between two
#' consecutive estimation steps is smaller than RelativeErrMax;
#' default = 1e-3.}
#' }
#'
#' \strong{methodR} Random numbers must be generated for each MC study. For
#' each distribution, different methods are available for this (partly also
#' depending on alpha). For more information, the documentation of the
#' respective \code{r...()} distribution can be called up. By default, the fastest
#' method is selected. Since the deviation error can amplify to the edges of
#' alpha depending on the method, it is recommended to check the generated
#' random numbers once for each distribution using the density function before
#' starting the simulation.
#'
#' \strong{Parallelization} Parallelization of the function is possible with
#' using [parallelizeMCsimulation()]. If someone wants to parallelize the
#' function manually, the parameter \code{MCparam} must be set to \code{1} and
#' the parameter \code{SeedOption} must be changed for each iteration.
#'
#' Since this package is structurally based on the \strong{"StableEstim"
#' package by Tarak Kharrat and Georgi N. Boshnakov}, more detailed
#' documentation can be found in their documentation.
#'
#' @seealso
#' \url{https://github.com/GeoBosh/StableEstim/blob/master/R/Simulation.R}
#'
#' @references
#' Massing, T. (2023), 'Parametric Estimation of Tempered Stable Laws'
#'
#' Kim, Y. s.; Rachev, S. T.; Bianchi, M. L. & Fabozzi, F. J. (2008), 'Financial
#' market models with lévy processes and time-varying volatility'
#' \doi{10.1016/j.jbankfin.2007.11.004}
#'
#' Hansen, L. P. (1982), 'Large sample properties of generalized method of
#' moments estimators' \doi{10.2307/1912775}
#'
#' Hansen, L. P.; Heaton, J. & Yaron, A. (1996), 'Finite-Sample Properties of
#' Some Alternative GMM Estimators' \doi{10.1080/07350015.1996.10524656}
#'
#' Feuerverger, A. & McDunnough, P. (1981), 'On the efficiency of empirical
#' characteristic function procedures'
#' \doi{10.1111/j.2517-6161.1981.tb01143.x}
#'
#' Carrasco, M. & Kotchoni, R. (2017), 'Efficient estimation using the
#' characteristic function' \doi{10.1017/S0266466616000025};
#'
#' Kuechler, U. & Tappe, S. (2013), 'Tempered stable distribution and processes'
#' \doi{10.1016/j.spa.2013.06.012}
#'
#' @param ParameterMatrix The matrix is to be composed of vectors, row by row.
#' Each vector must fit the pattern of theta of the \code{TemperedType}.
#' @param SampleSizes Sample sizes to be used to simulate the data. By default,
#' we use 200 (small sample size) and 1600 (large sample size);
#' vector of integer.
#' @param MCparam Number of Monte Carlo simulation for each couple of parameter,
#' default=100; integer
#' @param TemperedType A String. Either "CTS", "TSS", "NTS", "MTS", "GTS",
#' "KRTS", "RDTS".
#' @param Estimfct The estimation function to be used. A String.
#' Either "ML", "GMM", "Cgmm", or "GMC".
#' @param HandleError Logical flag: if set to TRUE, the simulation doesn't stop
#' when an error in the estimation function is encountered. A vector of
#' (size 4) NA is saved and the the simulation carries on. See details.
#' @param saveOutput Logical flag: if set to TRUE, a csv file (for each couple
#' of parameter) with the the estimation
#' information is saved in the current directory. See details.
#' @param SeedOptions List to control the seed generation. See details.
#' @param eps Numerical error tolerance. \code{1e-06} by default.
#' @param algo algorithm: For GMM: \code{"2SGMM"} is the two step GMM proposed
#' by Hansen (1982). \code{"CueGMM"} and \code{"ITGMM"} are respectively the
#' continuous updated and the iterative GMM proposed by Hansen, Eaton et Yaron
#' (1996) and adapted to the continuum case. For GMC: \code{"2SGMC", "CueGMC"}.
#' For Cgmm: \code{"2SCgmm", "CueCgmm", ...}.
#' @param regularization regularization scheme to be used for moment methods,
#' one of \code{"Tikhonov"} (Tikhonov), \code{"LF"} (Landweber-Fridmann) and
#' \code{"cut-off"} (spectral cut-off).
#' @param WeightingMatrix type of weighting matrix used to compute the
#' objective function for the GMM and GMC methods, one of \code{"OptAsym"} (the
#' optimal asymptotic), \code{"DataVar"} (the data driven, only for GMM) and
#' \code{"Id"} (the identity matrix).
#' @param t_scheme scheme used to select the points for the GMM method where the
#' moment conditions are evaluated, one of \code{"equally"} (equally placed),
#' \code{"NonOptAr"} (non optimal arithmetic placement), \code{"uniformOpt"}
#' (uniform optimal placement), \code{"ArithOpt"} (arithmetic optimal placement)
#' , \code{"Var Opt"} (optimal variance placement) and \code{"free"} (users need
#' to pass their own set of points in ...).
#' @param alphaReg value of the regularisation parameter; numeric. Example Value
#' could be ==0.01.
#' @param t_free sequence, if \code{t_scheme=="free"}.
#' @param subdivisions Number of subdivisions used to compute the different
#' integrals involved in the computation of the objective function for the Cgmm
#' method (to minimise); numeric.
#' @param IntegrationMethod Numerical integration method to be used to
#' approximate the (vectorial) integrals for the Cgmm method. Users can choose
#' between "Uniform" discretization or the "Simpson"'s rule (the 3-point
#' Newton-Cotes quadrature rule).
#' @param randomIntegrationLaw Probability measure associated to the Hilbert
#' space spanned by the moment conditions for the Cgmm method.
#' @param s_min,s_max Lower and Upper bounds of the interval where the moment
#' conditions are considered for the Cgmm method; numeric.
#' @param ncond Integer. Number of moment conditions (until order \code{ncond})
#' for the GMC method. Must not be less than 3 for TSS, 6 for CTS, 5 for NTS.
#' @param IterationControl only used if algo = "IT..." or algo = "Cue..."
#' to control the iterations. See Details.
#' @param nb_t integer, if you set \code{t_scheme <- "equally"}. nb_t could be
#' == 20 for example.
#' @param methodR A string. Method generates random variates of TS distribution.
#' "TM" by default. Switches automatically if the method is not applicable in
#' this way.
#' @param ... Other arguments to be passed to the estimation function.
#'
#' @return If \code{saveOutput == FALSE}, the return object is a list of 2.
#' Results of the simulation are listed in \code{$outputMat}. If \code{
#' saveOutput == TRUE}, only a csv file is saved and nothing is returned.
#'
#' @examples
#' \donttest{
#' TemperedEstim_Simulation(ParameterMatrix = rbind(c(1.5,1,1,1,1,0),
#' c(0.5,1,1,1,1,0)),
#' SampleSizes = c(4), MCparam = 4,
#' TemperedType = "CTS", Estimfct = "ML",
#' saveOutput = FALSE)
#'
#' TemperedEstim_Simulation(ParameterMatrix = rbind(c(1.5,1,1,1,1,0)),
#' SampleSizes = c(4), MCparam = 4,
#' TemperedType = "CTS", Estimfct = "GMM",
#' saveOutput = FALSE, algo = "2SGMM",
#' regularization = "cut-off",
#' WeightingMatrix = "OptAsym", t_scheme = "free",
#' alphaReg = 0.01,
#' t_free = seq(0.1,2,length.out=12))
#'
#' TemperedEstim_Simulation(ParameterMatrix = rbind(c(1.45,0.55,1,1,1,0)),
#' SampleSizes = c(4), MCparam = 4,
#' TemperedType = "CTS", Estimfct = "Cgmm",
#' saveOutput = FALSE, algo = "2SCgmm",
#' alphaReg = 0.01, subdivisions = 50,
#' IntegrationMethod = "Uniform",
#' randomIntegrationLaw = "unif",
#' s_min = 0, s_max= 1)
#'
#' TemperedEstim_Simulation(ParameterMatrix = rbind(c(1.45,0.55,1,1,1,0)),
#' SampleSizes = c(4), MCparam = 4,
#' TemperedType = "CTS", Estimfct = "GMC",
#' saveOutput = FALSE, algo = "2SGMC",
#' alphaReg = 0.01, WeightingMatrix = "OptAsym",
#' regularization = "cut-off", ncond = 8)
#' }
#'
#' @export
TemperedEstim_Simulation <- function(ParameterMatrix,
SampleSizes = c(200, 1600),
MCparam = 100,
TemperedType = c("CTS", "TSS",
"NTS", "MTS", "GTS",
"KRTS", "RDTS"),
Estimfct = c("ML", "GMM", "Cgmm", "GMC"),
HandleError = TRUE, saveOutput = FALSE,
SeedOptions = NULL, eps = 1e-06,
algo = NULL, regularization = NULL,
WeightingMatrix = NULL, t_scheme = NULL,
alphaReg = NULL, t_free = NULL,
nb_t = NULL, subdivisions = NULL,
IntegrationMethod = NULL,
randomIntegrationLaw = NULL, s_min = NULL,
s_max = NULL, ncond = NULL,
IterationControl = NULL,
methodR = "TM", ...) {
#seeAlso: https://github.com/GeoBosh/StableEstim/blob/master/R/Simulation.R
SeedVector <- getSeedVector(MCparam, SeedOptions)
Estimfct <- match.arg(arg = Estimfct,
choices = c("ML", "GMM", "Cgmm", "GMC"))
TemperedType <- match.arg(arg = TemperedType, choices =
c("CTS", "TSS","NTS", "MTS", "GTS",
"KRTS", "RDTS"))
nab <- nrow(ParameterMatrix)
npar <- ncol(ParameterMatrix)
lS <- length(SampleSizes)
nRowOutput <- nab * lS
OutputCollection <- empty_list <- vector(mode = "list", length = nab)
returnList <- empty_list <- vector(mode = "list")
#returnList <- matrix(data = NA, ncol = npar, nrow = nab*MCparam)
indexStatOutput <- 1
#Values that needs to be set in the beginning as checkpoint files are not
# allowed
ab <- 1
sample <- 1
mc <- 0
if(MCparam != 1){
CheckPointValues <- readCheckPoint(ParameterMatrix, TemperedType,
Estimfct, nab, npar, lS, MCparam,
eps = eps,
algo = algo,
regularization = regularization,
WeightingMatrix =
WeightingMatrix,
t_scheme = t_scheme,
alphaReg = alphaReg,
t_free = t_free,
nb_t = nb_t,
subdivisions = subdivisions,
IntegrationMethod =
IntegrationMethod,
randomIntegrationLaw =
randomIntegrationLaw,
s_min = s_min,
s_max = s_max,
ncond = ncond,
IterationControl = IterationControl,
ab = ab, sample = sample, mc = mc,
...)
}
else{
CheckPointValues <- list(ab = ab, nab = nab, npar = npar, sample = sample,
nSS = lS, mc = mc, MCparam = 1)
}
updatedCheckPointValues <- updateCheckPointValues(CheckPointValues, MCparam,
lS, nab)
# Ist im Zielordner bereits eine csv-Datei mit dem gleichen Namen (thetaT
# und MCparam sind gleich), wird die Datei aktuell um die weiteren
# Ergebnisse aktualisiert. Mention in Details
# Wird dieser Code aktiviert, würde das Überschreiben einer Datei verboten
# werden.
#
#if (updatedCheckPointValues$mc_start != 1){
# print("'Can't Compute Stat summary when the process doesn't
# start from the beginning!!")
#}
for (ab in updatedCheckPointValues$ab_start:nab) {
thetaT <- ParameterMatrix[ab, ]
outputString <- switch(TemperedType,
CTS = paste("Alpha=", thetaT[1] ,
" *** DeltaP=", thetaT[2],
" *** DeltaM=", thetaT[3],
" *** LambdaP=", thetaT[4],
" *** LambdaM=", thetaT[5],
" *** mu=", thetaT[6], sep = ""),
TSS = paste("Alpha=", thetaT[1] ,
" *** Delta=", thetaT[2],
" *** Lambda=", thetaT[3], sep = ""),
NTS = paste("Alpha=", thetaT[1] ,
" *** Beta=", thetaT[2],
" *** Delta=", thetaT[3],
" *** Lambda=", thetaT[4],
" *** mu=", thetaT[5], sep = ""),
MTS = paste("Alpha=", thetaT[1] ,
" *** Delta=", thetaT[2],
" *** LambdaP=", thetaT[3],
" *** LambdaM=", thetaT[4],
" *** mu=", thetaT[5], sep = ""),
GTS = paste("AlphaP=", thetaT[1] ,
" *** AlphaM=", thetaT[2],
" *** DeltaP=", thetaT[3],
" *** DeltaM=", thetaT[4],
" *** LambdaP=", thetaT[5],
" *** LambdaM=", thetaT[6],
" *** mu=", thetaT[7], sep = ""),
KRTS = paste("Alpha=", thetaT[1] ,
" *** kP=", thetaT[2],
" *** kM=", thetaT[3],
" *** rP=", thetaT[4],
" *** rM=", thetaT[5],
" *** pP=", thetaT[6],
" *** pM=", thetaT[7],
" *** mu=", thetaT[8], sep = ""),
RDTS = paste("Alpha=", thetaT[1] ,
" *** Delta=", thetaT[2],
" *** LambdaP=", thetaT[3],
" *** LambdaM=", thetaT[4],
" *** mu=", thetaT[5], sep = "")
# ,CGMY = paste("C=", thetaT[1] ,
# " *** G=", thetaT[2],
# " *** M=", thetaT[3],
# " *** Y=", thetaT[4], sep = "")
)
cat("---------------- ", outputString, " --------------- \n", sep = "")
if (saveOutput) initOutputFile(thetaT, MCparam, TemperedType,
Estimfct, ...)
EstimOutput <- ComputeMCSimForTempered(thetaT = thetaT,
MCparam = MCparam,
SampleSizes =
as.vector(SampleSizes),
SeedVector = SeedVector,
TemperedType = TemperedType,
Estimfct = Estimfct,
HandleError = HandleError,
ab_current = ab,
nab = nab,
npar = npar, ParameterMatrix,
CheckPointValues =
updatedCheckPointValues,
saveOutput = saveOutput,
eps = eps,
algo = algo,
regularization = regularization,
WeightingMatrix =
WeightingMatrix,
t_scheme = t_scheme,
alphaReg = alphaReg,
t_free = t_free,
nb_t = nb_t,
subdivisions = subdivisions,
IntegrationMethod =
IntegrationMethod,
randomIntegrationLaw =
randomIntegrationLaw,
s_min = s_min,
s_max = s_max,
ncond = ncond,
IterationControl =
IterationControl,
methodR = methodR,
...)
OutputCollection <- EstimOutput
if (saveOutput == FALSE){
if (length(returnList) == 0) returnList <- EstimOutput
else returnList <- Map(f = rbind, x = returnList, init = EstimOutput)
}
}
# if(MCparam != 1){
# deleteCheckPoint(ParameterMatrix, TemperedType, Estimfct, nab, npar, lS,
# MCparam,
# eps = eps,
# algo = algo,
# regularization = regularization,
# WeightingMatrix =
# WeightingMatrix,
# t_scheme = t_scheme,
# alphaReg = alphaReg,
# t_free = t_free,
# subdivisions = subdivisions,
# IntegrationMethod =
# IntegrationMethod,
# randomIntegrationLaw =
# randomIntegrationLaw,
# s_min = s_min,
# s_max = s_max,
# ncond = ncond,
# IterationControl = IterationControl,
# methodR = methodR,
# ...)
# }
if (saveOutput == FALSE){
return(returnList)
}
}
#' Function to parallelize the Monte Carlo Simulation
#'
#' Since the Monte Carlo Simulation is very computationally intensive, it may
#' be worthwhile to split it across all available processor cores. To do this,
#' simply pass all the parameters from the [TemperedEstim_Simulation()]
#' function to this function in the same way.
#'
#' In this function exactly the arguments must be passed, which are also needed
#' for the function [TemperedEstim_Simulation()]. However, a few functions of
#' [TemperedEstim_Simulation()] are not possible here. The restrictions are
#' described in more detail for the individual arguments.
#'
#' In addition to the arguments of function [TemperedEstim_Simulation()], the
#' argument "cores" can be assigned an integer value. This value determines how
#' many different processes are to be parallelized. If value is \code{NULL}, R
#' tries to read out how many cores the processor has and passes this
#' value to "cores".
#'
#' During the simulation, the progress of the simulation can be viewed in a
#' file in the workspace named "IterationControlForParallelization.txt".
#'
#' @param ParameterMatrix The matrix is to be composed of vectors, row by row.
#' Each vector must fit the pattern of theta of the \code{TemperedType}.
#' Compared to the function [TemperedEstim_Simulation()], the matrix here may
#' contain only one parameter vector.
#' @param MCparam Number of Monte Carlo simulation for each couple of parameter,
#' default=100; integer
#' @param SampleSizes Sample sizes to be used to simulate the data. By default,
#' we use 200 (small sample size). Vector of integer. Compared to the function
#' [TemperedEstim_Simulation()], the vector here may contain only one integer.
#' @param saveOutput Logical flag: In the function [TemperedEstim_Simulation()]
#' the argument can be true. Then an external csv file is created. Here the
#' argument must be false. The output of the values works in this function
#' exclusively via the return of the function.
#' @param SeedOptions is an argument what can be used in
#' [TemperedEstim_Simulation()] but must be NULL here.
#' @param cores size of cluster for parallelization. Positive Integer.
#' @param iterationDisplayToFileSystem creates a text file in your file system
#' that displays the current iteration of the simulation.
#' @param ... The function works only if all necessary arguments from the
#' function [TemperedEstim_Simulation()] are passed. See description and
#' details.
#'
#' @return The return object is a list of 2. Results of the simulation are
#' listed in \code{$outputMat}.
#'
#' @export
#' @importFrom foreach %dopar%
parallelizeMCsimulation <- function(
ParameterMatrix,
MCparam = 10000,
SampleSizes = c(200),
saveOutput = FALSE,
cores = 2,
SeedOptions = NULL,
iterationDisplayToFileSystem = FALSE,
...){
mc <- NULL
if (!is.null(SeedOptions)){
stop("SeedOptions is used by the function and cannot be passed as an
argument")
}
if(is.null(ParameterMatrix) || nrow(ParameterMatrix) != 1){
stop("Compared to the function TemperedEstim_Simulation(), the matrix here
may contain only one parameter vector.")
}
if(is.null(SampleSizes) || length(SampleSizes) != 1){
stop("Compared to the function TemperedEstim_Simulation(), the SampleSizes
vector here may contain only one integer.")
}
if(is.null(saveOutput) || saveOutput){
stop("Compared to the function TemperedEstim_Simulation(), saveOutput
must be FALSE.")
}
if (is.null(cores) || cores < 1){
cores <- parallel::detectCores()
}
if(MCparam < cores){
cores <- MCparam
}
R <- MCparam #MonteCarloRuns
cl <- parallel::makeCluster(cores)
doParallel::registerDoParallel(cl)
resultOfSimulation <- foreach::foreach(mc = 1:R, .combine = rbind,
.export = ls("package:TempStable",
all.names = TRUE)
)%dopar%{
returnValue <- TemperedEstim_Simulation(
ParameterMatrix = ParameterMatrix,
MCparam = 1,
SampleSizes = SampleSizes,
saveOutput = saveOutput,
SeedOptions = list(MCtot = R, seedStart = mc),
... = ...
)
if (iterationDisplayToFileSystem){
utils::write.table(
paste("Last Monte Carlo run: ", mc) ,
file = base::paste("IterationControlForParallelization.txt"),
sep = "\t", row.names = FALSE)
}
returnValue$outputMat
}
parallel::stopCluster(cl)
attr(resultOfSimulation, "rng") <- NULL
attr(resultOfSimulation, "doRNG_version") <- NULL
#Delete txt file
if (iterationDisplayToFileSystem){
base::unlink(x = base::paste("IterationControlForParallelization.txt"),
force = TRUE)
}
return(resultOfSimulation)
}
# No export.
getSeedVector <- function(Outputsize, SeedOptions = NULL) {
set.seed(345)
if (is.null(SeedOptions))
vec <- as.vector(sample.int(n = 3 * Outputsize, size = Outputsize))
else {
MCtot <- SeedOptions$MCtot
seedStart <- SeedOptions$seedStart
seedEnd <- seedStart + Outputsize
vec <-
as.vector(sample.int(n = 3 * MCtot, size = MCtot))[seedStart:seedEnd]
}
vec
}
# No export.
ComputeMCSimForTempered <- function(thetaT, MCparam, SampleSizes, SeedVector,
TemperedType, Estimfct, HandleError,
ab_current, nab, npar, ParameterMatrix,
CheckPointValues = NULL, saveOutput, eps,
algo, regularization, WeightingMatrix,
t_scheme, alphaReg, t_free,nb_t,
subdivisions, IntegrationMethod,
randomIntegrationLaw, s_min, s_max, ncond,
IterationControl, methodR, ...)
{
if (TemperedType == "CTS") {
Ncol <- 16
} else if (TemperedType == "TSS") {
Ncol <- 10
} else if (TemperedType == "NTS" || TemperedType == "MTS" ||
TemperedType == "RDTS") {
Ncol <- 14
} else if (TemperedType == "GTS") {
Ncol <- 18
} else if (TemperedType == "KRTS") {
Ncol <- 20
}
# else {
# Ncol <- 12
# }
nSS <- length(SampleSizes)
Nrow <- nSS * MCparam
Output <- matrix(data = NA, ncol = Ncol, nrow = Nrow)
if (TemperedType == "CTS") {
colnames(Output) <- c("alphaT", "delta+T", "delta-T", "lambda+T",
"lambda-T", "muT", "data size", "seed", "alphaE",
"delta+E", "delta-E", "lambda+E", "lambda-E",
"muE", "failure", "time")
} else if (TemperedType == "TSS") {
colnames(Output) <- c("alphaT", "deltaT", "lambdaT", "data size",
"seed", "alphaE", "deltaE", "lambdaE",
"failure", "time")
} else if (TemperedType == "NTS") {
colnames(Output) <- c("alphaT", "betaT", "deltaT", "lambdaT", "muT",
"data size", "seed", "alphaE", "betaE", "deltaE",
"lambdaE", "muE", "failure", "time")
} else if (TemperedType == "MTS") {
colnames(Output) <- c("alphaT", "delta", "lambda+T", "lambda-T", "muT",
"data size", "seed", "alphaE", "deltaE", "lambda+E",
"lambda-E", "muE", "failure", "time")
} else if (TemperedType == "GTS") {
colnames(Output) <- c("alpha+T", "alpha-T", "delta+T", "delta-T",
"lambda+T", "lambda-T", "muT", "data size", "seed",
"alpha+E", "alpha-E", "delta+E", "delta-E",
"lambda+E", "lambda-E", "muE", "failure", "time")
} else if (TemperedType == "KRTS") {
colnames(Output) <- c("alphaT", "k+T", "k-T", "r+T", "r-T",
"p+T", "p-T", "muT", "data size", "seed",
"alphaE", "k+E", "k-E", "r+E", "r-E",
"p+E", "p-E", "muE", "failure", "time")
} else if (TemperedType == "RDTS") {
colnames(Output) <- c("alphaT", "delta", "lambda+T", "lambda-T", "muT",
"data size", "seed", "alphaE", "deltaE", "lambda+E",
"lambda-E", "muE", "failure", "time")
}
# else {
# colnames(Output) <- c("C.T", "G.T", "M.T", "Y.T", "data size", "seed",
# "C.E", "G.E", "M.E", "Y.E", "failure", "time")
# }
if (ab_current == CheckPointValues$ab_start) {
sample_start = CheckPointValues$sample_start
mc_start = CheckPointValues$mc_start
} else {
sample_start = 1
mc_start = 1
}
for (sample in sample_start:nSS) {
size <- SampleSizes[sample]
if (sample != sample_start) mc_start = 1
for (mc in mc_start:MCparam) {
tIter <- getTime_()
iter <- mc + (sample - 1) * MCparam
set.seed(seed <- SeedVector[mc])
if (TemperedType == "CTS") {
x <- rCTS(n = size, alpha = thetaT[1], deltap = thetaT[2],
deltam = thetaT[3], lambdap = thetaT[4],
lambdam = thetaT[5], mu = thetaT[6], methodR = methodR, ...)
} else if (TemperedType == "TSS") {
x <- rTSS(n = size, alpha = thetaT[1], delta = thetaT[2],
lambda = thetaT[3], methodR = methodR, ...)
} else if (TemperedType == "NTS") {
x <- rNTS(n = size, alpha = thetaT[1], beta = thetaT[2],
delta = thetaT[3], lambda = thetaT[4], mu = thetaT[5],
methodR = methodR, ...)
} else if (TemperedType == "MTS") {
x <- rMTS(n = size, theta = thetaT, methodR = methodR, ...)
} else if (TemperedType == "GTS") {
x <- rGTS(n = size, theta = thetaT, methodR = methodR, ...)
} else if (TemperedType == "KRTS") {
x <- rKRTS(n = size, alpha = thetaT[1], kp = thetaT[2], km = thetaT[3],
rp = thetaT[4], rm = thetaT[5], pp = thetaT[6],
pm = thetaT[7], mu = thetaT[8], methodR = methodR, ...)
} else if (TemperedType == "RDTS") {
x <- rRDTS(n = size, theta = thetaT, methodR = methodR, ...)
}
# else {
# x <- rCGMY(n = size, C = thetaT[1], M = thetaT[2],
# G = thetaT[3], Y = thetaT[4])
# }
Estim <- getTempEstimation(thetaT = thetaT, x = x, seed = seed,
size = size, Ncol = Ncol,
TemperedType = TemperedType,
Estimfct = Estimfct,
HandleError = HandleError, eps = eps,
algo = algo, regularization = regularization,
WeightingMatrix =
WeightingMatrix,
t_scheme = t_scheme,
alphaReg = alphaReg,
t_free = t_free,
nb_t = nb_t,
subdivisions = subdivisions,
IntegrationMethod =
IntegrationMethod,
randomIntegrationLaw =
randomIntegrationLaw,
s_min = s_min,
s_max = s_max,
ncond = ncond,
IterationControl = IterationControl,
...)
Output[iter, ] <- Estim$outputMat
file <- Estim$file
#When checkpoint file should be availabe again
# if (!is.null(CheckPointValues) && MCparam != 1) {
# writeCheckPoint(ParameterMatrix, TemperedType, Estimfct, ab_current,
# nab, npar, sample, nSS, mc,
# MCparam,
# eps,
# algo,
# regularization,
# WeightingMatrix,
# t_scheme,
# alphaReg,
# t_free,
# subdivisions,
# IntegrationMethod,
# randomIntegrationLaw,
# s_min,
# s_max,
# ncond,
# IterationControl,
# ...)
# }
if (saveOutput) updateOutputFile(thetaT, MCparam, TemperedType,
Estim)
StableEstim::PrintEstimatedRemainingTime(iter, tIter, Nrow)
}
}
#End Sample
if (isFALSE(saveOutput)){
return(list(outputMat = Output))
}
else{
return(list(outputMat = Output, file = file))
}
}
# No export.
getTempEstimation <- function(thetaT, x, seed, size, Ncol, TemperedType,
Estimfct, HandleError, eps,
algo,
regularization,
WeightingMatrix,
t_scheme,
alphaReg,
t_free,
nb_t,
subdivisions,
IntegrationMethod,
randomIntegrationLaw,
s_min,
s_max,
ncond,
IterationControl,
...) {
output <- vector(length = Ncol)
if (TemperedType == "CTS") {
output[1:8] <- c(thetaT, size, seed)
} else if (TemperedType == "TSS") {
output[1:5] <- c(thetaT, size, seed)
} else if (TemperedType == "NTS") {
output[1:7] <- c(thetaT, size, seed)
} else if (TemperedType == "MTS") {
output[1:7] <- c(thetaT, size, seed)
} else if (TemperedType == "GTS") {
output[1:9] <- c(thetaT, size, seed)
} else if (TemperedType == "KRTS") {
output[1:10] <- c(thetaT, size, seed)
} else if (TemperedType == "RDTS") {
output[1:7] <- c(thetaT, size, seed)
}
# else { #CGMY
# output[1:6] <- c(thetaT, size, seed)
# }
theta0 <- thetaT - 0.1 #noise
EstimRes <- TemperedEstim(TemperedType = TemperedType,
EstimMethod = Estimfct,
data = x, theta0 = theta0, ComputeCov = FALSE,
HandleError = HandleError, eps = eps,
algo = algo, regularization = regularization,
WeightingMatrix =
WeightingMatrix,
t_scheme = t_scheme,
alphaReg = alphaReg,
t_free = t_free,
nb_t = nb_t,
subdivisions = subdivisions,
IntegrationMethod =
IntegrationMethod,
randomIntegrationLaw =
randomIntegrationLaw,
s_min = s_min,
s_max = s_max,
ncond = ncond,
IterationControl = IterationControl, ...)
if (TemperedType == "CTS") {
output[9:14] <- EstimRes@par
} else if (TemperedType == "TSS") {
output[6:8] <- EstimRes@par
} else if (TemperedType == "NTS") {
output[8:12] <- EstimRes@par
} else if (TemperedType == "MTS") {
output[8:12] <- EstimRes@par
} else if (TemperedType == "GTS") {
output[10:16] <- EstimRes@par
} else if (TemperedType == "KRTS") {
output[11:18] <- EstimRes@par
} else if (TemperedType == "RDTS") {
output[8:12] <- EstimRes@par
}
# else {
# output[7:10] <- EstimRes@par
# }
if (TemperedType == "CTS") {
output[15:16] <- c(EstimRes@failure, EstimRes@duration)
} else if (TemperedType == "TSS") {
output[9:10] <- c(EstimRes@failure, EstimRes@duration)
} else if (TemperedType == "NTS") {
output[13:14] <- c(EstimRes@failure, EstimRes@duration)
} else if (TemperedType == "MTS") {
output[13:14] <- c(EstimRes@failure, EstimRes@duration)
} else if (TemperedType == "GTS") {
output[17:18] <- c(EstimRes@failure, EstimRes@duration)
} else if (TemperedType == "KRTS") {
output[19:20] <- c(EstimRes@failure, EstimRes@duration)
} else if (TemperedType == "RDTS") {
output[13:14] <- c(EstimRes@failure, EstimRes@duration)
}
# else {
# output[11:12] <- c(EstimRes@failure, EstimRes@duration)
# }
list(outputMat = output, file = EstimRes@method)
}
# Function title
# Merge with EstimSimulation
#
# @examples
# ComputeMCSimForTempered_parallel(1,c(1.5, 1, 1, 1, 1, 0),10,"CTS","ML")
# ComputeMCSimForTempered_parallel(1,c(0.5, 1, 1),10,"TSS","Cgmm",
# IntegrationMethod = "Simpson",
# randomIntegrationLaw = "unif")
#
# No export.
# ComputeMCSimForTempered_parallel <- function(MCparam, thetaT, size,
# TemperedType = c("CTS",
# "TSS",
# "NTS","CGMY"),
# Estimfct = c("ML", "GMM", "Cgmm",
# "GMC"),
# HandleError = TRUE, eps, ...) {
# Estimfct <- match.arg(Estimfct)
# TemperedType <- match.arg(TemperedType)
# Ncol <- ifelse(TemperedType == "CTS", 15, 9)
# if (TemperedType == "CTS") {
# Ncol <- 15
# } else if (TemperedType == "TSS") {
# Ncol <- 9
# } else if (TemperedType == "NTS") {
# Ncol <- 13
# } else {
# Ncol <- 11
# }
# Output <- numeric(Ncol)
# if (TemperedType == "CTS") {
# names(Output) <- c("alphaT", "delta+T", "delta-T", "lambda+T",
# "lambda-T", "muT", "data size", "alphaE", "delta+E",
# "delta-E", "lambda+E", "lambda-E", "muE", "failure",
# "time")
# } else if (TemperedType == "TSS") {
# names(Output) <- c("alphaT", "deltaT", "lambdaT", "data size", "alphaE",
# "deltaE", "lambdaE", "failure", "time")
# } else if (TemperedType == "NTS") {
# names(Output) <- c("alphaT", "betaT", "deltaT", "lambdaT", "muT",
# "data size", "alphaE", "betaE", "deltaE", "lambdaE",
# "muE", "failure", "time")
# } else {
# names(Output) <- c("C.T", "G.T", "M.T", "Y.T", "data size", "C.E",
# "G.E", "M.E", "Y.E", "failure", "time")
# }
#
# if (TemperedType == "CTS") {
# x <- rCTS(n = size, alpha = thetaT[1], deltap = thetaT[2],
# deltam = thetaT[3], lambdap = thetaT[4], lambdam = thetaT[5],
# mu = thetaT[6])
# } else if (TemperedType == "TSS") {
# x <- rTSS(n = size, alpha = thetaT[1], delta = thetaT[2],
# lambda = thetaT[3])
# } else if (TemperedType == "NTS") {
# x <- rNTS(n = size, alpha = thetaT[1], beta = thetaT[2],
# delta = thetaT[3], lambda = thetaT[4], mu = thetaT[5])
# } else {
# x <- rCGMY(n = size, C = theta[1], G = theta[2], G = theta[3],
# Y = theta[4])
# }
# Estim <- getTempEstimation_parallel(thetaT = thetaT, x = x, size = size,
# Ncol = Ncol,
# TemperedType = TemperedType,
# Estimfct = Estimfct,
# HandleError = HandleError,
# eps = eps, ...)
# Output <- c(MCparam, Estim)
#
# return(Output)
# }
# No export.
# getTempEstimation_parallel <- function(thetaT, x, size, Ncol, TemperedType,
# Estimfct, HandleError, eps, ...) {
# output <- vector(length = Ncol)
# if (TemperedType == "CTS") {
# output[1:7] <- c(thetaT, size)
# } else if (TemperedType == "TSS") {
# output[1:4] <- c(thetaT, size)
# } else if (TemperedType == "NTS") {
# output[1:6] <- c(thetaT, size)
# } else {
# output[1:5] <- c(thetaT, size)
# }
# theta0 <- thetaT - 0.1 #noise
# EstimRes <- TemperedEstim_v2(TemperedType = TemperedType,
# EstimMethod = Estimfct, data = x,
# theta0 = theta0, ComputeCov = FALSE,
# HandleError = HandleError, eps = eps, ...)
# if (TemperedType == "CTS") {
# output[8:13] <- EstimRes$par
# } else if (TemperedType == "TSS") {
# output[5:7] <- EstimRes$par
# } else if (TemperedType == "NTS") {
# output[7:11] <- EstimRes$par
# } else {
# output[6:9] <- EstimRes$par
# }
# if (TemperedType == "CTS") {
# output[14:15] <- c(EstimRes$failure, EstimRes$duration)
# } else if (TemperedType == "TSS") {
# output[8:9] <- c(EstimRes$failure, EstimRes$duration)
# } else if (TemperedType == "NTS") {
# output[12:13] <- c(EstimRes$failure, EstimRes$duration)
# } else {
# output[10:11] <- c(EstimRes$failure, EstimRes$duration)
# }
# return(output)
# }
##### for statistical summary#####
# Function title
#
# Gap holder for description.
#
# Gap holder for details.
#
# @param EstimOutput A gap holder.
# @param FctsToApply A gap holder.
# @param SampleSizes A gap holder.
# @param CheckMat A gap holder.
# @param tolFailCheck A gap holder.
# @param MCparam A gap holder.
#
# @return Gap holder for return.
#
# @export
# ComputeStatOutput <- function(EstimOutput, FctsToApply, SampleSizes, CheckMat,
# tolFailCheck, MCparam, ...) {
# list(alpha = ComputeStatOutputPar(EstimOutput = EstimOutput,
# FctsToApply = FctsToApply, par = "alpha",
# SampleSizes = SampleSizes,
# CheckMat = CheckMat,
# tolFailCheck = tolFailCheck,
# MCparam = MCparam, ...),
# beta = ComputeStatOutputPar(EstimOutput = EstimOutput,
# FctsToApply = FctsToApply, par = "beta",
# SampleSizes = SampleSizes,
# CheckMat = CheckMat,
# tolFailCheck = tolFailCheck,
# MCparam = MCparam, ...),
# gamma = ComputeStatOutputPar(EstimOutput = EstimOutput,
# FctsToApply = FctsToApply, par = "gamma",
# SampleSizes = SampleSizes,
# CheckMat = CheckMat,
# tolFailCheck = tolFailCheck,
# MCparam = MCparam, ...),
# delta = ComputeStatOutputPar(EstimOutput = EstimOutput,
# FctsToApply = FctsToApply, par = "delta",
# SampleSizes = SampleSizes,
# CheckMat = CheckMat,
# tolFailCheck = tolFailCheck,
# MCparam = MCparam, ...))
# }
##### for Output File#####
# No export.
initOutputFile <- function(thetaT, MCparam, TemperedType, Estimfct, ...) {
method <- Estim_Des_Temp(TemperedType, Estimfct, ...)
fileName <- get_filename(thetaT, MCparam, TemperedType, method)
if (!file.exists(fileName)) {
x <- switch(TemperedType,
CTS = paste("alphaT", "delta+T", "delta-T", "lambda+T",
"lambda-T", "muT", "data size", "seed",
"alphaE", "delta+E", "delta-E", "lambda+E",
"lambda-E", "muE", "failure", "time",
sep = ","),
TSS = paste("alphaT", "deltaT", "lambdaT",
"data size", "seed", "alphaE", "deltaE",
"lambdaE", "failure", "time",
sep = ","),
NTS = paste("alphaT", "betaT", "deltaT", "lambdaT", "muT",
"data size", "seed", "alphaE", "betaE", "deltaE",
"lambdaE", "muE", "failure", "time",
sep = ","),
MTS = paste("alphaT", "delta", "lambda+T", "lambda-T", "muT",
"data size", "seed", "alphaE", "deltaE", "lambda+E",
"lambda-E", "muE", "failure", "time", sep = ""),
GTS = paste("alpha+T", "alpha-T", "delta+T", "delta-T",
"lambda+T", "lambda-T", "muT", "data size", "seed",
"alpha+E", "alpha-E", "delta+E", "delta-E",
"lambda+E", "lambda-E", "muE", "failure", "time",
sep = ""),
KRTS = paste("alphaT", "k+T", "k-T", "r+T", "r-T",
"p+T", "p-T", "muT", "data size", "seed",
"alphaE", "k+E", "k-E", "r+E", "r-E",
"p+E", "p-E", "muE", "failure", "time", sep = ""),
RDTS = paste("alphaT", "delta", "lambda+T", "lambda-T", "muT",
"data size", "seed", "alphaE", "deltaE",
"lambda+E", "lambda-E", "muE", "failure", "time",
sep = "")
# ,CGMY = paste("CT", "GT", "MT", "YT", "data size", "seed", "CE",
# "GE", "ME", "YE","failure", "time",
# sep = ",")
)
write(x, file = fileName, sep = "\n")
}
}
# No export.
Estim_Des_Temp <- function(TemperedType = c("CTS", "TSS", "NTS", "MTS", "GTS",
"KRTS", "RDTS"),
EstimMethod = c("ML", "GMM", "Cgmm", "GMC"),
eps,
algo,
regularization,
WeightingMatrix,
t_scheme,
alphaReg,
t_free,
nb_t,
subdivisions,
IntegrationMethod,
randomIntegrationLaw,
s_min,
s_max,
ncond,
IterationControl,
...) {
TemperedType <- match.arg(TemperedType)
EstimMethod <- match.arg(EstimMethod)
EstimFcts <- getTempEstimFcts(TemperedType, EstimMethod,
eps = eps,
algo = algo,
regularization = regularization,
WeightingMatrix =
WeightingMatrix,
t_scheme = t_scheme,
alphaReg = alphaReg,
t_free = t_free,
nb_t = nb_t,
subdivisions = subdivisions,
IntegrationMethod =
IntegrationMethod,
randomIntegrationLaw =
randomIntegrationLaw,
s_min = s_min,
s_max = s_max,
ncond = ncond,
IterationControl = IterationControl,
...)
EstimFcts$methodDes(eps = eps,
algo = algo,
regularization = regularization,
WeightingMatrix =
WeightingMatrix,
t_scheme = t_scheme,
alphaReg = alphaReg,
t_free = t_free,
nb_t = nb_t,
subdivisions = subdivisions,
IntegrationMethod =
IntegrationMethod,
randomIntegrationLaw =
randomIntegrationLaw,
s_min = s_min,
s_max = s_max,
ncond = ncond,
IterationControl = IterationControl,
...)
}
##### for Checkpoints#####
# No export.
readCheckPoint <- function(ParameterMatrix, TemperedType, Estimfct, nab, npar,
nSS, MCparam,
eps,
algo,
regularization,
WeightingMatrix,
t_scheme,
alphaReg,
t_free,
nb_t,
subdivisions,
IntegrationMethod,
randomIntegrationLaw,
s_min,
s_max,
ncond,
IterationControl,
ab,
sample,
mc,
...) {
method <- Estim_Des_Temp(TemperedType, Estimfct,
eps = eps,
algo = algo,
regularization = regularization,
WeightingMatrix =
WeightingMatrix,
t_scheme = t_scheme,
alphaReg = alphaReg,
t_free = t_free,
nb_t = nb_t,
subdivisions = subdivisions,
IntegrationMethod =
IntegrationMethod,
randomIntegrationLaw =
randomIntegrationLaw,
s_min = s_min,
s_max = s_max,
ncond = ncond,
IterationControl = IterationControl, ...)
#This code can be used to create a checkpoint text file.
# fileName <- get_filename_checkPoint_Temp(ParameterMatrix, nab, npar,
# MCparam, method)
# if (!file.exists(fileName)) {
# write(x = "## ab;nab;npar;sample;nSS;mc;MCparam", file = fileName,
# sep = "\n")
# ab <- 1
# sample <- 1
# mc <- 0
# write(x = paste("--", ab, nab, npar, sample, nSS, mc, MCparam,
# sep = ";"), file = fileName, sep = "\n", append = TRUE)
# } else {
# tab <- as.numeric(utils::read.table(file = fileName,
# header = F, sep = ";"))
# ab <- tab[2]
# sample <- tab[5]
# mc <- tab[7]
# n_ab <- tab[3]
# n_par <- tab[4]
# n_SS <- tab[6]
# mc_Param <- tab[8]
# stopifnot(nab == n_ab, npar == n_par, nSS == n_SS, mc_Param == MCparam)
# }
list(ab = ab, nab = nab, npar = npar, sample = sample, nSS = nSS, mc = mc,
MCparam = MCparam)
}
# No export.
# This function writes every value of the current checkpoint in a string as
# return.
get_filename_checkPoint_Temp <- function(ParameterMatrix, nab, npar, MCparam,
method) {
# This should be adapted for every Tempered Type.
# BUT: The filenames are getting too long and this results in errors.
# Version 0.1.0 will not feature checkpoints during calculation.
#
# case: StableEstim
# if (npar == 3) {
# MC <- paste(paste("alpha0=", ParameterMatrix[1, 1], sep = ""),
# paste("delta0=", ParameterMatrix[1, 2], sep = ""),
# paste("lambda0=", ParameterMatrix[1, 3], sep = ""),
# paste("alphan=", ParameterMatrix[nab, 1], sep = ""),
# paste("deltan=", ParameterMatrix[nab, 2], sep = ""),
# paste("lambdan=", ParameterMatrix[nab, 3], sep = ""),
# paste("MCparam", MCparam, sep = ""), sep = "_")
# } else {
# case: "CTS"
# MC <- paste(paste("alpha0=", ParameterMatrix[1, 1], sep = ""),
# paste("delta+0=", ParameterMatrix[1, 2], sep = ""),
# paste("delta-0=", ParameterMatrix[1, 3], sep = ""),
# paste("lambda+0=", ParameterMatrix[1, 4], sep = ""),
# paste("lambda-0=",ParameterMatrix[1, 5], sep = ""),
# paste("mu0=", ParameterMatrix[1, 6], sep = ""),
# paste("alphan=", ParameterMatrix[nab,1], sep = ""),
# paste("delta+n=", ParameterMatrix[nab, 2], sep = ""),
# paste("delta-n=", ParameterMatrix[nab,3], sep = ""),
# paste("lambda+n=", ParameterMatrix[nab, 4], sep = ""),
# paste("lambda-n=", ParameterMatrix[nab,5], sep = ""),
# paste("mun=", ParameterMatrix[nab, 5], sep = ""),
# paste("MCparam", MCparam, sep = ""), sep = "_")
# }
MC <- paste("Test", sep = "")
methodTrunc <- method
if (base::nchar(methodTrunc) > 30){
methodTrunc <- substring(methodTrunc,1,30)
}
fileName <- paste(MC, methodTrunc, "_CHECKPOINT.txt", sep = "")
fileName
}
# No export.
updateCheckPointValues <- function(CheckPointValues, MCparam, lS, nab) {
ab_start <- CheckPointValues$ab
sample_start <- CheckPointValues$sample
mc_start <- CheckPointValues$mc
if (CheckPointValues$mc == MCparam) {
mc_start = 1
if (CheckPointValues$sample == lS) {
sample_start = 1
if (CheckPointValues$ab == nab)
stop("Simulation finished already! check your output file")
else ab_start = CheckPointValues$ab + 1
} else sample_start = CheckPointValues$sample + 1
} else mc_start = mc_start + 1
list(ab_start = ab_start, sample_start = sample_start, mc_start = mc_start)
}
# No export.
deleteCheckPoint <- function(ParameterMatrix, TemperedType, Estimfct, nab, npar,
nSS, MCparam,
eps,
algo,
regularization,
WeightingMatrix,
t_scheme,
alphaReg,
t_free,
subdivisions,
IntegrationMethod,
randomIntegrationLaw,
s_min,
s_max,
ncond,
IterationControl,
...) {
method <- Estim_Des_Temp(TemperedType, Estimfct,
eps = eps,
algo = algo,
regularization = regularization,
WeightingMatrix =
WeightingMatrix,
t_scheme = t_scheme,
alphaReg = alphaReg,
t_free = t_free,
subdivisions = subdivisions,
IntegrationMethod =
IntegrationMethod,
randomIntegrationLaw =
randomIntegrationLaw,
s_min = s_min,
s_max = s_max,
ncond = ncond,
IterationControl = IterationControl,
...)
##This code can be used to create a checkpoint text file.
# fileName <- get_filename_checkPoint_Temp(ParameterMatrix, nab, npar,
# MCparam, method)
# unlink(x = fileName, force = TRUE)
}
# No export.
writeCheckPoint <- function(ParameterMatrix, TemperedType, Estimfct, ab, nab,
npar, sample, nSS, mc, MCparam, ...) {
method <- Estim_Des_Temp(TemperedType, Estimfct, ...)
fileName <- get_filename_checkPoint_Temp(ParameterMatrix, nab, npar,
MCparam, method)
line = readLines(fileName, -1)
line[2] = paste("--", ab, nab, npar, sample, nSS, mc, MCparam, sep = ";")
writeLines(line, fileName)
}
#Added by Cedric 20220726
# Currently not necessary. Not even adatpted
#NameStatOutput <- function(FctsToApply, StatOutput) {
# Names <- c("alpha", "beta", "n", names(FctsToApply), "failure", "time")
# lapply(StatOutput, function(x) {
# colnames(x) <- Names
# return(x)
# })
#}
#Added by Cedric 20220729
# No export.
get_filename <- function(thetaT, MCparam, TemperedType, method,
extension = ".csv") {
MC <- switch(TemperedType,
CTS = paste(paste("Alpha=", thetaT[1]),
paste("DeltaP=", thetaT[2]),
paste("DeltaM=", thetaT[3]),
paste("LambdaP=", thetaT[4]),
paste("LambdaM=", thetaT[5]),
paste("mu=", thetaT[6]),
"MCparam", MCparam, sep = "_"),
TSS = paste(paste("Alpha=", thetaT[1]),
paste("Delta=", thetaT[2]),
paste("Lambda=", thetaT[3]),
"MCparam", MCparam, sep = "_"),
NTS = paste(paste("Alpha=", thetaT[1]),
paste("Beta=", thetaT[2]),
paste("Delta=", thetaT[3]),
paste("Lambda=", thetaT[4]),
paste("mu=", thetaT[5]),
"MCparam", MCparam, sep = "_"),
MTS = paste(paste("Alpha=", thetaT[1]),
paste("Delta=", thetaT[2]),
paste("LambdaP=", thetaT[3]),
paste("LambdaM=", thetaT[4]),
paste("mu=", thetaT[5]),
"MCparam", MCparam, sep = "_"),
GTS = paste(paste("AlphaP=", thetaT[1]),
paste("AlphaM=", thetaT[2]),
paste("DeltaP=", thetaT[3]),
paste("DeltaM=", thetaT[4]),
paste("LambdaP=", thetaT[5]),
paste("LambdaM=", thetaT[6]),
paste("mu=", thetaT[7]),
"MCparam", MCparam, sep = "_"),
KRTS = paste(paste("Alpha=", thetaT[1]),
paste("kP=", thetaT[2]),
paste("kM=", thetaT[3]),
paste("rP=", thetaT[4]),
paste("rM=", thetaT[5]),
paste("pP=", thetaT[6]),
paste("pM=", thetaT[7]),
paste("mu=", thetaT[8]),
"MCparam", MCparam, sep = "_"),
RDTS = paste(paste("Alpha=", thetaT[1]),
paste("Delta=", thetaT[2]),
paste("LambdaP=", thetaT[3]),
paste("LambdaM=", thetaT[4]),
paste("mu=", thetaT[5]),
"MCparam", MCparam, sep = "_")
# ,CGMY = paste(paste("C=", thetaT[1]),
# paste("G=", thetaT[2]),
# paste("M=", thetaT[3]),
# paste("Y=", thetaT[4]),
# "MCparam", MCparam, sep = "_")
)
fileName <- paste(MC, method, extension, sep = "")
fileName
}
#Added by Cedric 20220805
# No export.
updateOutputFile <- function(thetaT, MCparam, TemperedType, Output){
method <- Output$file
fileName <- get_filename(thetaT, MCparam, TemperedType, method)
if (!file.exists(fileName)) {
x <- switch(TemperedType,
CTS = paste("alphaT", "delta+T", "delta-T", "lambda+T",
"lambda-T", "muT", "data size", "seed",
"alphaE", "delta+E", "delta-E", "lambda+E",
"lambda-E", "muE", "failure", "time",
sep = ","),
TSS = paste("alphaT", "deltaT", "lambdaT",
"data size", "seed", "alphaE", "deltaE",
"lambdaE", "failure", "time",
sep = ","),
NTS = paste("alphaT", "betaT", "deltaT", "lambdaT", "muT",
"data size", "seed", "alphaE", "betaE", "deltaE",
"lambdaE", "muE", "failure", "time",
sep = ","),
MTS = paste("alphaT", "delta", "lambda+T", "lambda-T", "muT",
"data size", "seed", "alphaE", "deltaE", "lambda+E",
"lambda-E", "muE", "failure", "time", sep = ","),
GTS = paste("alpha+T", "alpha-T", "delta+T", "delta-T",
"lambda+T", "lambda-T", "muT", "data size", "seed",
"alpha+E", "alpha-E", "delta+E", "delta-E",
"lambda+E", "lambda-E", "muE", "failure", "time",
sep = ","),
KRTS = paste("alphaT", "k+T", "k-T", "r+T", "r-T",
"p+T", "p-T", "muT", "data size", "seed",
"alphaE", "k+E", "k-E", "r+E", "r-E",
"p+E", "p-E", "muE", "failure", "time", sep = ","),
RDTS = paste("alphaT", "delta", "lambda+T", "lambda-T", "muT",
"data size", "seed", "alphaE", "deltaE",
"lambda+E", "lambda-E", "muE", "failure", "time",
sep = ",")
# ,CGMY = paste("CT", "GT", "MT", "YT", "data size", "seed", "CE",
# "GE", "ME", "YE","failure", "time",
# sep = ",")
)
write(x, file = fileName, sep = "\n")
}
write(x = paste(as.character(Output$outputMat), collapse=","),
file = fileName, sep="\n", append=TRUE)
}
#Added by Cedric 20221011
# No export.
getTime_ <- function() proc.time()[3]
Try the TempStable package in your browser
Any scripts or data that you put into this service are public.
R Package Documentation
Browse R Packages
We want your feedback!
Note that we can't provide technical support on individual packages. You should contact the package authors for that.
Embedding an R snippet on your website
Add the following code to your website.
For more information on customizing the embed code, read Embedding Snippets.