Estimation of functional data observations using fPCA, with applications to DTI tractography data
Sahil Seth, Dept. of Biostatistics, Johns Hopkins University
Estimation of functional data observations using fPCA, with applications to DTI tractography data
Sahil Seth, Dept. of Biostatistics, Johns Hopkins University
Help
Contact:
Sahil Seth
Dept of Biostatistics
Johns Hopkins Univ
sseth4 [at] jhu.edu
R script to run some examples:
# script to create a dataset
source("http://biostat.jhsph.edu/~saseth/R/mkDat.R")
# some plotting and other functions
source("http://biostat.jhsph.edu/~saseth/R/funcs.R")
## making a sample dataset
dat <- mkDat(I=25,noVis=5,D=30)
YFinal <- dat$Y
resFinal <- LFPCA(YFinal,dat$subject,dat$time)
################# Segment ##################
## remove some segments
fracMiss <- 0.5 #fraction of functions with missing segments
minLen <- 4 # min length of the segment
maxLen <- 10 # max length of the segment, should be less than D
rmSeg <- removeDat(YFinal,type="seg",fracMiss,minLen,maxLen,sim)
resMiss <- LFPCA(rmSeg$missDat, subject1, time1)
rmFuncs <- rmSeg$missFuncs
tsubs <- subject1[rmFuncs]
tvis <- visit1[rmFuncs]
ttim <- time1[rmFuncs]
resEst <- predict(resFinal,tsubs,ttim,tvis)
resPred <- predict(resMiss,tsubs,ttim,tvis)
## plot them in the working dir
plotFuncs(model=resFinal,est=resEst,
pred=resPred,file=”missSegs.pdf"
,type="seg",miss=rm,plotVis=TRUE)
## get rmpse for them
R2YEst <- getR2(YFinal[rmFuncs,],resEst$YHat)$R2
R2YPred <- getR2(YFinal[rmFuncs,],resPred$YHat)$R2
R2EstPred <- getR2(resEst$YHat,resPred$YHat)$R2
############ Tract ######################
## remove some tracts
rm <- removeDat(dat=YFinal,type="tract",fracMiss=fracMiss,
sim=sim,
subject=subject1,minVis=2)
## dat=YFinal;type="tract";fracMiss=fracMiss;sim=sim;
## subject=subject1;minVis=2
rmFuncs <- rm$missFuncs
tsubs <- subject1[rmFuncs]
tvis <- visit1[rmFuncs]
ttim <- time1[rmFuncs]
resMiss <- LFPCA( YFinal[-rmFuncs,], subject1[-rmFuncs],
time1[-rmFuncs])
## plot them in the working dir
resEst <- predict(resFinal,tsubs,ttim,tvis)
resPred <- predict(resMiss,tsubs,ttim)
plotFuncs(model=resFinal,est=resEst,
pred=resPred,file=”missTracts.pdf",
,type="tract",miss=rm,plotVis=FALSE)
## get rmpse for them
## both should be in the same order
R2YEst <- getR2(YFinal[rmFuncs,],resEst$YHat)$R2
## this gave spline error
R2YPred <- getR2(YFinal[rmFuncs,],resPred$YHat)$R2
R2EstPred <- getR2(resEst$YHat,resPred$YHat)$R2
Step 1: Install and load these dependencies:
require(mgcv) #gamm function
require(Matrix) #used in simulations
library(orthopolynom) #used in simulations
Step 2: Install the predictions package as follows:
Download it and cd to the directory
[http://biostat.jhsph.edu/~saseth/R/predLFPCA_1.0.tar.gz]:
R CMD INSTALL predLFPCA_1.0.tar.gz
Step 3:
Load the package:
library(predLFPCA)