library(fdaPDE)
## Load the hub data
data(hub2.5D)
hub2.5D.nodes = hub2.5D$hub2.5D.nodes
hub2.5D.triangles = hub2.5D$hub2.5D.triangles
mesh = create.mesh.2.5D(nodes = hub2.5D.nodes, triangles = hub2.5D.triangles)
## Create the Finite Element basis
FEMbasis = create.FEM.basis(mesh)
## Create a datamatrix
datamatrix = NULL
for(ii in 1:50){
a1 = rnorm(1, mean = 1, sd = 1)
a2 = rnorm(1, mean = 1, sd = 1)
a3 = rnorm(1, mean = 1, sd = 1)
func_evaluation = numeric(nrow(mesh$nodes))
for (i in 0:(nrow(mesh$nodes)-1)){
func_evaluation[i+1] = a1* sin(2*pi*mesh$nodes[i+1,1]) +
a2* sin(2*pi*mesh$nodes[i+1,2]) +
a3*sin(2*pi*mesh$nodes[i+1,3]) + 1
}
data = func_evaluation + rnorm(nrow(mesh$nodes), mean = 0, sd = 0.5)
datamatrix = rbind(datamatrix, data)
}
## Compute the mean of the datamatrix and subtract it to the data
data_bar = colMeans(datamatrix)
data_demean = matrix(rep(data_bar,50), nrow=50, byrow=TRUE)
datamatrix_demeaned = datamatrix - data_demean
## Set the smoothing parameter lambda
lambda = 0.00375
## Estimate the first 2 Principal Components
FPCA_solution = FPCA.FEM(datamatrix = datamatrix_demeaned,
FEMbasis = FEMbasis, lambda = lambda, nPC = 2)
## Plot the functional loadings of the estimated Principal Components
plot(FPCA_solution$loadings.FEM)
Run the code above in your browser using DataLab