R/denoise.R

In adlift: An Adaptive Lifting Scheme Algorithm

"denoise" <-
function(x,f,pred=LinearPred,neigh=1,int=TRUE,clo=FALSE,keep=2,rule="median",returnall=FALSE){

#denoises input using prediction scheme specified by other arguments
#artificial level thresholding, values generated by ebayesthresh using normalized detail coefficients

newcoeff<-NULL
ndetlist<-list()
tclist<-NULL


out<-fwtnp(x,f,LocalPred=pred,neighbours=neigh,intercept=int,closest=clo,nkeep=keep,varonly=TRUE)
indsd<-sqrt(out$v)
norcoeff<-out$coeff/indsd

lr<-out$lengthsremove    #vector deciding how to divide up coefficients into artificial levels
rem<-out$removelist      #used to convert output to original lr,rem)

al<-artlev(lr,rem)      #the list of indices of removelist separated into levels

levno<-length(al)
for (i in 1:levno){
	ndetlist[[i]]<-norcoeff[al[[i]]]
}

sd<-mad(ndetlist[[1]])    #uses the first (largest) level to estimate noise standard deviation
sd1<-mad(ndetlist[[1]],center=0)
sd2<-mad(norcoeff[rem]) #uses all the normalised details to estimate the noise std dev

for (i in 1:levno){
	tclist<-ebayesthresh(ndetlist[[i]],prior="cauchy",a=NA,sdev=sd,threshrule=rule)
	newcoeff[al[[i]]]<-tclist*indsd[al[[i]]]
}

newcoeff[out$pointsin]<-out$coeff[out$pointsin]

fhat<-invtnp(x,newcoeff,out$lengths,lr,out$pointsin,rem,out$neighbrs,out$schemehist,out$interhist,length(x)-keep,int,neigh,clo,pred)

if(returnall){
	return(list(fhat=fhat,w=out$W,indsd=indsd,al=al,sd=sd))
}
else{
	return(fhat$coeff)
}

}