.GRfiltreScan<-function(csc,idx=rep(1,length(csc)),alpha=5:10,perc=0.7){
cmpmiss<-function(i) sum(i)/length(i)
alpha=alpha[alpha<=min(max(alpha),length(csc))]
isFeat=rep(FALSE,length(csc))
for(add in alpha){
z <- embed(csc,add)
lrep=rowSums(z)/add
l=which(lrep >=perc)
l=l[which(idx[l]==idx[l+add-1])]
for(i in 0:(add-1)) isFeat[l+i]=TRUE
}
return(isFeat)
}
### collapse the eic attr
.GRgetEicInfos<-function(eicParams,eicloc=NULL){
if(!is.null(eicloc)) if(is.data.frame(eicloc)){
fileeics=paste(ifelse(is.null(eicParams$dirEic),"./",eicParams$dirEic),
unique(eicloc$GrpEic),
ifelse(is.null(eicParams$addEic),"./",eicParams$addEic),".rda",sep="")
names(fileeics)=unique(eicloc$GrpEic)
}
if(is.null(eicloc))
fileeics=list.files(ifelse(is.null(eicParams$dirEic),"./",eicParams$dirEic),pattern=paste(eicParams$addEic,".rda",sep=""),full.names=TRUE)
fileeics=fileeics[file.exists(fileeics)]
if(length(fileeics)==0) stop('No file found!')
allr=list()
for(k in 1:length(fileeics)){
cat(ifelse(k%%10==0,"X","."))
load(fileeics[k])
allr[[k]]=attr(dfeic,"eic")
}
newtabeic=do.call("rbind",allr)
# if(is.null(eicloc)){
grp=factor(newtabeic$GrpEic,names(sort(tapply(newtabeic$mzap,newtabeic$GrpEic,mean,na.rm=T))))
newtabeic=newtabeic[order(grp,newtabeic$mzap),]
# }
newtabeic
}
.GRestimBl<-function(tabeic,blSid,eicParams,ncl=1){
.infctBlfct<-function(ifile,blSid=NULL,minScan=11,quantBl=0.5){
load(ifile)
dfeic=dfeic[dfeic$samp%in%blSid,]
leics=attr(dfeic,"eic")[,"Id"]
matbl=matrix(NA,nrow=length(leics),ncol=length(blSid),dimnames=list(leics,blSid))
if("ineic" %in% names(dfeic)) dfeic=dfeic[which(dfeic$ineic==1),]
if(nrow(dfeic)==0) return(matbl)
eicbl=paste(dfeic$samp,dfeic$eic)
l2use=names(which(table(eicbl)>=minScan))
if(length(l2use)==0) return(matbl)
dfeic=dfeic[which(eicbl%in%l2use & dfeic$y>0 & !is.na(dfeic$y)),]
if(nrow(dfeic)==0) return(matbl)
blsamp=tapply(1:nrow(dfeic),dfeic$samp,function(x){
x=tapply(dfeic$y[x],dfeic$eic[x],quantile,quantBl)
x=x[match(rownames(matbl),names(x))]
names(x)=rownames(matbl)
x
})
for(i in names(blsamp)) matbl[names(blsamp[[i]]),i]=blsamp[[i]]
# print("OKKK")
rm(list=c("dfeic"))
return(matbl)
}
fileeics=paste(ifelse(is.null(eicParams$dirEic),"./",eicParams$dirEic),
unique(tabeic$GrpEic),
ifelse(is.null(eicParams$addEic),"./",eicParams$addEic),".rda",sep="")
names(fileeics)=unique(tabeic$GrpEic)
fileeics=fileeics[file.exists(fileeics)]
print(fileeics)
if(length(fileeics)==0) stop('No file found!')
d0=proc.time()[3]
cat("Started at ",date(),sep="")
if(ncl==1 | length(fileeics)==1){
cat(" on 1 processor\n",sep="")
allr=list()
for(k in 1:length(fileeics)){
cx=fileeics[k]
cat(cx," ")
allr[[k]]=.infctBlfct(cx,blSid,eicParams$minScanBl,eicParams$quantBl)
}
}
if(ncl>1 & length(fileeics)>1){
require("snowfall")
ncl=max(1,min(ncl,parallel:::detectCores()))
cat(" on ",ncl," processors\n",sep="")
sfInit(parallel=TRUE, cpus=ncl, type='SOCK')
sfExport( ".infctBlfct", local=TRUE )
allr=sfClusterApplyLB(fileeics,.infctBlfct,blSid,eicParams$minScanBl,eicParams$quantBl)
sfStop()
}
d1=proc.time()[3]
cat("\nCompleted at ",date()," -> took ",round(d1-d0,1)," secs \n",sep="")
allr=allr[which(!sapply(allr,is.null))]
cat("\n")
if(length(allr)==0) return(NULL)
blres=do.call("rbind",allr)
blres=blres[rownames(blres)%in%tabeic$Id,,drop=F]
blres=blres[match(tabeic$Id,rownames(blres)),,drop=F]
rownames(blres)=tabeic$Id
return(blres)
}
.GRcompEICst<-function(tabeic,eicParams,lSids,ncl=1){
.incompEICst<-function(cfi,blvect=NULL,eicParams){
load(cfi)
leics=attr(dfeic,"eic")$Id
if(is.null(blvect)) blvect=eicParams$Mint
if(length(blvect==1)){
iblvect=rep(blvect,length(leics))
names(iblvect)=leics
} else iblvect=blvect[names(blvect)%in%leics]
if(any(is.na(iblvect))) iblvect[is.na(iblvect)]=eicParams$Mint
dfeic$csc=as.vector(dfeic$y>=iblvect[dfeic$eic])
Ncons=Coda2=HMax=list()
lusamp=unique(dfeic$samp)
if("ineic" %in% names(dfeic)) unique(dfeic$samp[which(dfeic$ineic==1)])
for(isamp in lusamp){
idfeic=dfeic[which(dfeic$samp==isamp),]
if("ineic" %in% names(idfeic)) idfeic=idfeic[which(idfeic$ineic==1),]
re=.GRfiltreScan(idfeic$csc,idfeic$eic,alpha=eicParams$LSc,perc=eicParams$Perc)
Ncons[[isamp]]=tapply(re,idfeic$eic,sum)[leics]
Coda2[[isamp]]=tapply(idfeic$y,idfeic$eic,.GRcodadw2)[leics]
HMax[[isamp]]=tapply(idfeic$y,idfeic$eic,function(x) floor(max(x,na.rm=T)))[leics]
}
Ncons=do.call("rbind",Ncons)
Coda2=do.call("rbind",Coda2)
HMax=do.call("rbind",HMax)
colnames(HMax)=colnames(Coda2)=colnames(Ncons)=leics
return(list(Ncons=Ncons,Coda2=Coda2,HMax=HMax))
}
fileeics=paste(ifelse(is.null(eicParams$dirEic),"./",eicParams$dirEic),
unique(tabeic$GrpEic),
ifelse(is.null(eicParams$addEic),"./",eicParams$addEic),".rda",sep="")
names(fileeics)=unique(tabeic$GrpEic)
fileeics=fileeics[file.exists(fileeics)]
if(length(fileeics)==0) stop('No file found!')
if("Bl"%in%names(tabeic)) blvect=tabeic$Bl else blvect=rep(eicParams$Mint,nrow(tabeic))
names(blvect)=tabeic$Id
d0=proc.time()[3]
cat("Started at ",date(),sep="")
if(ncl==1 | length(fileeics)==1){
cat(" on 1 processor\n",sep="")
allr=list()
for(k in 1:length(fileeics)){
cx=fileeics[k]
cat(ifelse(k%%10==0,"X","."))
allr[[k]]=.incompEICst(cx,blvect,eicParams)
}
}
if(ncl>1 & length(fileeics)>1){
require("snowfall")
ncl=max(1,min(ncl,parallel:::detectCores()))
cat(" on ",ncl," processors\n",sep="")
sfInit(parallel=TRUE, cpus=ncl, type='SOCK')
sfExport(".incompEICst", local=TRUE )
sfLibrary(GRMeta)
# sfExport( ".GRfiltreScan", local=TRUE )
# sfExport( ".GRcodadw2", local=TRUE )
allr=sfClusterApplyLB(fileeics,.incompEICst,blvect,eicParams)
sfStop()
}
d1=proc.time()[3]
cat("\nCompleted at ",date()," -> took ",round(d1-d0,1)," secs \n",sep="")
allr=allr[which(!sapply(allr,is.null))]
cat("\n")
if(length(allr)==0) return(NULL)
if(missing(lSids)) lSids=sort(unique(unlist(lapply(allr,function(x) rownames(x[[1]])))))
allr=lapply(allr,function(x) lapply(x,function(y) y[match(lSids,rownames(y)),,drop=F]))
lmats=unique(unlist(lapply(allr,names)))
allmat=list()
for(imat in lmats){
tmpm=lapply(allr,function(x) x[[imat]])
tmpm=tmpm[which(!sapply(tmpm,is.null))]
tmpm=do.call("cbind",tmpm)
tmpm=tmpm[,match(tabeic$Id,colnames(tmpm)),drop=F]
dimnames(tmpm)=list(lSids,tabeic$Id)
allmat[[imat]]=tmpm
}
rm(list="allr")
return(allmat)
}
Add the following code to your website.
For more information on customizing the embed code, read Embedding Snippets.