Profile: Profile Estimation Functions
In CollocInfer: Collocation Inference for Dynamic Systems

Profiling Routines

R Documentation

Profile Estimation Functions

Description

These functions are wrappers that create lik and proc functions and run generalized profiling.

Usage

Profile.LS(fn,data,times,pars,coefs=NULL,basisvals=NULL,lambda,
                        fd.obj=NULL,more=NULL,weights=NULL,quadrature=NULL,
                        likfn = make.id(), likmore = NULL,
                        in.meth='nlminb',out.meth='nls',
                        control.in,control.out,eps=1e-6,active=NULL,posproc=FALSE,
                        poslik=FALSE,discrete=FALSE,names=NULL,sparse=FALSE)

Profile.multinorm(fn,data,times,pars,coefs=NULL,basisvals=NULL,var=c(1,0.01),
                        fd.obj=NULL,more=NULL,quadrature=NULL,
                        in.meth='nlminb',out.meth='optim',
                        control.in,control.out,eps=1e-6,active=NULL,
                        posproc=FALSE,poslik=FALSE,discrete=FALSE,names=NULL,sparse=FALSE)

Arguments

`fn`	A function giving the right hand side of a differential/difference equation. The function should have arguments times The times at which the RHS is being evaluated. x The state values at those times. p Parameters to be entered in the system. more An object containing additional inputs to `fn` It should return a matrix of the same dimension of `x` giving the right hand side values. If `fn` is given as a single function, its derivatives are estimated by finite-differencing with stepsize `eps`. Alternatively, a list can be supplied with elements: fn Function to calculate the right hand side should accept a matrix of state values at . dfdx Function to calculate the derivative with respect to `x` dfdp Function to calculate the derivative with respect to `p` d2fdx2 Function to calculate the second derivative with respect to `x` d2fdxdp Function to calculate the second derivative with respect to `x` and `p` These functions take the same arguments as `fn` and should output multidimensional arrays with the dimensions ordered according to time, state, deriv1, deriv2; here derivatives with respect to `x` always precede derivatives with respect to `p`.
`data`	Matrix of observed data values.
`times`	Vector observation times for the data.
`pars`	Initial values of parameters to be estimated processes.
`coefs`	Vector giving the current estimate of the coefficients in the spline.
`basisvals`	Values of the collocation basis to be used. This can either be a basis object from the `fda` package, or a list elements: bvals.obs A matrix giving the values of the basis at the observation times bvals A matrix giving the values of the basis at the quadrature times dbvals A matrix giving the derivative of the basis at the quadrature times
`lambda`	(`Profile.LS` only) Penalty value trading off fidelity to data with fidelity to differential equations.
`var`	(`profile.Cproc` or `profile.Dproc`) A vector of length 2, giving
`fd.obj`	(Optional) A functional data object; if this is non-null, `coefs` and `basisvals` is extracted from here.
`more`	An object specifying additional arguments to `fn`.
`weights`	(`Profile.LS` only)
`quadrature`	Quadrature points, should contain two elements (if not NULL) qpts Quadrature points; defaults to midpoints between knots qwts Quadrature weights; defaults to normalizing by the length of `qpts`.
`in.meth`	Inner optimization function to be used, currently one of 'nlminb', 'MaxNR', 'optim' or 'house'. The last calls `SplineEst.NewtRaph`. This is fast but has poor convergence.
`out.meth`	Outer optimization function to be used, depending on the type of method `Profile.LS` One of 'nls' or 'ProfileGN'; the latter calls `Profile.GausNewt` which is fast but may have poor convergence. `Profile.multinorm` One of 'optim' (defaults to BFGS routine in `optim` unless `control.out$meth` specifies otherwise), 'nlminb', or 'maxNR'.
`control.in`	Control object for inner optimization function.
`control.out`	Control object for outer optimization function.
`eps`	Finite differencing step size, if needed.
`active`	Incides indicating which parameters of `pars` should be estimated; defaults to all of them.
`posproc`	Should the state vector be constrained to be positive? If this is the case, the state is represented by an exponentiated basis expansion in the `proc` object.
`poslik`	Should the state be exponentiated before being compared to the data? When the state is represented on the log scale (`posproc=TRUE`), this is an alternative to taking the log of the data.
`discrete`	Is this a discrete-time or a continuous-time system? If discrete, the derivative is instead taken to be the value at the next time point.
`names`	The names of the state variables if not given by the column names of `coefs`.
`sparse`	Should sparse matrices be used for basis values? This option can save memory when `ProfileGausNewt` and `SplineEstNewtRaph` are called. Otherwise sparse matrices will be converted to full matrices and this can slow the code down.
`likfn`	Defines a map from the trajectory to the observations. This should be in the same form as `fn`. If a function is given, derivatives are estimated by finite differencing, otherwise a list is expected to provide the same derivatives as `fn`. If `poslik=TRUE`, the states are exponentiated before the `likfn` is evaluated and the derivatives are updated to account for this. Defaults to the identity transform.
`likmore`	A list containing additional inputs to `likfn` if needed, otherwise set to `NULL`

Details

These functional all carry out the profiled optimization method of Ramsay et al 2007. Profile.LS uses a sum of squared errors criteria for both fit to data and the fit of the derivatives to a differential equation. Profile.multinorm uses multivariate normal approximations. discrete changes the system to a discrete-time difference equation with the right hand side function giving the transition function.

Note that these all call outeropt, which creates temporary files 'curcoefs.tmp' and 'optcoefs.tmp' to update coefficients as pars evolves. These overwrite existing files of those names and are deleted before the function terminates.

Value

A list with elements

`pars`	Optimized parameters
`coefs`	Optimized coefficients at `pars`
`lik`	The `lik` object generated
`proc`	The `proc` item generated
`data`	The data used in doing the fitting.
`times`	The vector of times at which the observations were made

Examples


###############################
####   Data             #######
###############################

data(FhNdata)

###############################
####  Basis Object      #######
###############################

knots = seq(0,20,0.2)
norder = 3
nbasis = length(knots) + norder - 2
range = c(0,20)

bbasis = create.bspline.basis(range=range(FhNtimes),nbasis=nbasis,
	norder=norder,breaks=knots)


#### Start from pre-estimated values to speed up optimization

data(FhNest)

spars = FhNestPars
coefs = FhNestCoefs

lambda = 10000

res1 = Profile.LS(make.fhn(),data=FhNdata,times=FhNtimes,pars=spars,coefs=coefs,
  basisvals=bbasis,lambda=lambda,in.meth='nlminb',out.meth='nls')

Covar = Profile.covariance(pars=res1$pars,times=FhNtimes,data=FhNdata,
  coefs=res1$coefs,lik=res1$lik,proc=res1$proc) 


## Not run: 
## Alternative, starting from perturbed coefficients -- takes too long for 
# automatic checks in CRAN

# Initial values for coefficients will be obtained by smoothing

DEfd = smooth.basis(FhNtimes,FhNdata,fdPar(bbasis,1,0.5))   # Smooth to estimate
                                                            # coefficients first
coefs = DEfd$fd$coefs
colnames(coefs) = FhNvarnames


###############################
####     Optimization       ###
###############################

spars = c(0.25,0.15,2.5)          # Perturbed parameters
names(spars)=FhNparnames
lambda = 10000

res1 = Profile.LS(make.fhn(),data=FhNdata,times=FhNtimes,pars=spars,coefs=coefs,
  basisvals=bbasis,lambda=lambda,in.meth='nlminb',out.meth='nls')

par(mfrow=c(2,1))
plotfit.fd(FhNdata,FhNtimes,fd(res1$coefs,bbasis))

## End(Not run)  
  
 
 
  
## Not run: 
####################################################
###  An Explicitly Multivariate Normal Formation ### 
####################################################

var = c(1,0.0001)

res2 = Profile.multinorm(make.fhn(),FhNdata,FhNtimes,pars=res1$pars,
          res1$coefs,bbasis,var=var,out.meth='nlminb', in.meth='nlminb')

## End(Not run)

CollocInfer documentation built on April 4, 2025, 2:21 a.m.