Main function of this package : computes the Laplace deconvolution with noisy discrete non-equally spaced observations on a finite time interval (see reference).
numeric vector, the observed noisy observations of the Laplace convolution
g
numeric vector, the known kernel of the Laplace convolution
times
numeric vector, the observation times (default 1:length(Y))
sigma
numeric, the noise level
cpen
numeric, the penalization constant (default 2)
atab
numeric vector, an array of value for a (default -1, see Details)
Mmax
integer, the maximum degree (default 25)
ncores.max
numeric, number of used cores (default Inf)
verbose
boolean to control information output (default FALSE)
withplot
boolean to control plot output (default FALSE)
Details
atab defines the values of the scale parameter used to build the Laguerre functions basis. If atab is of length 1 and negative then atab<-seq(0.4,3,by=0.05/abs(atab))/sqrt(Tmax/10) where Tmax=max(times).
Value
a list containing:
f.hat, numeric vector, the estimate at observation times
q.hat, numeric vector, the reconstructed convolution of the kernel g and the estimate f.hat, computed at observation times
f.coef, numeric vector, the coefficients of f.hat in the selected Laguerre function basis
info, only for internal use — not documented
a.hat, numeric, the selected value of the parameter a
M0, numeric, the dimension of the selected model
sigma, numeric, the noise level
cpen, numeric, the penalization constant used in the penalty
Author(s)
Y. Rozenholc and M. Pensky
References
Laplace deconvolution on the basis of time domain data and its application to Dynamic Contrast Enhanced imaging by F. Comte, C-A. Cuenod, M. Pensky, Y. Rozenholc (ArXiv http://arxiv.org/abs/1405.7107)
Examples
## Not run:
#### AN ARTICIAL EXAMPLE ####
library(LaplaceDeconv)
par(mfrow=c(1,1))
set.seed(29102015)
sigma=0.02
a = 1
t = seq(0,5,l=100)
g = 20*t^2*exp(-5*t)
f.coef = c(0.4,0.02,0.01)
# compute the Laplace convolution from g, kernel computed at times t, and the function
# described by its decomposition in Laguerre function basis with scale a :
fg = LaguerreLaplaceConvolution(t,g,f.coef,a)
# the noisy observations :
Y = fg+sigma*rnorm(length(fg))
# estimation of f from the observation and the kernel :
L = LagLaplDeconv(Y,g,t,sigma)
matplot(t,cbind(g,MakeLaguerreMatrix(a,3)(t)%*%f.coef,fg,L$q.hat,L$f.hat,Y),lty=1,
type=c('b',rep('l',4),'p'),ylab='',pch='x')
# display results of estimation
legend('topright',lty=c(rep(1,5),0),pch=c('x',rep('',4),'x'),
legend=c(
'g: partially observed kernel',
'f: unknown',
'q=fxg: unknown convolution',
expression(hat(q)*': plug-in convolution'),
expression(hat(f)*': estimation of f'),
'Y: observations'),
col=1:6)
## End(Not run)
## Not run:
#### A REAL EXAMPLE USING DCE-MRI DATA FROM A TUMOR ####
library(LaplaceDeconv)
par(mfrow=c(1,2))
# load data from patient before the treatment
data(EX_DCEMRI_t0)
# display AIF and tumoral enhancements
matplot(ex_dcemri$times,
cbind(ex_dcemri$AIF,ex_dcemri$TUM_1,ex_dcemri$TUM_2,ex_dcemri$TUM_3),
ylab='',lty=1,type=c('b',rep('p',3)),pch='+',main='Observations')
legend('topright',pch='+',legend=c('AIF','TUM_1','TUM_2','TUM_3'),col=1:4)
# estimation of the contrast agent survival functions
L1 = LagLaplDeconv(ex_dcemri$TUM_1,ex_dcemri$AIF,ex_dcemri$times,ex_dcemri$sigma)
L2 = LagLaplDeconv(ex_dcemri$TUM_2,ex_dcemri$AIF,ex_dcemri$times,ex_dcemri$sigma)
L3 = LagLaplDeconv(ex_dcemri$TUM_3,ex_dcemri$AIF,ex_dcemri$times,ex_dcemri$sigma)
matlines(ex_dcemri$times,cbind(L1$q.hat,L2$q.hat,L3$q.hat),type='l',lty=1,col=2:4)
# display results of estimation
matplot(ex_dcemri$times,cbind(L1$f.hat,L2$f.hat,L3$f.hat),type='l',lty=1,col=2:4,
ylab='survival',main='Contrast agent survival fcts')
legend('topright',lty=1,col=2:4,
legend=c(
paste0('TUM_1 - a.hat=',round(L1$a.hat,digits=2)),
paste0('TUM_2 - a.hat=',round(L2$a.hat,digits=2)),
paste0('TUM_3 - a.hat=',round(L3$a.hat,digits=2))
)
)
## End(Not run)