source: trunk/grails-app/services/nl/tno/massSequencing/ClassificationService.groovy @ 58

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

Implemented importing of classifications

File size: 9.6 KB
Line 
1package nl.tno.massSequencing
2
3import java.io.BufferedWriter;
4import java.io.File;
5import java.io.Writer;
6import java.util.ArrayList;
7import java.util.List;
8import org.codehaus.groovy.grails.commons.ConfigurationHolder
9import org.hibernate.StatelessSession
10import java.util.zip.*
11import nl.tno.massSequencing.classification.*
12
13class ClassificationService implements nl.tno.massSequencing.imports.Importer {
14        def fileService
15        def sessionFactory
16
17        static transactional = false
18
19        /**
20         * Determines whether a file can be processed.
21         * @param filetype      Filetype of the file
22         * @see FileService.determineFileType()
23         * @return
24         */
25        public boolean canParseFileType( String filetype ) {
26                switch( filetype ) {
27                        case "groups":
28                        case "taxonomy":
29                                return true;
30                        default:
31                                return false;
32                }
33        }
34
35        /**
36         * Parses the given GROUPS or TAXONOMY file
37         *
38         * @param file                  File to parse
39         * @param filetype              Type of the given file
40         * @param onProgress    Closure to execute when progress indicators should be updated.
41         *                                              Has 2 parameters: numFilesProcessed and numBytesProcessed that indicate the number
42         *                                              of files and bytes that have been processed in this file (so the first parameter should
43         *                                              only be 1 when the file is finished)
44         *
45         * @return                              List structure. Examples:
46         *
47         *   [ success: true, filename: 'abc.fasta', type: 'fasta', numSequences: 200 ]
48         *   [ success: true, filename: 'abc.qual', type: 'qual', numSequences: 200, avgQuality: 36 ]
49         *   [ success: false, filename: 'abc.txt', type: 'txt', message: 'Filetype could not be parsed.' ]
50         */
51        public Map parseFile( File file, String filetype, Closure onProgress ) {
52                switch( filetype ) {
53                        case "groups":
54                        case "taxonomy":
55                                def fileInfo = countLines( file, onProgress );
56                                return [ success: true, filename: file.getName(), type: filetype, numLines: fileInfo.numLines, filesize: fileInfo.fileSize ]
57                        default:
58                                onProgress( file.length(), 0 );
59                                return [ success: false, type: filetype, message: 'Filetype could not be parsed.' ]
60                }
61        }
62
63        /**
64         * Parses the given GROUPS file
65         * @param file                  File to parse
66         * @param onProgress    Closure to execute when progress indicators should be updated.
67         *                                              Has 2 parameters: numFilesProcessed and numBytesProcessed that indicate the number
68         *                                              of files and bytes that have been processed in this file (so the first parameter should
69         *                                              only be 1 when the file is finished)
70         * @return                              List structure. Examples:
71         *
72         *   [ success: true, filename: 'abc.fasta', type: 'fasta', numSequences: 200 ]
73         *   [ success: false, filename: 'def.fasta', type: 'fasta', message: 'File is not a valid FASTA file' ]
74         */
75        protected Map countLines( File file, Closure onProgress ) {
76
77                long startTime = System.nanoTime();
78                log.trace "Start counting lines for " + file.getName()
79
80                // Count the number of lines. The real processing is done later on, because
81                // the GROUPS and TAXONOMY files must be handled together
82                long numLines = 0;
83                long bytesProcessed = 0;
84
85                file.eachLine { line ->
86                        numLines++;
87                       
88                        // Update progress every time 100kB is processed
89                        bytesProcessed += line.size();
90                        if( bytesProcessed > 100000 ) {
91                                onProgress( bytesProcessed, 0 );
92                                bytesProcessed = 0;
93                        }
94                }
95
96                // Update progress and say we're finished
97                onProgress( bytesProcessed, 0 );
98
99                log.trace "Finished counting lines for " + file.getName() + ": " + ( System.nanoTime() - startTime ) / 1000000L
100
101                return [ numLines: numLines, fileSize: file.size() ];
102        }
103       
104        /**
105         * Saves the classification of sequences into the database.
106         * @param inputFiles    List of maps with groups- and taxonomy files. 
107         * @return                              List of processed files
108         */
109        public List storeClassification( List inputFiles, Map excelMatches, Closure onProgress ) {
110                def filesStored = [];
111               
112                // In order to store classification, we need a .groups and a .taxonomy file with an
113                // equal number of lines
114                inputFiles.each { inputFile ->
115                        // Check whether this file has already been processed as a matching file for another file
116                        if( filesStored*.filename.contains( inputFile.filename ) )
117                                return;
118                               
119                        // Search for a matching file
120                        def searchForType = inputFile.type == 'groups' ? 'taxonomy' : 'groups';
121                        def matchingFile = findMatchingClassificationFile( inputFiles, searchForType, inputFile.numLines, filesStored*.filename );
122                       
123                        // If no matching file is found, give a message to the user
124                        if( !matchingFile ) {
125                                inputFile[ 'success' ] = false
126                                inputFile[ 'message' ] = "No matching " + searchForType + " file could be found for this file (or a matching file is already processed)."
127                               
128                                filesStored << inputFile;
129                                return;
130                        }
131                       
132                        println "Handling files: " + inputFile.filename + " - " + matchingFile.filename
133                        println excelMatches;
134                       
135                        // Now parse the two matching files
136                        def classificationInfo = []
137                        if( inputFile.type == 'groups' ) {
138                                classificationInfo = storeClassificationFromFiles( inputFile.filename, matchingFile.filename, excelMatches, onProgress );
139                        } else {
140                                classificationInfo = storeClassificationFromFiles( matchingFile.filename, inputFile.filename, excelMatches, onProgress );
141                        }
142                        filesStored += classificationInfo;
143                }
144               
145                return filesStored;
146        }
147       
148        /**
149         * Store the classification from the given files
150         * @param groupsFile    Filename
151         * @param taxonomyFile
152         * @return
153         */
154        protected List storeClassificationFromFiles( String groupsFilename, String taxonomyFilename, Map excelMatches, Closure onProgress ) {
155                def returnList = []
156               
157                File groupsFile = fileService.get( groupsFilename );
158                File taxonomyFile = fileService.get( taxonomyFilename );
159               
160                if( !groupsFile ) {
161                        returnList << [ success: false, type: 'groups', filename: groupsFilename, message: 'Groups file doesn\'t exist.' ]
162                }
163                if( !taxonomyFile ) {
164                        returnList << [ success: false, type: 'taxonomy', filename: groupsFilename, message: 'Taxonomy file doesn\'t exist.' ]
165                }
166               
167                if( !groupsFile || !taxonomyFile )
168                        return returnList;
169                       
170                // Open files for reading and create temporary variables
171                def taxonomyReader = taxonomyFile.newReader();
172                def groupsReader = groupsFile.newReader();
173                       
174                def parts, sequenceName, classification, sampleName, taxa, taxon, line, groupsLine
175               
176                def classificationScoreRegex = /(\(\d+\))/
177                def startLevel = 0    // The first level that is present in the file. Sometimes, the root element (level 0) is not mentioned in the file, in that case, this values should be set to 1
178               
179                def i = 0;
180                def start = System.currentTimeMillis();
181                def lapTime = start;
182               
183                // Create a stateless session in order to speed up inserts
184                StatelessSession statelessSession = sessionFactory.openStatelessSession();
185                statelessSession.beginTransaction();
186               
187                println "Starting at " + start
188                       
189                // For each line in the taxonomy file, also read the groups file. If the sequenceNames don't match,
190                // discard the sequence
191                def discardedSequences = 0;
192                def importedSequences = 0;
193                long bytesProcessed = 0;
194                while( ( line = taxonomyReader.readLine() ) ) {
195                        // Also read a line from the groupsFile
196                        groupsLine = groupsReader.readLine();
197                       
198                        bytesProcessed += line.size() + groupsLine.size();
199                       
200                        // Find the taxon for this line
201                        parts = line.tokenize( "\t" );
202                       
203                        if( parts.size() != 2 ) {
204                                // Skip this line because it is incorrect
205                                continue;
206                        }
207                       
208                        sequenceName = parts[ 0 ];
209                        classification = parts[ 1 ];
210                       
211                        if( !groupsLine ) {
212                                // Skip this line because the groups file doesn't contain a new line
213                                discardedSequences++;
214                                continue
215                        }
216
217                        // Find the sequence name and sample name in the groups file
218                        parts = groupsLine.tokenize( "\t" );
219                        if( parts.size() != 2 ) {
220                                // Skip this line because the groups file line is incorrect
221                                discardedSequences++;
222                                continue;
223                        }
224                       
225                        if( parts[ 0 ] != sequenceName ) {
226                                log.trace( "Sequence names in taxonomy and groups files don't match: " + sequenceName + " " + parts[ 0 ] )
227                                discardedSequences++;
228                                continue;
229                        }
230                        sampleName = parts[ 1 ];
231                       
232                        // Split classification to check whether the taxon already exists
233                        taxa = classification.replaceAll( classificationScoreRegex, '' ).tokenize( ';' ).findAll { it }
234                       
235                        // Find taxon or create one if needed
236                        taxon = Taxon.findOrCreateTaxonByPath( taxa, startLevel );
237                       
238                        // Determine assaySample
239                        def assaySample = excelMatches[ sampleName ]
240                       
241                        if( !assaySample ) {
242                                log.trace( "Sample name from GROUPS file is not found in excel file: " + sampleName )
243                                discardedSequences++;
244                                continue;
245                        }
246                       
247                        // Create a new sequence record
248                        def s = new Sequence()
249                        s.name = sequenceName
250                        s.classification = taxon
251                        s.assaySample = assaySample
252                       
253                        statelessSession.insert(s);
254                       
255                        //println "Sequence " + sequenceName + " / classification " + classification + " / sample: " + assaySample
256
257                        if( i % 100 == 0 ) {
258                                onProgress( bytesProcessed );
259                                bytesProcessed = 0;
260                        }
261                       
262                        importedSequences++;
263                       
264                        i++;
265                }
266               
267                onProgress( bytesProcessed );
268               
269                taxonomyReader.close();
270                groupsReader.close();
271               
272                returnList << [ filename: taxonomyFilename, type: 'taxonomy', success: true, importedSequences: importedSequences, discardedSequences: discardedSequences ]
273                returnList << [ filename: groupsFilename, type: 'groups', success: true, importedSequences: importedSequences, discardedSequences: discardedSequences ]
274               
275                return returnList;
276        }
277       
278        /**
279         * Find a file that matches the parameters and is not already stored
280         * @param inputFiles
281         * @param type
282         * @param numLines
283         * @param alreadyStored
284         * @return
285         */
286        protected Map findMatchingClassificationFile( def inputFiles, String type, long numLines, def alreadyStored ) {
287                return inputFiles.find { it.type == type && it.numLines == numLines && !alreadyStored?.contains( it.filename ) };
288        }
289}
Note: See TracBrowser for help on using the repository browser.