bootstrap(kkk, y_no, ktype, K, ode_par, intp_data, www = NULL)

## Not run: 
require(mvtnorm)
noise = 0.1  ## set the variance of noise
SEED = 19537
set.seed(SEED)
## Define ode function, we use lotka-volterra model in this example.
## we have two ode states x[1], x[2] and four ode parameters alpha, beta, gamma and delta.
LV_fun = function(t,x,par_ode){
  alpha=par_ode[1]
  beta=par_ode[2]
  gamma=par_ode[3]
  delta=par_ode[4]
  as.matrix( c( alpha*x[1]-beta*x[2]*x[1] , -gamma*x[2]+delta*x[1]*x[2] ) )
}
## Define the gradient of ode function against ode parameters
## df/dalpha,  df/dbeta, df/dgamma, df/ddelta where f is the differential equation.
LV_grlNODE= function(par,grad_ode,y_p,z_p) {
alpha = par[1]; beta= par[2]; gamma = par[3]; delta = par[4]
dres= c(0)
dres[1] = sum( -2*( z_p[1,]-grad_ode[1,])*y_p[1,]*alpha )
dres[2] = sum( 2*( z_p[1,]-grad_ode[1,])*y_p[2,]*y_p[1,]*beta)
dres[3] = sum( 2*( z_p[2,]-grad_ode[2,])*gamma*y_p[2,] )
dres[4] = sum( -2*( z_p[2,]-grad_ode[2,])*y_p[2,]*y_p[1,]*delta)
dres
}

## create a ode class object
kkk0 = ode$new(2,fun=LV_fun,grfun=LV_grlNODE)
## set the initial values for each state at time zero.
xinit = as.matrix(c(0.5,1))
## set the time interval for the ode numerical solver.
tinterv = c(0,6)
## solve the ode numerically using predefined ode parameters. alpha=1, beta=1, gamma=4, delta=1.
kkk0$solve_ode(c(1,1,4,1),xinit,tinterv)

## Add noise to the numerical solution of the ode model and use it as the noisy observation.
n_o = max( dim( kkk0$y_ode) )
t_no = kkk0$t
y_no =  t(kkk0$y_ode) + rmvnorm(n_o,c(0,0),noise*diag(2))

## Create a ode class object by using the simulation data we created from the ode numerical solver.
## If users have experiment data, they can replace the simulation data with the experiment data.
## Set initial value of ode parameters.
init_par = rep(c(0.1),4)
init_yode = t(y_no)
init_t = t_no
kkk = ode$new(1,fun=LV_fun,grfun=LV_grlNODE,t=init_t,ode_par= init_par, y_ode=init_yode )

## The following examples with CPU or elapsed time > 10s

##Use function 'rkg' to estimate the ode parameters. The standard gradient matching method is coded
##in the the 'rkg' function. The parameter estimations are stored in the returned vector of 'rkg'.
## Choose a kernel type for 'rkhs' interpolation. Two options are provided 'rbf' and 'mlp'.
ktype ='rbf'
rkgres = rkg(kkk,y_no,ktype)
## show the results of ode parameter estimation using the standard gradient matching
kkk$ode_par

## Perform bootstrap procedure to estimate the median absolute deviations of ode parameters
# here we get the resulting interpolation from gradient matching using 'rkg' for each ode state
bbb = rkgres$bbb
nst = length(bbb)
intp_data = list()
for( i in 1:nst) {
    intp_data[[i]] = bbb[[i]]$predictT(bbb[[i]]$t)$pred
}
K = 12 # the number of bootstrap replicates
mads = bootstrap(kkk, y_no, ktype, K, ode_par, intp_data)

## show the results of ode parameter estimation and its uncertainty
## using the standard gradient matching
ode_par
mads

############# gradient matching + ODE regularisation
crtype='i'
lam=c(10,1,1e-1,1e-2,1e-4)
lamil1 = crossv(lam,kkk,bbb,crtype,y_no)
lambdai1=lamil1[[1]]
res = third(lambdai1,kkk,bbb,crtype)
oppar = res$oppar

### do bootstrap here for gradient matching + ODE regularisation
ode_par = oppar
K = 12
intp_data = list()
for( i in 1:nst) {
    intp_data[[i]] = res$rk3$rk[[i]]$predictT(bbb[[i]]$t)$pred
}
mads = bootstrap(kkk, y_no, ktype, K, ode_par, intp_data)
ode_par
mads

############# gradient matching + ODE regularisation + warping
###### warp state
peod = c(6,5.3) #8#9.7     ## the guessing period
eps= 1          ## the standard deviation of period
fixlens=warpInitLen(peod,eps,rkgres)
kkkrkg = kkk$clone()
www = warpfun(kkkrkg,bbb,peod,eps,fixlens,y_no,kkkrkg$t)

### do bootstrap here for gradient matching + ODE regularisation + warping
nst = length(bbb)
K = 12
ode_par = www$wkkk$ode_par
intp_data = list()
for( i in 1:nst) {
    intp_data[[i]] = www$bbbw[[i]]$predictT(www$wtime[i, ])$pred
}
mads = bootstrap(kkk, y_no, ktype, K, ode_par, intp_data,www)
ode_par
mads


## End(Not run)