source: trunk/grails-app/services/nl/tno/metagenomics/FastaService.groovy @ 3

Last change on this file since 3 was 3, checked in by robert@…, 9 years ago

Externalized configuration; improved assay view (detail views of runs and samples); implemented uploading and parsing of FASTA and QUAL files

File size: 13.4 KB
Line 
1package nl.tno.metagenomics
2
3import java.io.File;
4import java.util.ArrayList;
5import org.codehaus.groovy.grails.commons.ConfigurationHolder
6
7class FastaService {
8        def fileService
9        def fuzzySearchService
10
11        static transactional = true
12
13        /**
14         * Parses uploaded files and checks them for FASTA and QUAL files
15         * @param filenames             List of filenames currently existing in the upload directory
16         * @param onProgress    Closure to execute when progress indicators should be updated.
17         *                                              Has 2 parameters: numFilesProcessed and numBytesProcessed that indicate the number
18         *                                              of files and bytes that have been processed in total
19         * @param directory             Directory to move the files to
20         * @return                              Structure with information about the parsed files. The 'success' files are
21         *                                              moved to the given directory
22         *
23         * [
24         *              success: [
25         *                      [filename: 'abc.fasta', type: FASTA, numSequences: 190]
26         *                      [filename: 'cde.fasta', type: FASTA, numSequences: 140]
27         *                      [filename: 'abc.qual', type: QUAL, numSequences: 190, avgQuality: 38]
28         *                      [filename: 'cde.qual', type: QUAL, numSequences: 140, avgQuality: 29]
29         *              ],
30         *              failure: [
31         *                      [filename: 'testing.xls', type: 'unknown', message: 'Type not recognized']
32         *              ]
33         * ]
34         *
35         */
36        def parseFiles( ArrayList filenames, Closure onProgress, File directory = null ) {
37                if( filenames.size() == 0 ) {
38                        return [ success: [], failure: [] ];
39                }
40
41                if( !directory ) {
42                        directory = fileService.absolutePath( ConfigurationHolder.config.metagenomics.fileDir )
43                }
44
45                def success = [];
46                def failure = [];
47
48                long filesProcessed = 0;
49                long bytesProcessed = 0;
50
51                // Loop through all filenames
52                filenames.each { filename ->
53                        // TODO: Handle zipped files
54                        def file = fileService.get( filename );
55                        String filetype = determineFileType( file );
56
57                        if( !fileTypeValid( filetype ) ) {
58                                fileService.delete(filename);
59                                failure << [ filename: filename, originalfilename: fileService.originalFilename( filename ), type: filetype, message: 'File type not accepted' ];
60                        } else {
61                                try {
62                                        def contents = parseFile( file, filetype, { files, bytes ->
63                                                filesProcessed += files;
64                                                bytesProcessed += bytes;
65
66                                                onProgress( filesProcessed, bytesProcessed );
67                                        } );
68
69                                        contents.filename = file.getName();
70                                        contents.originalfilename = fileService.originalFilename( contents.filename )
71
72                                        if( contents.success ) {
73                                                success << contents;
74                                        } else {
75                                                fileService.delete(filename);
76                                                failure << contents;
77                                        }
78                                } catch( Exception e ) {
79                                        // If anything fails during parsing, return an error message
80                                        fileService.delete(filename);
81                                        failure << [ filename: filename, originalfilename: fileService.originalFilename( filename ), type: filetype, message: e.getMessage() ];
82                                }
83
84                        }
85                }
86
87                return [ success: success, failure: failure ];
88        }
89
90        /**
91         * Matches uploaded fasta and qual files and combines them with the samples they probably belong to
92         * @param parsedFiles   Parsed files
93         * [
94         *              [filename: 'abc.fasta', type: FASTA, numSequences: 190]
95         *              [filename: 'cde.fasta', type: FASTA, numSequences: 140]
96         *              [filename: 'abc.qual', type: QUAL, numSequences: 190, avgQuality: 38]
97         *              [filename: 'cde.qual', type: QUAL, numSequences: 140, avgQuality: 29]
98         * ]
99         * @param samples               AssaySample objects to which the files should be matched.
100         * @return                              Structure with information about the matching.
101         * [
102         *              [
103         *                      fasta:  [filename: 'abc.fasta', type: FASTA, numSequences: 190],
104         *                      qual:   [filename: 'abc.qual', type: QUAL, numSequences: 190, avgQuality: 38],
105         *                      feasibleQuals: [
106         *                              [filename: 'abc.qual', type: QUAL, numSequences: 190, avgQuality: 38],
107         *                              [filename: 'def.qual', type: QUAL, numSequences: 190, avgQuality: 21]
108         *                      ]
109         *                      sample: AssaySample object     
110         */
111        def matchFiles( def parsedFiles, def samples ) {
112                def fastas = parsedFiles.findAll { it.type == "fasta" }
113                def quals = parsedFiles.findAll { it.type == "qual" }
114                samples = samples.toList()
115
116                def files = [];
117
118                fastas.each { fastaFile ->
119                        // Remove extension
120                        def matchWith = fastaFile.originalfilename.substring( 0, fastaFile.originalfilename.lastIndexOf( '.' ) )
121
122                        // Determine feasible quals (based on number of sequences )
123                        def feasibleQuals = quals.findAll { it.numSequences == fastaFile.numSequences }
124
125                        // Best matching qual file
126                        def qualIdx = fuzzySearchService.mostSimilarWithIndex( matchWith + '.qual', feasibleQuals.originalfilename );
127
128                        def qual = null
129                        if( qualIdx != null )
130                                qual = feasibleQuals[ qualIdx ];
131
132                        // Best matching sample
133                        def sampleIdx = fuzzySearchService.mostSimilarWithIndex( matchWith, samples.sample.name );
134                        def assaySample = null
135                        if( sampleIdx != null ) {
136                                assaySample = samples[ sampleIdx ];
137                        }
138
139                        files << [
140                                                fasta: fastaFile,
141                                                feasibleQuals: feasibleQuals,
142                                                qual: qual,
143                                                assaySample: assaySample
144                                        ]
145                }
146
147                return files;
148
149
150        }
151
152
153        /**
154         * Determines the file type of a given file, based on the extension.
155         * @param file  File to assess
156         * @return
157         */
158        protected String determineFileType( File file ) {
159                if( !file )
160                        return null
161
162                // Determine the extension of the file
163                String name = file.getName();
164                if( name.lastIndexOf( '.' ) == -1 ) {
165                        // Without an extension the file type is unknown
166                        return ""
167                }
168
169                String extension = name.substring( name.lastIndexOf( '.' ) + 1 ).toLowerCase();
170
171                switch( extension ) {
172                        case "fasta":
173                        case "fna":
174                                return "fasta";
175                        case "qual":
176                        case "fqa":
177                                return "qual";
178                        case "xls":
179                                return "excel";
180                        case "zip":
181                                return "zip";
182                        case "gz":
183                                return "gzip";
184                        default:
185                                return "";              // meaning 'unknown'
186                }
187        }
188
189        /**
190         * Determines whether a file can be processed.
191         * @param filetype      Filetype of the file
192         * @see determineFileType()
193         * @return
194         */
195        protected boolean fileTypeValid( String filetype ) {
196                switch( filetype ) {
197                        case "fasta":
198                        case "qual":
199                                return true;
200                        default:
201                                return false;
202                }
203        }
204
205        /**
206         * Parses the given file
207         * @param file                  File to parse
208         * @param filetype              Type of the given file
209         * @param onProgress    Closure to execute when progress indicators should be updated.
210         *                                              Has 2 parameters: numFilesProcessed and numBytesProcessed that indicate the number
211         *                                              of files and bytes that have been processed in this file (so the first parameter should
212         *                                              only be 1 when the file is finished)
213         *
214         * @return                              List structure. Examples:
215         *
216         *   [ success: true, filename: 'abc.fasta', type: 'fasta', numSequences: 200 ]
217         *   [ success: true, filename: 'abc.qual', type: 'qual', numSequences: 200, avgQuality: 36 ]
218         *   [ success: false, filename: 'abc.txt', type: 'txt', message: 'Filetype could not be parsed.' ]
219         */
220        protected def parseFile( File file, String filetype, Closure onProgress ) {
221                switch( filetype ) {
222                        case "fasta":
223                                return parseFasta( file, onProgress );
224                        case "qual":
225                                return parseQual( file, onProgress );
226                        default:
227                                onProgress( 1, file.length() );
228                                return [ success: false, type: filetype, message: 'Filetype could not be parsed.' ]
229                }
230        }
231
232        /**
233         * Parses the given FASTA file
234         * @param file                  File to parse
235         * @param onProgress    Closure to execute when progress indicators should be updated.
236         *                                              Has 2 parameters: numFilesProcessed and numBytesProcessed that indicate the number
237         *                                              of files and bytes that have been processed in this file (so the first parameter should
238         *                                              only be 1 when the file is finished)
239         * @return                              List structure. Examples:
240         *
241         *   [ success: true, filename: 'abc.fasta', type: 'fasta', numSequences: 200 ]
242         *   [ success: false, filename: 'def.fasta', type: 'fasta', message: 'File is not a valid FASTA file' ]
243         */
244        protected def parseFasta( File file, Closure onProgress ) {
245
246                long startTime = System.nanoTime();
247                log.trace "Start parsing FASTA " + file.getName()
248
249                // Count the number of lines, starting with '>' (and where the following line contains a character other than '>')
250                long numSequences = 0;
251                long bytesProcessed = 0;
252                boolean lookingForSequence = false;
253
254                file.eachLine { line ->
255                        if( line ) {
256                                if( !lookingForSequence && line[0] == '>' ) {
257                                        lookingForSequence = true;
258                                } else if( lookingForSequence ) {
259                                        if( line[0] != '>' ) {
260                                                numSequences++;
261                                                lookingForSequence = false;
262                                        }
263                                }
264
265
266                                // Update progress every time 1MB is processed
267                                bytesProcessed += line.size();
268                                if( bytesProcessed > 1000000 ) {
269                                        onProgress( 0, bytesProcessed );
270                                        bytesProcessed = 0;
271                                }
272                        }
273
274
275
276                }
277
278                // Update progress and say we're finished
279                onProgress( 1, bytesProcessed );
280
281                log.trace "Finished parsing FASTA " + file.getName() + ": " + ( System.nanoTime() - startTime ) / 1000000L
282
283                return [ success: true, type: "fasta", numSequences: numSequences ];
284        }
285
286        /**
287         * Parses the given QUAL file
288         * @param file                  File to parse
289         * @param onProgress    Closure to execute when progress indicators should be updated. 
290         *                                              Has 2 parameters: numFilesProcessed and numBytesProcessed that indicate the number
291         *                                              of files and bytes that have been processed in this file (so the first parameter should
292         *                                              only be 1 when the file is finished)
293         * @return                              List structure. Examples:
294         *
295         *   [ success: true, filename: 'abc.qual', type: 'qual', numSequences: 200, avgQuality: 31 ]
296         *   [ success: false, filename: 'def.qual', type: 'qual', message: 'File is not a valid QUAL file' ]
297         */
298        protected def parseQual( File file, Closure onProgress ) {
299                long startTime = System.nanoTime();
300                log.trace "Start parsing QUAL " + file.getName()
301
302                // Count the number of lines, starting with '>'. After we've found such a character, we continue looking for
303                // quality scores
304                long numSequences = 0;
305                long bytesProcessed = 0;
306                def quality = [ quality: 0.0, number: 0L ]
307
308                boolean lookingForFirstQualityScores = false;
309                file.eachLine { line ->
310                        if( line ) {
311                                if( !lookingForFirstQualityScores && line[0] == '>' ) {
312                                        lookingForFirstQualityScores = true;
313                                } else if( lookingForFirstQualityScores ) {
314                                        if( line[0] != '>' ) {
315                                                numSequences++;
316                                                lookingForFirstQualityScores = false;
317
318                                                quality = updateQuality( quality, line );
319                                        }
320                                } else {
321                                        quality = updateQuality( quality, line );
322                                }
323                               
324                                // Update progress every time 1MB is processed
325                                bytesProcessed += line.size();
326                                if( bytesProcessed > 1000000 ) {
327                                        onProgress( 0, bytesProcessed );
328                                        bytesProcessed = 0;
329                                }
330                        }
331
332                }
333
334                // Update progress and say we're finished
335                onProgress( 1, bytesProcessed );
336
337                log.trace "Finished parsing QUAL " + file.getName() + ": " + ( System.nanoTime() - startTime ) / 1000000L
338
339                return [ success: true, type: "qual", numSequences: numSequences, avgQuality: quality.quality ];
340        }
341
342        /**
343         * Parses the given line and updates the average quality based on the scores
344         * @param quality       [quality: 0.0f, number: 0L]
345         * @param line          String of integer quality scores, separated by a whitespace
346         * @return                      [quality: 0.0f, number: 0L]
347         */
348        protected def updateQuality( def quality, String line ) {
349                // Determine current average
350                List tokens = line.tokenize();
351                Long total = 0;
352
353                tokens.each {
354                        total += Integer.parseInt( it );
355                }
356
357                int numTokens = tokens.size();
358
359                // Update the given average
360                if( numTokens > 0 ) {
361                        quality.number += numTokens;
362                        quality.quality = quality.quality + ( ( total / numTokens as double ) - quality.quality ) / quality.number * numTokens;
363                }
364
365                return quality
366        }
367
368        /**
369         * Moves a fasta and qual file to their permanent location, and returns information about these files
370         *
371         * @param fastaFile                     Filename of the fasta file in the temporary directory
372         * @param qualFile                      Filename of the fasta file in the temporary directory
373         * @param processedFiles        Structure with data about uploaded files and matches
374         * @return      [ fasta: <fasta filename>, qual: <qual filename>, numSequences: <number of sequences>, avgQuality: <average quality> ]
375         */
376        public def savePermanent( String fastaFile, String qualFile, def processedFiles ) {
377                File permanentDirectory = fileService.absolutePath( ConfigurationHolder.config.metagenomics.fileDir );
378                def returnStructure = [:];
379
380                if( fileService.fileExists( fastaFile ) ) {
381                        // Lookup the original filename
382                        def fastaData = processedFiles.parsed.success.find { it.filename == fastaFile };
383                        if( fastaData ) {
384                                returnStructure.fasta = fileService.moveFileToUploadDir( fileService.get( fastaFile ), fastaData.originalfilename, permanentDirectory );
385                                returnStructure.numSequences = fastaData.numSequences;
386                        } else {
387                                throw new Exception( "Fasta file wasn't uploaded the right way. Maybe the session has expired/" );
388                        }
389                } else {
390                        // If the file doesn't exist, we can't save anything to the database.
391                        throw new Exception( "Fasta file to save doesn't exist on disk" );
392                }
393
394                if( qualFile && fileService.fileExists( qualFile ) ) {
395                        // Lookup the original filename
396                        def qualData = processedFiles.parsed.success.find { it.filename == qualFile };
397
398                        if( qualData ) {
399                                returnStructure.qual = fileService.moveFileToUploadDir( fileService.get(qualFile ), qualData.originalfilename, permanentDirectory );
400                                returnStructure.avgQuality = qualData.avgQuality
401                        } else {
402                                // Delete the uploaded fasta file, since this is a serious error
403                                fileService.delete( returnStructure.fasta, permanentDirectory );
404                                throw new Exception( "Qual file wasn't uploaded the right way. Maybe the session has expired" );
405                        }
406                } else {
407                        // If the file doesn't exist, we don't save any quality information
408                        returnStructure.qual = null
409                        returnStructure.avgQuality = 0;
410                }
411
412                return returnStructure;
413        }
414}
Note: See TracBrowser for help on using the repository browser.