Nothing
#' @title Double sigmoidal fit function with h0.
#'
#' @param dataInput A dataframe or a list containing the dataframe. The dataframe should be composed of at least two columns. One represents time, and the other represents intensity. The data should be normalized with the normalize data function sicegar::normalizeData() before being imported into this function.
#' @param tryCounter A counter that shows the number of times the data was fit via maximum likelihood function.
#' @param startList The vector containing the initial set of parameters that algorithm tries for the first fit attempt for the relevant parameters. The vector composes of seven elements; 'finalAsymptoteIntensityRatio', 'maximum', 'slope1Param', 'midPoint1Param' , 'slope2Param', 'midPointDistanceParam', and 'h0'. Detailed explanations of those parameters can be found in vignettes. Defaults are finalAsymptoteIntensityRatio = 0, maximum = 1, slope1Param = 1, midPoint1Param = 0.33, slope2Param = 1, midPointDistanceParam=0.29, and h0 = 0. The numbers are in normalized time intensity scale.
#' @param lowerBounds The lower bounds for the randomly generated start parameters. The vector is composed of seven elements; 'finalAsymptoteIntensityRatio', 'maximum', 'slope1Param', 'midPoint1Param' , 'slope2Param', 'midPointDistanceParam', and 'h0'. Detailed explanations of those parameters can be found in vignettes. Defaults are finalAsymptoteIntensityRatio = 0, maximum = 0.3, slope1Param = .01, midPoint1Param = -0.52, slope2Param = .01, midPointDistanceParam = 0.04, and h0 = -0.1. The numbers are in normalized time intensity scale.
#' @param upperBounds The upper bounds for the randomly generated start parameters. The vector is composed of seven elements; 'finalAsymptoteIntensityRatio', 'maximum', 'slope1Param', 'midPoint1Param' , 'slope2Param', 'midPointDistanceParam', and 'h0'. Detailed explanations of those parameters can be found in vignettes. Defaults are finalAsymptoteIntensityRatio = 1, maximum = 1.5, slope1Param = 180, midPoint1Param = 1.15, slope2Param = 180, midPointDistanceParam = 0.63, and h0 = 0.3. The numbers are in normalized time intensity scale.
#' @param min_Factor Defines the minimum step size used by the fitting algorithm. Default is 1/2^20.
#' @param n_iterations Defines maximum number of iterations used by the fitting algorithm. Default is 1000.
#'
#' @description The function fits a double sigmoidal curve to given data by using likelihood maximization (LM) algorithm and provides the parameters (maximum, final asymptote intensity, slope1Param, midpoint1Param, slope2Param, and midpoint distance) describing the double sigmoidal fit as output. It also contains information about the goodness of fits such as AIC, BIC, residual sum of squares, and log likelihood.
#' @return Returns the fitted parameters and goodness of fit metrics.
#' @export
#' @examples
#' # runif() is used here for consistency with previous versions of the sicegar package. However,
#' # rnorm() will generate symmetric errors, producing less biased numerical parameter estimates.
#' # We recommend errors generated with rnorm() for any simulation studies on sicegar.
#'time=seq(3, 24, 0.1)
#'
#'#simulate intensity data and add noise
#'noise_parameter <- 0.2
#'intensity_noise <- stats::runif(n = length(time), min = 0, max = 1) * noise_parameter
#'intensity <- doublesigmoidalFitFormula_h0(time,
#' finalAsymptoteIntensityRatio = .3,
#' maximum = 4,
#' slope1Param = 1,
#' midPoint1Param = 7,
#' slope2Param = 1,
#' midPointDistanceParam = 8,
#' h0 = 1)
#'intensity <- intensity + intensity_noise
#'
#'dataInput <- data.frame(intensity = intensity, time = time)
#'normalizedInput <- normalizeData(dataInput)
#'parameterVector <- doublesigmoidalFitFunction_h0(normalizedInput, tryCounter = 1)
#'
#'
#'#Check the results
#'# doublesigmoidalFitFunction_h0() is run on the startList param values (because 'tryCounter = 1')
#'# use multipleFitFunction() for multiple random starts in order to optimize
#'if(parameterVector$isThisaFit){
#' intensityTheoretical <-
#' doublesigmoidalFitFormula_h0(
#' time,
#' finalAsymptoteIntensityRatio = parameterVector$finalAsymptoteIntensityRatio_Estimate,
#' maximum = parameterVector$maximum_Estimate,
#' slope1Param = parameterVector$slope1Param_Estimate,
#' midPoint1Param = parameterVector$midPoint1Param_Estimate,
#' slope2Param = parameterVector$slope2Param_Estimate,
#' midPointDistanceParam = parameterVector$midPointDistanceParam_Estimate,
#' h0 = parameterVector$h0_Estimate)
#'
#' comparisonData <- cbind(dataInput, intensityTheoretical)
#'
#' require(ggplot2)
#' ggplot(comparisonData) +
#' geom_point(aes(x = time, y = intensity)) +
#' geom_line(aes(x = time, y = intensityTheoretical)) +
#' expand_limits(x = 0, y = 0)}
#'
#'if(!parameterVector$isThisaFit) {print(parameterVector)}
#'
doublesigmoidalFitFunction_h0 <- function(dataInput,
tryCounter,
startList=list(finalAsymptoteIntensityRatio = 0,
maximum = 1,
slope1Param = 1,
midPoint1Param = 0.33,
slope2Param = 1,
midPointDistanceParam = 0.29,
h0 = 0.1),
lowerBounds=c(finalAsymptoteIntensityRatio = 0,
maximum = 0.3,
slope1Param = .01,
midPoint1Param = -0.52,
slope2Param = .01,
midPointDistanceParam = 0.04,
h0 = 0),
upperBounds=c(finalAsymptoteIntensityRatio = 1,
maximum = 1.5,
slope1Param = 180,
midPoint1Param = 1.15,
slope2Param = 180,
midPointDistanceParam = 0.63,
h0 = 0.3),
min_Factor = 1/2^20,
n_iterations = 1000)
{
isalist = (is.list(dataInput) & !is.data.frame(dataInput))
if(isalist){
dataFrameInput <- dataInput$timeIntensityData
}
isadataframe <- (is.data.frame(dataInput))
if(isadataframe){
dataFrameInput <- dataInput
}
if(tryCounter==1){
counterDependentStartList <- startList
}
if(tryCounter!=1){
randomVector <- stats::runif(length(startList), 0, 1)
names(randomVector) <- c("finalAsymptoteIntensityRatio",
"maximum",
"slope1Param",
"midPoint1Param",
"slope2Param",
"midPointDistanceParam",
"h0")
counterDependentStartVector <- randomVector * (upperBounds - lowerBounds) + lowerBounds
counterDependentStartList <- as.list(counterDependentStartVector)
}
theFitResult <- try(minpack.lm::nlsLM(intensity ~ doublesigmoidalFitFormula_h0(time,
finalAsymptoteIntensityRatio,
maximum,
slope1Param,
midPoint1Param,
slope2Param,
midPointDistanceParam,
h0),
dataFrameInput,
start = counterDependentStartList,
control = list(maxiter = n_iterations, minFactor = min_Factor),
lower = lowerBounds,
upper = upperBounds,
trace=F), silent = TRUE)
if(!inherits(theFitResult, "try-error"))
{
parameterMatrix <- summary(theFitResult)$parameters
colnames(parameterMatrix) <- c("Estimate", "Std_Error", "t_value", "Pr_t")
parameterVector <- c(t(parameterMatrix))
names(parameterVector) <- c("finalAsymptoteIntensityRatio_N_Estimate", "finalAsymptoteIntensityRatio_Std_Error", "finalAsymptoteIntensityRatio_t_value", "finalAsymptoteIntensityRatio_Pr_t",
"maximum_N_Estimate", "maximum_Std_Error", "maximum_t_value", "maximum_Pr_t",
"slope1Param_N_Estimate", "slope1Param_Std_Error", "slope1Param_t_value", "slope1Param_Pr_t",
"midPoint1Param_N_Estimate", "midPoint1Param_Std_Error", "midPoint1Param_t_value", "midPoint1Param_Pr_t",
"slope2Param_N_Estimate", "slope2Param_Std_Error", "slope2Param_t_value", "slope2Param_Pr_t",
"midPointDistanceParam_N_Estimate", "midPointDistanceParam_Std_Error", "midPointDistanceParam_t_value", "midPointDistanceParam_Pr_t",
"h0_N_Estimate","h0_Std_Error","h0_t_value","h0_Pr_t" # ← appended
)
parameterVector <- c(parameterVector,
residual_Sum_of_Squares = sum((as.vector(stats::resid(theFitResult)))^2)[1],
log_likelihood = as.vector(stats::logLik(theFitResult))[1],
AIC_value = as.vector(stats::AIC(theFitResult))[1],
BIC_value = as.vector(stats::BIC(theFitResult))[1])
parameterList <- as.list(parameterVector)
parameterList$isThisaFit <- TRUE
parameterList$startVector <- counterDependentStartList
if(isalist){
parameterList$dataScalingParameters <- as.list(dataInput$dataScalingParameters)
}
parameterList$model <- as.character("doublesigmoidal")
parameterList$additionalParameters <- FALSE
parameterDf <- as.data.frame(parameterList)
#Renormalize Parameters
parameterDf <- doublesigmoidalRenormalizeParameters_h0(parameterDf, isalist)
}
if(inherits(theFitResult, "try-error"))
{
parameterVector <- rep(NA, 28)
names(parameterVector) <- c("finalAsymptoteIntensityRatio_N_Estimate", "finalAsymptoteIntensityRatio_Std_Error", "finalAsymptoteIntensityRatio_t_value", "finalAsymptoteIntensityRatio_Pr_t",
"maximum_N_Estimate", "maximum_Std_Error", "maximum_t_value", "maximum_Pr_t",
"slope1Param_N_Estimate", "slope1Param_Std_Error", "slope1Param_t_value", "slope1Param_Pr_t",
"midPoint1Param_N_Estimate", "midPoint1Param_Std_Error", "midPoint1Param_t_value", "midPoint1Param_Pr_t",
"slope2Param_N_Estimate", "slope2Param_Std_Error", "slope2Param_t_value", "slope2Param_Pr_t",
"midPointDistanceParam_N_Estimate", "midPointDistanceParam_Std_Error", "midPointDistanceParam_t_value", "midPointDistanceParam_Pr_t",
"h0_N_Estimate","h0_Std_Error","h0_t_value","h0_Pr_t" # ← appended
)
parameterVector <- c(parameterVector,
residual_Sum_of_Squares = Inf,
log_likelihood = NA,
AIC_value = NA,
BIC_value = NA)
parameterList <- as.list(parameterVector)
parameterList$isThisaFit <- FALSE
parameterList$startVector <- counterDependentStartList
if(isalist) {
parameterList$dataScalingParameters <- as.list(dataInput$dataScalingParameters)
}
parameterList$model <- "doublesigmoidal"
parameterDf <- as.data.frame(parameterList)
#Renormalize Parameters
parameterDf <- doublesigmoidalRenormalizeParameters_h0(parameterDf, isalist)
}
return(parameterDf)
}
#**************************************
#**************************************
#' @title Double Sigmoidal Formula with h0
#'
#' @param x the "time" column of the dataframe
#' @param finalAsymptoteIntensityRatio this is the ratio between the final asymptote intensity and maximum intensity of the fitted curve.
#' @param maximum the maximum intensity that the double sigmoidal function reaches.
#' @param slope1Param the slope parameter of the sigmoidal function at the steepest point in the exponential phase of the viral production.
#' @param midPoint1Param the x axis value of the steepest point in the function.
#' @param slope2Param the slope parameter of the sigmoidal function at the steepest point in the lysis phase. i.e when the intensity is decreasing.
#' @param midPointDistanceParam the distance between the time of steepest increase and steepest decrease in the intensity data. In other words the distance between the x axis values of arguments of slope1Param and slope2Param.
#' @param h0 the lower asymptote (baseline) intensity
#'
#' @description Calculates intensities using the double sigmoidal model fit and the parameters (maximum, final asymptote intensity, slope1Param, midpoint1Param, slope2Param, and midpoint distance).
#' @return Returns the predicted intensities for the given time points with the double sigmoidal fitted parameters for the double sigmoidal fit.
#' @export
#' @examples
#' # runif() is used here for consistency with previous versions of the sicegar package. However,
#' # rnorm() will generate symmetric errors, producing less biased numerical parameter estimates.
#' # We recommend errors generated with rnorm() for any simulation studies on sicegar.
#'time <- seq(3, 24, 0.1)
#'
#'#simulate intensity data and add noise
#'noise_parameter <- 0.2
#'intensity_noise <- stats::runif(n = length(time), min = 0, max = 1) * noise_parameter
#'intensity <- doublesigmoidalFitFormula_h0(time,
#' finalAsymptoteIntensityRatio = .3,
#' maximum = 4,
#' slope1Param = 1,
#' midPoint1Param = 7,
#' slope2Param = 1,
#' midPointDistanceParam = 8,
#' h0 = 1)
#'intensity <- intensity + intensity_noise
#'
#'dataInput <- data.frame(intensity = intensity, time = time)
#'normalizedInput <- normalizeData(dataInput)
#'parameterVector <- doublesigmoidalFitFunction_h0(normalizedInput, tryCounter = 1)
#'
#'
#'#Check the results
#'# doublesigmoidalFitFunction_h0() is run on the startList param values (because 'tryCounter = 1')
#'# use multipleFitFunction() for multiple random starts in order to optimize
#'if(parameterVector$isThisaFit){
#' intensityTheoretical <-
#' doublesigmoidalFitFormula_h0(
#' time,
#' finalAsymptoteIntensityRatio = parameterVector$finalAsymptoteIntensityRatio_Estimate,
#' maximum = parameterVector$maximum_Estimate,
#' slope1Param = parameterVector$slope1Param_Estimate,
#' midPoint1Param = parameterVector$midPoint1Param_Estimate,
#' slope2Param = parameterVector$slope2Param_Estimate,
#' midPointDistanceParam = parameterVector$midPointDistanceParam_Estimate,
#' h0 = parameterVector$h0_Estimate)
#'
#' comparisonData <- cbind(dataInput, intensityTheoretical)
#'
#' require(ggplot2)
#' ggplot(comparisonData) +
#' geom_point(aes(x = time, y = intensity)) +
#' geom_line(aes(x = time, y = intensityTheoretical)) +
#' expand_limits(x = 0, y = 0)
#' }
#'
#'if(!parameterVector$isThisaFit){
#' print(parameterVector)
#' }
#'
#'
doublesigmoidalFitFormula_h0 <- function(x,
finalAsymptoteIntensityRatio,
maximum,
slope1Param,
midPoint1Param,
slope2Param,
midPointDistanceParam,
h0)
{
if(slope1Param < 0){
stop("slope1Param should be a positive number")
}
if(slope2Param < 0){
stop("slope2Param should be a positive number. It is the absolute value of the second slopeParam")
}
if(midPointDistanceParam < 0){
stop("midPointDistanceParam should be a positive number. It is the distance between two steepest points of exponential phase and lysis")
}
if(finalAsymptoteIntensityRatio < 0 | finalAsymptoteIntensityRatio > 1){
stop("finalAsymptoteIntensityRatio should be a number between 0 and 1")
}
if(maximum < 0){
stop("maximum should be a positive number")
}
optimizeIntervalMin <- midPoint1Param - 2 * midPointDistanceParam
optimizeIntervalMax <- midPoint1Param + 3 * midPointDistanceParam
xmax <- stats::optimize(f1,
c(optimizeIntervalMin, optimizeIntervalMax),
tol = 0.0001,
B1 = slope1Param, M1 = midPoint1Param, B2 = slope2Param, L = midPointDistanceParam, maximum = TRUE);
argumentt <- xmax$maximum;
constt <- f0(argumentt, slope1Param, midPoint1Param, slope2Param, midPointDistanceParam);
y <- f2_h0(x,
finalAsymptoteIntensityRatio, maximum,
slope1Param, midPoint1Param,
slope2Param, midPointDistanceParam,
constt, argumentt, h0);
return(y)
}
#**************************************
#**************************************
# Those four functions are necessary internal steps for double sigmoidal model
# All four are internal functions
# Name convention is different than the rest of the package
# A2 correspond to finalAsymptoteIntensityRatio
# Ka correspond to maximum
# B1 correspond to slope1Param
# B2 correspond to slope2Param
# M1 correspond to midpoint1Param
# L correspond to midPointDistanceParam
# argument correspond to maximum_x (i.e the time when intensity reaches maximum with
# respect to the double sigmoidal model)
# const correspond to the maximum value the function reaches when the input data intensity and
# is normalized between 0,1
# f0, f1, and f0mid are unchanged and are available in doublesigmoidlFitFunctions.R
# f0 <- function (x, B1, M1, B2, L) {
# (1/( ( 1+exp(-B1 * (x - M1)) ) * ( 1 + exp(B2 * (x - (M1 + L))) ) ))
# }
# f1 <- function (x, B1, M1, B2, L) {
# log((1 / ( ( 1 + exp(-B1 * (x - M1)) ) * ( 1 + exp(B2 * (x - (M1 + L))) ) )))
# }
f2_h0 <- function (x, A2, Ka, B1, M1, B2, L, const, argument, h0){
fBasics::Heaviside(x - argument) * (f0(x, B1, M1, B2, L) * ((Ka - A2 * Ka)/(const)) + A2 * Ka) +
(1 - fBasics::Heaviside(x - argument)) * (f0(x, B1, M1, B2, L) * ((Ka - h0)/(const)) + h0)
}
# f0mid <- function (x, B1, M1, B2, L,const) {
# -const / 2 + 1 / ( ( 1 + exp(-B1 * (x-M1)) ) * ( 1 + exp(B2*(x-(M1+L))) ) )
# }
#**************************************
#**************************************
# @title f_argmax_doublesigmoidal_h0
#
# @param parameterVector is a dataFrame that includes columns with names:
#
# *"dataScalingParameters.timeRange" that represents maximum time of intensity measurements,
#
# parameters of double sigmoidal model:
# *"slope1Param_Estimate",
# *"midPoint1Param_Estimate",
# *"slope2Param_Estimate",
# *"midPointDistanceParam_Estimate".
# (Note: additional parameters in Df does not cause problems)
#
# finally a column with the used model information model="doublesigmoidal".
# @description This function is designed to compensate for the small discrepancy in numerical parameters in the double sigmoidal model. This function is called by numericalReCalculation.
# @return return the x value of the maximum in between time 0 and maximum time measured in sample intensities. In other words the time that the double sigmoidal function reaches its maximum.
#
# @examples
# # Related example
#
## Initial Command to Reset the System
#rm(list = ls())
#if (is.integer(dev.list())){dev.off()}
#cat("\014")
#
## Generate a parameters DF
#parameterDf = data.frame(dataScalingParameters.timeRange = 24,
# slope1Param_Estimate = 1,
# midPoint1Param_Estimate = 7,
# slope2Param_Estimate = 1,
# midPointDistanceParam_Estimate = 8,
# model = "doublesigmoidal")
#
# # Find the x value of the maximum for this model
# # with in the range of [0, dataScalingParameters.timeRange]
# maximum_x = f_argmax_doublesigmoidal_h0(parameterDf)
#
f_argmax_doublesigmoidal_h0 <- function(parameterVector){
slope1Param <- parameterVector$slope1Param_N_Estimate
slope2Param <- parameterVector$slope2Param_N_Estimate
midPointDistanceParam <- parameterVector$midPointDistanceParam_N_Estimate
finalAsymptoteIntensityRatio <- parameterVector$finalAsymptoteIntensityRatio_N_Estimate
maximum <- parameterVector$maximum_N_Estimate
midPoint1Param <- parameterVector$midPoint1Param_N_Estimate
timeRange <- parameterVector$dataScalingParameters.timeRange
if(slope1Param <0 ){
stop("slope1Param should be a positive number")
}
if(slope2Param < 0){
stop("slope2Param should be a positive number. It is the absolute value of the slope2Param")
}
if(midPointDistanceParam < 0){
stop("midPointDistanceParam should be a positive number. It is the distance between two steepest points of exponential phase and lysis")
}
if(finalAsymptoteIntensityRatio < 0 | finalAsymptoteIntensityRatio > 1){
stop("finalAsymptoteIntensityRatio should be a number between 0 and 1")
}
if(maximum < 0){
stop("maximum should be a positive number")
}
optimizeIntervalMin <- midPoint1Param - 2 * midPointDistanceParam
optimizeIntervalMax <- midPoint1Param + 3 * midPointDistanceParam
xmax <- stats::optimize(f1,
c(optimizeIntervalMin, optimizeIntervalMax),
tol = 0.0001,
B1 = slope1Param, M1 = midPoint1Param, B2 = slope2Param, L = midPointDistanceParam, maximum = TRUE);
argumentt <- xmax$maximum * timeRange;
return(argumentt)
}
#**************************************
#**************************************
# f_mid1_doublesigmoidal is unchanged and is available in doublesigmoidalFitFunctions.R
# # @title f_mid1_doublesigmoidal
# #
# # @param parameterDf is a dataFrame that includes columns with names:
# #
# # *"dataScalingParameters.timeRange" that represents maximum time of intensity measurements,
# #
# # parameters of double sigmoidal model:
# # *"slope1Param_Estimate",
# # *"midPoint1Param_Estimate",
# # *"slope2Param_Estimate",
# # *"midPointDistanceParam_Estimate".
# # (Note: additional parameters in Df does not cause problems)
# #
# # finally a column with the used model information model="doublesigmoidal".
# #
# # @description This function is designed to compensate for the small discrepancy the in numerical parameters in double sigmoidal model. This function is called by numericalReCalculation.
# # @return returns the x values of midpoint1Param, by assuming the y of midpoint1Param is the half of the maximum value it reaches.
# # @export
# #
# # @examples
# # # Related example
# #
# ## Initial Command to Reset the System
# #rm(list = ls())
# #if (is.integer(dev.list())){dev.off()}
# #cat("\014")
# #
# ## Generate a parameters DF
# #parameterDf <- data.frame(dataScalingParameters.timeRange = 24,
# # slope1Param_Estimate = 1,
# # midPoint1Param_Estimate = 7,
# # slope2Param_Estimate = 1,
# # midPointDistanceParam_Estimate = 8,
# # model = "doublesigmoidal")
# #
# # # Find the x value of the maximum for this model
# # # with in the range of [0, dataScalingParameters.timeRange]
# # midPoint1Param_x <- f_mid1_doublesigmoidal(parameterDf)
# #
# f_mid1_doublesigmoidal <- function(parameterDf){
#
# max_x <- parameterDf$dataScalingParameters.timeRange
# xmax <- stats::optimize(f1,
# interval = c(-1.125 * max_x, max_x * 3),
# tol = 0.0001,
# B1 = parameterDf$slope1Param_Estimate,
# M1 = parameterDf$midPoint1Param_Estimate,
# B2 = parameterDf$slope2Param_Estimate,
# L = parameterDf$midPointDistanceParam_Estimate,
# maximum = TRUE);
# argumentt <- xmax$maximum;
# constt <- f0(argumentt,
# B1 = parameterDf$slope1Param_Estimate,
# M1 = parameterDf$midPoint1Param_Estimate,
# B2 = parameterDf$slope2Param_Estimate,
# L = parameterDf$midPointDistanceParam_Estimate);
#
# mid1x <- try(expr = stats::uniroot(f0mid,
# interval = c(-1.125 * max_x, argumentt),
# tol = 0.0001,
# B1 = parameterDf$slope1Param_Estimate,
# M1 = parameterDf$midPoint1Param_Estimate,
# B2 = parameterDf$slope2Param_Estimate,
# L = parameterDf$midPointDistanceParam_Estimate, const = constt), silent = TRUE);
# if(inherits(mid1x, "try-error"))
# {
# rangeList <- seq(-1,30)
# rangeListCounter <- -1
# while(inherits(mid1x, "try-error") && rangeListCounter <= 30)
# {
# #print(rangeListCounter)
# rangeListSelection <- rangeList[rangeListCounter + 2];
# mid1x <- try(expr = stats::uniroot(f0mid,
# interval = c((-1) * (2^rangeListSelection) * max_x, argumentt),
# tol = 0.0001,
# B1 = parameterDf$slope1Param_Estimate,
# M1 = parameterDf$midPoint1Param_Estimate,
# B2 = parameterDf$slope2Param_Estimate,
# L = parameterDf$midPointDistanceParam_Estimate, const = constt), silent = TRUE);
# rangeListCounter <- rangeListCounter + 1
# }
# }
#
# mid1x <- mid1x$root
#
# return(mid1x)
# }
#**************************************
#**************************************
# f_mid2_doublesigmoidal is unchanged and is available in doublesigmoidalFitFunctions.R
# # @title f_mid2_doublesigmoidal
# #
# # @param parameterDf is a dataFrame that includes columns with names:
# #
# # *"dataScalingParameters.timeRange" that represents maximum time of intensity measurements,
# #
# # parameters of double sigmoidal model:
# # *"slope1Param_Estimate",
# # *"midPoint1Param_Estimate",
# # *"slope2Param_Estimate",
# # *"midPointDistanceParam_Estimate".
# # (Note: additional parameters in Df does not cause problems)
# #
# # finally a column with the used model information model="doublesigmoidal".
# #
# # @description This function is designed to compensate for the small discrepancy in numerical parameters in double sigmoidal model. This function is called by numericalReCalculation.
# # @return returns the x values of midpoint2, assuming the y value of midpoint2 is half of the distance between the maximum value it reaches and the final asymptote.
# # @export
# #
# # @examples
# # # Related example
# #
# ## Initial Command to Reset the System
# #rm(list = ls())
# #if (is.integer(dev.list())){dev.off()}
# #cat("\014")
# #
# ## Generate a parameters DF
# #parameterDf <- data.frame(dataScalingParameters.timeRange = 24,
# # slope1Param_Estimate = 1,
# # midPoint1Param_Estimate = 7,
# # slope2Param_Estimate = 1,
# # midPointDistanceParam_Estimate = 8,
# # model = "doublesigmoidal")
# #
# # # Find the x value of the maximum for this model
# # # with in the range of [0, dataScalingParameters.timeRange]
# # midPoint2_x <- f_mid2_doublesigmoidal(parameterDf)
# #
# f_mid2_doublesigmoidal <- function(parameterDf){ # also untouched for now, but also a place we added try-error stop...
#
# max_x <- parameterDf$dataScalingParameters.timeRange
# xmax <- stats::optimize(f1,
# interval = c(-1.125 * max_x, max_x * 3),
# tol = 0.0001,
# B1 = parameterDf$slope1Param_Estimate,
# M1 = parameterDf$midPoint1Param_Estimate,
# B2 = parameterDf$slope2Param_Estimate,
# L = parameterDf$midPointDistanceParam_Estimate,
# maximum = TRUE);
#
# argumentt <- xmax$maximum;
# constt <- f0(argumentt,
# B1 = parameterDf$slope1Param_Estimate,
# M1 = parameterDf$midPoint1Param_Estimate,
# B2 = parameterDf$slope2Param_Estimate ,
# L = parameterDf$midPointDistanceParam_Estimate);
#
# mid2x <- try(stats::uniroot(f0mid,
# interval = c(argumentt, max_x * (1.25)),
# tol = 0.0001,
# B1 = parameterDf$slope1Param_Estimate,
# M1 = parameterDf$midPoint1Param_Estimate,
# B2 = parameterDf$slope2Param_Estimate,
# L = parameterDf$midPointDistanceParam_Estimate, const=constt), silent = TRUE);
#
# if(inherits(mid2x, "try-error")){
# rangeList <- seq(-1,30)
# rangeListCounter <- -1
# while(inherits(mid2x, "try-error") && rangeListCounter <= 30){
# rangeListSelection <- rangeList[rangeListCounter + 2];
# mid2x <- try(stats::uniroot(f0mid,
# interval = c(argumentt, max_x * (2^rangeListSelection)),
# tol = 0.0001, B1 = parameterDf$slope1Param_Estimate,
# M1 = parameterDf$midPoint1Param_Estimate,
# B2 = parameterDf$slope2Param_Estimate,
# L = parameterDf$midPointDistanceParam_Estimate, const = constt), silent = TRUE);
# rangeListCounter <- rangeListCounter + 1
# }
# }
#
# mid2x <- mid2x$root
#
# return(mid2x)
# }
#**************************************
#**************************************
# @title f_slope_doublesigmoidal_h0
# @param x as a single time point or a sequence of time points that one wants to calculate the derivative of.
# @param parameterDf is a dataFrame that includes columns with names:
#
# *"dataScalingParameters.timeRange" which represents the maximum time of intensity measurements,
#
# parameters of double sigmoidal model:
# *"finalAsymptoteIntensityRatio_Estimate",
# *"maximum_Estimate",
# *"slope1Param_Estimate",
# *"midPoint1Param_Estimate",
# *"slope2Param_Estimate",
# *"midPointDistanceParam_Estimate".
#
# finally a column with the used model information model="doublesigmoidal".
# @param timeStep is the time step for the derivative calculation increment "h". Default value is 0.00001
#
#
# @description The function calculates the numerical slope of double sigmoidal function with given parameters by using 5 points derivative. It is designed to compensate for the small discrepancy in numerical parameters in the double sigmoidal model. The function is called by numericalReCalculation.
# @return return the numerical slope of double sigmoidal for the given x value.
# @export
#
# @examples
# # Related examples
#
## Example 1
## Initial Command to Reset the System
#rm(list = ls())
#if (is.integer(dev.list())){dev.off()}
#cat("\014")
#
## Generate a parameters DF
#parameterDataframe <- data.frame(dataScalingParameters.timeRange = 24,
# finalAsymptoteIntensityRatio_Estimate = .3,
# maximum_Estimate = 4,
# slope1Param_Estimate = 1,
# midPoint1Param_Estimate = 7,
# slope2Param_Estimate = 1,
# midPointDistanceParam_Estimate = 8,
# model = "doublesigmoidal")
#
# # Find the x value of the maximum for this model
# # with in the range of [0, dataScalingParameters.timeRange]
# slopeParam_Estimate <- f_slope_doublesigmoidal_h0(x = 6.5,
# parameterDf = parameterDataframe,
# timeStep = 0.00001)
#
## Example 2
## x might be vector as well
## Initial Command to Reset the System
#rm(list = ls())
#if (is.integer(dev.list())){dev.off()}
#cat("\014")
#
## Generate a parameters DF
#parameterDataframe <- data.frame(dataScalingParameters.timeRange = 24,
# finalAsymptoteIntensityRatio_Estimate = .3,
# maximum_Estimate = 4,
# slope1Param_Estimate = 1,
# midPoint1Param_Estimate = 7,
# slope2Param_Estimate = 1,
# midPointDistanceParam_Estimate = 8,
# model = "doublesigmoidal")
#
# # Find the x value of the maximum for this model
# # with in the range of [0, dataScalingParameters.timeRange]
# slopeParam_Estimate <- f_slope_doublesigmoidal_h0(x = seq(-6.5, 50, 0.5),
# parameterDf = parameterDataframe,
# timeStep = 0.00001)
## Note that some points in the x sequence are out of the range of
## time <- [0, dataScalingParameters.timeRange]
#
f_slope_doublesigmoidal_h0 <- function(x, parameterDf, timeStep = 0.00001){
fxp2h <- doublesigmoidalFitFormula_h0(x = x + 2 * timeStep,
finalAsymptoteIntensityRatio = parameterDf$finalAsymptoteIntensityRatio_Estimate,
maximum = parameterDf$maximum_Estimate,
slope1Param = parameterDf$slope1Param_Estimate,
midPoint1Param = parameterDf$midPoint1Param_Estimate,
slope2Param = parameterDf$slope2Param_Estimate,
midPointDistanceParam = parameterDf$midPointDistanceParam_Estimate,
h0 = parameterDf$h0_Estimate)
fxph <- doublesigmoidalFitFormula_h0(x = x + timeStep,
finalAsymptoteIntensityRatio = parameterDf$finalAsymptoteIntensityRatio_Estimate,
maximum = parameterDf$maximum_Estimate,
slope1Param = parameterDf$slope1Param_Estimate,
midPoint1Param = parameterDf$midPoint1Param_Estimate,
slope2Param = parameterDf$slope2Param_Estimate,
midPointDistanceParam = parameterDf$midPointDistanceParam_Estimate,
h0 = parameterDf$h0_Estimate)
fxmh <- doublesigmoidalFitFormula_h0(x = x - timeStep,
finalAsymptoteIntensityRatio = parameterDf$finalAsymptoteIntensityRatio_Estimate,
maximum = parameterDf$maximum_Estimate,
slope1Param = parameterDf$slope1Param_Estimate,
midPoint1Param = parameterDf$midPoint1Param_Estimate,
slope2Param = parameterDf$slope2Param_Estimate,
midPointDistanceParam = parameterDf$midPointDistanceParam_Estimate,
h0 = parameterDf$h0_Estimate)
fxm2h <- doublesigmoidalFitFormula_h0(x = x - 2 * timeStep,
finalAsymptoteIntensityRatio = parameterDf$finalAsymptoteIntensityRatio_Estimate,
maximum = parameterDf$maximum_Estimate,
slope1Param = parameterDf$slope1Param_Estimate,
midPoint1Param = parameterDf$midPoint1Param_Estimate,
slope2Param = parameterDf$slope2Param_Estimate,
midPointDistanceParam = parameterDf$midPointDistanceParam_Estimate,
h0 = parameterDf$h0_Estimate)
der <- (-1 * fxp2h + 8 * fxph - 8 * fxmh + 1 * fxm2h) / (12 * timeStep)
return(der)
}
#**************************************
#**************************************
# @title sigmoidalRenormalizeParameters_h0 (This is an internal function)
# @param parametersDf is the parameter data frame generated by sigmoidal fit function and
# includes the parameters named
# *maximum_N_Estimate (normalized Maximum Estimate)
# *finalAsymptoteIntensityRatio_N_Estimate (normalized final asymptote intensity estimate)
# *slope1Param_N_Estimate (normalized Slope1Param Estimate)
# *midPoint1Param_N_Estimate (normalized Midpoint1Param Estimate)
# *slope2Param_N_Estimate (normalized Slope2Param Estimate)
# *midPointDistanceParam_N_Estimate (normalized distance between midpoints estimate)
#. *h0_N_Estimate (normalized estimate of lower asymptote)
# *dataScalingParameters.intensityRange the range of initial unnormalized intensity. Provided if the data is normalized
# *parameterDF$dataScalingParameters.intensityMin the minimum of unnormalized intensity. Provided if the data is normalized
# *parameterDF$dataScalingParameters.timeRange the maximum time that the experiment reach. Provided if the data is normalized
# @param isalist checks if the input provided is a list (i.e normalized) or a data frame (i.e not normalized)
# @details If the fit was done on normalized data, parameter estimates will be on normalized scale.
# This function unnormalizes those parameter estimations.
doublesigmoidalRenormalizeParameters_h0 <- function(parameterDF, isalist)
{
if(isalist){
parameterDF$finalAsymptoteIntensityRatio_Estimate <- (parameterDF$finalAsymptoteIntensityRatio_N_Estimate * parameterDF$maximum_N_Estimate * parameterDF$dataScalingParameters.intensityRange + parameterDF$dataScalingParameters.intensityMin) /
(parameterDF$maximum_N_Estimate * parameterDF$dataScalingParameters.intensityRange + parameterDF$dataScalingParameters.intensityMin)
parameterDF$maximum_Estimate <- parameterDF$maximum_N_Estimate * parameterDF$dataScalingParameters.intensityRange + parameterDF$dataScalingParameters.intensityMin
parameterDF$slope1Param_Estimate <- parameterDF$slope1Param_N_Estimate / parameterDF$dataScalingParameters.timeRange
parameterDF$midPoint1Param_Estimate <- parameterDF$midPoint1Param_N_Estimate * parameterDF$dataScalingParameters.timeRange
parameterDF$slope2Param_Estimate <- parameterDF$slope2Param_N_Estimate / parameterDF$dataScalingParameters.timeRange
parameterDF$midPointDistanceParam_Estimate <- parameterDF$midPointDistanceParam_N_Estimate * parameterDF$dataScalingParameters.timeRange
parameterDF$h0_Estimate <- parameterDF$h0_N_Estimate * parameterDF$dataScalingParameters.intensityRange + parameterDF$dataScalingParameters.intensityMin
}
if(!isalist){
parameterDF$finalAsymptoteIntensityRatio_Estimate <- parameterDF$finalAsymptoteIntensityRatio_N_Estimate
parameterDF$maximum_Estimate <- parameterDF$maximum_N_Estimate
parameterDF$slope1Param_Estimate <- parameterDF$slope1Param_N_Estimate
parameterDF$midPoint1Param_Estimate <- parameterDF$midPoint1Param_N_Estimate
parameterDF$slope2Param_Estimate <- parameterDF$slope2Param_N_Estimate
parameterDF$midPointDistanceParam_Estimate <- parameterDF$midPointDistanceParam_N_Estimate
parameterDF$h0_Estimate <- parameterDF$h0_N_Estimate
}
return(parameterDF)
}
#**************************************
Any scripts or data that you put into this service are public.
Add the following code to your website.
For more information on customizing the embed code, read Embedding Snippets.