
Implements functional generalized additive models for functional and scalar covariates and scalar responses.
Additionally implements functional linear models. This function is a wrapper for mgcv's gam
and its siblings to fit models of the general form
fgam(formula, fitter = NA, tensortype = c("te", "t2"), ...)
a fitted fgam-object, which is a gam
-object with some additional information
in a fgam-entry. If fitter is "gamm" or "gamm4", only the $gam part of the
returned list is modified in this way.
a formula with special terms as for gam, with additional special terms
af
(), lf
(), re
().
the name of the function used to estimate the model. Defaults to gam
if the matrix of functional responses has less than 2e5 data points and to
bam
if not. "gamm" (see gamm
) and "gamm4"
(see gamm4
) are valid options as well.
additional arguments that are valid for gam
or bam
; for example,
specify a gamma
> 1 to increase amount of smoothing when using GCV to choose smoothing
parameters or method="REML"
to change to REML for estimation of smoothing parameters
(default is GCV).
Binomial responses should be specified as a numeric vector rather than as a matrix or a factor.
Mathew W. McLean mathew.w.mclean@gmail.com and Fabian Scheipl
McLean, M. W., Hooker, G., Staicu, A.-M., Scheipl, F., and Ruppert, D. (2014). Functional generalized additive models. Journal of Computational and Graphical Statistics, 23 (1), pp. 249-269. Available at https://www.ncbi.nlm.nih.gov/pmc/articles/PMC3982924/.
af
, lf
, predict.fgam
, vis.fgam
data(DTI)
## only consider first visit and cases (no PASAT scores for controls)
y <- DTI$pasat[DTI$visit==1 & DTI$case==1]
X <- DTI$cca[DTI$visit==1 & DTI$case==1, ]
X_2 <- DTI$rcst[DTI$visit==1 & DTI$case==1, ]
## remove samples containing missing data
ind <- rowSums(is.na(X)) > 0
ind2 <- rowSums(is.na(X_2)) > 0
y <- y[!(ind | ind2)]
X <- X[!(ind | ind2), ]
X_2 <- X_2[!(ind | ind2), ]
N <- length(y)
## fit fgam using FA measurements along corpus callosum
## as functional predictor with PASAT as response
## using 8 cubic B-splines for marginal bases with third
## order marginal difference penalties
## specifying gamma > 1 enforces more smoothing when using
## GCV to choose smoothing parameters
#fit <- fgam(y ~ af(X, k = c(8, 8), m = list(c(2, 3), c(2, 3))), gamma = 1.2)
## fgam term for the cca measurements plus an flm term for the rcst measurements
## leave out 10 samples for prediction
test <- sample(N, 10)
#fit <- fgam(y ~ af(X, k = c(7, 7), m = list(c(2, 2), c(2, 2))) +
# lf(X_2, k=7, m = c(2, 2)), subset=(1:N)[-test])
#plot(fit)
## predict the ten left outs samples
#pred <- predict(fit, newdata = list(X=X[test, ], X_2 = X_2[test, ]), type='response',
# PredOutOfRange = TRUE)
#sqrt(mean((y[test] - pred)^2))
## Try to predict the binary response disease status (case or control)
## using the quantile transformed measurements from the rcst tract
## with a smooth component for a scalar covariate that is pure noise
y <- DTI$case[DTI$visit==1]
X <- DTI$cca[DTI$visit==1, ]
X_2 <- DTI$rcst[DTI$visit==1, ]
ind <- rowSums(is.na(X)) > 0
ind2 <- rowSums(is.na(X_2)) > 0
y <- y[!(ind | ind2)]
X <- X[!(ind | ind2), ]
X_2 <- X_2[!(ind | ind2), ]
z1 <- rnorm(length(y))
## select=TRUE allows terms to be zeroed out of model completely
#fit <- fgam(y ~ s(z1, k = 10) + af(X_2, k=c(7,7), m = list(c(2, 1), c(2, 1)),
# Qtransform=TRUE), family=binomial(), select=TRUE)
#plot(fit)
Run the code above in your browser using DataLab