From 0326bf3e8a7b2ad89fd8db4d9dbb89e630e02dd1 Mon Sep 17 00:00:00 2001 From: clarkliming Date: Wed, 1 Nov 2023 08:41:35 +0000 Subject: [PATCH] deploy: 38e4aa3d26f0be0cacaec180893274737d7d38bc --- main/coverage-report/index.html | 12935 +++++++++++++++--------------- 1 file changed, 6436 insertions(+), 6499 deletions(-) diff --git a/main/coverage-report/index.html b/main/coverage-report/index.html index bae77291b..db3d2c0f0 100644 --- a/main/coverage-report/index.html +++ b/main/coverage-report/index.html @@ -94,7 +94,7 @@ font-size: 11px; }
-

mmrm coverage - 96.71%

+

mmrm coverage - 96.63%

  • @@ -106,8 +106,8 @@

    mmrm coverage - 96.71%

-
- +
+
@@ -11642,3470 +11642,3442 @@

mmrm coverage - 96.71%

69 - - 
#' - `visits_zero_inds`: length `n` `integer` containing zero-based visits indices.
- - - - 70 - 
#' - `n_visits`: `int` with the number of visits, which is the dimension of the
- 71 + 70 
#'      covariance matrix.
- 72 + 71 
#' - `n_subjects`: `int` with the number of subjects.
- 73 + 72 
#' - `subject_zero_inds`: length `n_subjects` `integer` containing the zero-based start
- 74 + 73 
#'     indices for each subject.
- 75 + 74 
#' - `subject_n_visits`: length `n_subjects` `integer` containing the number of
- 76 + 75 
#'     observed visits for each subjects. So the sum of this vector equals `n`.
- 77 + 76 
#' - `cov_type`: `string` value specifying the covariance type.
- 78 + 77 
#' - `is_spatial_int`: `int` specifying whether the covariance structure is spatial(1) or not(0).
- 79 + 78 
#' - `reml`: `int` specifying whether REML estimation is used (1), otherwise ML (0).
- 80 + 79 
#' - `subject_groups`: `factor` specifying the grouping for each subject.
- 81 + 80 
#' - `n_groups`: `int` with the number of total groups
- 82 + 81 
#'
- 83 + 82 
#' @details Note that the `subject_var` must not be factor but can also be character.
- 84 + 83 
#'   If it is character, then it will be converted to factor internally. Here
- 85 + 84 
#'   the levels will be the unique values, sorted alphabetically and numerically if there
- 86 + 85 
#'   is a common string prefix of numbers in the character elements. For full control
- 87 + 86 
#'   on the order please use a factor.
- 88 + 87 
#'
- 89 + 88 
#' @keywords internal
- 90 + 89 
h_mmrm_tmb_data <- function(formula_parts,
- 91 + 90 
                            data,
- 92 + 91 
                            weights,
- 93 + 92 
                            reml,
- 94 + 93 
                            singular = c("drop", "error", "keep"),
- 95 + 94 
                            drop_visit_levels,
- 96 + 95 
                            allow_na_response = FALSE,
- 97 + 96 
                            drop_levels = TRUE) {
- 98 + 97 297x 
  assert_class(formula_parts, "mmrm_tmb_formula_parts")
- 99 + 98 297x 
  assert_data_frame(data)
- 100 + 99 297x 
  varname <- formula_parts[grepl("_var", names(formula_parts))]
- 101 + 100 297x 
  assert_names(
- 102 + 101 297x 
    names(data),
- 103 + 102 297x 
    must.include = unlist(varname, use.names = FALSE)
- 104 + 103 
  )
- 105 + 104 297x 
  assert_true(is.factor(data[[formula_parts$subject_var]]) || is.character(data[[formula_parts$subject_var]]))
- 106 + 105 297x 
  assert_numeric(weights, len = nrow(data))
- 107 + 106 297x 
  assert_flag(reml)
- 108 + 107 297x 
  singular <- match.arg(singular)
- 109 + 108 297x 
  assert_flag(drop_visit_levels)
- 110 + 109 

                     
                   
                   
-                    111
+                    110
                     297x
                     
                       
  if (is.character(data[[formula_parts$subject_var]])) {

                     
                   
                   
-                    112
+                    111
                     4x
                     
                       
    data[[formula_parts$subject_var]] <- factor(

                     
                   
                   
-                    113
+                    112
                     4x
                     
                       
      data[[formula_parts$subject_var]],

                     
                   
                   
-                    114
+                    113
                     4x
                     
                       
      levels = stringr::str_sort(unique(data[[formula_parts$subject_var]]), numeric = TRUE)

                     
                   
                   
-                    115
+                    114
                     
                     
                       
    )

                     
                   
                   
-                    116
+                    115
                     
                     
                       
  }

                     
                   
                   
-                    117
+                    116
                     297x
                     
                       
  data_order <- if (formula_parts$is_spatial) {

                     
                   
                   
-                    118
+                    117
                     15x
                     
                       
    order(data[[formula_parts$subject_var]])

                     
                   
                   
-                    119
+                    118
                     
                     
                       
  } else {

                     
                   
                   
-                    120
+                    119
                     282x
                     
                       
    subject_visit_data <- data[, c(formula_parts$subject_var, formula_parts$visit_var)]

                     
                   
                   
-                    121
+                    120
                     282x
                     
                       
    is_duplicated <- duplicated(subject_visit_data)

                     
                   
                   
-                    122
+                    121
                     282x
                     
                       
    if (any(is_duplicated)) {

                     
                   
                   
-                    123
+                    122
                     1x
                     
                       
      stop(

                     
                   
                   
-                    124
+                    123
                     1x
                     
                       
        "time points have to be unique for each subject, detected following duplicates in data:\n",

                     
                   
                   
-                    125
+                    124
                     1x
                     
                       
        paste(utils::capture.output(print(subject_visit_data[is_duplicated, ])), collapse = "\n")

                     
                   
                   
-                    126
+                    125
                     
                     
                       
      )

                     
                   
                   
-                    127
+                    126
                     
                     
                       
    }

                     
                   
                   
-                    128
+                    127
                     281x
                     
                       
    order(data[[formula_parts$subject_var]], data[[formula_parts$visit_var]])

                     
                   
                   
-                    129
+                    128
                     
                     
                       
  }

                     
                   
                   
-                    130
+                    129
                     296x
                     
                       
  if (identical(formula_parts$is_spatial, FALSE)) {

                     
                   
                   
-                    131
+                    130
                     281x
                     
                       
    h_confirm_large_levels(length(levels(data[[formula_parts$visit_var]])))

                     
                   
                   
-                    132
+                    131
                     
                     
                       
  }

                     
                   
                   
-                    133
+                    132
                     295x
                     
                       
  data <- data[data_order, ]

                     
                   
                   
-                    134
+                    133
                     295x
                     
                       
  weights <- weights[data_order]

                     
                   
                   
-                    135
+                    134
                     295x
                     
                       
  data <- data.frame(data, weights)

                     
                   
                   
-                    136
+                    135
                     
                     
                       
  # Weights is always the last column.

                     
                   
                   
-                    137
+                    136
                     295x
                     
                       
  weights_name <- colnames(data)[ncol(data)]

                     
                   
                   
-                    138
+                    137
                     
                     
                       
  # If `y` is allowed to be NA, then first replace y with 1:n, then replace it with original y.

                     
                   
                   
-                    139
+                    138
                     295x
                     
                       
  if (allow_na_response) {

                     
                   
                   
-                    140
+                    139
                     69x
                     
                       
    y_original <- eval(formula_parts$full_formula[[2]], envir = data)

                     
                   
                   
-                    141
+                    140
                     69x
                     
                       
    vn <- deparse(formula_parts$full_formula[[2]])

                     
                   
                   
-                    142
+                    141
                     69x
                     
                       
    data[[vn]] <- seq_len(nrow(data))

                     
                   
                   
-                    143
+                    142
                     
                     
                       
  } else {

                     
                   
                   
-                    144
+                    143
                     226x
                     
                       
    h_warn_na_action()

                     
                   
                   
-                    145
+                    144
                     
                     
                       
  }

                     
                   
                   
-                    146
+                    145
                     295x
                     
                       
  full_frame <- eval(

                     
                   
                   
-                    147
+                    146
                     295x
                     
                       
    bquote(stats::model.frame(

                     
                   
                   
-                    148
+                    147
                     295x
                     
                       
      formula_parts$full_formula,

                     
                   
                   
-                    149
+                    148
                     295x
                     
                       
      data = data,

                     
                   
                   
-                    150
+                    149
                     295x
                     
                       
      weights = .(as.symbol(weights_name)),

                     
                   
                   
-                    151
+                    150
                     295x
                     
                       
      na.action = stats::na.omit

                     
                   
                   
-                    152
+                    151
                     
                     
                       
    ))

                     
                   
                   
-                    153
+                    152
                     
                     
                       
  )

                     
                   
                   
-                    154
+                    153
                     295x
                     
                       
  if (drop_levels) {

                     
                   
                   
-                    155
+                    154
                     228x
                     
                       
    full_frame <- droplevels(full_frame, except = formula_parts$visit_var)

                     
                   
                   
-                    156
+                    155
                     
                     
                       
  }

                     
                   
                   
-                    157
+                    156
                     
                     
                       
  # If `y` is allowed to be NA, replace it with original y.

                     
                   
                   
-                    158
+                    157
                     295x
                     
                       
  if (allow_na_response) {

                     
                   
                   
-                    159
+                    158
                     69x
                     
                       
    full_frame[[vn]] <- y_original[full_frame[[vn]]]

                     
                   
                   
-                    160
+                    159
                     
                     
                       
  }

                     
                   
                   
-                    161
+                    160
                     295x
                     
                       
  if (drop_visit_levels && !formula_parts$is_spatial && is.factor(full_frame[[formula_parts$visit_var]])) {

                     
                   
                   
-                    162
+                    161
                     212x
                     
                       
    old_levels <- levels(full_frame[[formula_parts$visit_var]])

                     
                   
                   
-                    163
+                    162
                     212x
                     
                       
    full_frame[[formula_parts$visit_var]] <- droplevels(full_frame[[formula_parts$visit_var]])

                     
                   
                   
-                    164
+                    163
                     212x
                     
                       
    new_levels <- levels(full_frame[[formula_parts$visit_var]])

                     
                   
                   
-                    165
+                    164
                     212x
                     
                       
    dropped <- setdiff(old_levels, new_levels)

                     
                   
                   
-                    166
+                    165
                     212x
                     
                       
    if (length(dropped) > 0) {

                     
                   
                   
-                    167
+                    166
                     3x
                     
                       
      message("In ", formula_parts$visit_var, " there are dropped visits: ", toString(dropped))

                     
                   
                   
-                    168
+                    167
                     
                     
                       
    }

                     
                   
                   
-                    169
+                    168
                     
                     
                       
  }

                     
                   
                   
-                    170
+                    169
                     295x
                     
                       
  x_matrix <- stats::model.matrix(formula_parts$model_formula, data = full_frame)

                     
                   
                   
-                    171
+                    170
                     294x
                     
                       
  x_cols_aliased <- stats::setNames(rep(FALSE, ncol(x_matrix)), nm = colnames(x_matrix))

                     
                   
                   
-                    172
+                    171
                     294x
                     
                       
  qr_x_mat <- qr(x_matrix)

                     
                   
                   
-                    173
+                    172
                     294x
                     
                       
  if (qr_x_mat$rank < ncol(x_matrix)) {

                     
                   
                   
-                    174
+                    173
                     26x
                     
                       
    cols_to_drop <- utils::tail(qr_x_mat$pivot, ncol(x_matrix) - qr_x_mat$rank)

                     
                   
                   
-                    175
+                    174
                     26x
                     
                       
    if (identical(singular, "error")) {

                     
                   
                   
-                    176
+                    175
                     1x
                     
                       
      stop(

                     
                   
                   
-                    177
+                    176
                     1x
                     
                       
        "design matrix only has rank ", qr_x_mat$rank, " and ", length(cols_to_drop),

                     
                   
                   
-                    178
+                    177
                     1x
                     
                       
        " columns (", toString(colnames(x_matrix)[cols_to_drop]), ") could be dropped",

                     
                   
                   
-                    179
+                    178
                     1x
                     
                       
        " to achieve full rank ", ncol(x_matrix), " by using `accept_singular = TRUE`"

                     
                   
                   
-                    180
+                    179
                     
                     
                       
      )

                     
                   
                   
-                    181
+                    180
                     25x
                     
                       
    } else if (identical(singular, "drop")) {

                     
                   
                   
-                    182
+                    181
                     11x
                     
                       
      assign_attr <- attr(x_matrix, "assign")

                     
                   
                   
-                    183
+                    182
                     11x
                     
                       
      contrasts_attr <- attr(x_matrix, "contrasts")

                     
                   
                   
-                    184
+                    183
                     11x
                     
                       
      x_matrix <- x_matrix[, -cols_to_drop, drop = FALSE]

                     
                   
                   
-                    185
+                    184
                     11x
                     
                       
      x_cols_aliased[cols_to_drop] <- TRUE

                     
                   
                   
-                    186
+                    185
                     11x
                     
                       
      attr(x_matrix, "assign") <- assign_attr[-cols_to_drop]

                     
                   
                   
-                    187
+                    186
                     11x
                     
                       
      attr(x_matrix, "contrasts") <- contrasts_attr

                     
                   
                   
-                    188
+                    187
                     
                     
                       
    }

                     
                   
                   
-                    189
+                    188
                     
                     
                       
  }

                     
                   
                   
-                    190
+                    189
                     
                     
                       

                     
                   
                   
-                    191
+                    190
                     293x
                     
                       
  y_vector <- as.numeric(stats::model.response(full_frame))

                     
                   
                   
-                    192
+                    191
                     293x
                     
                       
  weights_vector <- as.numeric(stats::model.weights(full_frame))

                     
                   
                   
-                    193
+                    192
                     293x
                     
                       
  n_subjects <- length(unique(full_frame[[formula_parts$subject_var]]))

                     
                   
                   
-                    194
+                    193
                     293x
                     
                       
  subject_zero_inds <- which(!duplicated(full_frame[[formula_parts$subject_var]])) - 1L

                     
                   
                   
-                    195
+                    194
                     293x
                     
                       
  subject_n_visits <- c(utils::tail(subject_zero_inds, -1L), nrow(full_frame)) - subject_zero_inds

                     
                   
                   
-                    196
+                    195
                     
                     
                       
  # It is possible that `subject_var` is factor with more levels (and this does not affect fit)

                     
                   
                   
-                    197
+                    196
                     
                     
                       
  # so no check is needed for `subject_visits`.

                     
                   
                   
-                    198
+                    197
                     293x
                     
                       
  assert_true(all(subject_n_visits > 0))

                     
                   
                   
-                    199
+                    198
                     293x
                     
                       
  if (!is.null(formula_parts$group_var)) {

                     
                   
                   
-                    200
+                    199
                     39x
                     
                       
    assert_factor(data[[formula_parts$group_var]])

                     
                   
                   
-                    201
+                    200
                     39x
                     
                       
    subject_groups <- full_frame[[formula_parts$group_var]][subject_zero_inds + 1L]

                     
                   
                   
-                    202
+                    201
                     39x
                     
                       
    n_groups <- nlevels(subject_groups)

                     
                   
                   
-                    203
+                    202
                     
                     
                       
  } else {

                     
                   
                   
-                    204
+                    203
                     254x
                     
                       
    subject_groups <- factor(rep(0L, n_subjects))

                     
                   
                   
-                    205
+                    204
                     254x
                     
                       
    n_groups <- 1L

                     
                   
                   
-                    206
+                    205
                     
                     
                       
  }

                     
                   
                   
-                    207
+                    206
                     293x
                     
                       
  coordinates <- full_frame[, formula_parts$visit_var, drop = FALSE]

                     
                   
                   
-                    208
+                    207
                     293x
                     
                       
  if (formula_parts$is_spatial) {

                     
                   
                   
-                    209
+                    208
                     15x
                     
                       
    lapply(coordinates, assert_numeric)

                     
                   
                   
-                    210
+                    209
                     15x
                     
                       
    coordinates_matrix <- as.matrix(coordinates)

                     
                   
                   
-                    211
-                    15x
-                    
-                      
    visits_zero_inds <- 0L

-                    
-                  
-                  
-                    212
+                    210
                     15x
                     
                       
    n_visits <- max(subject_n_visits)

                     
                   
                   
-                    213
+                    211
                     
                     
                       
  } else {

                     
                   
                   
-                    214
+                    212
                     278x
                     
                       
    assert(identical(ncol(coordinates), 1L))

                     
                   
                   
-                    215
+                    213
                     278x
                     
                       
    assert_factor(coordinates[[1L]])

                     
                   
                   
-                    216
-                    278x
-                    
-                      
    visits_zero_inds <- as.integer(coordinates[[1L]]) - 1L

-                    
-                  
-                  
-                    217
+                    214
                     278x
                     
-                      
    coordinates_matrix <- as.matrix(visits_zero_inds, ncol = 1)

+                      
    coordinates_matrix <- as.matrix(as.integer(coordinates[[1L]]) - 1, ncol = 1)

                     
                   
                   
-                    218
+                    215
                     278x
                     
                       
    n_visits <- nlevels(coordinates[[1L]])

                     
                   
                   
-                    219
+                    216
                     278x
                     
                       
    assert_true(all(subject_n_visits <= n_visits))

                     
                   
                   
-                    220
+                    217
                     
                     
                       
  }

                     
                   
                   
-                    221
+                    218
                     293x
                     
                       
  structure(

                     
                   
                   
-                    222
+                    219
                     293x
                     
                       
    list(

                     
                   
                   
-                    223
+                    220
                     293x
                     
                       
      full_frame = full_frame,

                     
                   
                   
-                    224
+                    221
                     293x
                     
                       
      data = data,

                     
                   
                   
-                    225
+                    222
                     293x
                     
                       
      x_matrix = x_matrix,

                     
                   
                   
-                    226
+                    223
                     293x
                     
                       
      x_cols_aliased = x_cols_aliased,

                     
                   
                   
-                    227
+                    224
                     293x
                     
                       
      coordinates = coordinates_matrix,

                     
                   
                   
-                    228
+                    225
                     293x
                     
                       
      y_vector = y_vector,

                     
                   
                   
-                    229
+                    226
                     293x
                     
                       
      weights_vector = weights_vector,

                     
                   
                   
-                    230
-                    293x
-                    
-                      
      visits_zero_inds = visits_zero_inds,

-                    
-                  
-                  
-                    231
+                    227
                     293x
                     
                       
      n_visits = n_visits,

                     
                   
                   
-                    232
+                    228
                     293x
                     
                       
      n_subjects = n_subjects,

                     
                   
                   
-                    233
+                    229
                     293x
                     
                       
      subject_zero_inds = subject_zero_inds,

                     
                   
                   
-                    234
+                    230
                     293x
                     
                       
      subject_n_visits = subject_n_visits,

                     
                   
                   
-                    235
+                    231
                     293x
                     
                       
      cov_type = formula_parts$cov_type,

                     
                   
                   
-                    236
+                    232
                     293x
                     
                       
      is_spatial_int = as.integer(formula_parts$is_spatial),

                     
                   
                   
-                    237
+                    233
                     293x
                     
                       
      reml = as.integer(reml),

                     
                   
                   
-                    238
+                    234
                     293x
                     
                       
      subject_groups = subject_groups,

                     
                   
                   
-                    239
+                    235
                     293x
                     
                       
      n_groups = n_groups

                     
                   
                   
-                    240
+                    236
                     
                     
                       
    ),

                     
                   
                   
-                    241
+                    237
                     293x
                     
                       
    class = "mmrm_tmb_data"

                     
                   
                   
-                    242
+                    238
                     
                     
                       
  )

                     
                   
                   
-                    243
+                    239
                     
                     
                       
}

                     
                   
                   
-                    244
+                    240
                     
                     
                       

                     
                   
                   
-                    245
+                    241
                     
                     
                       
#' Start Parameters for `TMB` Fit

                     
                   
                   
-                    246
+                    242
                     
                     
                       
#'

                     
                   
                   
-                    247
+                    243
                     
                     
                       
#' @param formula_parts (`mmrm_tmb_formula_parts`)\cr produced by

                     
                   
                   
-                    248
+                    244
                     
                     
                       
#'  [h_mmrm_tmb_formula_parts()].

                     
                   
                   
-                    249
+                    245
                     
                     
                       
#' @param tmb_data (`mmrm_tmb_data`)\cr produced by [h_mmrm_tmb_data()].

                     
                   
                   
-                    250
+                    246
                     
                     
                       
#' @param start (`numeric` or `NULL`)\cr optional start values for variance

                     
                   
                   
-                    251
+                    247
                     
                     
                       
#'   parameters.

                     
                   
                   
-                    252
+                    248
                     
                     
                       
#' @param n_groups (`int`)\cr number of groups.

                     
                   
                   
-                    253
+                    249
                     
                     
                       
#' @return List with element `theta` containing the start values for the variance

                     
                   
                   
-                    254
+                    250
                     
                     
                       
#'   parameters.

                     
                   
                   
-                    255
+                    251
                     
                     
                       
#'

                     
                   
                   
-                    256
+                    252
                     
                     
                       
#' @keywords internal

                     
                   
                   
-                    257
+                    253
                     
                     
                       
h_mmrm_tmb_parameters <- function(formula_parts,

                     
                   
                   
-                    258
+                    254
                     
                     
                       
                                  tmb_data,

                     
                   
                   
-                    259
+                    255
                     
                     
                       
                                  start,

                     
                   
                   
-                    260
+                    256
                     
                     
                       
                                  n_groups = 1L) {

                     
                   
                   
-                    261
+                    257
                     217x
                     
                       
  assert_class(formula_parts, "mmrm_tmb_formula_parts")

                     
                   
                   
-                    262
+                    258
                     217x
                     
                       
  assert_class(tmb_data, "mmrm_tmb_data")

                     
                   
                   
-                    263
+                    259
                     
                     
                       

                     
                   
                   
-                    264
+                    260
                     217x
                     
                       
  m <- tmb_data$n_visits

                     
                   
                   
-                    265
+                    261
                     217x
                     
                       
  start_value <- switch(formula_parts$cov_type,

                     
                   
                   
-                    266
+                    262
                     217x
                     
                       
    us = rep(0, m * (m + 1) / 2),

                     
                   
                   
-                    267
+                    263
                     217x
                     
                       
    toep = rep(0, m),

                     
                   
                   
-                    268
+                    264
                     217x
                     
                       
    toeph = rep(0, 2 * m - 1),

                     
                   
                   
-                    269
+                    265
                     217x
                     
                       
    ar1 = c(0, 0.5),

                     
                   
                   
-                    270
+                    266
                     217x
                     
                       
    ar1h = c(rep(0, m), 0.5),

                     
                   
                   
-                    271
+                    267
                     217x
                     
                       
    ad = rep(0, m),

                     
                   
                   
-                    272
+                    268
                     217x
                     
                       
    adh = rep(0, 2 * m - 1),

                     
                   
                   
-                    273
+                    269
                     217x
                     
                       
    cs = rep(0, 2),

                     
                   
                   
-                    274
+                    270
                     217x
                     
                       
    csh = rep(0, m + 1),

                     
                   
                   
-                    275
+                    271
                     217x
                     
                       
    sp_exp = rep(0, 2)

                     
                   
                   
-                    276
+                    272
                     
                     
                       
  )

                     
                   
                   
-                    277
+                    273
                     217x
                     
                       
  theta_dim <- length(start_value) * n_groups

                     
                   
                   
-                    278
+                    274
                     217x
                     
                       
  if (!is.null(start)) {

                     
                   
                   
-                    279
+                    275
                     4x
                     
                       
    assert_numeric(start, len = theta_dim, any.missing = FALSE, finite = TRUE)

                     
                   
                   
-                    280
+                    276
                     
                     
                       
  } else {

                     
                   
                   
-                    281
+                    277
                     213x
                     
                       
    start <- rep(start_value, n_groups)

                     
                   
                   
-                    282
+                    278
                     
                     
                       
  }

                     
                   
                   
-                    283
+                    279
                     217x
                     
                       
  list(theta = start)

                     
                   
                   
-                    284
+                    280
                     
                     
                       
}

                     
                   
                   
-                    285
+                    281
                     
                     
                       

                     
                   
                   
-                    286
+                    282
                     
                     
                       
#' Asserting Sane Start Values for `TMB` Fit

                     
                   
                   
-                    287
+                    283
                     
                     
                       
#'

                     
                   
                   
-                    288
+                    284
                     
                     
                       
#' @param tmb_object (`list`)\cr created with [TMB::MakeADFun()].

                     
                   
                   
-                    289
+                    285
                     
                     
                       
#'

                     
                   
                   
-                    290
+                    286
                     
                     
                       
#' @return Nothing, only used for assertions.

                     
                   
                   
-                    291
+                    287
                     
                     
                       
#'

                     
                   
                   
-                    292
+                    288
                     
                     
                       
#' @keywords internal

                     
                   
                   
-                    293
+                    289
                     
                     
                       
h_mmrm_tmb_assert_start <- function(tmb_object) {

                     
                   
                   
-                    294
+                    290
                     204x
                     
                       
  assert_list(tmb_object)

                     
                   
                   
-                    295
+                    291
                     204x
                     
                       
  assert_subset(c("fn", "gr", "par"), names(tmb_object))

                     
                   
                   
-                    296
+                    292
                     
                     
                       

                     
                   
                   
-                    297
+                    293
                     204x
                     
                       
  if (is.na(tmb_object$fn(tmb_object$par))) {

                     
                   
                   
-                    298
+                    294
                     1x
                     
                       
    stop("negative log-likelihood is NaN at starting parameter values")

                     
                   
                   
-                    299
+                    295
                     
                     
                       
  }

                     
                   
                   
-                    300
+                    296
                     203x
                     
                       
  if (any(is.na(tmb_object$gr(tmb_object$par)))) {

                     
                   
                   
-                    301
+                    297
                     1x
                     
                       
    stop("some elements of gradient are NaN at starting parameter values")

                     
                   
                   
-                    302
+                    298
                     
                     
                       
  }

                     
                   
                   
-                    303
+                    299
                     
                     
                       
}

                     
                   
                   
-                    304
+                    300
                     
                     
                       

                     
                   
                   
-                    305
+                    301
                     
                     
                       
#' Checking the `TMB` Optimization Result

                     
                   
                   
-                    306
+                    302
                     
                     
                       
#'

                     
                   
                   
-                    307
+                    303
                     
                     
                       
#' @param tmb_opt (`list`)\cr optimization result.

                     
                   
                   
-                    308
+                    304
                     
                     
                       
#' @param mmrm_tmb (`mmrm_tmb`)\cr result from [h_mmrm_tmb_fit()].

                     
                   
                   
-                    309
+                    305
                     
                     
                       
#'

                     
                   
                   
-                    310
+                    306
                     
                     
                       
#' @return Nothing, only used to generate warnings in case that the model

                     
                   
                   
-                    311
+                    307
                     
                     
                       
#' did not converge.

                     
                   
                   
-                    312
+                    308
                     
                     
                       
#'

                     
                   
                   
-                    313
+                    309
                     
                     
                       
#' @keywords internal

                     
                   
                   
-                    314
+                    310
                     
                     
                       
h_mmrm_tmb_check_conv <- function(tmb_opt, mmrm_tmb) {

                     
                   
                   
-                    315
+                    311
                     199x
                     
                       
  assert_list(tmb_opt)

                     
                   
                   
-                    316
+                    312
                     199x
                     
                       
  assert_subset(c("par", "objective", "convergence", "message"), names(tmb_opt))

                     
                   
                   
-                    317
+                    313
                     199x
                     
                       
  assert_class(mmrm_tmb, "mmrm_tmb")

                     
                   
                   
-                    318
+                    314
                     
                     
                       

                     
                   
                   
-                    319
+                    315
                     199x
                     
                       
  if (!is.null(tmb_opt$convergence) && tmb_opt$convergence != 0) {

                     
                   
                   
-                    320
+                    316
                     1x
                     
                       
    warning("Model convergence problem: ", tmb_opt$message, ".")

                     
                   
                   
-                    321
+                    317
                     1x
                     
                       
    return()

                     
                   
                   
-                    322
+                    318
                     
                     
                       
  }

                     
                   
                   
-                    323
+                    319
                     198x
                     
                       
  theta_vcov <- mmrm_tmb$theta_vcov

                     
                   
                   
-                    324
+                    320
                     198x
                     
                       
  if (is(theta_vcov, "try-error")) {

                     
                   
                   
-                    325
+                    321
                     3x
                     
                       
    warning("Model convergence problem: hessian is singular, theta_vcov not available.")

                     
                   
                   
-                    326
+                    322
                     3x
                     
                       
    return()

                     
                   
                   
-                    327
+                    323
                     
                     
                       
  }

                     
                   
                   
-                    328
+                    324
                     195x
                     
                       
  if (!all(is.finite(theta_vcov))) {

                     
                   
                   
-                    329
+                    325
                     1x
                     
                       
    warning("Model convergence problem: theta_vcov contains non-finite values.")

                     
                   
                   
-                    330
+                    326
                     1x
                     
                       
    return()

                     
                   
                   
-                    331
+                    327
                     
                     
                       
  }

                     
                   
                   
-                    332
+                    328
                     194x
                     
                       
  qr_rank <- qr(theta_vcov)$rank

                     
                   
                   
-                    333
+                    329
                     194x
                     
                       
  if (qr_rank < ncol(theta_vcov)) {

                     
                   
                   
-                    334
+                    330
                     1x
                     
                       
    warning("Model convergence problem: theta_vcov is numerically singular.")

                     
                   
                   
-                    335
+                    331
                     
                     
                       
  }

                     
                   
                   
-                    336
+                    332
                     
                     
                       
}

                     
                   
                   
-                    337
+                    333
                     
                     
                       

                     
                   
                   
-                    338
+                    334
                     
                     
                       
#' Extract covariance matrix from `TMB` report and input data

                     
                   
                   
-                    339
+                    335
                     
                     
                       
#'

                     
                   
                   
-                    340
+                    336
                     
                     
                       
#' This helper does some simple post-processing to extract covariance matrix or named

                     
                   
                   
-                    341
+                    337
                     
                     
                       
#' list of covariance matrices if the fitting is using grouped covariance matrices.

                     
                   
                   
-                    342
+                    338
                     
                     
                       
#'

                     
                   
                   
-                    343
+                    339
                     
                     
                       
#' @param tmb_report (`list`)\cr report created with [TMB::MakeADFun()] report function.

                     
                   
                   
-                    344
+                    340
                     
                     
                       
#' @param tmb_data (`mmrm_tmb_data`)\cr produced by [h_mmrm_tmb_data()].

                     
                   
                   
-                    345
+                    341
                     
                     
                       
#' @param visit_var (`character`)\cr character vector of the visit variable

                     
                   
                   
-                    346
+                    342
                     
                     
                       
#' @param is_spatial (`flag`)\cr indicator whether the covariance structure is spatial.

                     
                   
                   
-                    347
+                    343
                     
                     
                       
#' @return Return a simple covariance matrix if there is no grouping, or a named

                     
                   
                   
-                    348
+                    344
                     
                     
                       
#' list of estimated grouped covariance matrices,

                     
                   
                   
-                    349
+                    345
                     
                     
                       
#' with its name equal to the group levels.

                     
                   
                   
-                    350
+                    346
                     
                     
                       
#'

                     
                   
                   
-                    351
+                    347
                     
                     
                       
#' @keywords internal

                     
                   
                   
-                    352
+                    348
                     
                     
                       
h_mmrm_tmb_extract_cov <- function(tmb_report, tmb_data, visit_var, is_spatial) {

                     
                   
                   
-                    353
+                    349
                     198x
                     
                       
  d <- dim(tmb_report$covariance_lower_chol)

                     
                   
                   
-                    354
+                    350
                     198x
                     
                       
  visit_names <- if (!is_spatial) {

                     
                   
                   
-                    355
+                    351
                     186x
                     
                       
    levels(tmb_data$full_frame[[visit_var]])

                     
                   
                   
-                    356
+                    352
                     
                     
                       
  } else {

                     
                   
                   
-                    357
+                    353
                     12x
                     
                       
    c(0, 1)

                     
                   
                   
-                    358
+                    354
                     
                     
                       
  }

                     
                   
                   
-                    359
+                    355
                     198x
                     
                       
  cov <- lapply(

                     
                   
                   
-                    360
+                    356
                     198x
                     
                       
    seq_len(d[1] / d[2]),

                     
                   
                   
-                    361
+                    357
                     198x
                     
                       
    function(i) {

                     
                   
                   
-                    362
+                    358
                     231x
                     
                       
      ret <- tcrossprod(tmb_report$covariance_lower_chol[seq(1 + (i - 1) * d[2], i * d[2]), ])

                     
                   
                   
-                    363
+                    359
                     231x
                     
                       
      dimnames(ret) <- list(visit_names, visit_names)

                     
                   
                   
-                    364
+                    360
                     231x
                     
                       
      return(ret)

                     
                   
                   
-                    365
+                    361
                     
                     
                       
    }

                     
                   
                   
-                    366
+                    362
                     
                     
                       
  )

                     
                   
                   
-                    367
+                    363
                     198x
                     
                       
  if (identical(tmb_data$n_groups, 1L)) {

                     
                   
                   
-                    368
+                    364
                     165x
                     
                       
    cov <- cov[[1]]

                     
                   
                   
-                    369
+                    365
                     
                     
                       
  } else {

                     
                   
                   
-                    370
+                    366
                     33x
                     
                       
    names(cov) <- levels(tmb_data$subject_groups)

                     
                   
                   
-                    371
+                    367
                     
                     
                       
  }

                     
                   
                   
-                    372
+                    368
                     198x
                     
                       
  return(cov)

                     
                   
                   
-                    373
+                    369
                     
                     
                       
}

                     
                   
                   
-                    374
+                    370
                     
                     
                       

                     
                   
                   
-                    375
+                    371
                     
                     
                       
#' Build `TMB` Fit Result List

                     
                   
                   
-                    376
+                    372
                     
                     
                       
#'

                     
                   
                   
-                    377
+                    373
                     
                     
                       
#' This helper does some simple post-processing of the `TMB` object and

                     
                   
                   
-                    378
+                    374
                     
                     
                       
#' optimization results, including setting names, inverting matrices etc.

                     
                   
                   
-                    379
+                    375
                     
                     
                       
#'

                     
                   
                   
-                    380
+                    376
                     
                     
                       
#' @param tmb_object (`list`)\cr created with [TMB::MakeADFun()].

                     
                   
                   
-                    381
+                    377
                     
                     
                       
#' @param tmb_opt (`list`)\cr optimization result.

                     
                   
                   
-                    382
+                    378
                     
                     
                       
#' @param formula_parts (`mmrm_tmb_formula_parts`)\cr produced by

                     
                   
                   
-                    383
+                    379
                     
                     
                       
#'  [h_mmrm_tmb_formula_parts()].

                     
                   
                   
-                    384
+                    380
                     
                     
                       
#' @param tmb_data (`mmrm_tmb_data`)\cr produced by [h_mmrm_tmb_data()].

                     
                   
                   
-                    385
+                    381
                     
                     
                       
#'

                     
                   
                   
-                    386
+                    382
                     
                     
                       
#' @return List of class `mmrm_tmb` with:

                     
                   
                   
-                    387
+                    383
                     
                     
                       
#'   - `cov`: estimated covariance matrix, or named list of estimated group specific covariance matrices.

                     
                   
                   
-                    388
+                    384
                     
                     
                       
#'   - `beta_est`: vector of coefficient estimates.

                     
                   
                   
-                    389
+                    385
                     
                     
                       
#'   - `beta_vcov`: Variance-covariance matrix for coefficient estimates.

                     
                   
                   
-                    390
+                    386
                     
                     
                       
#'   - `beta_vcov_inv_L`: Lower triangular matrix `L` of the inverse variance-covariance matrix decomposition.

                     
                   
                   
-                    391
+                    387
                     
                     
                       
#'   - `beta_vcov_inv_D`: vector of diagonal matrix `D` of the inverse variance-covariance matrix decomposition.

                     
                   
                   
-                    392
+                    388
                     
                     
                       
#'   - `theta_est`: vector of variance parameter estimates.

                     
                   
                   
-                    393
+                    389
                     
                     
                       
#'   - `theta_vcov`: variance-covariance matrix for variance parameter estimates.

                     
                   
                   
-                    394
+                    390
                     
                     
                       
#'   - `neg_log_lik`: obtained negative log-likelihood.

                     
                   
                   
-                    395
+                    391
                     
                     
                       
#'   - `formula_parts`: input.

                     
                   
                   
-                    396
+                    392
                     
                     
                       
#'   - `data`: input.

                     
                   
                   
-                    397
+                    393
                     
                     
                       
#'   - `weights`: input.

                     
                   
                   
-                    398
+                    394
                     
                     
                       
#'   - `reml`: input as a flag.

                     
                   
                   
-                    399
+                    395
                     
                     
                       
#'   - `opt_details`: list with optimization details including convergence code.

                     
                   
                   
-                    400
+                    396
                     
                     
                       
#'   - `tmb_object`: original `TMB` object created with [TMB::MakeADFun()].

                     
                   
                   
-                    401
+                    397
                     
                     
                       
#'   - `tmb_data`: input.

                     
                   
                   
-                    402
+                    398
                     
                     
                       
#'

                     
                   
                   
-                    403
+                    399
                     
                     
                       
#' @details Instead of inverting or decomposing `beta_vcov`, it can be more efficient to use its robust

                     
                   
                   
-                    404
+                    400
                     
                     
                       
#'   Cholesky decomposition `LDL^T`, therefore we return the corresponding two components `L` and `D`

                     
                   
                   
-                    405
+                    401
                     
                     
                       
#'   as well since they have been available on the `C++` side already.

                     
                   
                   
-                    406
+                    402
                     
                     
                       
#'

                     
                   
                   
-                    407
+                    403
                     
                     
                       
#' @keywords internal

                     
                   
                   
-                    408
+                    404
                     
                     
                       
h_mmrm_tmb_fit <- function(tmb_object,

                     
                   
                   
-                    409
+                    405
                     
                     
                       
                           tmb_opt,

                     
                   
                   
-                    410
+                    406
                     
                     
                       
                           formula_parts,

                     
                   
                   
-                    411
+                    407
                     
                     
                       
                           tmb_data) {

                     
                   
                   
-                    412
+                    408
                     196x
                     
                       
  assert_list(tmb_object)

                     
                   
                   
-                    413
+                    409
                     196x
                     
                       
  assert_subset(c("fn", "gr", "par", "he"), names(tmb_object))

                     
                   
                   
-                    414
+                    410
                     196x
                     
                       
  assert_list(tmb_opt)

                     
                   
                   
-                    415
+                    411
                     196x
                     
                       
  assert_subset(c("par", "objective", "convergence", "message"), names(tmb_opt))

                     
                   
                   
-                    416
+                    412
                     196x
                     
                       
  assert_class(formula_parts, "mmrm_tmb_formula_parts")

                     
                   
                   
-                    417
+                    413
                     196x
                     
                       
  assert_class(tmb_data, "mmrm_tmb_data")

                     
                   
                   
-                    418
+                    414
                     
                     
                       

                     
                   
                   
-                    419
+                    415
                     196x
                     
                       
  tmb_report <- tmb_object$report(par = tmb_opt$par)

                     
                   
                   
-                    420
+                    416
                     196x
                     
                       
  x_matrix_cols <- colnames(tmb_data$x_matrix)

                     
                   
                   
-                    421
+                    417
                     196x
                     
                       
  cov <- h_mmrm_tmb_extract_cov(tmb_report, tmb_data, formula_parts$visit_var, formula_parts$is_spatial)

                     
                   
                   
-                    422
+                    418
                     196x
                     
                       
  beta_est <- tmb_report$beta

                     
                   
                   
-                    423
+                    419
                     196x
                     
                       
  names(beta_est) <- x_matrix_cols

                     
                   
                   
-                    424
+                    420
                     196x
                     
                       
  beta_vcov <- tmb_report$beta_vcov

                     
                   
                   
-                    425
+                    421
                     196x
                     
                       
  dimnames(beta_vcov) <- list(x_matrix_cols, x_matrix_cols)

                     
                   
                   
-                    426
+                    422
                     196x
                     
                       
  beta_vcov_inv_L <- tmb_report$XtWX_L # nolint

                     
                   
                   
-                    427
+                    423
                     196x
                     
                       
  beta_vcov_inv_D <- tmb_report$XtWX_D # nolint

                     
                   
                   
-                    428
+                    424
                     196x
                     
                       
  theta_est <- tmb_opt$par

                     
                   
                   
-                    429
+                    425
                     196x
                     
                       
  names(theta_est) <- NULL

                     
                   
                   
-                    430
+                    426
                     196x
                     
                       
  theta_vcov <- try(solve(tmb_object$he(tmb_opt$par)), silent = TRUE)

                     
                   
                   
-                    431
+                    427
                     196x
                     
                       
  opt_details_names <- setdiff(

                     
                   
                   
-                    432
+                    428
                     196x
                     
                       
    names(tmb_opt),

                     
                   
                   
-                    433
+                    429
                     196x
                     
                       
    c("par", "objective")

                     
                   
                   
-                    434
+                    430
                     
                     
                       
  )

                     
                   
                   
-                    435
+                    431
                     196x
                     
                       
  structure(

                     
                   
                   
-                    436
+                    432
                     196x
                     
                       
    list(

                     
                   
                   
-                    437
+                    433
                     196x
                     
                       
      cov = cov,

                     
                   
                   
-                    438
+                    434
                     196x
                     
                       
      beta_est = beta_est,

                     
                   
                   
-                    439
+                    435
                     196x
                     
                       
      beta_vcov = beta_vcov,

                     
                   
                   
-                    440
+                    436
                     196x
                     
                       
      beta_vcov_inv_L = beta_vcov_inv_L,

                     
                   
                   
-                    441
+                    437
                     196x
                     
                       
      beta_vcov_inv_D = beta_vcov_inv_D,

                     
                   
                   
-                    442
+                    438
                     196x
                     
                       
      theta_est = theta_est,

                     
                   
                   
-                    443
+                    439
                     196x
                     
                       
      theta_vcov = theta_vcov,

                     
                   
                   
-                    444
+                    440
                     196x
                     
                       
      neg_log_lik = tmb_opt$objective,

                     
                   
                   
-                    445
+                    441
                     196x
                     
                       
      formula_parts = formula_parts,

                     
                   
                   
-                    446
+                    442
                     196x
                     
                       
      data = tmb_data$data,

                     
                   
                   
-                    447
+                    443
                     196x
                     
                       
      weights = tmb_data$weights_vector,

                     
                   
                   
-                    448
+                    444
                     196x
                     
                       
      reml = as.logical(tmb_data$reml),

                     
                   
                   
-                    449
+                    445
                     196x
                     
                       
      opt_details = tmb_opt[opt_details_names],

                     
                   
                   
-                    450
+                    446
                     196x
                     
                       
      tmb_object = tmb_object,

                     
                   
                   
-                    451
+                    447
                     196x
                     
                       
      tmb_data = tmb_data

                     
                   
                   
-                    452
+                    448
                     
                     
                       
    ),

                     
                   
                   
-                    453
+                    449
                     196x
                     
                       
    class = "mmrm_tmb"

                     
                   
                   
-                    454
+                    450
                     
                     
                       
  )

                     
                   
                   
-                    455
+                    451
                     
                     
                       
}

                     
                   
                   
-                    456
+                    452
                     
                     
                       

                     
                   
                   
-                    457
+                    453
                     
                     
                       
#' Low-Level Fitting Function for MMRM

                     
                   
                   
-                    458
+                    454
                     
                     
                       
#'

                     
                   
                   
-                    459
+                    455
                     
                     
                       
#' @description `r lifecycle::badge("experimental")`

                     
                   
                   
-                    460
+                    456
                     
                     
                       
#'

                     
                   
                   
-                    461
+                    457
                     
                     
                       
#' This is the low-level function to fit an MMRM. Note that this does not

                     
                   
                   
-                    462
+                    458
                     
                     
                       
#' try different optimizers or adds Jacobian information etc. in contrast to

                     
                   
                   
-                    463
+                    459
                     
                     
                       
#' [mmrm()].

                     
                   
                   
-                    464
+                    460
                     
                     
                       
#'

                     
                   
                   
-                    465
+                    461
                     
                     
                       
#' @param formula (`formula`)\cr model formula with exactly one special term

                     
                   
                   
-                    466
+                    462
                     
                     
                       
#'   specifying the visits within subjects, see details.

                     
                   
                   
-                    467
+                    463
                     
                     
                       
#' @param data (`data.frame`)\cr input data containing the variables used in

                     
                   
                   
-                    468
+                    464
                     
                     
                       
#'   `formula`.

                     
                   
                   
-                    469
+                    465
                     
                     
                       
#' @param weights (`vector`)\cr input vector containing the weights.

                     
                   
                   
-                    470
+                    466
                     
                     
                       
#' @inheritParams h_mmrm_tmb_data

                     
                   
                   
-                    471
+                    467
                     
                     
                       
#' @param covariance (`cov_struct`)\cr A covariance structure type definition,

                     
                   
                   
-                    472
+                    468
                     
                     
                       
#'   or value that can be coerced to a covariance structure using

                     
                   
                   
-                    473
+                    469
                     
                     
                       
#'   [as.cov_struct()]. If no value is provided, a structure is derived from

                     
                   
                   
-                    474
+                    470
                     
                     
                       
#'   the provided formula.

                     
                   
                   
-                    475
+                    471
                     
                     
                       
#' @param control (`mmrm_control`)\cr list of control options produced by

                     
                   
                   
-                    476
+                    472
                     
                     
                       
#'   [mmrm_control()].

                     
                   
                   
-                    477
+                    473
                     
                     
                       
#' @inheritParams fit_single_optimizer

                     
                   
                   
-                    478
+                    474
                     
                     
                       
#'

                     
                   
                   
-                    479
+                    475
                     
                     
                       
#' @return List of class `mmrm_tmb`, see [h_mmrm_tmb_fit()] for details.

                     
                   
                   
-                    480
+                    476
                     
                     
                       
#'   In addition, it contains elements `call` and `optimizer`.

                     
                   
                   
-                    481
+                    477
                     
                     
                       
#'

                     
                   
                   
-                    482
+                    478
                     
                     
                       
#' @details

                     
                   
                   
-                    483
+                    479
                     
                     
                       
#' The `formula` typically looks like:

                     
                   
                   
-                    484
+                    480
                     
                     
                       
#'

                     
                   
                   
-                    485
+                    481
                     
                     
                       
#' `FEV1 ~ RACE + SEX + ARMCD * AVISIT + us(AVISIT | USUBJID)`

                     
                   
                   
-                    486
+                    482
                     
                     
                       
#'

                     
                   
                   
-                    487
+                    483
                     
                     
                       
#' which specifies response and covariates as usual, and exactly one special term

                     
                   
                   
-                    488
+                    484
                     
                     
                       
#' defines which covariance structure is used and what are the visit and

                     
                   
                   
-                    489
+                    485
                     
                     
                       
#' subject variables.

                     
                   
                   
-                    490
+                    486
                     
                     
                       
#'

                     
                   
                   
-                    491
+                    487
                     
                     
                       
#' Always use only the first optimizer if multiple optimizers are provided.

                     
                   
                   
-                    492
+                    488
                     
                     
                       
#'

                     
                   
                   
-                    493
+                    489
                     
                     
                       
#' @export

                     
                   
                   
-                    494
+                    490
                     
                     
                       
#'

                     
                   
                   
-                    495
+                    491
                     
                     
                       
#' @examples

                     
                   
                   
-                    496
+                    492
                     
                     
                       
#' formula <- FEV1 ~ RACE + SEX + ARMCD * AVISIT + us(AVISIT | USUBJID)

                     
                   
                   
-                    497
+                    493
                     
                     
                       
#' data <- fev_data

                     
                   
                   
-                    498
+                    494
                     
                     
                       
#' system.time(result <- fit_mmrm(formula, data, rep(1, nrow(fev_data))))

                     
                   
                   
-                    499
+                    495
                     
                     
                       
fit_mmrm <- function(formula,

                     
                   
                   
-                    500
+                    496
                     
                     
                       
                     data,

                     
                   
                   
-                    501
+                    497
                     
                     
                       
                     weights,

                     
                   
                   
-                    502
+                    498
                     
                     
                       
                     reml = TRUE,

                     
                   
                   
-                    503
+                    499
                     
                     
                       
                     covariance = NULL,

                     
                   
                   
-                    504
+                    500
                     
                     
                       
                     tmb_data,

                     
                   
                   
-                    505
+                    501
                     
                     
                       
                     formula_parts,

                     
                   
                   
-                    506
+                    502
                     
                     
                       
                     control = mmrm_control()) {

                     
                   
                   
-                    507
+                    503
                     204x
                     
                       
  if (missing(formula_parts) || missing(tmb_data)) {

                     
                   
                   
-                    508
+                    504
                     67x
                     
                       
    covariance <- h_reconcile_cov_struct(formula, covariance)

                     
                   
                   
-                    509
+                    505
                     65x
                     
                       
    formula_parts <- h_mmrm_tmb_formula_parts(formula, covariance)

                     
                   
                   
-                    510
+                    506
                     
                     
                       

                     
                   
                   
-                    511
+                    507
                     65x
                     
                       
    if (!formula_parts$is_spatial && !is.factor(data[[formula_parts$visit_var]])) {

                     
                   
                   
-                    512
+                    508
                     1x
                     
                       
      stop("Time variable must be a factor for non-spatial covariance structures")

                     
                   
                   
-                    513
+                    509
                     
                     
                       
    }

                     
                   
                   
-                    514
+                    510
                     
                     
                       

                     
                   
                   
-                    515
+                    511
                     64x
                     
                       
    assert_class(control, "mmrm_control")

                     
                   
                   
-                    516
+                    512
                     64x
                     
                       
    assert_list(control$optimizers, min.len = 1)

                     
                   
                   
-                    517
+                    513
                     64x
                     
                       
    assert_numeric(weights, any.missing = FALSE)

                     
                   
                   
-                    518
+                    514
                     64x
                     
                       
    assert_true(all(weights > 0))

                     
                   
                   
-                    519
+                    515
                     64x
                     
                       
    tmb_data <- h_mmrm_tmb_data(

                     
                   
                   
-                    520
+                    516
                     64x
                     
                       
      formula_parts, data, weights, reml,

                     
                   
                   
-                    521
+                    517
                     64x
                     
                       
      singular = if (control$accept_singular) "drop" else "error", drop_visit_levels = control$drop_visit_levels

                     
                   
                   
-                    522
+                    518
                     
                     
                       
    )

                     
                   
                   
-                    523
+                    519
                     
                     
                       
  } else {

                     
                   
                   
-                    524
+                    520
                     137x
                     
                       
    assert_class(tmb_data, "mmrm_tmb_data")

                     
                   
                   
-                    525
+                    521
                     137x
                     
                       
    assert_class(formula_parts, "mmrm_tmb_formula_parts")

                     
                   
                   
-                    526
+                    522
                     
                     
                       
  }

                     
                   
                   
-                    527
+                    523
                     201x
                     
                       
  tmb_parameters <- h_mmrm_tmb_parameters(formula_parts, tmb_data, start = control$start, n_groups = tmb_data$n_groups)

                     
                   
                   
-                    528
+                    524
                     
                     
                       

                     
                   
                   
-                    529
+                    525
                     201x
                     
                       
  tmb_object <- TMB::MakeADFun(

                     
                   
                   
-                    530
+                    526
                     201x
                     
                       
    data = tmb_data,

                     
                   
                   
-                    531
+                    527
                     201x
                     
                       
    parameters = tmb_parameters,

                     
                   
                   
-                    532
+                    528
                     201x
                     
                       
    hessian = TRUE,

                     
                   
                   
-                    533
+                    529
                     201x
                     
                       
    DLL = "mmrm",

                     
                   
                   
-                    534
+                    530
                     201x
                     
                       
    silent = TRUE

                     
                   
                   
-                    535
+                    531
                     
                     
                       
  )

                     
                   
                   
-                    536
+                    532
                     201x
                     
                       
  h_mmrm_tmb_assert_start(tmb_object)

                     
                   
                   
-                    537
+                    533
                     201x
                     
                       
  used_optimizer <- control$optimizers[[1L]]

                     
                   
                   
-                    538
+                    534
                     201x
                     
                       
  used_optimizer_name <- names(control$optimizers)[1L]

                     
                   
                   
-                    539
+                    535
                     201x
                     
                       
  args <- with(

                     
                   
                   
-                    540
+                    536
                     201x
                     
                       
    tmb_object,

                     
                   
                   
-                    541
+                    537
                     201x
                     
                       
    c(

                     
                   
                   
-                    542
+                    538
                     201x
                     
                       
      list(par, fn, gr),

                     
                   
                   
-                    543
+                    539
                     201x
                     
                       
      attr(used_optimizer, "args")

                     
                   
                   
-                    544
+                    540
                     
                     
                       
    )

                     
                   
                   
-                    545
+                    541
                     
                     
                       
  )

                     
                   
                   
-                    546
+                    542
                     201x
                     
                       
  if (identical(attr(used_optimizer, "use_hessian"), TRUE)) {

                     
                   
                   
-                    547
+                    543
                     4x
                     
                       
    args$hessian <- tmb_object$he

                     
                   
                   
-                    548
+                    544
                     
                     
                       
  }

                     
                   
                   
-                    549
+                    545
                     201x
                     
                       
  tmb_opt <- do.call(

                     
                   
                   
-                    550
+                    546
                     201x
                     
                       
    what = used_optimizer,

                     
                   
                   
-                    551
+                    547
                     201x
                     
                       
    args = args

                     
                   
                   
-                    552
+                    548
                     
                     
                       
  )

                     
                   
                   
-                    553
+                    549
                     
                     
                       
  # Ensure negative log likelihood is stored in `objective` element of list.

                     
                   
                   
-                    554
+                    550
                     194x
                     
                       
  if ("value" %in% names(tmb_opt)) {

                     
                   
                   
-                    555
+                    551
                     189x
                     
                       
    tmb_opt$objective <- tmb_opt$value

                     
                   
                   
-                    556
+                    552
                     189x
                     
                       
    tmb_opt$value <- NULL

                     
                   
                   
-                    557
+                    553
                     
                     
                       
  }

                     
                   
                   
-                    558
+                    554
                     194x
                     
                       
  fit <- h_mmrm_tmb_fit(tmb_object, tmb_opt, formula_parts, tmb_data)

                     
                   
                   
-                    559
+                    555
                     194x
                     
                       
  h_mmrm_tmb_check_conv(tmb_opt, fit)

                     
                   
                   
-                    560
+                    556
                     194x
                     
                       
  fit$call <- match.call()

                     
                   
                   
-                    561
+                    557
                     194x
                     
                       
  fit$call$formula <- formula_parts$formula

                     
                   
                   
-                    562
+                    558
                     194x
                     
                       
  fit$optimizer <- used_optimizer_name

                     
                   
                   
-                    563
+                    559
                     194x
                     
                       
  fit

                     
                   
                   
-                    564
+                    560
                     
                     
                       
}

@@ -15114,14 +15086,14 @@ 
mmrm coverage - 96.71%

                 
               
             

-