730 covariates with different extrapolation method

NEWS.md

@@ -4,13 +4,37 @@
 
 - Error when specifying `wd` without `modName`
 
+- With Linear and midpoint of a time between two points, how `rxode2`
+  handles missing values has changed.  When the missing value is lower
+  than the requested time, it will look backward until it finds the
+  first non-missing value (or if all are missing start looking
+  forward).  When the missing value is higher than the requested time,
+  the algorithm will look forward until it finds the first non-missing
+  value (or if all are missing, start looking backward).
+
 ## Possible breaking changes (though unlikely)
 
 - `iCov` is no longer merged to the event dataset.  This makes solving
   with `iCov` slightly faster (#743)
 
 ## New features
 
+- Now you can specify the interpolation method per covariate in the model:
+
+  - `linear(var1, var2)` says both `var1` and `var2` would use linear
+    interpolation when they are a time-varying covariate. You could
+    also use `linear(var1)`
+
+  - `locf()` declares variables using last observation carried forward
+
+  - `nocb()` declares variables using next observation carried backward
+
+  - `midpoint()` declares variables using midpoint interpolation
+
+- `linear()`, `locf()`, `locb()`, `midpoint()`, `params()`, `cmt()`
+  and `dvid()` declarations are now ignored when loading a `rxode2`
+  model with `rxS()`
+
 - Empty arguments to `rxRename()` give a warning (#688)
 
 - Promoting from covariates to parameters with model piping (via `ini()`) now
@@ -38,6 +62,9 @@
 
 - When assigning reserved variables, the parser will error. See #744
 
+- Linear interpolation will now adjust the times as well as the values
+  when `NA` values are observed.
+
 ## Big change
 
 - At the request of CRAN, combine `rxode2parse`, `rxode2random`, and

R/RcppExports.R

@@ -187,7 +187,7 @@ etTransEvidIsObs <- function(isObsSexp) {
 #'   steady concentration at the actual time of dose, otherwise when
 #'   `FALSE` the doses are shifted
 #'
-#' @inheritParams rxSolve::rxSolve
+#' @inheritParams rxSolve
 #'
 #' @return Object for solving in rxode2
 #'

R/build.R

@@ -62,6 +62,22 @@ d/dt(blood)     = a*intestine - b*blood
 
 .rxodeBuildCode <- function() {
   # This builds the code needed for rxode2
+  message("Generate grammar include file")
+  dparser::mkdparse(devtools::package_file("inst/tran.g"),
+                    devtools::package_file("src/"),
+                    grammar_ident="rxode2parse")
+  l <- readLines(devtools::package_file("src/tran.g.d_parser.c"))
+  .w <- which(grepl("#line ", l))
+  if (.w > 1L) {
+    .w <- .w[1L]
+    l[.w] <- sub("[#]line([^\"]*\").*(src.*)", "#line\\1\\2",l[.w])
+  }
+  tran.g.h <- file(devtools::package_file("src/tran.g.d_parser.h"), "wb")
+  writeLines(l, tran.g.h)
+  close(tran.g.h)
+  unlink(devtools::package_file("src/tran.g.d_parser.c"))
+
+
   # generate control
   try({
     message("generate defines")

R/rxsolve.R

@@ -229,14 +229,38 @@
 #'     by solving the line between the observed covariates and extrapolating the new
 #'     covariate value.
 #'
-#' * `"constant"` -- Last observation carried forward (the default).
+#' * `"locf"` -- Last observation carried forward (the default).
 #'
-#' * `"NOCB"` -- Next Observation Carried Backward.  This is the same method
+#' * `"nocb"` -- Next Observation Carried Backward.  This is the same method
 #'       that NONMEM uses.
 #'
 #' * `"midpoint"` Last observation carried forward to midpoint; Next observation
 #'   carried backward to midpoint.
 #'
+#'   For time-varying covariates where a missing value is present, the
+#'   interpolation method will use either "locf" or "nocb" throughout
+#'   if they are the type of covariate interpolation that is selected.
+#'
+#'   When using the linear or midpoint interpolation, the lower point
+#'   in the interpolation will use locf to interpolate missing
+#'   covariates and the upper point will use the nocb to interpolate
+#'   missing covariates.
+#'
+#' @param naInterpolation specifies the interpolation method for
+#'   time-varying covariates when the instantaneous value is `NA` (not
+#'   during an explicit interpolation) and the `covsInterpolation` is
+#'   either `"midpoint"` or `"linear"`. This can be:
+#'
+#'   * `"locf"` -- last observation carried forward (default)
+#'
+#'   * `"nocb"` -- next observation carried backward.
+#'
+#'   This will look for the prior value (backwards/locf) when
+#'   instantaneously missing, or the next value when instantaneously
+#'   missing.  If all the covariates are missing and you find the
+#'   end/beginning of the individual record, switch direction.  If all
+#'   are really missing, then return missing.
+#'
 #' @param addCov A boolean indicating if covariates should be added
 #'     to the output matrix or data frame. By default this is
 #'     disabled.
@@ -661,6 +685,7 @@ rxSolve <- function(object, params = NULL, events = NULL, inits = NULL,
                     hmaxSd = 0, hini = 0, maxordn = 12L, maxords = 5L, ...,
                     cores,
                     covsInterpolation = c("locf", "linear", "nocb", "midpoint"),
+                    naInterpolation = c("locf", "nocb"),
                     addCov = TRUE, sigma = NULL, sigmaDf = NULL,
                     sigmaLower = -Inf, sigmaUpper = Inf,
                     nCoresRV = 1L, sigmaIsChol = FALSE,
@@ -794,6 +819,11 @@ rxSolve <- function(object, params = NULL, events = NULL, inits = NULL,
     } else {
       covsInterpolation <- c("linear"=0L, "locf"=1L, "nocb"=2L, "midpoint"=3L)[match.arg(covsInterpolation)]
     }
+    if (checkmate::testIntegerish(naInterpolation, len=1, lower=0, upper=2, any.missing=FALSE)) {
+      naInterpolation <- as.integer(naInterpolation)
+    } else {
+      naInterpolation <- c("locf"=1L, "nocb"=0L)[match.arg(naInterpolation)]
+    }
     if (missing(naTimeHandle) && !is.null(getOption("rxode2.naTimeHandle", NULL))) {
       naTimeHandle <- getOption("rxode2.naTimeHandle")
     }
@@ -1127,8 +1157,8 @@ rxSolve <- function(object, params = NULL, events = NULL, inits = NULL,
       addlDropSs=addlDropSs,
       ssAtDoseTime=ssAtDoseTime,
       ss2cancelAllPending=ss2cancelAllPending,
+      naInterpolation=naInterpolation,
       .zeros=unique(.zeros)
-
     )
     class(.ret) <- "rxControl"
     return(.ret)
@@ -1744,7 +1774,6 @@ rxSolve.default <- function(object, params = NULL, events = NULL, inits = NULL,
     }
     .minfo(sprintf("omega/sigma items treated as zero: '%s'", paste(.ctl$.zeros, collapse="', '")))
   }
-
   if (rxode2.debug) {
     .envReset$ret <- .collectWarnings(rxSolveSEXP(object, .ctl, .nms, .xtra,
                                                   params, events, inits,

R/symengine.R

@@ -1441,6 +1441,8 @@ rxToSE <- function(x, envir = NULL, progress = FALSE,
                )
         }
       } else {
+        if (.fun %in% c("param", "dvid", "cmt", "locf", "nocb",
+                        "midpoint", "linear")) return(NULL)
         .udf <- try(get(.fun, envir = .rxToSE.envir$parent, mode="function"), silent =TRUE)
         if (inherits(.udf, "try-error")) {
           .udf <- try(get(.fun, envir = rxode2::.udfEnvSet(NULL), mode="function"), silent =TRUE)

inst/include/rxode2_control.h

@@ -96,5 +96,6 @@
 #define Rxc_addlDropSs 91
 #define Rxc_ssAtDoseTime 92
 #define Rxc_ss2cancelAllPending 93
-#define Rxc__zeros 94
+#define Rxc_naInterpolation 94
+#define Rxc__zeros 95
 #endif // __rxode2_control_H__

inst/include/rxode2parseStruct.h

@@ -40,14 +40,14 @@ typedef struct {
   int ncov;
   char modNamePtr[1000];
   int *par_cov;
+  int *par_cov_interp;
+
   double *inits;
   double *scale;
   bool do_par_cov;
   // approx fun options
-  double f1;
-  double f2;
-  int kind;
   int is_locf;
+  int instant_backward;
   int cores;
   int doesRandom;
   int extraCmt;
@@ -214,6 +214,8 @@ typedef struct {
   //double *extraDoseIi; // ii doses unsupported
   bool lastIsSs2;
   double *timeThread;
+  int idxLow;
+  int idxHi;
 } rx_solving_options_ind;
 
 typedef struct {

inst/include/rxode2parse_control.h

@@ -95,7 +95,8 @@
 #define Rxc_addlDropSs 91
 #define Rxc_ssAtDoseTime 92
 #define Rxc_ss2cancelAllPending 93
-#define Rxc__zeros 94
+#define Rxc_naInterpolation 94
+#define Rxc__zeros 95
 #define RxMv_params 0
 #define RxMv_lhs 1
 #define RxMv_state 2
@@ -117,8 +118,9 @@
 #define RxMv_slhs 18
 #define RxMv_alag 19
 #define RxMv_udf 20
-#define RxMv_timeId 21
-#define RxMv_md5 22
+#define RxMv_interp 21
+#define RxMv_timeId 22
+#define RxMv_md5 23
 #define RxMvFlag_ncmt 0
 #define RxMvFlag_ka 1
 #define RxMvFlag_linB 2

inst/tran.g

@@ -1,22 +1,23 @@
 //loop
 statement_list : (statement)+ ;
 
-statement 
+statement
   : assignment end_statement
   | ini        end_statement
   | ini0       end_statement
   | ini0f      end_statement
   | fbio       end_statement
   | alag       end_statement
   | rate       end_statement
-  | dur        end_statement 
+  | dur        end_statement
   | derivative end_statement
   | dfdy       end_statement
   | mtime      end_statement
   | mat0       end_statement
   | matF       end_statement
   | printf_statement end_statement
   | param_statement end_statement
+  | interp_statement end_statement
   | cmt_statement end_statement
   | dvid_statementI end_statement
   | break_statement end_statement
@@ -46,6 +47,8 @@ cmt_statement
 param_statement
     : "params?" '(' (identifier_r | theta0 | theta | eta) (',' (identifier_r | theta0 | theta | eta) )*  ')';
 
+interp_statement: ('locf' | 'linear' | 'nocb' | 'midpoint') '(' (identifier_r | theta0 | theta | eta) (',' (identifier_r | theta0 | theta | eta) )* ')';
+
 printf_statement
   : printf_command '(' string (',' logical_or_expression )* ')';
 
@@ -86,10 +89,10 @@ matF: '_rxF' '=' logical_or_expression;
 
 mtime     : 'mtime' '(' identifier_r_no_output ')' ('=' | '<-' | '~') logical_or_expression;
 
-logical_or_expression : logical_and_expression 
+logical_or_expression : logical_and_expression
     (('||' | '|')  logical_and_expression)* ;
 
-logical_and_expression : equality_expression0 
+logical_and_expression : equality_expression0
     (('&&' | '&') equality_expression0)* ;
 
 equality_expression0 : equality_expression |
@@ -99,23 +102,25 @@ equality_expression0 : equality_expression |
     '(' equality_str ')' |
     '!' '(' equality_str ')' |
     '(' '!' identifier_r ')' |
-    '!' identifier_r | 
+    '!' identifier_r |
     '!' function ;
 
 equality_str : equality_str1 | equality_str2;
 equality_str1 : string ('!=' | '==' ) identifier_r;
 equality_str2 : identifier_r ('!=' | '==' ) string;
 
-equality_expression : relational_expression 
+equality_expression : relational_expression
     (('!=' | '==' ) relational_expression)* ;
 
+relational_op: '<' | '>' | '<=' | '>=' | '<-' | '->';
+
 relational_expression : additive_expression
-    (('<' | '>' | '<=' | '>=') additive_expression)* ;
+    (relational_op additive_expression)* ;
 
 additive_expression : multiplicative_expression
     (('+' | '-') multiplicative_expression)* ;
 
-multiplicative_expression : unary_expression 
+multiplicative_expression : unary_expression
     (mult_part)* ;
 
 mult_part : ('*' | '/') unary_expression ;
@@ -128,7 +133,7 @@ power_expression : primary_expression power_operator exponent_expression;
 
 power_operator   : ('^' | '**');
 
-primary_expression 
+primary_expression
   : constant
   | identifier_r
   | theta0
@@ -179,4 +184,3 @@ identifier_r_no_output_2: "[.]+[a-zA-Z_][a-zA-Z0-9_.]*" $term -4;
 identifier: "[a-zA-Z][a-zA-Z0-9_.]*" $term -4;
 whitespace: ( "[ \t\r\n]+" | singleLineComment )*;
 singleLineComment: '#' "[^\n]*";
-