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
Run the code above in your browser using DataLab