FLM: Functional Linear Models

View source: R/FLM.R

FLMR Documentation

Functional Linear Models

Description

Functional linear models for scalar or functional responses and functional predictors.

Usage

FLM(Y, X, XTest = NULL, optnsListY = NULL, optnsListX = NULL, nPerm = NULL)

Arguments

Y

Either an n-dimensional vector whose elements consist of scalar responses, or a list which contains functional responses in the form of a list LY and the time points LT at which they are observed (i.e., list(Ly = LY,Lt = LT)).

X

A list of lists which contains the observed functional predictors list Lxj and the time points list Ltj at which they are observed. It needs to be of the form list(list(Ly = Lx1,Lt = Lxt1),list(Ly = Lx2,Lt = Lxt2),...)

XTest

A list which contains the values of functional predictors for a held-out testing set.

optnsListY

A list of options control parameters for the response specified by list(name=value). See ‘Details’ in FPCA.

optnsListX

A list of options control parameters for the predictors specified by list(name=value). See ‘Details’ in FPCA.

nPerm

If this argument is specified, perform a permutation test to obtain the (global) p-value for the test of regression relationship between X and Y. Recommend to set to 1000 or larger if specified.

Value

A list of the following:

alpha

A length-one numeric if the response Y is scalar. Or a vector of length(workGridY) of the fitted constant alpha(t) in the linear model if Y is functional.

betaList

A list of fitted beta(s) vectors, one per predictor, if Y is scalar. Each of dimension length(workGridX[[j]]).

Or a list of fitted beta(s,t) matrices, one per predictor, if Y is functional. Each of dimension length(workGridX[[j]]) times length(workGridY).

R2

The functional R2

pv

Permutation p-value based on the functional R2

yHat

A length n vector if Y is scalar.

Or an n by length(workGridY) matrix of fitted Y's from the model if Y is functional.

yPred

Same as YHat if XTest is not provided.

Or a length length(XTest[[1]]$Ly) vector of predicted Y's if Y is scalar.

Or a length(XTest[[1]]$Ly) by length(workGridY) matrix of predicted Y's if Y is functional.

estXi

A list of n by k_j matrices of estimated functional principal component scores of predictors, where k_j is the number of eigenfunctions selected for each predictor.

testXi

A list of n by k_j matrices of estimated functional principal component scores of predictors in XTest, with eigenfunctions fitted only with X.

lambdaX

A length sum_j k_j vector of estimated eigenvalues for predictors.

workGridX

A list of vectors, each is a working grid for a predictor.

phiY

A length(workGridY) by k_y the estimated eigenfunctions of Y's, where k_y is number of eigenfunctions selected for Y. NULL if Y is scalar.

workGridY

A vector of working grid of the response Y's. NULL if Y is scalar

References

Yao, F., Müller, H.G., Wang, J.L. (2005). Functional linear regression analysis for longitudinal data. Annals of Statistics 33, 2873–2903. Hall, P., Horowitz, J.L. (2007). Methodology and convergence rates for functional linear regression. The Annals of Statistics, 35(1), 70–91.

Examples

set.seed(1000)

library(MASS)

### functional covariate
phi1 <- function(t,k) sqrt(2)*sin(2*pi*k*t)
phi2 <- function(t,k) sqrt(2)*cos(2*pi*k*t)

lambdaX <- c(1,0.7)

# training set
n <- 50
Xi <- matrix(rnorm(2*n),nrow=n,ncol=2)

denseLt <- list(); denseLy <- list()
sparseLt <- list(); sparseLy <- list()

t0 <- seq(0,1,length.out=51)
for (i in 1:n) {
 denseLt[[i]] <- t0
  denseLy[[i]] <- lambdaX[1]*Xi[i,1]*phi1(t0,1) + lambdaX[2]*Xi[i,2]*phi1(t0,2)
  
  ind <- sort(sample(1:length(t0),3))
  sparseLt[[i]] <- t0[ind]
  sparseLy[[i]] <- denseLy[[i]][ind]
}

denseX <- list(Ly=denseLy,Lt=denseLt)
sparseX <- list(Ly=sparseLy,Lt=sparseLt)

denseX <- list(X=denseX)
sparseX <- list(X=sparseX)

# test set
N <- 30

XiTest <- matrix(rnorm(2*N),nrow=N,ncol=2)

denseLtTest <- list(); denseLyTest <- list()

sparseLtTest <- list(); sparseLyTest <- list()

t0 <- seq(0,1,length.out=51)
for (i in 1:N) {
  denseLtTest[[i]] <- t0
  denseLyTest[[i]] <- lambdaX[1]*XiTest[i,1]*phi1(t0,1) + lambdaX[2]*XiTest[i,2]*phi1(t0,2)
  
  ind <- sort(sample(1:length(t0),5))
  sparseLtTest[[i]] <- t0[ind]
  sparseLyTest[[i]] <- denseLyTest[[i]][ind]
}

denseXTest <- list(Ly=denseLyTest,Lt=denseLtTest)
sparseXTest <- list(Ly=sparseLyTest,Lt=sparseLtTest)

denseXTest <- list(X=denseXTest)
sparseXTest <- list(X=sparseXTest)


### scalar response
beta <- c(1, -1)
Y <- c(Xi%*%diag(lambdaX)%*%beta) + rnorm(n,0,0.5)
YTest <- c(XiTest%*%diag(lambdaX)%*%beta) + rnorm(N,0,0.5)

## dense
denseFLM <- FLM(Y=Y,X=denseX,XTest=denseXTest,optnsListX=list(FVEthreshold=0.95))

trueBetaList <- list()
trueBetaList[[1]] <- cbind(phi1(denseFLM$workGridX[[1]],1),phi1(denseFLM$workGridX[[1]],2))%*%beta

# coefficient function estimation error (L2-norm)
plot(denseFLM$workGridX[[1]],denseFLM$betaList[[1]],type='l',xlab='t',ylab=paste('beta',1,sep=''))
points(denseFLM$workGridX[[1]],trueBetaList[[1]],type='l',col=2)

denseEstErr <-
  sqrt(trapzRcpp(denseFLM$workGridX[[1]],(denseFLM$betaList[[1]] - trueBetaList[[1]])^2))
denseEstErr

op <- par(mfrow=c(1,2))
plot(denseFLM$yHat,Y,xlab='fitted Y', ylab='observed Y')
abline(coef=c(0,1),col=8)
plot(denseFLM$yPred,YTest,xlab='predicted Y', ylab='observed Y')
abline(coef=c(0,1),col=8)
par(op)

# prediction error
densePredErr <- sqrt(mean((YTest - denseFLM$yPred)^2))
densePredErr


functionaldata/tPACE documentation built on Aug. 16, 2022, 8:27 a.m.