Package 'SensIAT'

Plot for Estimated Treatment Effect for SensIAT_fulldata_jackknife_results Objects

Description

The horizontal and vertical axes represent the sensitivity parameter alpha for the control and treatment groups, respectively. The plot shows at each combination of alpha values zero if the level confidence interval contains zero, otherwise the bound of the confidence interval that is closest to zero.

Usage

## S3 method for class 'SensIAT_fulldata_jackknife_results'
autoplot(object, level = 0.05, ..., include.rugs = NA)

Arguments

object	A SensIAT_fulldata_jackknife_results object.
level	Significance level of the confidence interval.
...	Additional arguments passed to predict.
include.rugs	If TRUE, adds rugs to the plot. If FALSE, no rugs are added. When NA, rugs are added only if the number of distinct values of alpha_control and alpha_treatment is less than or equal to 10.

Examples

# Note: fitting the jackknife is computationally expensive,
#       so this example is here for reference.
## Not run: 
full.object <-
    fit_SensIAT_fulldata_model(
        data = SensIAT_example_fulldata,
        trt = Treatment_group == "treatment",
        outcome_modeler = fit_SensIAT_single_index_fixed_coef_model,
        id = Subject_ID,
        outcome = Outcome,
        time = Time,
        knots = c(60, 260, 460),
        alpha = c(-0.6, -0.3, 0, 0.3, 0.6)
    )
jk.full.model <- jackknife(full.object, time = 180)
ggplot2::autoplot(jk.full.model)

## End(Not run)

---

Plot for Estimated Treatment Effect for SensIAT_fulldata_model Objects

Description

The horizontal and vertical axes represent the sensitivity parameter alpha for the control and treatment groups, respectively. The contour plot shows the estimated treatment effect at each combination of alpha values.

Usage

## S3 method for class 'SensIAT_fulldata_model'
autoplot(object, time, include.rugs = NA, ...)

Arguments

object	A SensIAT_fulldata_model object.
time	Time at which to plot the estimates.
include.rugs	If TRUE, adds rugs indicating the locations where the sensitivity was evaluated to the plot. If FALSE, no rugs are added. If NA, rugs are added only if the number of distinct values of alpha_control and alpha_treatment is less than or equal to 10.
...	Additional arguments passed to predict.

Value

A ggplot2 object.

Examples

full.object <-
    fit_SensIAT_fulldata_model(
        data = SensIAT_example_fulldata,
        trt = Treatment_group == "treatment",
        outcome_modeler = fit_SensIAT_single_index_fixed_coef_model,
        id = Subject_ID,
        outcome = Outcome,
        time = Time,
        knots = c(60, 260, 460),
        alpha = c(-0.6, -0.3, 0, 0.3, 0.6)
    )
ggplot2::autoplot(full.object, time = 180)

---

Plot a SensIAT_within_group_model Object

Description

This creates a line plot for a SensIAT_within_group_model object. The horizontal axis represents time, and the vertical axis represents the expected marginal outcome given the sensitivity parameter alpha.

Usage

## S3 method for class 'SensIAT_within_group_model'
autoplot(object, ...)

Arguments

object	A SensIAT_within_group_model object.
...	currently ignored

Value

A ggplot2 object.

Examples

# Note: example takes a few seconds to run.

object <-
    fit_SensIAT_within_group_model(
        group.data = SensIAT_example_data,
        outcome_modeler = fit_SensIAT_single_index_fixed_bandwidth_model,
        id = Subject_ID,
        outcome = Outcome,
        time = Time,
        knots = c(60, 260, 460),
        End = 830,
        alpha = c(-0.6, -0.3, 0, 0.3, 0.6),
        intensity.args = list(bandwidth = 30)
    )
ggplot2::autoplot(object) +
    # Title not included
    ggplot2::ggtitle("SensIAT within group model") +
    # Nor are bounds on reasonable values of alpha
    ggplot2::geom_hline(yintercept = c(1.2, 3), linetype = "dotted", linewidth = 1.5)

---

Plot Estimates at Given Times for SensIAT_withingroup_jackknife_results Objects

Description

Horizontal axis represents time, and the vertical axis represents the outcome from the model. Point plotted is the mean estimate, and the error bars show the level confidence interval using the variance estimated from the jackknife.

Usage

## S3 method for class 'SensIAT_withingroup_jackknife_results'
autoplot(object, level = 0.95, width = NULL, ...)

Arguments

object	A SensIAT_withingroup_jackknife_results object produced from SensIAT_jackknife.
level	Significance level of the confidence interval.
width	Width of the dodge for position, default is half the minimum distance between time evaluation points.
...	Ignored.

Value

A ggplot2 object.

Examples

# Note: fitting the jackknife is computationally expensive,
#       so this example is here for reference.
## Not run: 
fitted <-
    fit_SensIAT_within_group_model(
        group.data = SensIAT_example_data,
        outcome_modeler = fit_SensIAT_single_index_fixed_coef_model,
        alpha = c(-0.6, -0.3, 0, 0.3, 0.6),
        id = Subject_ID,
        outcome = Outcome,
        time = Time,
        intensity.args = list(bandwidth = 30),
        knots = c(60, 260, 460),
        End = 830
    )
jackknife.estimates <- SensIAT_jackknife(fitted, time = c(90, 180, 270, 360, 450))
ggplot2::autoplot(jackknife.estimates)

## End(Not run)

---

Benchmark Term2 Integration Methods

Description

Efficiently compares different term2 integration methods by isolating just the term2 computation from the full model fitting. This is useful for performance analysis and method selection without the overhead of repeated model fitting.

Usage

benchmark_term2_methods(
  data,
  id,
  time,
  outcome.model,
  knots,
  spline_degree = 3L,
  alpha = 0,
  impute_data,
  link = c("identity", "log"),
  methods = c("fast", "original", "fixed_grid", "seeded_adaptive"),
  grid_sizes = c(50, 100, 200),
  n_patients = NULL,
  n_iterations = 3,
  reference_method = "fast",
  seed = 42
)

Arguments

data	Data frame with longitudinal observations
id	Unquoted column name for subject identifier
time	Numeric vector of observation times
outcome.model	A fitted outcome model (e.g., from fit_SensIAT_single_index_fixed_coef_model())
knots	Knots for marginal mean model spline basis
spline_degree	Degree of B-spline basis (default: 3)
alpha	Numeric vector of sensitivity parameters to test
impute_data	Function to impute data at arbitrary times: ⁠function(t, patient_data) -> data.frame⁠
link	Link function: "identity" or "log"
methods	Character vector of methods to benchmark. Options: "fast", "original", "fixed_grid", "seeded_adaptive", "gauss_legendre"
grid_sizes	For grid-based methods, vector of grid sizes to test (default: c(50, 100, 200))
n_patients	Number of patients to use for benchmark (NULL = all patients)
n_iterations	Number of timing iterations per method
reference_method	Method to use as accuracy reference (default: "fast")
seed	Random seed for reproducibility

Value

A list with components:

timing

Data frame with timing results per method

accuracy

Data frame with accuracy metrics vs reference

reference_results

List of reference results per patient

setup_info

List with setup parameters and dimensions

Examples

## Not run: 
# Setup data with lag variables
data_with_lags <- SensIAT_example_data |>
  dplyr::group_by(Subject_ID) |>
  dplyr::mutate(
    ..prev_outcome.. = dplyr::lag(Outcome, order_by = Time),
    ..prev_time.. = dplyr::lag(Time, default = 0, order_by = Time),
    ..delta_time.. = Time - dplyr::lag(Time, order_by = Time)
  ) |>
  dplyr::ungroup()

# Fit outcome model
outcome.model <- fit_SensIAT_single_index_fixed_coef_model(
  Outcome ~ splines::ns(..prev_outcome.., df = 3) + ..delta_time.. - 1,
  id = Subject_ID,
  data = data_with_lags |> dplyr::filter(Time > 0)
)

# Imputation function
impute_fn <- function(t, df) {
  extrapolate_from_last_observation(t, df, "Time",
    slopes = c("..delta_time.." = 1))
}

# Run benchmark
results <- benchmark_term2_methods(
  data = data_with_lags,
  id = Subject_ID,
  time = data_with_lags$Time,
  outcome.model = outcome.model,
  knots = c(60, 260, 460),
  alpha = c(-0.1, 0, 0.1),
  impute_data = impute_fn,
  methods = c("fast", "fixed_grid", "seeded_adaptive"),
  grid_sizes = c(50, 100),
  n_patients = 10
)

# View timing results
print(results$timing)

## End(Not run)

---

Compute Influence Terms

Description

This function computes the influence terms for the marginal outcome model sensitivity analysis. It is a generic function that can handle different types of outcome models.

Usage

compute_influence_terms(outcome.model, intensity.model, alpha, data, ...)

## Default S3 method:
compute_influence_terms(
  outcome.model,
  intensity.model,
  alpha,
  data,
  id,
  base,
  ...
)

## S3 method for class 'glm'
compute_influence_terms(
  outcome.model,
  intensity.model,
  alpha,
  data,
  base,
  tolerance = .Machine$double.eps^(1/3),
  na.action = na.fail,
  id = NULL,
  time = NULL,
  ...
)

## S3 method for class ''SensIAT::Single-index-outcome-model''
compute_influence_terms(
  outcome.model,
  intensity.model,
  alpha,
  data,
  base,
  tolerance = .Machine$double.eps^(1/3),
  na.action = na.fail,
  id = NULL,
  time = NULL,
  ...
)

Arguments

outcome.model	The outcome model fitted to the data.
intensity.model	The intensity model fitted to the data.
alpha	A numeric vector representing the sensitivity parameter.
data	A data frame containing the observations.
...	Additional arguments passed to the method.
id	A variable representing the patient identifier.
base	A spline basis object.
tolerance	Numeric value indicating the tolerance for integration, default is .Machine$double.eps^(1/3).
na.action	Function to handle missing values, default is na.fail.
time	Variable indicating the time variable in the data, by Default will be extracted from the intensity model response.

Methods (by class)

• compute_influence_terms(default): Generic method, which throws a not implemented error.

• compute_influence_terms(glm): Method for Generalized Linear Models (GLM).

• compute_influence_terms(`SensIAT::Single-index-outcome-model`): Optimized method for the single index model.

---

Compute Conditional Expected Values based on Outcome Model

Description

Compute Conditional Expected Values based on Outcome Model

Usage

compute_SensIAT_expected_values(
  model,
  alpha = 0,
  new.data = model.frame(model),
  ...
)

## S3 method for class 'lm'
compute_SensIAT_expected_values(model, alpha, new.data, ...)

## S3 method for class 'glm'
compute_SensIAT_expected_values(
  model,
  alpha,
  new.data,
  ...,
  y.max = NULL,
  eps = .Machine$double.eps
)

## S3 method for class 'negbin'
compute_SensIAT_expected_values(
  model,
  alpha,
  new.data,
  ...,
  y.max = NULL,
  eps = .Machine$double.eps^(1/4)
)

Arguments

model	An object representing the output of the outcome model.
alpha	The sensitivity parameter
new.data	Data to compute conditional means for, defaults to the model frame for the fitted model.
...	passed onto methods.
y.max	The maximum value of the outcome variable for the Poisson and Negative Binomial models. If omitted it is chosen from the quantile function for the distribution at 1-eps.
eps	The tolerance for the quantile function used to estimate y.max, default is .Machine$double.eps.

Details

Compute the conditional expectations needed for predictions in the models. Two additional values/expectations are computed:

• ⁠$E \big[ Y(t) \exp \{ \alpha Y(t) \} | A(t)=1, \bar{O}(t) \big]$⁠, returned as E_Yexp_alphaY, and

• ⁠$E \big[ \exp \{ \alpha Y(t) \} \ | A(t)=1, \bar{O}(t) \big]$⁠, returned as E_exp_alphaY.

For the methods shown here

Value

The new.data frame with additional columns alpha, E_Yexp_alphaY, and E_exp_alphaY appended.

Methods (by class)

• compute_SensIAT_expected_values(lm): (Gaussian) Linear Model method Uses closed-form analytical solutions for the conditional expectations. For Y distributed as Normal(mu, sigma^2): E[exp(alpha*Y)] = exp(alpha*mu + alpha^2*sigma^2/2) and E[Y*exp(alpha*Y)] = (mu + alpha*sigma^2) * exp(alpha*mu + alpha^2*sigma^2/2)

• compute_SensIAT_expected_values(glm): Generalized Linear Model method

• compute_SensIAT_expected_values(negbin): Negative Binomial Model method

Examples

model <- lm(mpg ~ as.factor(cyl) + disp + wt, data = mtcars)
compute_SensIAT_expected_values(model, alpha = c(-0.3, 0, 0.3), new.data = mtcars[1:5, ])
model <- glm(cyl ~ mpg + disp + wt, data = mtcars, family = poisson())
compute_SensIAT_expected_values(model, alpha = c(-0.3, 0, 0.3), new.data = mtcars[1:5, ]) |>
    dplyr::mutate("E(y|alpha)" = .data$E_Yexp_alphaY / .data$E_exp_alphaY)

---

Fit SensIAT Marginal Mean Model (Unified)

Description

This function fits a marginal mean model for sensitivity analysis, supporting both single-index (linear) and generalized linear outcome models.

Usage

fit_marginal_model(
  data,
  outcome.model,
  intensity.model,
  id,
  time,
  alpha = 0,
  knots,
  ...
)

Arguments

data	Data frame containing all required variables.
outcome.model	Outcome model object or formula (GLM, single-index, or LM).
intensity.model	Intensity model object (e.g., coxph).
id	Subject identifier variable.
time	Time variable.
alpha	Sensitivity parameter(s).
knots	Spline knot locations.
...	Additional arguments passed to underlying methods.

Value

A fitted SensIAT marginal mean model object.

---

Produce fitted model for group (treatment or control)

Description

Produces a fitted model that may be used to produce estimates of mean and variance for the given group.

Usage

fit_SensIAT_fulldata_model(data, trt, ...)

fit_SensIAT_within_group_model(
  group.data,
  outcome_modeler,
  id,
  outcome,
  time,
  knots,
  alpha = 0,
  End = NULL,
  intensity.args = list(),
  outcome.args = list(),
  influence.args = list(),
  spline.degree = 3,
  add.terminal.observations = TRUE,
  loss = c("lp_mse", "quasi-likelihood"),
  link = c("identity", "log", "logit"),
  term2_method = c("fast", "original", "fixed_grid", "seeded_adaptive", "gauss_legendre"),
  impute_data = NULL
)

Arguments

data	the full data set.
trt	an expression that determine what is treated as the treatment. Everything not treatment is considered control.
...	common arguments passed to fit_SensIAT_within_group_model.
group.data	The data for the group that is being analyzed. Preferably passed in as a single tibble that internally is subsetted/filtered as needed.
outcome_modeler	function for fitting the outcome model. Called with a formula, data argument and outcome.args list.
id	The variable that identifies the patient.
outcome	The variable that contains the outcome.
time	The variable that contains the time.
knots	knot locations for defining the spline basis.
alpha	The sensitivity parameter.
End	The end time for this data analysis, we need to set the default value as the max value of the time.
intensity.args	A list of optional arguments for intensity model. See the Intensity Arguments section.
outcome.args	parameters as needed passed into the outcome_modeler. One special element may be 'model.modifications' which, if present, should be a formula that will be used to modify the outcome model per, update.formula.
influence.args	A list of optional arguments used when computing the influence. See the Influence Arguments section.
spline.degree	The degree of the spline basis.
add.terminal.observations	Logical indicating whether to add terminal observations to the data. If TRUE, data may not contain any NAs. if FALSE, data will be assumed to already include the terminal observations
loss	The loss function to use. Options are "lp_mse" (default) and "quasi-likelihood". Only used when link is not "identity".
link	The link function to use. Options are "identity" (default), "log", and "logit". When "identity", the original analytic formula is used; otherwise, the generalized iterative solver is used.
term2_method	Method for computing term2 influence components. Options are: "fast" (default), "original", "fixed_grid", "seeded_adaptive", "gauss_legendre". Only used when link is not "identity".
impute_data	A function that takes ⁠(t, df)⁠ and returns the imputed data at time t. If NULL (default), a standard imputation function is generated automatically.

Details

This function should be agnostic to whether it is being provided a treatment or control group.

Value

a list with class SensIAT-fulldata-fitted-model with two components, control and treatment, each of which is an independently fitted SensIAT-within-group-fitted-model fit with the fit_within_group_model function.

A SensIAT_within_group_model object containing:

• models: List with intensity and outcome sub-models

• data: The transformed data used for fitting

• variables: List of variable names (id, outcome, time)

• End: The end time for the analysis

• influence: Influence function values by patient

• alpha: Sensitivity parameter values

• coefficients: Estimated spline coefficients for each alpha

• coefficient.variance: Variance of coefficients for each alpha

• base: The SplineBasis object

• V_inverse: Inverse of the Gram matrix

• link: The link function used

• loss: The loss function used

• term2_method: The term2 integration method used

Functions

• fit_SensIAT_fulldata_model(): Fit the Marginal Mean Sensitivity Analysis with Tilting Assumption for Both Treatment and Control Groups.

Intensity Arguments

The intensity.args list may contain the following elements:

• model.modifications A formula that will be used to modify the intensity model from it's default, per update.formula.

• kernel The kernel function for the intensity model. Default is the Epanechnikov kernel.

• bandwidth The bandwidth for the intensity model kernel.

Influence Arguments

The influence.args list may contain the following elements:

• method The method for integrating, adaptive or fixed quadrature. Default is 'adaptive'.

• tolerance The tolerance when using adaptive quadrature.

• delta The bin width for fixed quadrature.

• resolution alternative to delta by specifying the number of bins.

• fix_discontinuity Whether to account for the discontinuity in the influence at observation times.

Examples

# Example 1: Default identity link (original analytic solution)
model <-
    fit_SensIAT_within_group_model(
        group.data = SensIAT_example_data,
        outcome_modeler = fit_SensIAT_single_index_fixed_coef_model,
        alpha = c(-0.6, -0.3, 0, 0.3, 0.6),
        id = Subject_ID,
        outcome = Outcome,
        time = Time,
        End = 830,
        knots = c(60, 260, 460),
    )

# Example 2: Log link with generalized iterative solver
model_log <-
    fit_SensIAT_within_group_model(
        group.data = SensIAT_example_data,
        outcome_modeler = fit_SensIAT_single_index_fixed_coef_model,
        alpha = 0,
        id = Subject_ID,
        outcome = Outcome,
        time = Time,
        End = 830,
        knots = c(60, 260, 460),
        link = "log",
        loss = "lp_mse"
    )

---

Fit the Marginal Means Model

Description

Fit the Marginal Means Model

Usage

fit_SensIAT_marginal_mean_model(
  data,
  id,
  alpha,
  knots,
  outcome.model,
  intensity.model,
  spline.degree = 3L,
  ...
)

Arguments

data	Data for evaluation of the model. Should match the data used to fit the intensity and outcome models.
id	The subject identifier variable in the data. Lazy evaluation is used, so it can be a symbol or a string.
alpha	Sensitivity parameter, a vector of values.
knots	Location of spline knots. If a SplineBasis object is provided, it is used directly.
outcome.model	The observed effects model.
intensity.model	The assessment time intensity model.
spline.degree	The degree of the spline basis, default is 3 (cubic splines).
...	Additional arguments passed to compute_influence_terms.

Value

a list with the fitted model, including the coefficients and their variances for each alpha value.

Examples

# Note: example takes approximately 30 seconds to run.

library(survival)
library(dplyr)
library(splines)
# Create followup data with lags
# added variables `..prev_time..`, `..delta_time..` and `..prev_outcome..`
# have special interpretations when computing the influence.
data_with_lags <- SensIAT_example_data |>
    dplyr::group_by(Subject_ID) |>
    dplyr::mutate(
        ..prev_outcome.. = dplyr::lag(Outcome, default = NA_real_, order_by = Time),
        ..prev_time.. = dplyr::lag(Time, default = 0, order_by = Time),
        ..delta_time.. = Time - dplyr::lag(.data$Time, default = NA_real_, order_by = Time)
    )

# Create the observation time intensity model
intensity.model <-
    coxph(Surv(..prev_time.., Time, !is.na(Outcome)) ~ ..prev_outcome.. + strata(Visit),
        data = data_with_lags |> dplyr::filter(.data$Time > 0)
    )

# Create the observed outcome model
outcome.model <-
    fit_SensIAT_single_index_fixed_coef_model(
        Outcome ~ ns(..prev_outcome.., df = 3) + ..delta_time.. - 1,
        id = Subject_ID,
        data = data_with_lags |> filter(Time > 0)
    )

# Fit the marginal outcome model
mm <- fit_SensIAT_marginal_mean_model(
    data = data_with_lags,
    id = Subject_ID,
    alpha = c(-0.6, -0.3, 0, 0.3, 0.6),
    knots = c(60, 260, 460),
    intensity.model = intensity.model,
    time.vars = c("..delta_time.."),
    outcome.model = outcome.model
)

---

Fit the marginal mean model for generalize outcomes.

Description

This function supports multiple integration methods for term2 computation, including adaptive and fixed-grid approaches. The implementation includes numerical stability improvements (exp(-μ) multiplication vs division) and extensive caching optimizations for repeated expected value computations.

Usage

fit_SensIAT_marginal_mean_model_generalized(
  data,
  time,
  id,
  alpha,
  knots,
  outcome.model,
  intensity.model,
  impute_data,
  loss = c("lp_mse", "quasi-likelihood"),
  link = c("identity", "log", "logit"),
  spline.degree = 3L,
  ...,
  BBsolve.control = list(maxit = 1000, tol = 1e-06, trace = FALSE),
  term2_method = c("fast", "original", "fixed_grid", "seeded_adaptive", "gauss_legendre"),
  term2_grid_n = 100,
  use_expected_cache = TRUE
)

Arguments

data	Data for evaluation of the model. Should match the data used to fit the intensity and outcome models.
time	The time variable in the data. Can be provided as a column name or vector.
id	The subject identifier variable in the data. Lazy evaluation is used, so it can be a symbol or a string.
alpha	Sensitivity parameter, a vector of values.
knots	Location of spline knots. If a SplineBasis object is provided, it is used directly.
outcome.model	The observed effects model.
intensity.model	The assessment time intensity model.
impute_data	A function that takes (t, df) and returns the imputed data at time t. Should handle extrapolation from the last observed time point.
loss	The loss function to use. Options are "lp_mse", "mean_mse", and "quasi-likelihood".
link	The link function to use. Options are "identity", "log", and "logit".
spline.degree	The degree of the spline basis, default is 3 (cubic splines).
...	Additional arguments passed to compute_influence_terms.
BBsolve.control	Control parameters for the BB::sane optimizer, including maxit and tol.
term2_method	Method for computing term2 influence components. Options are: "fast": Optimized closure-based integrand with adaptive Simpson's (default) "original": Standard implementation with adaptive Simpson's "fixed_grid": Pre-computed expected values on fixed grid with composite Simpson's rule "seeded_adaptive": Adaptive Simpson's seeded with pre-computed grid points "gauss_legendre": Gauss-Legendre quadrature (requires statmod package)
term2_grid_n	Number of grid points/nodes for fixed_grid, seeded_adaptive, and gauss_legendre methods (default 100)
use_expected_cache	Logical; whether to cache expected values for performance (default TRUE)

Details

Integration Methods for Term2

The function offers four integration methods with different performance/accuracy tradeoffs:

Adaptive Methods (fast, original):

• Use adaptive Simpson's quadrature with automatic subdivision

• Best accuracy for irregular integrands

• fast method uses optimized closure-based integrand construction

Fixed-Grid Method (fixed_grid):

• Pre-computes expected values at fixed grid points (once per alpha)

• Uses composite Simpson's rule for integration

• 2-5x faster when optimizing over beta (multiple iterations)

• Best for smooth integrands with sufficient grid density

Seeded Adaptive Method (seeded_adaptive):

• Combines pre-computation with adaptive refinement

• Starts with pre-computed grid, subdivides where needed

• Good balance of speed and accuracy

Gauss-Legendre Method (gauss_legendre):

• Uses Gauss-Legendre quadrature via statmod::gauss.quad()

• Highly accurate for smooth integrands with fewer evaluation points

• Exact for polynomials up to degree 2n-1 using n points

• Requires the statmod package

Outcome Model Compatibility

Unlike simulation code that assumes specific single-index model formulas, this function supports any outcome model with a compute_SensIAT_expected_values method, including:

• Single-index models (⁠fit_SensIAT_single_index_*_model⁠)

• Generalized linear models (GLM)

• Linear models (LM)

• Negative binomial models

Performance Optimizations

Term1 Optimizations (alpha-independent):

• Y-scaled observations pre-computed once

• Patient index mappings cached for O(1) lookups

• Identity link: weights pre-computed (don't depend on beta or alpha)

• Single-index models: global PMF constants extracted once

Term2 Optimizations:

• Integration grids pre-computed (alpha-independent)

• Basis evaluations at grid points pre-computed

• Expected values computed once per alpha (for grid methods)

• Per-patient caching of expected values

• Weight functions use numerically stable exp(-μ) multiplication

• Weight functions use numerically stable exp(-μ) multiplication

• Fast method uses closure-based integration with reduced allocations

Examples

library(survival)
library(splines)

data_with_lags <- SensIAT_example_data |>
    dplyr::group_by(Subject_ID) |>
    dplyr::mutate(
        ..prev_outcome.. = dplyr::lag(Outcome, default = NA_real_, order_by = Time),
        ..prev_time.. = dplyr::lag(Time, default = 0, order_by = Time),
        ..delta_time.. = Time - dplyr::lag(.data$Time, default = NA_real_, order_by = Time)
    )

# Create the observation time intensity model
intensity.model <-
    coxph(Surv(..prev_time.., Time, !is.na(Outcome)) ~ ..prev_outcome.. + strata(Visit),
        data = data_with_lags |> dplyr::filter(.data$Time > 0)
    )

# Create the observed outcome model
outcome.model <-
    fit_SensIAT_single_index_fixed_coef_model(
        Outcome ~ ns(..prev_outcome.., df = 3) + ..delta_time.. - 1,
        id = Subject_ID,
        data = data_with_lags |> dplyr::filter(.data$Time > 0)
    )
fit_SensIAT_marginal_mean_model_generalized(
    data = data_with_lags,
    time = data_with_lags$Time,
    id = data_with_lags$Subject_ID,
    alpha = 0,
    knots = c(60, 260, 460),
    outcome.model = outcome.model,
    intensity.model = intensity.model,
    loss = "lp_mse",
    link = "log",
    impute_data = \(t, df){
        data_wl <- df |>
            dplyr::mutate(
                ..prev_time.. = Time,
                ..prev_outcome.. = Outcome,
                ..delta_time.. = 0
            )
        extrapolate_from_last_observation(t, data_wl, "Time", slopes = c("..delta_time.." = 1))
    }
)

---

Outcome Modeler for SensIAT Single Index Model.

Description

Outcome Modeler for SensIAT Single Index Model.

Usage

fit_SensIAT_single_index_fixed_coef_model(
  formula,
  data,
  kernel = "K2_Biweight",
  method = "nmk",
  id = ..id..,
  initial = NULL,
  ...
)

fit_SensIAT_single_index_fixed_bandwidth_model(
  formula,
  data,
  kernel = "K2_Biweight",
  method = "nmk",
  id = ..id..,
  initial = NULL,
  ...
)

Arguments

formula	The outcome model formula
data	The data to fit the outcome model to. Should only include follow-up data, i.e. time > 0.
kernel	The kernel to use for the outcome model.
method	The optimization method to use for the outcome model, either "optim", "nlminb", or "nmk".
id	The patient identifier variable for the data.
initial	Either a vector of initial values or a function to estimate initial values. If NULL (default), the initial values are estimated using the MAVE::mave.compute function.
...	Currently ignored, included for future compatibility.

Value

Object of class SensIAT::Single-index-outcome-model which contains the outcome model portion.

Functions

• fit_SensIAT_single_index_fixed_bandwidth_model(): for fitting with a fixed bandwidth

Examples

# A basic example using fixed intensity bandwidth.
object <-
    fit_SensIAT_within_group_model(
        group.data = SensIAT_example_data,
        outcome_modeler = fit_SensIAT_single_index_fixed_bandwidth_model,
        id = Subject_ID,
        outcome = Outcome,
        time = Time,
        knots = c(60, 260, 460),
        End = 830,
        intensity.args = list(bandwidth = 30)
    )

# A basic example using variable bandwidth but with fixed first coefficient.
object.bw <-
    fit_SensIAT_within_group_model(
        group.data = SensIAT_example_data,
        outcome_modeler = fit_SensIAT_single_index_fixed_coef_model,
        id = Subject_ID,
        outcome = Outcome,
        time = Time,
        knots = c(60, 260, 460),
        End = 830,
        intensity.args = list(bandwidth = 30)
    )

---

Single Index Model using MAVE and Optimizing Bandwidth.

Description

Single index model estimation using minimum average variance estimation (MAVE). A direction is estimated using MAVE, and then the bandwidth is selected by minimization of the cross-validated pseudo-integrated squared error. Optionally, the initial coefficients of the outcome model can be re-estimated by optimization on a spherical manifold. This option requires the ManifoldOptim package (see ManifoldOptim::manifold.optim()).

Usage

fit_SensIAT_single_index_norm1coef_model(
  formula,
  data,
  kernel = "K2_Biweight",
  mave.method = "meanMAVE",
  id = ..id..,
  bw.selection = c("ise", "mse"),
  bw.method = c("optim", "grid", "optimize"),
  bw.range = c(0.01, 1.5),
  reestimate.coef = 0,
  ...
)

Arguments

formula	The outcome model formula
data	The data to fit the outcome model to. Should only include follow-up data, i.e. time > 0.
kernel	The kernel to use for the outcome model.
mave.method	The method to use for the MAVE estimation.
id	The patient identifier variable for the data.
bw.selection	The criteria for bandwidth selection, either 'ise' for Integrated Squared Error or 'mse' for Mean Squared Error.
bw.method	The method for bandwidth selection, either 'optim' for using optimization or 'grid' for grid search.
bw.range	A numeric vector of length 2 indicating the range of bandwidths to consider for selection as a multiple of the standard deviation of the single index predictor.
reestimate.coef	number of iterations to go through.
...	Additional arguments to be passed to optim.

Value

Object of class SensIAT::Single-index-outcome-model which contains the outcome model portion.

---

Perform Jackknife Resampling on an Object

Description

Perform Jackknife Resampling on an Object

Usage

jackknife(object, ...)

## S3 method for class 'SensIAT_within_group_model'
jackknife(object, time, ...)

## S3 method for class 'SensIAT_fulldata_model'
jackknife(object, time, ...)

Arguments

object	An object to cross validate on.
...	Additional arguments passed to the method.
time	Time points for which to estimate the response.

Value

A data frame of the jackknife resampling results. For SensIAT objects, a tibble with columns alpha, time, jackknife_mean, and jackknife_var, where jackknife_mean is the mean of the jackknife estimates and jackknife_var is the estimated variances of the response at the given time points for the specified alpha values.

Methods (by class)

• jackknife(SensIAT_within_group_model): Perform jackknife resampling on a SensIAT_within_group_model object.

• jackknife(SensIAT_fulldata_model): Perform jackknife resampling on a SensIAT_fulldata_model object.

Examples

## Not run: 
object <-
    fit_SensIAT_within_group_model(
        group.data = SensIAT_example_data,
        outcome_modeler = fit_SensIAT_single_index_fixed_coef_model,
        alpha = c(-0.6, -0.3, 0, 0.3, 0.6),
        id = Subject_ID,
        outcome = Outcome,
        time = Time,
        intensity.args = list(bandwidth = 30),
        knots = c(60, 260, 460),
        End = 830
    )
jackknife.estimates <- jackknife(object, time = c(90, 180, 270, 360, 450))

## End(Not run)

---

Give the Marginal Mean Estimate and its Estimated Asymptotic Variance

Description

Give the marginal mean model estimate

Usage

## S3 method for class 'SensIAT_fulldata_model'
predict(object, time, ...)

## S3 method for class 'SensIAT_within_group_model'
predict(object, time, include.var = TRUE, ..., base = object$base)

Arguments

object	SensIAT_within_group_model object
time	Time points of interest
...	Currently ignored.
include.var	Logical. If TRUE, the variance of the outcome is also returned
base	A SplineBasis object used to evaluate the basis functions.

Value

If include.var is TRUE, a tibble with columns time, mean, and var is returned. otherwise if include.var is FALSE, only the mean vector is returned.

Functions

• predict(SensIAT_fulldata_model): For each combination of time and alpha estimate the mean response and variance for each group as well as estimate the mean treatment effect and variance.

Examples

model <-
    fit_SensIAT_within_group_model(
        group.data = SensIAT_example_data,
        outcome_modeler = fit_SensIAT_single_index_fixed_coef_model,
        alpha = c(-0.6, -0.3, 0, 0.3, 0.6),
        id = Subject_ID,
        outcome = Outcome,
        time = Time,
        End = 830,
        knots = c(60, 260, 460),
    )
predict(model, time = c(90, 180))

---

Prepare Data for Sensitivity Analysis with Irregular Assessment Times

Description

This function prepares the data for SensIAT analysis by transforming it into a format suitable for the SensIAT models.

Usage

prepare_SensIAT_data(
  data,
  id.var,
  time.var,
  outcome.var,
  End,
  add.terminal.observations = TRUE
)

Arguments

data	A data frame containing the data to be prepared.
id.var	The variable in data that identifies the subject.
time.var	The variable in data that identifies the time of the observation.
outcome.var	The variable in data that contains the outcome of interest.
End	The end time for the analysis. Observations with time greater than End will be filtered out.
add.terminal.observations	Logical indicating whether to add terminal observations to the data (TRUE), or terminal observations have already been added (FALSE).

Value

A data frame with the following transformations:

• Data filtered to time less than or equal to End.

• Observations are arranged by id.var and time.var.

• Terminal observations added if add.terminal.observations is TRUE, with ..time.. set to End and ..outcome.. set to NA, if the subject has less observations than the maximum number of observations.

• New variables created:

  • ..id.. aliases id.var,

  • ..time.. aliases time.var,

  • ..outcome.. aliases outcome.var,

  • ..visit_number.. is the visit number within each subject derived from time.var,

  • ..prev_outcome.., i.e. lag-outcome, the outcome from the previous visit,

  • ..prev_time.., i.e. lag-time, the time from the previous visit,

  • ..delta_time.., the difference in time between the current and previous visit.

Examples

prepare_SensIAT_data(SensIAT_example_data, Subject_ID, Time, Outcome, 830)

exdata <- tibble::tibble(
    ID = rep(1:2, c(3, 5)),
    Time = c(
        0, 30, 60,
        0, 30, 60, 90, 120
    ),
    Outcome = floor(runif(8, 1, 100))
)

prepare_SensIAT_data(exdata, ID, Time, Outcome, 120)

---

SensIAT Example Data

Description

A simulated dataset for use in the SensIAT tutorial, testing and documentation.

Usage

SensIAT_example_data

SensIAT_example_fulldata

Format

SensIAT_example_data is a data frame with 779 rows and 4 variables consisting of 200 simulated patients. Each row in the data represents a visit for the patient. The columns are:

Subject_ID

A unique identifier for each patient.

Visit

The ordinal number of the visit for the patient. Baseline observation is 0.

Time

The time of the visit in days, since baseline.

Outcome

The outcome of interest.

SensIAT_example_fulldata is a data frame with 1614 rows and 5 variables consisting of 400 simulated patients, 200 for each treatment arm. Each row in the data represents a visit for the patient. The columns are:

Subject_ID

A unique identifier for each patient.

Visit

The ordinal number of the visit for the patient. Baseline observation is 0.

Time

The time of the visit in days, since baseline.

Outcome

The outcome of interest.

Treatment_group

Treatment or control group.

Functions

• SensIAT_example_fulldata: A simulated dataset with both treatment and control groups.

---

Simulate SensIAT Data

Description

Generates simulated longitudinal data following the SensIAT model structure, alternating between observation time generation (Cox proportional hazards model) and outcome generation.

Usage

simulate_SensIAT_data(
  n_subjects,
  End,
  intensity_coef = -0.05,
  outcome_coef = list(prev_outcome = c(0.7, -0.1, 0.05), time = -0.001, delta_time =
    -0.002, intercept = 2),
  baseline_hazard = 0.005,
  outcome_sd = 1.5,
  initial_outcome_mean = 5,
  initial_outcome_sd = 2,
  max_visits = 50,
  seed = NULL,
  link = "identity"
)

Arguments

n_subjects	Number of subjects to simulate.
End	Maximum follow-up time.
intensity_coef	Coefficient for the effect of previous outcome on observation intensity. Can be a scalar or vector (one per visit number stratum).
outcome_coef	Named list of coefficients for outcome model including: prev_outcome - coefficients for natural spline of previous outcome (length 3) time - coefficient for time since baseline delta_time - coefficient for time since last observation intercept - intercept term
baseline_hazard	Baseline hazard function. Either a function of time and visit number, or a numeric value for constant baseline hazard.
outcome_sd	Standard deviation of the outcome residuals.
initial_outcome_mean	Mean of the initial (baseline) outcome.
initial_outcome_sd	Standard deviation of the initial outcome.
max_visits	Maximum number of visits per subject (to prevent infinite loops).
seed	Random seed for reproducibility.
link	Link function for outcome model. One of "identity", "log", or "logit". Determines the scale on which the outcome model operates.

Value

A tibble with columns:

• Subject_ID - Subject identifier

• Time - Observation time

• Outcome - Observed outcome value

• Additional columns may be added for internal use

Examples

# Simulate data with default parameters
sim_data <- simulate_SensIAT_data(
    n_subjects = 100,
    End = 830,
    intensity_coef = -0.05,
    outcome_coef = list(
        prev_outcome = c(0.7, -0.1, 0.05),
        time = -0.001,
        delta_time = -0.002,
        intercept = 2
    ),
    baseline_hazard = 0.005,
    outcome_sd = 1.5
)

---

Simulate Treatment and Control Groups

Description

Generate simulated data for both treatment and control groups with potentially different parameters.

Usage

simulate_SensIAT_two_groups(
  n_subjects,
  End,
  intensity_coef = -0.05,
  outcome_coef = list(prev_outcome = c(0.7, -0.1, 0.05), time = -0.001, delta_time =
    -0.002, intercept = 2),
  baseline_hazard = 0.005,
  outcome_sd = 1.5,
  initial_outcome_mean = 5,
  initial_outcome_sd = 2,
  max_visits = 50,
  treatment_effect = 0,
  treatment_intensity_effect = 1,
  seed = NULL,
  link = "identity"
)

Arguments

n_subjects	Number of subjects to simulate.
End	Maximum follow-up time.
intensity_coef	Coefficient for the effect of previous outcome on observation intensity. Can be a scalar or vector (one per visit number stratum).
outcome_coef	Named list of coefficients for outcome model including: prev_outcome - coefficients for natural spline of previous outcome (length 3) time - coefficient for time since baseline delta_time - coefficient for time since last observation intercept - intercept term
baseline_hazard	Baseline hazard function. Either a function of time and visit number, or a numeric value for constant baseline hazard.
outcome_sd	Standard deviation of the outcome residuals.
initial_outcome_mean	Mean of the initial (baseline) outcome.
initial_outcome_sd	Standard deviation of the initial outcome.
max_visits	Maximum number of visits per subject (to prevent infinite loops).
treatment_effect	Additive treatment effect on outcomes (added to intercept on link scale).
treatment_intensity_effect	Multiplicative effect on observation intensity (values < 1 mean fewer observations in treatment group).
seed	Random seed for reproducibility.
link	Link function for outcome model. One of "identity", "log", or "logit".

Value

A tibble with an additional Treatment column indicating group assignment.

Examples

# Simulate data with treatment effect
sim_data <- simulate_SensIAT_two_groups(
    n_subjects = 100,
    End = 830,
    treatment_effect = 1.5,
    treatment_intensity_effect = 0.9
)

---