Introduction to BLiSS method"

knitr::opts_chunk$set(
  collapse = TRUE,
  comment = "#>"
)
  library(bliss)

This vignette describes step by step how to use the BLiSS method, by:

One single functional covariate

Simulation of the data set

To obtain data, several characteristics must be specified to simulate the data, such as n (number of observations), p (number of measurement moments) curves $x_{i}(.)$), beta$_$types (the shape of the coefficient function), and b$_$inf and b$_$sup (to define the domain of curves $x_{i}(.)$). Based on these parameters, we can use the sim function in the following code to simulate $x_{i}(.)$ curves, and real values $y_{i}$, from the functional linear regression model.

  param <- list(                        # define the "param" to simulate data
                Q=1,                    # the number of functional covariate
                n=100,                  # n is the sample size and p is the
                p=c(50),                # number of time observations of the curves
                beta_types=c("smooth"), # define the shape of the "true" coefficient function
                grids_lim=list(c(0,1))) # Give the beginning and the end of the observation's domain of the functions.

  data <- sim(param) # Simulate the data

Apply the Bliss method

To obtain an a posteriori sample, we use the Gibbs algorithm. We use the main function fit$_$Bliss which calls sub-functions that allow us

This main function requires a param list containing

  param <- list(            # define the required values of the Bliss method.
                iter=1e3,   # The number of iteration of the main numerical algorithm of Bliss.
                burnin=2e2, # The number of burnin iteration for the Gibbs Sampler
                K=c(3))     # The number of intervals of the beta


  res_bliss<-fit_Bliss(data=data,param=param,verbose=TRUE)

  # Structure of a Bliss object
  str(res_bliss)

Plot the result

We give here the code to obtain representations of the a posteriori distribution. First, we give the code to obtain a posteriori probabilities $\alpha(t|D)$, relative to the support. Then, the image$_$Bliss function is used to represent the subsequent distribution of the coefficient function.

  param$ylim <- range(range(res_bliss$beta_posterior_density[[1]]$grid_beta_t),
                      c(-5,5))
  param$cols <- rev(heat.colors(100))
  image_Bliss(res_bliss$beta_posterior_density,param,q=1)
  lines(res_bliss$data$grids[[1]],res_bliss$Bliss_estimate[[1]],type="s",lwd=2)
  lines(res_bliss$data$grids[[1]],res_bliss$data$betas[[1]],col=2,lwd=2,type="s")
  ylim <- range(range(res_bliss$Bliss_estimate[[1]]),
                range(res_bliss$Smooth_estimate[[1]]))
  plot_bliss(res_bliss$data$grids[[1]],
             res_bliss$Bliss_estimate[[1]],lwd=2)
  lines_bliss(res_bliss$data$grids[[1]],
             res_bliss$Smooth_estimate[[1]],lty=2)

Several functional covariates

To avoid execution lengths, this section is not executed. Please, give it a try.

Simulate the dataset

   param <- list(Q=2,
                 n=300,
                 p=c(40,60),
                 beta_shapes=c("simple","smooth"),
                 grids_lim=list(c(0,1),c(0,2)))

  data <- sim(param)

Apply the Bliss method

  param <- list(       # define the required values of the Bliss method.
     iter=1e3,         # The number of iteration of the main numerical algorithm of Bliss.
     burnin=2e2,       # The number of burnin iteration for the Gibbs Sampler
     K=c(3,3))         # The number of intervals of the beta

  res_Bliss_mult <- fit_Bliss(data=data,param=param)

Plot the result

   q <- 1
   param$ylim <- range(range(res_Bliss_mult$beta_posterior_density[[q]]$grid_beta_t),
                       c(-5,5))
   param$cols <- rev(heat.colors(100))
   image_Bliss(res_Bliss_mult$beta_posterior_density,param,q=q)
   lines(res_Bliss_mult$data$grids[[q]],res_Bliss_mult$Bliss_estimate[[q]],type="s",lwd=2)
   lines(res_Bliss_mult$data$grids[[q]],res_Bliss_mult$data$betas[[q]],col=2,lwd=2,type="s")

  ylim <- range(range(res_Bliss_mult$Bliss_estimate[[q]]),
                 range(res_Bliss_mult$Smooth_estimate[[q]]))
   plot_bliss(res_Bliss_mult$data$grids[[q]],
              res_Bliss_mult$Bliss_estimate[[q]],lwd=2,ylim=ylim)
   lines(res_Bliss_mult$data$grids[[q]],
         res_Bliss_mult$Smooth_estimate[[q]],lty=2)


   q <- 2
   param$ylim <- range(range(res_Bliss_mult$beta_posterior_density[[q]]$grid_beta_t),
                       c(-5,5))
   param$cols <- rev(heat.colors(100))
   image_Bliss(res_Bliss_mult$beta_posterior_density,param,q=q)
   lines(res_Bliss_mult$data$grids[[q]],res_Bliss_mult$Bliss_estimate[[q]],type="s",lwd=2)
   lines(res_Bliss_mult$data$grids[[q]],res_Bliss_mult$data$betas[[q]],col=2,lwd=2,type="l")

   ylim <- range(range(res_Bliss_mult$Bliss_estimate[[q]]),
                 range(res_Bliss_mult$Smooth_estimate[[q]]))
   plot_bliss(res_Bliss_mult$data$grids[[q]],
              res_Bliss_mult$Bliss_estimate[[q]],lwd=2,ylim=ylim)
   lines(res_Bliss_mult$data$grids[[q]],
         res_Bliss_mult$Smooth_estimate[[q]],lty=2)

Session informations

  sessionInfo()


Try the bliss package in your browser

Any scripts or data that you put into this service are public.

bliss documentation built on July 13, 2020, 9:06 a.m.