Nothing
#' @export
rxUiGet.saemMuRefCovariateDataFrame <- function(x, ...) {
.ui <- x[[1]]
if (!rxode2::rxGetControl(x[[1]], "muRefCov", getOption("nlmixr2.saemMuRefCov", TRUE))) {
return(data.frame(theta = character(0), covariate = character(0), covariateParameter = character(0)))
}
if (exists("muRefFinal", .ui)) {
.cov <- .ui$muRefFinal
} else {
.cov <- .ui$muRefCovariateDataFrame
}
.iniDf <- .ui$iniDf
.rm <- NULL
for (.i in seq_along(.cov$covariateParameter)) {
.cp <- .cov$covariateParameter[.i]
.w <- which(.iniDf$name == .cp)
if (.iniDf$fix[.w]) {
.rm <- c(.rm, -.i)
}
}
if (!is.null(.rm)) {
.cov <- .cov[.rm, ]
}
if (.saemGeneralLik(.ui) && length(.cov$theta) > 0L) {
# refinePhi0Lik only refines an intercept-only phi0, so a covariate on a
# theta without an eta (a <- add.sd + WT*cov.sd) stays in the model instead
.cov <- .cov[.cov$theta %in% .ui$muRefDataFrame$theta, , drop = FALSE]
}
.cov
}
attr(rxUiGet.saemMuRefCovariateDataFrame, "rstudio") <- NA
#' @export
rxUiGet.saemInParsAndMuRefCovariates <- function(x, ...) {
.ui <- x[[1]]
# mu ref final removes time varying covariates
.muRefFinal <- rxUiGet.saemMuRefCovariateDataFrame(x, ...)
.cov <- .ui$covariates
.muCov <- unique(.muRefFinal$covariate)
.cov <- .cov[!(.cov %in% .muCov)]
if (length(.ui$predDf$cond) > 1) {
.cov <- unique(c("CMT", .cov))
}
# A general-likelihood endpoint's own saem_mod formula references DV
# directly (rx_pred_=llikNorm(DV, ...) or a literal ll() expression), the
# same way predNoLhs/innerHess2 do -- but DV is never part of .ui$covariates
# (rxode2's etTran.cpp deliberately excludes any "dv"-named column from
# covariate detection, see CLAUDE.md), so it would otherwise never appear
# here. .configsaem's own DV-append step (R/saem_fit.R) only re-appends DV
# to the solve data when inPars already lists it -- omitting it here means
# DV silently never reaches saem_mod's par_ptr for a general-lik fit.
if (.saemGeneralLik(.ui)) {
.cov <- unique(c(.cov, "DV"))
}
list(inPars = .cov, covars = .muCov)
}
#attr(rxUiGet.saemInParsAndMuRefCovariates, "desc") <- "Get inPars and covars for saem"
#' @export
rxUiGet.saemInPars <- function(x, ...) {
.ret <- rxUiGet.saemInParsAndMuRefCovariates(x, ...)
.ret$inPars
}
#attr(rxUiGet.saemInPars, "desc") <- "get inPars"
attr(rxUiGet.saemInPars, "rstudio") <- "char"
#' @export
rxUiGet.saemCovars <- function(x, ...) {
.ret <- rxUiGet.saemInParsAndMuRefCovariates(x, ...)
.ret$covars
}
#attr(rxUiGet.saemInPars, "desc") <- "get saemn mu-referenced non-time varying covariates"
attr(rxUiGet.saemCovars, "rstudio") <- "char"
#' @export
rxUiGet.saemFunctionModPredQuote <- function(x, ...) {
.mod <- .nlmixr2estRxode2(rxUiGet.saemModel(x, ...), "rxSaem")
.fnPred <- bquote(function(a, b, c) {
rxode2::rxLoad(.(.mod))
rxode2::rxLock(.(.mod))
rxode2::rxAllowUnload(FALSE)
on.exit({
rxode2::rxUnlock(.(.mod))
rxode2::rxAllowUnload(TRUE)
rxode2::rxSolveFree()
})
.Call(`_nlmixr2est_saem_do_pred`, a, b, c)
})
list(.mod, .fnPred)
}
#' @export
rxUiGet.saemFunction <- function(x, ...) {
# This function depends on the number of time varying covariates in the data
.ui <- x[[1]]
.mod <- rxUiGet.saemFunctionModPredQuote(x, ...)
.fnPred <- .mod[[2]]
.mod <- .mod[[1]]
.fn <- bquote(function(a, b, c) {
rxode2::rxLoad(.(.mod))
rxode2::rxLock(.(.mod))
on.exit({
rxode2::rxUnlock(.(.mod))
rxode2::rxAllowUnload(TRUE)
rxode2::rxSolveFree()
})
if (missing(b) && missing(c)) {
.ret <- .Call(`_nlmixr2est_saem_fit`, a, PACKAGE = "nlmixr2est")
attr(.ret, "dopred") <- .(.fnPred)
return(.ret)
} else {
.curFn <- .(.fnPred)
return(.curFn(a, b, c))
}
})
.inPars <- rxUiGet.saemInPars(x, ...)
.param <- rxode2::rxParam(.mod)
.estParam <- rxUiGet.saemParamsToEstimate(x, ...)
.cov <- rxUiGet.saemMuRefCovariateDataFrame(x, ...)
if (length(.cov$covariateParameter) > 0) {
.param <- .param[!(.param %in% .cov$covariateParameter)]
.estParam <- .estParam[!(.estParam %in% .cov$covariateParameter)]
}
.parmUpdate <- vapply(
.param,
function(x) {
if (x %in% .estParam) {
return(1L)
} else {
return(0L)
}
},
integer(1),
USE.NAMES = FALSE
)
.nendpnt <- length(.ui$predDf$cond)
.fn <- eval(.fn)
attr(.fn, "form") <- "ode" ## Not sure this is necessary any more
attr(.fn, "neq") <- length(rxode2::rxState(.mod))
attr(.fn, "nlhs") <- length(rxode2::rxLhs(.mod))
attr(.fn, "nrhs") <- sum(.parmUpdate)
attr(.fn, "paramUpdate") <- .parmUpdate
attr(.fn, "rx") <- .mod
attr(.fn, "inPars") <- .inPars
attr(.fn, "nendpnt") <- .nendpnt
.fn
}
attr(rxUiGet.saemFunction, "rstudio") <- function() {}
#' @export
rxUiGet.saemFixed <- function(x, ...) {
.ui <- x[[1]]
.df <- .ui$iniDf
.dft <- .df[!is.na(.df$ntheta), ]
.dft <- .dft[.saemIsEstimableThetaRow(.ui, .dft), ]
if (length(.ui$mixProbs) > 0) {
.dft <- .dft[!(.dft$name %in% .ui$mixProbs), ]
}
.dft <- setNames(.dft$fix, paste(.dft$name))
.cov <- rxUiGet.saemMuRefCovariateDataFrame(x, ...)
if (length(.cov$theta) > 0) {
.theta <- .dft
.theta <- .theta[!(names(.theta) %in% .cov$covariateParameter)]
.allCovs <- rxUiGet.saemCovars(x, ...)
.lc <- length(.allCovs)
.m <- matrix(rep(NA_character_, .lc * length(.theta)), ncol = .lc)
dimnames(.m) <- list(names(.theta), .allCovs)
for (.c in seq_along(.cov$covariateParameter)) {
.curTheta <- .cov[.c, "theta"]
.curCov <- .cov[.c, "covariate"]
.curPar <- .cov[.c, "covariateParameter"]
.m[.curTheta, .curCov] <- .curPar
}
.m <- cbind(matrix(names(.theta), ncol = 1), .m)
.m <- as.vector(t(.m))
.dft <- .dft[.m[!is.na(.m)]]
}
.extra <- .ui$nonMuEtas
.extra <- setNames(rep(TRUE, length(.extra)), .ui$nonMuEtas)
c(.dft, .extra)
}
#attr(rxUiGet.saemFixed, "desc") <- "Get the saem fixed parameters"
attr(rxUiGet.saemFixed, "rstudio") <- function() {}
.saemEtaTrans <- function(x, ..., nonMu = FALSE) {
.ui <- x[[1]]
.etas <- .ui$iniDf[!is.na(.ui$iniDf$neta1), ]
.etas <- .etas$name[.etas$neta1 == .etas$neta2]
.thetas <- rxUiGet.saemParamsToEstimateCov(x, ...)
if (nonMu) {
.thetas <- .thetas[!(.thetas %in% .ui$nonMuEtas)]
}
.muRefDataFrame <- .ui$muRefDataFrame
vapply(
.etas,
function(eta) {
.w <- which(eta == .muRefDataFrame$eta)
if (length(.w) == 1L) {
.muTheta <- .muRefDataFrame$theta[.w]
.w <- which(.muTheta == .thetas)
if (length(.w) == 1L) return(.w)
}
if (nonMu && eta %in% .ui$nonMuEtas) {
.w <- which(eta == .etas)
if (length(.w) == 1L) return(-.w)
}
.w <- which(eta == .thetas)
if (length(.w) == 1L) {
return(.w)
}
return(NA_integer_)
},
integer(1),
USE.NAMES = FALSE
)
}
#' @export
rxUiGet.saemEtaTrans <- function(x, ...) {
.saemEtaTrans(x, ...)
}
attr(rxUiGet.saemEtaTrans, "rstudio") <- c(1L, 3L)
#' @export
rxUiGet.saemEtaTransPred <- function(x, ...) {
.saemEtaTrans(x, ..., nonMu = TRUE)
}
#attr(rxUiGet.saemEtaTrans, "desc") <- "Get the saem eta to theta translation"
attr(rxUiGet.saemEtaTransPred, "rstudio") <- c(1L, 3L)
#' @export
rxUiGet.saemOmegaTrans <- function(x, ...) {
.etaTrans <- rxUiGet.saemEtaTrans(x, ...)
.o <- order(.etaTrans)
.etaTrans2 <- .etaTrans
.c <- 1
for (i in .o) {
# An eta with no phi parameter (NA) owns no Gamma2_phi1 column, so it gets
# no slot either. Handing it a dense rank is what used to index past the
# matrix in .getSaemOmega() (#1047).
if (is.na(.etaTrans[i])) {
next
}
.etaTrans2[i] <- .c
.c <- .c + 1
}
.etaTrans2
}
#attr(rxUiGet.saemOmegaTrans, "desc") <- "Get the saem omega to UI omega translation"
attr(rxUiGet.saemOmegaTrans, "rstudio") <- c(1L, 3L)
#' @export
rxUiGet.saemOmegaShare <- function(x, ...) {
.ui <- x[[1]]
.etaNames <- rxUiGet.saemEtaNames(x, ...)
.ret <- rep(0L, length(.etaNames))
if (length(.ui$mixProbs) == 0L) {
return(.ret)
}
.allEtas <- .ui$iniDf[!is.na(.ui$iniDf$neta1), ]
.allEtas <- .allEtas[.allEtas$neta1 == .allEtas$neta2, "name"]
.mixCalls <- do.call(c, lapply(.ui$lstExpr, .findMixCalls))
if (length(.mixCalls) == 0L) {
return(.ret)
}
.groupId <- 1L
for (.mc in .mixCalls) {
.args <- as.list(.mc)[-1]
.comps <- .args[seq(1, length(.args), by = 2)]
.grpEtas <- unique(unlist(lapply(.comps, .extractEtas, etas = .allEtas)))
if (length(.grpEtas) > 1L) {
for (.eta in .grpEtas) {
.w <- which(.eta == .etaNames)
if (length(.w) == 1L) {
.ret[.w] <- .groupId
}
}
.groupId <- .groupId + 1L
}
}
.ret
}
#' @export
rxUiGet.saemOmegaShareSubpop <- function(x, ...) {
.ui <- x[[1]]
.etaNames <- rxUiGet.saemEtaNames(x, ...)
.ret <- rep(0L, length(.etaNames))
if (length(.ui$mixProbs) == 0L) {
return(.ret)
}
.allEtas <- .ui$iniDf[!is.na(.ui$iniDf$neta1), ]
.allEtas <- .allEtas[.allEtas$neta1 == .allEtas$neta2, "name"]
.mixCalls <- do.call(c, lapply(.ui$lstExpr, .findMixCalls))
if (length(.mixCalls) == 0L) {
return(.ret)
}
# Collect, for each eta, every component index it is referenced from, across
# ALL of the model's mix() calls. An eta seen from more than one component --
# or used anywhere outside a component, where it applies to every component --
# is SHARED, not owned: marking it as a component (the assignment below used
# to take whichever mix() call mentioned it last) sends a shared-eta mixture
# down the split-ETA code paths, which weight that eta's theta/omega update by
# a single component's responsibilities.
.outside <- unique(unlist(lapply(.ui$lstExpr, .extractEtasOutsideMix, etas = .allEtas)))
.compsOf <- list()
for (.mc in .mixCalls) {
.args <- as.list(.mc)[-1]
.comps <- .args[seq(1, length(.args), by = 2)]
for (.j in seq_along(.comps)) {
for (.eta in .extractEtas(.comps[[.j]], etas = .allEtas)) {
.compsOf[[.eta]] <- unique(c(.compsOf[[.eta]], .j))
}
}
}
for (.eta in names(.compsOf)) {
if (length(.compsOf[[.eta]]) != 1L || .eta %in% .outside) {
next
}
.w <- which(.eta == .etaNames)
if (length(.w) == 1L) {
.ret[.w] <- .compsOf[[.eta]]
}
}
.ret
}
#' @export
rxUiGet.saemModelOmega <- function(x, ...) {
.ui <- x[[1]]
.thetas <- rxUiGet.saemParamsToEstimateCov(x, ...)
.etaTrans <- rxUiGet.saemEtaTrans(x, ...)
.dm <- length(.thetas)
.mat <- matrix(rep(0, .dm * .dm), .dm)
.iniDf <- .ui$iniDf
.etd <- .iniDf[which(!is.na(.iniDf$neta1)), ]
for (i in seq_along(.etd$neta1)) {
.mat[.etaTrans[.etd$neta1[i]], .etaTrans[.etd$neta2[i]]] <-
.mat[.etaTrans[.etd$neta2[i]], .etaTrans[.etd$neta1[i]]] <- 1
}
.mat
}
#attr(rxUiGet.saemModelOmega, "desc") <- "Get the saem model omega"
attr(rxUiGet.saemModelOmega, "rstudio") <- lotri::lotri(a + b ~ c(1, 0.1, 1))
#' @export
rxUiGet.saemModelOmegaFixed <- function(x, ...) {
.ui <- x[[1]]
.thetas <- rxUiGet.saemParamsToEstimateCov(x, ...)
.etaTrans <- rxUiGet.saemEtaTrans(x, ...)
.dm <- length(.thetas)
.mat <- matrix(rep(0, .dm * .dm), .dm)
.iniDf <- .ui$iniDf
.etd <- .iniDf[which(!is.na(.iniDf$neta1)), ]
for (i in seq_along(.etd$neta1)) {
.mat[.etaTrans[.etd$neta1[i]], .etaTrans[.etd$neta2[i]]] <-
.mat[.etaTrans[.etd$neta2[i]], .etaTrans[.etd$neta1[i]]] <- as.integer(.etd$fix[i])
}
.mat
}
#attr(rxUiGet.saemModelOmegaFixed, "desc") <- "Get the indicator for saem model omega fixed components"
attr(rxUiGet.saemModelOmegaFixed, "rstudio") <- lotri::lotri(a + b ~ c(1, 0.1, 1))
#' @export
rxUiGet.saemModelOmegaFixedValues <- function(x, ...) {
.ui <- x[[1]]
.thetas <- rxUiGet.saemParamsToEstimateCov(x, ...)
.etaTrans <- rxUiGet.saemEtaTrans(x, ...)
.dm <- length(.thetas)
.mat <- matrix(rep(0, .dm * .dm), .dm)
.iniDf <- .ui$iniDf
.etd <- .iniDf[which(!is.na(.iniDf$neta1)), ]
for (i in seq_along(.etd$neta1)) {
.mat[.etaTrans[.etd$neta1[i]], .etaTrans[.etd$neta2[i]]] <-
.mat[.etaTrans[.etd$neta2[i]], .etaTrans[.etd$neta1[i]]] <- .etd$est[i]
}
.mat
}
#attr(rxUiGet.saemModelOmegaFixedValues, "desc") <- "Get the omega values may be fixed"
attr(rxUiGet.saemModelOmegaFixedValues, "rstudio") <- lotri::lotri(a + b ~ c(1, 0.1, 1))
#' @export
rxUiGet.saemLow <- function(x, ...) {
.ui <- x[[1]]
.ui$predDf$trLow
}
#attr(rxUiGet.saemLow, "desc") <- "Get the saem error transformation lower boundary"
attr(rxUiGet.saemLow, "rstudio") <- -Inf
#' @export
rxUiGet.saemHi <- function(x, ...) {
.ui <- x[[1]]
.ui$predDf$trHi
}
#attr(rxUiGet.saemHi, "desc") <- "Get the saem error transformation higher boundary"
attr(rxUiGet.saemHi, "rstudio") <- Inf
#' @export
rxUiGet.saemPropT <- function(x, ...) {
.ui <- x[[1]]
as.integer((.ui$predDf$errTypeF == "transformed") * 1L)
}
#attr(rxUiGet.saemPropT, "desc") <- "Get the saem transformation type for the function"
attr(rxUiGet.saemPropT, "rstudio") <- 1L
#' @export
rxUiGet.saemYj <- function(x, ...) {
.ui <- x[[1]]
as.integer(.ui$predDf$transform) - 1
}
#attr(rxUiGet.saemYj, "desc") <- "Get the saem transformation type"
attr(rxUiGet.saemYj, "rstudio") <- 1L
#' @export
rxUiGet.saemResMod <- function(x, ...) {
.ui <- x[[1]]
.predDf <- .ui$predDf
# A "norm" condition mixed with a genuine general-likelihood condition is
# promoted to "dnorm" at model-generation time (.createSaemLineObject,
# R/saemRxUiGetModel.R) and compiled with rx_r_ ~ 0 like every other
# general-likelihood endpoint -- so once the WHOLE fit is general-lik
# (.saemGeneralLik(.ui)), every condition has no residual error parameter,
# not just the ones whose own raw distribution already said so.
.genLik <- .saemGeneralLik(.ui)
vapply(
seq_along(.predDf$errType),
function(i) {
# general log-likelihood endpoint (ll() ~ expr, or any non-normal
# distribution() family): no residual error parameter -- the inner
# supplies the likelihood (distribution=4 path)
if (.genLik || .predDf$distribution[i] != "norm") {
return(0L)
}
.errType <- as.integer(.predDf$errType[i])
.hasLambda <- .predDf$transform[i] %in%
c("boxCox", "yeoJohnson", "logit + yeoJohnson", "probit + yeoJohnson", "logit + boxCox", "probit + boxCox")
if (.hasLambda) {
return(.errType + 5L)
} else {
return(.errType)
}
},
integer(1),
USE.NAMES = FALSE
)
}
#attr(rxUiGet.saemResMod, "desc") <- "saem res.mod component"
attr(rxUiGet.saemResMod, "rstudio") <- c(1L, 2L)
#' @export
rxUiGet.saemArActive <- function(x, ...) {
.ui <- x[[1]]
.predDf <- .ui$predDf
.iniDf <- .ui$iniDf
# Some rxode2 versions omit the predDf `ar` column entirely, expressing AR
# only through iniDf (err == "ar"); index it only when it is actually present
# for this condition, otherwise fall back to iniDf.
.ar <- .predDf$ar
vapply(
seq_along(.predDf$cond),
function(i) {
if (length(.ar) >= i && !is.na(.ar[i])) {
return(1L)
}
if (isTRUE(any(.iniDf$err == "ar" & .iniDf$condition == .predDf$cond[i]))) {
return(1L)
}
0L
},
integer(1),
USE.NAMES = FALSE
)
}
attr(rxUiGet.saemArActive, "rstudio") <- 0L
#' @export
rxUiGet.saemArCor <- function(x, ...) {
.ui <- x[[1]]
.predDf <- .ui$predDf
.iniDf <- .ui$iniDf
.ar <- .predDf$ar
vapply(
seq_along(.predDf$cond),
function(i) {
if (length(.ar) >= i && !is.na(.ar[i])) {
.v <- suppressWarnings(as.numeric(.ar[i]))
if (!is.na(.v)) return(.v)
}
.w <- which(.iniDf$err == "ar" & .iniDf$condition == .predDf$cond[i])
if (length(.w) == 1L) {
return(.iniDf$est[.w])
}
0.0
},
numeric(1),
USE.NAMES = FALSE
)
}
attr(rxUiGet.saemArCor, "rstudio") <- 0.0
#' @export
rxUiGet.saemModNumEst <- function(x, ...) {
.resMod <- rxUiGet.saemResMod(x, ...)
vapply(
.resMod,
function(i) {
if (i == 0L) {
return(0L)
} # general log-likelihood endpoint: no residual params
switch(
i,
1L, # add = 1
1L, # prop = 2
2L, # pow = 3
2L, # add + prop = 4
3L, # add + pow = 5
2L, # add + lambda = 6
2L, # prop + lambda = 7
3L, # pow + lambda = 8
3L, # add + prop + lambda = 9
4L # add + pow + lambda = 10
)
},
integer(1),
USE.NAMES = TRUE
)
}
#attr(rxUiGet.saemModNumEst, "desc") <- "saem number of parameters that can be estimated for each component"
attr(rxUiGet.saemModNumEst, "rstudio") <- c(1L, 2L)
#' @export
rxUiGet.saemModResOffset <- function(x, ...) {
# res_offset
cumsum(c(0, rxUiGet.saemModNumEst(x, ...)))
}
#attr(rxUiGet.saemModResOffset, "desc") <- "saem residual parameters offset"
attr(rxUiGet.saemModResOffset, "rstudio") <- c(1, 2)
#' @export
rxUiGet.saemModResTotalResiduals <- function(x, ...) {
# res_offset
sum(rxUiGet.saemModNumEst(x, ...))
}
#attr(rxUiGet.saemModResTotalResiduals, "desc") <- "saem total number of residuals"
attr(rxUiGet.saemModResTotalResiduals, "rstudio") <- c(1, 2)
#' Get the ares name based on the condition
#'
#' @param iniDf Ini data frame
#' @param cond Condition of the parameter
#' @param types The types of errors that support this distribution
#' @return Name of the ares for saem
#' @author Matthew L. Fidler
#' @noRd
.saemGetIniDfResNameFromType <- function(iniDf, cond, types, column = "name") {
.ini <- iniDf[iniDf$condition == cond, ]
.ini[.ini$err %in% types, column]
}
.saemGetIniDfAResName <- function(iniDf, cond, column = "name") {
.saemGetIniDfResNameFromType(iniDf, cond, c("add", "lnorm", "probitNorm", "logitNorm"), column = column)
}
.saemGetIniDfBResName <- function(iniDf, cond, column = "name") {
.saemGetIniDfResNameFromType(iniDf, cond, c("prop", "propT", "propF", "pow", "powT", "powF"), column = column)
}
.saemGetIniDfCResName <- function(iniDf, cond, column = "name") {
.saemGetIniDfResNameFromType(iniDf, cond, c("pow2", "powF2", "powT2"), column = column)
}
.saemGetIniDfLResName <- function(iniDf, cond, column = "name") {
.saemGetIniDfResNameFromType(iniDf, cond, c("boxCox", "yeoJohnson"), column = column)
}
#' Get saem residual item
#' @param ui rxode2 UI
#' @param column column name
#' @return Residual names
#' @author Matthew L. Fidler
#' @noRd
.saemGetResItem <- function(ui, column = "name") {
.predDf <- ui$predDf
.iniDf <- ui$iniDf
.numEst <- ui$saemModNumEst
.resMod <- ui$saemResMod
do.call(
"c",
lapply(seq_along(.numEst), function(i) {
.num <- .numEst[i]
.cond <- .predDf$cond[i]
.ret <- switch(
.resMod[i],
.saemGetIniDfAResName(.iniDf, .cond, column), # add = 1
.saemGetIniDfBResName(.iniDf, .cond, column), # prop = 2
c(.saemGetIniDfBResName(.iniDf, .cond, column), .saemGetIniDfCResName(.iniDf, .cond, column)), # pow = 3
c(.saemGetIniDfAResName(.iniDf, .cond, column), .saemGetIniDfBResName(.iniDf, .cond, column)), # add + prop = 4
c(
.saemGetIniDfAResName(.iniDf, .cond, column),
.saemGetIniDfBResName(.iniDf, .cond, column),
.saemGetIniDfCResName(.iniDf, .cond, column)
), # add + pow = 5
# add + lambda = 6
c(.saemGetIniDfAResName(.iniDf, .cond, column), .saemGetIniDfLResName(.iniDf, .cond, column)),
# prop + lambda = 7
c(.saemGetIniDfBResName(.iniDf, .cond, column), .saemGetIniDfLResName(.iniDf, .cond, column)),
c(
.saemGetIniDfBResName(.iniDf, .cond, column),
.saemGetIniDfCResName(.iniDf, .cond, column),
.saemGetIniDfLResName(.iniDf, .cond, column)
), # pow + lambda = 8
c(
.saemGetIniDfAResName(.iniDf, .cond, column),
.saemGetIniDfBResName(.iniDf, .cond, column),
.saemGetIniDfLResName(.iniDf, .cond, column)
), # add + prop + lambda = 9
c(
.saemGetIniDfAResName(.iniDf, .cond, column),
.saemGetIniDfBResName(.iniDf, .cond, column),
.saemGetIniDfCResName(.iniDf, .cond, column),
.saemGetIniDfLResName(.iniDf, .cond, column)
) # add + pow + lambda = 10
)
if (.num != length(.ret)) {
stop("endpoint '", .cond, "' for saem cannot locate the residual error(s) correctly", call. = FALSE)
}
.ret
})
)
}
#' @export
rxUiGet.saemResNames <- function(x, ...) {
.ui <- x[[1]]
.saemGetResItem(.ui, column = "name")
}
#attr(rxUiGet.saemResNames, "desc") <- "Get error names for SAEM"
attr(rxUiGet.saemResNames, "rstudio") <- "add.err"
#' @export
rxUiGet.saemResFixed <- function(x, ...) {
.ui <- x[[1]]
as.integer(.saemGetResItem(.ui, column = "fix"))
}
#attr(rxUiGet.saemResFixed, "desc") <- "Integer vector of residual fixed components"
attr(rxUiGet.saemResFixed, "rstudio") <- c(1L, 2L)
#' @export
rxUiGet.saemParHistResNames <- function(x, ...) {
.fix <- (rxUiGet.saemResFixed(x, ...) == 0L)
rxUiGet.saemResNames(x, ...)[.fix]
}
#attr(rxUiGet.saemParHistResNames, "desc") <-"Get the SAEM parameter history residual names"
attr(rxUiGet.saemParHistResNames, "rstudio") <- "add.sd"
#' @export
rxUiGet.saemResValue <- function(x, ...) {
.ui <- x[[1]]
.saemGetResItem(.ui, column = "est")
}
attr(rxUiGet.saemResValue, "rstudio") <- 0.7
#' @export
rxUiGet.saemEtaNames <- function(x, ...) {
.ui <- x[[1]]
.etaNames <- .ui$iniDf[!is.na(.ui$iniDf$neta1), ]
.etaNames <- .etaNames[.etaNames$neta1 == .etaNames$neta2, "name"]
## .etaTrans <- rxUiGet.saemOmegaTrans(x, ...)
.etaTrans <- rxUiGet.saemEtaTrans(x, ...)
.names <- rxUiGet.saemParamsToEstimateCov(x, ...)
.names <- rep("", length(.names))
for (.i in seq_along(.etaTrans)) {
.names[.etaTrans[.i]] <- .etaNames[.i]
}
.names <- .names[.names != ""]
.names
}
#attr(rxUiGet.saemParHistEtaNames, "desc") <- "Get ETA names for SAEM based on theta order"
attr(rxUiGet.saemEtaNames, "rstudio") <- "eta.names"
#' @export
rxUiGet.saemParHistOmegaKeep <- function(x, ...) {
.ui <- x[[1]]
.etaNames <- .ui$iniDf[!is.na(.ui$iniDf$neta1), ]
.etaNames <- .etaNames[.etaNames$neta1 == .etaNames$neta2, ]
.names <- rxUiGet.saemEtaNames(x, ...)
vapply(
.names,
function(etaName) {
.w <- which(.etaNames$name == etaName)
if (length(.w) == 1) {
return(1L - as.integer(.etaNames$fix[.w]))
} else {
stop("cannot figure out saemParHistOmegaKeep", call. = FALSE)
}
},
integer(1)
)
}
#attr(rxUiGet.saemOmegaKeep, "desc") <- "Get the etas that are kept for SAEM based on theta order"
attr(rxUiGet.saemParHistOmegaKeep, "rstudio") <- c("eta.ka" = 1)
#' @export
rxUiGet.saemParHistEtaNames <- function(x, ...) {
.names <- rxUiGet.saemParHistOmegaKeep(x, ...)
.names <- .names[.names == 1L]
if (length(.names) == 0) {
return(NULL)
}
paste0("V(", names(.names), ")")
}
#attr(rxUiGet.saemParHistEtaNames, "desc") <- "Get the parameter history eta names"
attr(rxUiGet.saemParHistEtaNames, "rstudio") <- "V(ka)"
#' Off-diagonal Omega covariances recorded in the SAEM iteration history
#'
#' Returns the estimated (non-fixed) off-diagonal Omega elements, mapped to the eta
#' order used by `Gamma2_phi1` (0-indexed row/col pairs for the C++ layer), and their
#' `cov.<eta>.<eta>` names (matching the covariance-matrix naming).
#' @return list(pairs = integer matrix (row, col), names = character)
#' @noRd
#' @export
rxUiGet.saemParHistOmegaOffInfo <- function(x, ...) {
.ui <- x[[1]]
.idf <- .ui$iniDf
.etaN <- rxUiGet.saemEtaNames(x, ...) # phi1 / Gamma2_phi1 order
.diag <- .idf[!is.na(.idf$neta1) & .idf$neta1 == .idf$neta2, , drop = FALSE]
.num2name <- stats::setNames(.diag$name, .diag$neta1)
.off <- .idf[!is.na(.idf$neta1) & .idf$neta1 != .idf$neta2 & !.idf$fix, , drop = FALSE]
if (nrow(.off) == 0L) {
return(list(pairs = matrix(integer(0), ncol = 2L), names = character(0)))
}
.pi <- match(.num2name[as.character(.off$neta1)], .etaN)
.pj <- match(.num2name[as.character(.off$neta2)], .etaN)
.hi <- pmax(.pi, .pj)
.lo <- pmin(.pi, .pj)
.ord <- order(.hi, .lo)
.hi <- .hi[.ord]
.lo <- .lo[.ord]
list(pairs = cbind(.hi - 1L, .lo - 1L), names = paste0("cov.", .etaN[.hi], ".", .etaN[.lo]))
}
attr(rxUiGet.saemParHistOmegaOffInfo, "rstudio") <- "off-diagonal omega history"
#' @export
rxUiGet.saemParHistNames <- function(x, ...) {
#join_cols(join_cols(Plambda, Gamma2_phi1.diag()), vcsig2).t();
.plambda <- rxUiGet.saemParamsToEstimate(x, ...)
.plambda <- .plambda[!rxUiGet.saemFixed(x, ...)]
.ui <- x[[1]]
c(
.plambda,
rxUiGet.saemParHistEtaNames(x, ...),
rxUiGet.saemParHistOmegaOffInfo(x, ...)$names,
rxUiGet.saemParHistResNames(x, ...),
.ui$mixProbs
)
}
attr(rxUiGet.saemParHistNames, "rstudio") <- c("ka", "add.sd")
#' @export
rxUiGet.saemMixProb <- function(x, ...) {
.ui <- x[[1]]
if (length(.ui$mixProbs) > 0) {
.probs <- .ui$theta[.ui$mixProbs]
# Defense-in-depth: upstream checks (rxode2::mix(), rxUiGet.thetaIniMix)
# don't catch an individually out-of-range probability whose sum still
# looks valid, so re-validate before deriving the implicit last component.
if (any(.probs < 0) || any(.probs > 1) || sum(.probs) > 1) {
stop(
"initial mixture probabilities must each be in [0, 1] and sum to no more than 1 (got: ",
paste(signif(.probs, 3), collapse = ", "),
")",
call. = FALSE
)
}
# To get the full simplex, we append 1 - sum(probs)
return(c(.probs, 1.0 - sum(.probs)))
} else {
return(numeric(0))
}
}
#' @export
rxUiGet.saemNMix <- function(x, ...) {
.ui <- x[[1]]
if (length(.ui$mixProbs) > 0) {
return(length(.ui$mixProbs) + 1L)
} else {
return(1L)
}
}
#' @export
rxUiGet.saemAres <- function(x, ...) {
.ui <- x[[1]]
.predDf <- .ui$predDf
.ini <- .ui$iniDf
.ini <- .ini[!is.na(.ini$err), ]
vapply(
.predDf$cond,
function(x) {
.tmp <- .ini[which(.ini$condition == x), ]
.w <- which(vapply(
.tmp$err,
function(x) {
x %in%
c(
"add",
"norm",
"dnorm",
"dpois",
"pois",
"dbinom",
"binom",
"dbern",
"bern",
"lnorm",
"dlnorm",
"logn",
"dlogn",
"logitNorm",
"probitNorm"
)
},
logical(1),
USE.NAMES = FALSE
))
if (length(.w) == 1) {
return(.tmp$est[.w])
} else {
return(10)
}
},
numeric(1),
USE.NAMES = FALSE
)
}
#attr(rxUiGet.saemAres, "desc") <- "ares initial estimates for saem"
attr(rxUiGet.saemAres, "rstudio") <- 0.7
#' @export
rxUiGet.saemBres <- function(x, ...) {
.ui <- x[[1]]
.predDf <- .ui$predDf
.ini <- .ui$iniDf
.ini <- .ini[!is.na(.ini$err), ]
vapply(
.predDf$cond,
function(x) {
.tmp <- .ini[which(.ini$condition == x), ]
.w <- which(vapply(
.tmp$err,
function(x) (any(x == "prop") || any(x == "propT") || any(x == "propF")),
logical(1),
USE.NAMES = FALSE
))
if (length(.w) == 1) {
return(.tmp$est[.w])
} else {
.w <- which(vapply(
.tmp$err,
function(x) (any(x == "pow") || any(x == "powT") || any(x == "powF")),
logical(1),
USE.NAMES = FALSE
))
if (length(.w) == 1) {
return(.tmp$est[.w])
} else {
return(1)
}
}
},
numeric(1),
USE.NAMES = FALSE
)
}
#attr(rxUiGet.saemBres, "desc") <- "bres initial estimates for saem"
attr(rxUiGet.saemBres, "rstudio") <- 0.7
#' @export
rxUiGet.saemCres <- function(x, ...) {
.ui <- x[[1]]
.predDf <- .ui$predDf
.ini <- .ui$iniDf
.ini <- .ini[!is.na(.ini$err), ]
vapply(
.predDf$cond,
function(x) {
.tmp <- .ini[which(.ini$condition == x), ]
.w <- which(vapply(
.tmp$err,
function(x) (any(x == "pow2") || any(x == "powT2") || any(x == "powF2")),
logical(1),
USE.NAMES = FALSE
))
if (length(.w) == 1) {
return(.tmp$est[.w])
} else {
return(1)
}
},
numeric(1),
USE.NAMES = FALSE
)
}
#attr(rxUiGet.saemCres, "desc") <- "cres initial estimates for saem"
attr(rxUiGet.saemCres, "rstudio") <- 0.7
#' @export
rxUiGet.saemLres <- function(x, ...) {
.ui <- x[[1]]
.predDf <- .ui$predDf
.ini <- .ui$iniDf
.ini <- .ini[!is.na(.ini$err), ]
vapply(
.predDf$cond,
function(x) {
.tmp <- .ini[which(.ini$condition == x), ]
.boxCox <- which(.tmp$err == "boxCox")
if (length(.boxCox) == 1L) {
return(.tmp$est[.boxCox])
}
.yeoJohnson <- which(.tmp$err == "yeoJohnson")
if (length(.yeoJohnson) == 1L) {
return(.tmp$est[.yeoJohnson])
}
return(1.0)
},
numeric(1),
USE.NAMES = FALSE
)
}
#attr(rxUiGet.saemLres, "desc") <- "lres (lambda) initial estimates for saem"
attr(rxUiGet.saemLres, "rstudio") <- 0.7
#' @export
rxUiGet.saemLogEta <- function(x, ...) {
.ui <- x[[1]]
.thetas <- rxUiGet.saemParamsToEstimate(x, ...)
.ce <- .ui$muRefCurEval
.cov <- rxUiGet.saemMuRefCovariateDataFrame(x, ...)
.thetas <- .thetas[!(.thetas %in% .cov$covariateParameter)]
vapply(
.thetas,
function(x) {
.w <- which(.ce$parameter == x)
if (length(.w) == 1L) {
return(.ce$curEval[.w] == "exp")
}
FALSE
},
logical(1)
)
}
#attr(rxUiGet.saemLogEta, "desc") <- "saem's log.eta for saem"
attr(rxUiGet.saemLogEta, "rstudio") <- c(tka = TRUE)
#' @export
rxUiGet.saemModelList <- function(x, ...) {
.ui <- x[[1]]
.mod <- list(saem_mod = .ui$saemFunction)
.covars <- rxUiGet.saemCovars(x, ...)
if (length(.covars) > 0) {
.mod$covars <- .covars
}
.mod$res.mod <- rxUiGet.saemResMod(x, ...)
.mod$arActive <- rxUiGet.saemArActive(x, ...)
.mod$arCor <- rxUiGet.saemArCor(x, ...)
.mod$log.eta <- rxUiGet.saemLogEta(x, ...)
.mod$ares <- rxUiGet.saemAres(x, ...)
.mod$bres <- rxUiGet.saemBres(x, ...)
.mod$omega <- rxUiGet.saemModelOmega(x, ...)
.mod
}
#attr(rxUiGet.saemModelList "desc") <- "saem's log.eta for saem"
#' @export
rxUiGet.saemInitTheta <- function(x, ...) {
.logEta <- rxUiGet.saemLogEta(x, ...)
.names <- names(.logEta)
.ui <- x[[1]]
.iniDf <- .ui$iniDf
.w <- .saemIsEstimableThetaRow(.ui, .iniDf)
.est <- setNames(.iniDf[.w, "est"], .iniDf[.w, "name"])
.cov <- rxUiGet.saemMuRefCovariateDataFrame(x, ...)
.est <- .est[!(names(.est) %in% .cov$covariateParameter)]
if (length(.ui$mixProbs) > 0) {
.est <- .est[!(names(.est) %in% .ui$mixProbs)]
}
.etaNames <- .iniDf[is.na(.iniDf$ntheta), ]
.etaNames <- .iniDf[.iniDf$neta1 == .iniDf$neta2, "name"]
.fixed <- rxUiGet.saemFixed(x, ...)
.theta <- .fixed
.theta <- .theta[!(names(.theta) %in% .cov$covariateParameter)]
.logEta <- .logEta[!(names(.logEta) %in% .cov$covariateParameter)]
.n <- vapply(.theta, function(x) ifelse(x, "FIXED", ""), character(1), USE.NAMES = FALSE)
.ret <- vapply(
seq_along(.logEta),
function(i) {
.isEta <- any(.names[i] %in% .etaNames)
.curName <- .names[i]
if (.logEta[i]) {
if (.isEta) {
return(1)
} else {
return(exp(.est[.curName]))
}
} else {
if (.isEta) {
return(0)
} else {
return(.est[.curName])
}
}
},
numeric(1),
USE.NAMES = FALSE
)
if (length(.cov$theta) > 0) {
.allCovs <- rxUiGet.saemCovars(x, ...)
.lc <- length(.allCovs)
.m <- matrix(rep(NA, .lc * length(.ret)), ncol = .lc)
dimnames(.m) <- list(names(.theta), .allCovs)
for (.c in seq_along(.cov$theta)) {
.curTheta <- .cov[.c, "theta"]
.curCov <- .cov[.c, "covariate"]
.curPar <- .cov[.c, "covariateParameter"]
.w <- which(.iniDf$name == .curPar)
.est <- .iniDf$est[.w]
.m[.curTheta, .curCov] <- .est
}
.ret <- c(.ret, as.vector(.m))
.ret <- setNames(.ret, c(.n, rep("", length(.ret) - length(.n))))
#.ret <- setNames(, c(.n, rep("", .lc + 1)))
} else {
.ret <- setNames(.ret, .n)
}
.ret
}
#attr(rxUiGet.saemInitTheta, "desc") <- "initialization for saem's theta"
attr(rxUiGet.saemInitTheta, "rstudio") <- c(" " = 1, tcl = 1)
#' @export
rxUiGet.saemInitOmega <- function(x, ...) {
.ui <- x[[1]]
.iniDf <- .ui$iniDf
.eta <- .iniDf[is.na(.iniDf$ntheta), ]
.eta <- .eta[.eta$neta1 == .eta$neta2, ]
.eta <- setNames(.eta$est, .eta$name)
.pars <- rxUiGet.saemParamsToEstimateCov(x, ...)
.ret <- rep(1.0, length(.pars))
.etaTrans <- rxUiGet.saemEtaTrans(x, ...)
for (i in seq_along(.etaTrans)) {
.ret[.etaTrans[i]] <- .eta[i]
}
.ret <- setNames(.ret, .pars)
.cov <- rxUiGet.saemMuRefCovariateDataFrame(x, ...)
if (length(.cov$covariateParameter) > 0) {
.ret <- .ret[!(names(.ret) %in% .cov$covariateParameter)]
}
.ret
}
#attr(rxUiGet.saemInitOmega, "desc") <- "initialization for saem's omega"
attr(rxUiGet.saemInitOmega, "rstudio") <- c(tka = 0.6)
#' @export
rxUiGet.saemInit <- function(x, ...) {
list(theta = rxUiGet.saemInitTheta(x, ...), omega = rxUiGet.saemInitOmega(x, ...))
}
#attr(rxUiGet.saemInit, "desc") <- "initialization for saem's theta and omega"
#' @export
rxUiGet.saemThetaDataFrame <- function(x, ...) {
.ui <- x[[1]]
.theta <- .ui$theta
.fixed <- .ui$iniDf[!is.na(.ui$iniDf$ntheta), "fix"]
data.frame(lower = -Inf, theta = .theta, fixed = .fixed, upper = Inf, row.names = names(.theta))
}
#attr(rxUiGet.saemThetaDataFrame, "desc") <- "Get theta data frame"
attr(rxUiGet.saemThetaDataFrame, "rstudio") <- NA
#' @export
rxUiGet.saemParHistThetaKeep <- function(x, ...) {
1L - as.integer(rxUiGet.saemFixed(x, ...))
}
#attr(rxUiGet.saemParHistThetaKeep, "desc") <- "The thetas that are kept in the parameter history"
attr(rxUiGet.saemParHistThetaKeep, "rstudio") <- c(1L, 2L)
#' @export
rxUiGet.saemAddProp <- function(x, ...) {
.ui <- x[[1]]
.parDf <- .ui$predDf
.addProp <- as.integer(.parDf$addProp)
.w <- which(.addProp == 3L)
if (length(.w) > 0) {
.default <- c(combined1 = 1, combined2 = 2)[rxode2::rxGetControl(.ui, "addProp", "combined2")]
.addProp[.w] <- .default
}
.addProp
}
#attr(rxUiGet.saemParHistThetaKeep, "desc") <- "Get the saem addProp integer vector"
attr(rxUiGet.saemAddProp, "rstudio") <- 2
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.