Merge branch 'issue-#68': fixes the main issue of #68

NAMESPACE

@@ -46,6 +46,11 @@ S3method(inv_link,gam)
 S3method(inv_link,gamm)
 S3method(inv_link,glm)
 S3method(inv_link,list)
+S3method(is_factor_term,bam)
+S3method(is_factor_term,gam)
+S3method(is_factor_term,gamm)
+S3method(is_factor_term,list)
+S3method(is_factor_term,terms)
 S3method(link,bam)
 S3method(link,family)
 S3method(link,gam)
@@ -72,6 +77,9 @@ S3method(smooth_dim,gamm)
 S3method(smooth_dim,mgcv.smooth)
 S3method(smooth_samples,default)
 S3method(smooth_samples,gam)
+S3method(term_variables,bam)
+S3method(term_variables,gam)
+S3method(term_variables,terms)
 S3method(vcov,scam)
 S3method(which_smooths,bam)
 S3method(which_smooths,default)
@@ -99,6 +107,7 @@ export(inv_link)
 export(is_by_smooth)
 export(is_continuous_by_smooth)
 export(is_factor_by_smooth)
+export(is_factor_term)
 export(is_mgcv_smooth)
 export(is_mrf_smooth)
 export(is_offset)

R/datasets.R

@@ -77,3 +77,22 @@ NULL
 ##' @keywords data
 ##' @docType data
 NULL
+
+##' Data from the General Social Survey (GSS) from the National Opinion Research
+##' Center of the University of Chicago
+##'
+##' A subset of the data from the [`carData::GSSvocab`] dataset from the
+##' `carData` package, containing observations from 2016 only.
+##'
+##' @format A data frame with 1858 rows and 3 variables:
+##'   * `vocab`: numeric; the number of words out of 10 correct on a vocabulary
+##'     test.
+##'   * `nativeBorn`: factor; Was the respondent born in the US? A factor with
+##'     levels `no` and `yes`.
+##'   * `ageGroup`: factor; grouped age of the respondent with levels `18-29`
+##'       `30-39`, `40-49`, `50-59`, and `60+`.##' 
+##'
+##' @keywords data
+##' @name gss_vocab
+##' @docType data
+NULL

R/evaluate_smooth.R

@@ -569,33 +569,55 @@
     }
 
     mf <- model.frame(object)  # data used to fit model
-    is_fac <- is.factor(mf[[term]]) # is term a factor?
+
+    ## is_fac <- is.factor(mf[[term]]) # is term a factor?
+    is_fac <- is_factor_term(tt, term)
 
     ## match the specific term, with term names mgcv actually uses
     ## for example in a model with multiple linear predictors, terms in
     ## nth linear predictor (for n > 1) get appended .{n-1}  
     ind <- match(term, vars)
-
-    ## take the actual mgcv version of the names for the `terms` argument
-    evaluated <- as.data.frame(predict(object, newdata = mf, type = 'terms',
-                                       terms = mgcv_names[ind], se = TRUE,
-                                       unconditional = unconditional))
-    evaluated <- setNames(evaluated, c("partial", "se"))
-    evaluated <- as_tibble(evaluated)
-
+
     if (is_fac) {
-        levs <- levels(mf[, term])
-        newd <- setNames(data.frame(fac = factor(levs, levels = levs)), "value")
-        spl <- lapply(split(evaluated, mf[, term]), `[`, i = 1, j = )
-        evaluated <- bind_rows(spl)
+        ## check order of term; if > 1 interaction and not handled
+        ord <- attr(tt, "order")[match(term, attr(tt, "term.labels"))]
+        if (ord > 1) {
+            stop("Interaction terms are not currently supported.")
+        }
+        ## facs <- attr(tt, 'factors')[, term]
+        newd <- unique(mf[, term, drop = FALSE])
+        ## ##fac_vars <- rownames(facs)
+        ## fac_vars <- names(facs)[as.logical(facs)]
+        ## facs <- attr(tt, 'factors')[, term]
+        ## newd <- unique(mf[, names(facs)[as.logical(facs)], drop = FALSE])
+        ## ##fac_vars <- rownames(facs)
+        ## fac_vars <- names(facs)[as.logical(facs)]
+        ## ##newd <- unique(mf[, fac_vars, drop = FALSE])
+        other_vars <- setdiff(names(mf), term)
+        other_data <- as_tibble(lapply(mf[other_vars], value_closest_to_median))
+        pred_data <- exec(expand_grid, !!!list(newd, other_data))
+        evaluated <- as.data.frame(predict(object, newdata = pred_data,
+                             type = 'terms',
+                             terms = term, se = TRUE,
+                             unconditional = unconditional,
+                             newdata.guaranteed = FALSE))
+        evaluated <- setNames(evaluated, c("partial", "se"))
+        evaluated <- as_tibble(evaluated)
         nr <- NROW(evaluated)
+        newd <- setNames(newd, "value")
         evaluated <- bind_cols(term = rep(term, nr),
-                               type = rep(ifelse(is_fac, "factor", "numeric"), nr),
+                               type = rep("factor", nr),
                                newd, evaluated)
     } else {
+        ## take the actual mgcv version of the names for the `terms` argument
+        evaluated <- as.data.frame(predict(object, newdata = mf, type = 'terms',
+                                           terms = mgcv_names[ind], se = TRUE,
+                                           unconditional = unconditional))
+        evaluated <- setNames(evaluated, c("partial", "se"))
+        evaluated <- as_tibble(evaluated)
         nr <- NROW(evaluated)
         evaluated <- bind_cols(term = rep(term, nr),
-                               type = rep(ifelse(is_fac, "factor", "numeric"), nr),
+                               type = rep("numeric", nr),
                                value = mf[[term]],
                                evaluated)
     }

R/utililties.R

@@ -805,3 +805,117 @@
     out <- TRUE
     out
 }
+
+##' Is a model term a factor (categorical)?
+##'
+##' Given the name (a term label) of a term in a model, identify if the term is a
+##' factor term or numeric. This is useful when considering interactions, where
+##' terms like `fac1:fac2` or `num1:fac1` may be requested by the user. Only for
+##' terms of the type `fac1:fac2` will this function return `TRUE`.
+##'
+##' @param object an R object on which method dispatch is performed
+##' @param term character; the name of a model term, in the sense of
+##'   `attr(terms(object), "term.labels")`. Currently not checked to see if the
+##'   term exists in the model.
+##' @param ... arguments passed to other methods.
+##'
+##' @return A logical: `TRUE` if and only if all variables involved in the term
+##'   are factors, otherwise `FALSE`.
+##'
+##' @export
+`is_factor_term` <- function(object, term, ...) {
+    UseMethod("is_factor_term", object)
+}
+
+##' @rdname is_factor_term
+##' @export
+`is_factor_term.terms` <- function(object, term, ...) {
+    if (missing(term)) {
+        stop("Argument 'term' must be provided.")
+    }
+    facs <- attr(object, "factors")
+    out <- if (term %in% colnames(facs)) {
+        facs <- facs[, term, drop = FALSE]
+        take <- rownames(facs)[as.logical(facs)]
+        data_types <- attr(object, 'dataClasses')[take]
+        all(data_types == "factor")
+    } else {
+        NULL
+    }
+    out
+}
+
+##' @rdname is_factor_term
+##' @export
+`is_factor_term.gam` <- function(object, term, ...) {
+    object <- terms(object)
+    is_factor_term(object, term, ...)
+}
+
+##' @rdname is_factor_term
+##' @export
+`is_factor_term.bam` <- function(object, term, ...) {
+    object <- terms(object)
+    is_factor_term(object, term, ...)
+}
+
+##' @rdname is_factor_term
+##' @export
+`is_factor_term.gamm` <- function(object, term, ...) {
+    object <- terms(object$gam)
+    is_factor_term(object, term, ...)
+}
+
+##' @rdname is_factor_term
+##' @export
+`is_factor_term.list` <- function(object, term, ...) {
+    if (!is_gamm4(object)) {
+        if (all(vapply(object, inherits, logical(1), "terms"))) {
+            out <- any(unlist(lapply(object, is_factor_term, term)))
+        } else {
+            stop("Don't know how to handle generic list objects.")
+        }
+    } else {
+        object <- terms(object$gam)
+        out <- is_factor_term(object, term, ...)
+    }
+    out
+}
+
+##' Names of variables involved in a specified model term
+##'
+##' Given the name (a term label) of a term in a model, returns the names
+##' of the variables involved in ther term.
+##'
+##' @param object an R object on which method dispatch is performed
+##' @param term character; the name of a model term, in the sense of
+##'   `attr(terms(object), "term.labels")`. Currently not checked to see if the
+##'   term exists in the model.
+##' @param ... arguments passed to other methods.
+##'
+##' @return A character vector of variable names.
+##'
+`term_variables` <- function(object, term, ...) {
+    UseMethod("terms_variables", object)
+}
+
+##' @rdname term_variables
+##' @export
+`term_variables.terms` <- function(object, term, ...) {
+    facs <- attr(object, "factors")[ , term]
+    names(facs)[as.logical(facs)]
+}
+
+##' @rdname term_variables
+##' @export
+`term_variables.gam` <- function(object, term, ...) {
+    object <- terms(object)
+    term_variables(object, term, ...)
+}
+
+##' @rdname term_variables
+##' @export
+`term_variables.bam` <- function(object, term, ...) {
+    object <- terms(object)
+    term_variables(object, term, ...)
+}