view export_to_fastq/export_to_fastq_past.R @ 0:97792524cc9c default tip

Migrated tool version 0.1 from old tool shed archive to new tool shed repository
author louise
date Tue, 07 Jun 2011 17:21:49 -0400
parents
children
line wrap: on
line source

my.args <- commandArgs(trailingOnly = TRUE)
# ----- Check input and output directories -----
if(!file.exists(my.args[1])){ 
   stop("The provided project directory does not exist!")
}
inputDirectory=my.args[1] #Directory where input data are. (e.g export or fastq ...)

output_file = my.args[3]

threshold=as.numeric(my.args[2]) #threshold for nFilter

#print(my.args)

library(ShortRead)
#source('/home/galaxy/galaxy_dev/tools/EMBL_tools/HTS_helper_src_for_export_to_fastq.R')

####Solution temporaire, chastityFilter sera dans le RNASeq package 
chastityFilter <- function(.name="Illumina Chastity Filter") 
{
  srFilter(function(x){
    if(any(rownames(varMetadata(alignData(x))) == "filtering")){
      keep<-alignData(x)$filtering=="Y"
    } else {
      warning(paste("The '",.name,"' filter is only valid for Illumina reads.",sep=""))
      keep<-rep(TRUE,length(x))
    }
    return(keep)
  },name=.name)
}

"summarize.by.transcripts" <- function(sample,annotation){
  
  transcripts <- do.call(rbind,lapply(names(sample),function(chr,sample,annotation){
    counts<-stats:::aggregate(sample[[chr]],list(transcript=annotation[chr]$transcript),sum)
  },sample,annotation))

  colnames(transcripts)[2] <- "counts"
  
  return(transcripts)
}
###

#----- FILTER ----
filter<- compose(chastityFilter(),nFilter(threshold=threshold))

#----- ALIGN ----
# call the readAligned function with this filter
aln<-readAligned(inputDirectory, type='SolexaExport',filter=filter, withAll=TRUE )
writeFastq(aln,file=output_file,mode='a')