test/Simulation_REG-T.R

In tulliapadellini/energytree: Conditional Inference trees for mixed type data

##########################
#Load the libraries

library(fda.usc)
library(roahd)
library(energy)
library(entropy)
library(partykit)

##Load the classification tree function
#source("functionsREG.R")


################



############
#Load the function that create the multivariate functional object

# source("gen_data_REG_CATS.R")
# #Build the  functional object ----- define parameters
# size <- c(50,50,50,50)
#
# P <- 50
# n.var <- 1

#Generate data
# data=list()
#
# M=100
#
#
#
# for(i in 1:M){
#   alpha <- matrix( round(runif(4,0.1,1),2),
#                    nrow=4, ncol=1)
#   beta <- matrix( round(runif(4,0.1,1),2),
#                   nrow=4, ncol=1)
#   data[[i]]=gen_data_reg(size=size,P=P,n.var=n.var,alpha=alpha,beta=beta,
#                  sigma=1)
#
#
# }
# save(data, file = "Regression_Sim_Dataset.RData")



load("test/Regression_Sim_Dataset.RData")

MEPmy <- c()

MEP_b <- c()

MEP_pc <- c()

simulations = function(M){
for(i in 1:M){
  print(i)

  #select train e test set
  #id.train <- sample(1:(3*size[1]), size=0.7*3*size[1])
  #list.train <- lapply(data[[i]], function(x) {x[id.train]})
  #list.test <- lapply(data[[i]], function(x) {x[-id.train]})

  nb <- 15

  resp.list <- lapply(data, function(x) x$Y)
  cov.list <- lapply(data, function(x) x$V1)

  #REGRESSION FUNCTIONAL TREE
  myREG <- mytree(response = resp.list[[i]],
                  covariates = cov.list[i],
                  case.weights = NULL,
                  minbucket = 195,
                  alpha = 0.05,
                  R = 1000,
                  rnd.sel = T,
                  rnd.splt = TRUE,
                  nb = nb)

  #FUNCTIONAL LINEAR MODEL BASIS
  x <- cov.list[[i]]
  y <- resp.list[[i]]
  tt=x[["argvals"]]

  dataf=as.data.frame(y)

  nbasis.x=50

  basis1=create.bspline.basis(rangeval=range(tt),nbasis=nbasis.x)


  basis.x=list("x"=basis1)


  res=fregre.basis(x,y,basis.x=basis1)



  #########FUNCTIONAL LINEAR MODEL --- PC


  res2 <- fregre.pc(x, y, kmax = 7)




  ###################
  ##PREDICTION MYREG##


  foo <- min.basis(data[[i]]$V1, numbasis = 50)

  fd3 <- fdata2fd(foo$fdata.est, type.basis = "bspline", nbasis = 50)

  # TULLIA: changed here
  #foo <- min.basis(data[[i]]$V1, numbasis = nb)
  #fd3 <- fdata2fd(foo$fdata.est, type.basis = "bspline", nbasis = foo$numbasis.opt)

  m.coef <- data.frame(t(fd3$coefs))
  for(j in 1: dim(m.coef)[2]){
    names(m.coef)[j] <- j
  }
  y_pred=predict.party(myREG, newdata = m.coef)
  MEPmy[i] <- (sum((y-y_pred)^2)/length(y))/(var(y))

  #################
  ##PREDICTION FLM

  y_pred1 <- predict.fregre.fd(res, new.fdataobj=data[[i]]$V1)
  MEP_b[i] <- (sum((y-y_pred1)^2)/length(y))/(var(y))


  ###PREDICTION PC
  y_pred2 <- predict.fregre.fd(res2,data[[i]]$V1)
  MEP_pc[i] <- (sum((y-y_pred2)^2)/length(y))/(var(y))



}
  return(list(MEP_b, MEP_pc, MEPmy))

}


simulations(3)


#save(MEPmy,MEP_b,MEP_pc, file="results.RData")