Package 'fdarep'

Title: Two-dimensional FPCA, Marginal FPCA, and Product FPCA for Repeated Functional Data
Description: Provides an implementation of two-dimensional functional principal component analysis (FPCA), Marginal FPCA, and Product FPCA for repeated functional data. Marginal and Product FPCA implementations are done for both dense and sparsely observed functional data. References: Chen, K., Delicado, P., & Müller, H. G. (2017) <doi:10.1111/rssb.12160>. Chen, K., & Müller, H. G. (2012) <doi:10.1080/01621459.2012.734196>. Hall, P., Müller, H.G. and Wang, J.L. (2006) <doi:10.1214/009053606000000272>. Yao, F., Müller, H. G., & Wang, J. L. (2005) <doi:10.1198/016214504000001745>.
Authors: Poorbita Kundu [aut, cre], Changbo Zhou [aut], Kehui Chen [aut], Pedro Delicado [aut], Su I Iao [aut], Hang Chen [aut], Han Chen [aut], Muqing Cui [aut], Hans-Georg Müller [cph, ths, aut]
Maintainer: Poorbita Kundu <[email protected]>
License: BSD_3_clause + file LICENSE
Version: 0.1.0
Built: 2025-02-23 03:47:18 UTC
Source: https://github.com/functionaldata/tfdarep

Help Index


Two-dimensional Functional Principal Component Analysis for dense functional data.

Description

Note: The code works for dense functional data observed on a regular grid, missing values are allowed.

Usage

Dense2dFPCA(
  X.age.year,
  n,
  num.years,
  num.ages,
  fpca.op = list(),
  pc.num = NULL
)

Arguments

X.age.year

An n by (num.years*num.ages) input data matrix, such that the ith row of the matrix gives the observed values for the ith individual. The values in each row are sorted first by years (dimension 1) and then by ages (dimension 2) within each year.

n

The sample size.

num.years

Dimension 1

num.ages

Dimension 2

fpca.op

A list of options control parameters specified by list(name=value) for the two-dimesnional FPCA; check fdapace::FPCA for available control options and default settings.

pc.num

A scalar denoting the maximum number of components to consider for the two-dimensional FPCA; default: chosen by FVE if NULL.

Details

The code works for dense functional data (with missing values), with density in both the direction of (age) dimension 2 and (year) dimension 1.

Value

A list containing the following fields:

mu

An num.ages by num.years matrix containing the bivariate mean function estimate.

pc.num

A scalar denoting the selected number of components for the two-dimensional FPCA.

res.2D.FPCA

A list containing the FPCA output for the fitted two-dimensional FPCA.

scores

An n by pc.num matrix of the estimated scores, such that the ith row of the matrix comprises estimated scores for the ith individual.

eig

An (num.years*num.ages) by pc.num matrix of the estimated product eigen functions. The estimated eigenfunctions in the otput eig are in the form of a vector rather than a matrix. For example, the first column in eig gives the first estimated eigenfunction such that gamma(s,t) -> eig[ ( (s-1)*num.ages + t ), 1] where LHS is the bivariate function in the usual form and RHS gives the corresponding element in the output vector. The rows are sorted first by years (dimension 1) and then by ages (dimension 2) within each year.

FVE

A vector of length pc.num, indicating the fraction of total variance explained by each product function, with corresponding 'FVEthreshold'.

References

  • Chen, K., Delicado, P., & Müller, H. G. (2017). Modelling function-valued stochastic processes, with applications to fertility dynamics. Journal of the Royal Statistical Society Series B: Statistical Methodology, 79(1), 177-196.

  • Chen, K., & Müller, H. G. (2012). Modeling repeated functional observations. Journal of the American Statistical Association, 107(500), 1599-1609.

  • Hall, P., Müller, H.G. and Wang, J.L. (2006). Properties of principal component methods for functional and longitudinal data analysis. Annals of Statistics, 34(3), 1493-1517.

  • Yao, F., Müller, H. G., & Wang, J. L. (2005). Functional data analysis for sparse longitudinal data. Journal of the American statistical association, 100(470), 577-590.

Examples

n <- 100 ### sample size
N <- 10
num.ages <- 20 ### dimension 2
num.years <- 15 ### dimension 1
dense_grid <- seq(0,1,length=N)
Lt <- list()
Ly <- list()
for (i in 1:n) {
  Lt[[i]] <- dense_grid ### dense time grid
  y_temp <- matrix(0,num.ages,num.years)
  for (s in 1:num.ages) {
    for (t in 1:num.years) {
      y_temp[s,t] <- y_temp[s,t]+cos(Lt[[i]][t])+rnorm(1,0,0.5)
    }
  }
  Ly[[i]] <- y_temp ### dense functional data
}
X.age.year <- matrix(0,n,num.years*num.ages)
for (i in 1:n) {
  X.age.year[i,] <- as.vector(Ly[[i]]) ### data matrix
}
res <- Dense2dFPCA(X.age.year, n , 15, 20, fpca.op=NULL,pc.num=2)
# Basic output
res$mu
res$pc.num
res$res.2D.FPCA
res$eig
res$FVE
res$pc.num
cumsum(res$FVE)

Marginal Functional Principal Component Analysis for dense functional data.

Description

Note: The code works for dense functional data observed on a regular grid, missing values are allowed, written by Kehui Chen 10/09/2015, based on the original code by Pedro Delicado.

Usage

DenseMarginalFPCA(
  X.age.year,
  n,
  num.years,
  num.ages,
  fpca.op1 = list(),
  fpca.op2 = list(),
  pc.j = NULL,
  pc.k = NULL
)

Arguments

X.age.year

An n by (num.years*num.ages) input data matrix, such that the ith row of the matrix gives the observed values for the ith individual. The values in each row are sorted first by years (dimension 1) and then by ages (dimension 2) within each year.

n

The sample size.

num.years

Dimension 1

num.ages

Dimension 2

fpca.op1

A list of options control parameters specified by list(name=value) for FPCA in direction of (age) dimension 2; check fdapace::FPCA for available control options and default settings.

fpca.op2

A list of options control parameters specified by list(name=value) for FPCA in direction of (year) dimension 1; check fdapace::FPCA for available control options and default settings.

pc.j

A scalar denoting the maximum number of components to consider for marginal FPCA in direction of (age) dimension 2; default: chosen by FVE if NULL

pc.k

A vector of length pc.j denoting the maximum number of components to consider for nested FPCA in direction of (year) dimension 1; default: chosen by FVE if NULL

Details

The code works for dense functional data (with missing values), with density in both the direction of (age) dimension 2 and (year) dimension 1.

Value

A list containing the following fields:

Xest

An n by (num.years*num.ages) estimated matrix, based on the fitted Marginal FPCA model. The ith row of the matrix gives the observed values for the ith individual. The values in each row are sorted first by years (dimension 1) and then by ages (dimension 2) within each year.

mu

An num.ages by num.years matrix containing the bivariate mean function estimate.

pc.j

A scalar denoting the selected number of components for FPCA in direction of (age) dimension 2.

pc.k

A vector of length pc.j, denoting the selected number of components for FPCA in direction of (year) dimension 1.

scores

An n by sum(pc.k) matrix of the estimated scores, such that the ith row of the matrix comprises estimated scores chi_1,1,chi_1,2,...chi_1,pc.k[1],chi_2,1,chi_2,2,...,chi_2,pc.k[2],...,chi_pc.j,1,chi_pc.j,2,...,chi_pc.j,pc.k[j] for the ith individual.

res.psi

A list containing the FPCA output for FPCA in direction of (age) dimension 2.

res.phi

A list containing the FPCA output for FPCA in direction of (year) dimension 1.

eig

An (num.years*num.ages) by sum(pc.k) matrix of the estimated product eigen functions. The rows are sorted first by years (dimension 1) and then by ages (dimension 2) within each year. The columns are sorted as in the scores.

psi

An num.ages by pc.j matrix containing the estimated eigenfunctions for FPCA in direction of (age) dimension 2.

phi

An num.years by sum(pc.k) matrix, containing the estimated eigenfunctions for FPCA in direction of (year) dimension 1.

FVE

A vector of length sum(pc.k), indicating the fraction of total variance explained by each product function, with corresponding 'FVEthreshold'.

VarOrdered

Variance explained by each term. The terms are ordered by var(chi_jk). One can select the best model by truncating at a desired level of FVE, and use names(VarOrdered) to see the corresponding model terms.

@examples n <- 100 ### sample size num.ages <- 20 ### dimension 2 num.years <- 15 ### dimension 1 dense_grid <- seq(0,1,length=N) Lt <- list() Ly <- list() for (i in 1:n) Lt[[i]] <- dense_grid ### dense time grid y_temp <- matrix(0,num.ages,num.years) for (s in 1:num.ages) for (t in 1:num.years) y_temp[s,t] <- y_temp[s,t]+cos(Lt[[i]][t])+rnorm(1,0,0.5) Ly[[i]] <- y_temp ### dense functional data X.age.year <- matrix(0,n,num.years*num.ages) for (i in 1:n) X.age.year[i,] <- as.vector(Ly[[i]]) ### data matrix res<-DenseMarginalFPCA(X.age.year, n, 15, 20, fpca.op1=NULL, fpca.op2=NULL, pc.j = NULL, pc.k = NULL) # Basic output res$Xest res$mu res$pc.j res$pc.k res$scores res$res.psi res$psi res$FVE res$VarOrdered #Additional scores fpca_psi <- res$res.psi xi_i_j <- fpca_psi$xiEst str(fpca_psi$xiEst) #xi1,(t) scores for 1st individual xi_i_j[1:num.years,] #xi1,(t) scores for 2nd individual xi_i_j[(num.years+1):(2*num.years),]

References

  • Chen, K., Delicado, P., & Müller, H. G. (2017). Modelling function-valued stochastic processes, with applications to fertility dynamics. Journal of the Royal Statistical Society Series B: Statistical Methodology, 79(1), 177-196.

  • Chen, K., & Müller, H. G. (2012). Modeling repeated functional observations. Journal of the American Statistical Association, 107(500), 1599-1609.

  • Hall, P., Müller, H.G. and Wang, J.L. (2006). Properties of principal component methods for functional and longitudinal data analysis. Annals of Statistics, 34(3), 1493-1517.

  • Yao, F., Müller, H. G., & Wang, J. L. (2005). Functional data analysis for sparse longitudinal data. Journal of the American statistical association, 100(470), 577-590.


Product Functional Principal Component Analysis for dense functional data.

Description

Note: The code works for dense functional data observed on a regular grid, missing values are allowed, written by Kehui Chen 10/09/2015, based on the original code by Pedro Delicado.

Usage

DenseProductFPCA(
  X.age.year,
  n,
  num.years,
  num.ages,
  fpca.op1 = list(),
  fpca.op2 = list(),
  pc.j = NULL,
  pc.k = NULL
)

Arguments

X.age.year

An n by (num.years*num.ages) input data matrix, such that the ith row of the matrix gives the observed values for the ith individual. The values in each row are sorted first by years (dimension 1) and then by ages (dimension 2) within each year.

n

The sample size.

num.years

Dimension 1

num.ages

Dimension 2

fpca.op1

A list of options control parameters specified by list(name=value) for FPCA in direction of (age) dimension 2; check fdapace::FPCA for available control options and default settings.

fpca.op2

A list of options control parameters specified by list(name=value) for FPCA in direction of (year) dimension 1; check fdapace::FPCA for available control options and default settings.

pc.j

A scalar denoting the maximum number of components to consider for FPCA in direction of (age) dimension 2; default: chosen by FVE if NULL.

pc.k

A scalar denoting the maximum number of components to consider for FPCA in direction of (year) dimension 1; default: chosen by FVE if NULL.

Details

The code works for dense functional data (with missing values), with density in both the direction of (age) dimension 2 and (year) dimension 1.

Value

A list containing the following fields:

Xest

An n by (num.years*num.ages) estimated matrix, based on the fitted Product FPCA model. The ith row of the matrix gives the observed values for the ith individual. The values in each row are sorted first by years (dimension 1) and then by ages (dimension 2) within each year.

mu

An num.ages by num.years matrix containing the bivariate mean function estimate.

pc.j

A scalar denoting the selected number of components for FPCA in direction of (age) dimension 2.

pc.k

A scalar denoting the selected number of components for FPCA in direction of (year) dimension 1.

scores

An n by (pc.k*pc.j) matrix of the estimated scores, such that the ith row of the matrix comprises estimated scores chi_1,1,chi_1,2,...chi_1,pc.k,chi_2,1,chi_2,2,...,chi_2,pc.k,...,chi_pc.j,1,chi_pc.j,2,...,chi_pc.j,pc.k for the ith individual.

res.psi

A list containing the FPCA output for FPCA in direction of (age) dimension 2.

res.phi

A list containing the FPCA output for FPCA in direction of (year) dimension 1.

eig

An (num.years*num.ages) by (pc.k*pc.j) matrix of the estimated product eigen functions. The rows are sorted first by years (dimension 1) and then by ages (dimension 2) within each year. The columns are sorted as in the scores.

psi

An num.ages by pc.j matrix containing the estimated eigenfunctions from FPCA in direction of (age) dimension 2.

phi

An num.years by pc.k matrix containing the estimated eigenfunctions from FPCA in direction of (year) dimension 1.

FVE

A vector of length (pc.k*pc.j), indicating the fraction of total variance explained by each product function, with corresponding 'FVEthreshold'.

VarOrdered

Variance explained by each term. The terms are ordered by var(chi_jk). One can select the best model by truncating at a desired level of FVE, and use names(VarOrdered) to see the corresponding model terms.

@examples n <- 100 ### sample size num.ages <- 20 ### dimension 2 num.years <- 15 ### dimension 1 dense_grid <- seq(0,1,length=N) Lt <- list() Ly <- list() for (i in 1:n) Lt[[i]] <- dense_grid ### dense time grid y_temp <- matrix(0,num.ages,num.years) for (s in 1:num.ages) for (t in 1:num.years) y_temp[s,t] <- y_temp[s,t]+cos(Lt[[i]][t])+rnorm(1,0,0.5) Ly[[i]] <- y_temp ### dense functional data X.age.year <- matrix(0,n,num.years*num.ages) for (i in 1:n) X.age.year[i,] <- as.vector(Ly[[i]]) ### data matrix res<-DenseProductFPCA(X.age.year, n, 15, 20, fpca.op1=NULL, fpca.op2=NULL, pc.j = NULL, pc.k = NULL) # Basic output res$Xest res$mu res$pc.j res$pc.k res$scores res$res.psi res$psi res$FVE res$VarOrdered

References

  • Chen, K., Delicado, P., & Müller, H. G. (2017). Modelling function-valued stochastic processes, with applications to fertility dynamics. Journal of the Royal Statistical Society Series B: Statistical Methodology, 79(1), 177-196.

  • Chen, K., & Müller, H. G. (2012). Modeling repeated functional observations. Journal of the American Statistical Association, 107(500), 1599-1609.

  • Hall, P., Müller, H.G. and Wang, J.L. (2006). Properties of principal component methods for functional and longitudinal data analysis. Annals of Statistics, 34(3), 1493-1517.

  • Yao, F., Müller, H. G., & Wang, J. L. (2005). Functional data analysis for sparse longitudinal data. Journal of the American statistical association, 100(470), 577-590.


Marginal Functional Principal Component Analysis for sparse functional data.

Description

Note: The code works for sparse functional data, written by Changbo Zhou 10/09/2015, based on the original code by Pedro Delicado.

Usage

SparseMarginalFPCA(
  sSup,
  Lt,
  Ly,
  fpca.op1 = list(),
  fpca.op2 = list(),
  pc.j = NULL,
  pc.k = NULL,
  bw_mu_min = NULL,
  bw_mu_max = NULL
)

Arguments

sSup

A vector of length num.ages representing the common grid for every individual in the direction of (age) dimension 2.

Lt

A list of n vectors containing the observation time points for every individual, such that the ith element of the list comprises the num.years.i points in the direction of (year) dimension 1 at which the functional-valued stochastic process is observed for the ith individual.

Ly

A list of n matrices containing the observed values for every individual, such that the ith element is an num.ages by num.years.i matrix of observed values for the ith individual.

fpca.op1

A list of options control parameters specified by list(name=value) for FPCA in direction of (age) dimension 2; check fdapace::FPCA for available control options and default settings.

fpca.op2

A list of options control parameters specified by list(name=value) for FPCA in direction of (year) dimension 1; check fdapace::FPCA for available control options and default settings.

pc.j

A scalar denoting the maximum number of components to consider for FPCA in direction of (age) dimension 2; default: chosen by FVE if NULL.

pc.k

A vector of length pc.j denoting the maximum number of components to consider for nested FPCA in direction of (year) dimension 1; default: chosen by FVE if NULL.

bw_mu_min

The minimum bandwidth value considered for bandwidth selection in mean function estimation, such that the final bandwidth chosen by 5-fold cross validation is above this minimum value; default:NULL, bandwidth chosen by 5-fold cross validation over the default range.

bw_mu_max

The maximum bandwidth value considered for bandwidth selection in mean function estimation, such that the final bandwidth chosen by 5-fold cross validation is below this value; default:NULL, bandwidth chosen by 5-fold cross validation over the default range.

Details

This code works for sparse functional data, with the notion of sparsity defined as follows. Sparsity in the year direction (dimension 1) means that the years at which the data are observed for a country (or individual unit) are sparsely distributed. However for the ith county (or individual unit), if the data are available for a particular year (dimension 1), then it is available for all the ages (dimension 2) in sSup corresponding to that specific year. Thus along (age) dimension 2, data type is dense. The 'usergrid' control option in FPCA indicates whether to use observation grid for fitting, if false FPCA will use equidistant grid. logical - default:FALSE. Along (age) dimension 2, FPCA is done for only for sSup as observation grid. Depending on the choice of usergrid for 'fpca.op2', FPCA in (year) dimension 1 is either fitted on the observed (pooled) grid or on the internal regular grid of default length 51.

Value

A list containing the following fields:

age.grid

A vector of length num.ages, representing the grid used for fitting FPCA in the direction of (age) dimension 2, same as sSup.

year.grid

A vector of length nWorkGrid, representing the grid used for fitting FPCA in the direction of (year) dimension 1.

mu

An num.ages by nWorkGrid matrix containing bivariate mean function estimate.

bwMu

The selected bandwidth for mean function estimation.

pc.j

A scalar denoting the selected number of components for FPCA in direction of (age) dimension 2.

pc.k

A vector of length pc.j, denoting the selected number of components for FPCA in direction of (year) dimension 1.

res.psi

A list containing the FPCA output for FPCA in direction of (age) dimension 2.

res.phi

A list containing the FPCA output for FPCA in direction of (year) dimension 1.

scores

A list of pc.j matrices containing the estimated scores, such that the jth element of the list is an n by pc.k[j] matrix with its ith row comprising the estimated scores chi_j,1,...chi_j,pc.k[j] for the ith individual.

psi

An num.ages by pc.j matrix containing the estimated eigenfunctions from FPCA in direction of (age) dimension 2.

phi

A list of pc.j matrices containing the estimated eigen functions from FPCA in direction of (year) dimension 1, such that the jth element of the list is an nWorkGrid by pc.k[j] matrix.

VarOrdered

A list of pc.j vectors, containing the variance explained by each term. The terms are ordered by var(chi_jk). One can select the best model by truncating at a desired level of FVE, and use names(VarOrdered) to see the corresponding model terms.

References

  • Chen, K., Delicado, P., & Müller, H. G. (2017). Modelling function-valued stochastic processes, with applications to fertility dynamics. Journal of the Royal Statistical Society Series B: Statistical Methodology, 79(1), 177-196.

  • Chen, K., & Müller, H. G. (2012). Modeling repeated functional observations. Journal of the American Statistical Association, 107(500), 1599-1609.

  • Hall, P., Müller, H.G. and Wang, J.L. (2006). Properties of principal component methods for functional and longitudinal data analysis. Annals of Statistics, 34(3), 1493-1517.

  • Yao, F., Müller, H. G., & Wang, J. L. (2005). Functional data analysis for sparse longitudinal data. Journal of the American statistical association, 100(470), 577-590.

Examples

Ly <- lapply(1:20, function(i){matrix(rnorm(13*(i)), 13, i)})
Lt <- lapply(1:20, function(i){1:(i)})
sSup <- c(1:13)
pc.j <- 2
pc.k <- c(2,3)
fpca.op1 <- NULL
fpca.op2 <- NULL
bw_mu_max <- NULL
bw_mu_min <- NULL
res <- SparseMarginalFPCA(sSup, Lt, Ly, fpca.op1, fpca.op2, pc.j, pc.k, bw_mu_min, bw_mu_max)

Product Functional Principal Component Analysis for sparse functional data.

Description

Note: The code works for sparse functional data, written by Changbo Zhou 10/09/2015, based on the original code by Pedro Delicado.

Usage

SparseProductFPCA(
  sSup,
  Lt,
  Ly,
  fpca.op1 = list(),
  fpca.op2 = list(),
  pc.j = NULL,
  pc.k = NULL,
  bw_mu_min = NULL,
  bw_mu_max = NULL
)

Arguments

sSup

A vector of length num.ages representing the common grid for every individual in the direction of (age) dimension 2.

Lt

A list of n vectors containing the observation time points for every individual, such that the ith element of the list gives the num.years.i points in the direction of (year) dimension 1 at which the functional-valued stochastic process is observed for the ith individual.

Ly

A list of n matrices containing the observed values for every individual, such that the ith element is an num.ages by num.years.i matrix of observed values for the i-th individual.

fpca.op1

A list of options control parameters specified by list(name=value) for FPCA in direction of (age) dimension 2; check fdapace::FPCA for available control options and default settings.

fpca.op2

A list of options control parameters specified by list(name=value) for FPCA in direction of (year) dimension 1; check fdapace::FPCA for available control options and default settings.

pc.j

A scalar denoting the maximum number of components to consider for FPCA in direction of (age) dimension 2; default: chosen by FVE if NULL

pc.k

A scalar denoting the maximum number of components to consider for FPCA in direction of (year) dimension 1; default: chosen by FVE if NULL.

bw_mu_min

The minimum bandwidth value considered for bandwidth selection in mean function estimation, such that the final bandwidth chosen by 5-fold cross validation is above this minimum value; default:NULL, bandwidth chosen by 5-fold cross validation over the default range.

bw_mu_max

The maximum bandwidth value considered for bandwidth selection in mean function estimation, such that the final bandwidth chosen by 5-fold cross validation is below this value; default:NULL, bandwidth chosen by 5-fold cross validation over the default range.

Details

This code works for sparse functional data, with the notion of sparsity defined as follows. Sparsity in the year direction (dimension 1) means that the years at which the data are observed for a country (or individual unit) are sparsely distributed. However for the ith county (or individual unit), if the data are available for a particular year (dimension 1), then it is available for all the ages (dimension 2) in sSup corresponding to that specific year. Thus along (age) dimension 2, data type is dense. The 'usergrid' control option in FPCA indicates whether to use observation grid for fitting, if false FPCA will use equidistant grid. logical - default:FALSE. Along (age) dimension 2, FPCA is done for only for sSup as observation grid. Depending on the choice of usergrid for 'fpca.op2', FPCA in (year) dimension 1 is either fitted on the observed (pooled) grid or on the internal regular grid of default length 51.

Value

A list containing the following fields:

age.grid

A vector of length num.ages, representing the grid used for fitting FPCA in the direction of (age) dimension 2, same as sSup.

year.grid

A vector of length nWorkGrid, representing the grid used for fitting FPCA in the direction of (year) dimension 1.

mu

An num.ages by nWorkGrid matrix containing the bivariate mean function estimate.

bwMu

The selected bandwidth for mean function estimation.

pc.j

A scalar denoting the selected number of components for FPCA in direction of (age) dimension 2.

pc.k

A scalar denoting the selected number of components for FPCA in direction of (year) dimension 1.

res.psi

A list containing the FPCA output for FPCA in direction of (age) dimension 2.

res.phi

A list containing the FPCA output for FPCA in direction of (year) dimension 1.

scores

A list of pc.j matrices containing the estimated scores, such that the jth element of the list is an n by pc.k matrix with its ith row comprising the estimated scores chi_j,1,...chi_j,pc.k for the ith individual.

psi

An num.ages by pc.j matrix containing the estimated eigenfunctions from FPCA in direction of (age) dimension 2.

phi

An nWorkGrid by pc.k matrix, containing the estimated eigenfunctions from FPCA in direction of (year) dimension 1.

VarOrdered

A list of pc.j vectors each of length pc.k, containing the variance explained by each term. The terms are ordered by var(chi_jk). One can select the best model by truncating at a desired level of FVE, and use names(VarOrdered) to see the corresponding model terms.

References

  • Chen, K., Delicado, P., & Müller, H. G. (2017). Modelling function-valued stochastic processes, with applications to fertility dynamics. Journal of the Royal Statistical Society Series B: Statistical Methodology, 79(1), 177-196.

  • Chen, K., & Müller, H. G. (2012). Modeling repeated functional observations. Journal of the American Statistical Association, 107(500), 1599-1609.

  • Hall, P., Müller, H.G. and Wang, J.L. (2006). Properties of principal component methods for functional and longitudinal data analysis. Annals of Statistics, 34(3), 1493-1517.

  • Yao, F., Müller, H. G., & Wang, J. L. (2005). Functional data analysis for sparse longitudinal data. Journal of the American statistical association, 100(470), 577-590.

Examples

Ly <- lapply(1:20, function(i){matrix(rnorm(13*(i)), 13, i)})
Lt <- lapply(1:20, function(i){1:(i)})
sSup <- c(1:13)
pc.j <- 2
pc.k <- 3
fpca.op1 <- NULL
fpca.op2 <- NULL
bw_mu_max <- 5.625000/2
bw_mu_min <- NULL
res <- SparseProductFPCA(sSup, Lt, Ly, fpca.op1, fpca.op2, pc.j, pc.k, bw_mu_min, bw_mu_max)