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

Last change on this file since 59 was 59, checked in by robert@…, 10 years ago
  • Improved speed by using ajax calls in pagination
  • Importing taxonomy files now works by adding it to sequenceData objects
File size: 9.2 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 org.hibernate.Transaction
11import java.util.zip.*
12import nl.tno.massSequencing.classification.*
13
14class ClassificationService implements nl.tno.massSequencing.imports.Importer {
15        def fileService
16        def sessionFactory
17
18        static transactional = false
19
20        /**
21         * Determines whether a file can be processed.
22         * @param filetype      Filetype of the file
23         * @see FileService.determineFileType()
24         * @return
25         */
26        public boolean canParseFileType( String filetype ) {
27                switch( filetype ) {
28                        case "groups":
29                        case "taxonomy":
30                                return true;
31                        default:
32                                return false;
33                }
34        }
35
36        /**
37         * Parses the given GROUPS or TAXONOMY file
38         *
39         * @param file                  File to parse
40         * @param filetype              Type of the given file
41         * @param onProgress    Closure to execute when progress indicators should be updated.
42         *                                              Has 2 parameters: numFilesProcessed and numBytesProcessed that indicate the number
43         *                                              of files and bytes that have been processed in this file (so the first parameter should
44         *                                              only be 1 when the file is finished)
45         *
46         * @return                              List structure. Examples:
47         *
48         *   [ success: true, filename: 'abc.fasta', type: 'fasta', numSequences: 200 ]
49         *   [ success: true, filename: 'abc.qual', type: 'qual', numSequences: 200, avgQuality: 36 ]
50         *   [ success: false, filename: 'abc.txt', type: 'txt', message: 'Filetype could not be parsed.' ]
51         */
52        public Map parseFile( File file, String filetype, Closure onProgress ) {
53                switch( filetype ) {
54                        case "groups":
55                        case "taxonomy":
56                                def fileInfo = countLines( file, onProgress );
57                                return [ success: true, filename: file.getName(), type: filetype, numLines: fileInfo.numLines, filesize: fileInfo.fileSize ]
58                        default:
59                                onProgress( file.length(), 0 );
60                                return [ success: false, type: filetype, message: 'Filetype could not be parsed.' ]
61                }
62        }
63
64        /**
65         * Parses the given GROUPS file
66         * @param file                  File to parse
67         * @param onProgress    Closure to execute when progress indicators should be updated.
68         *                                              Has 2 parameters: numFilesProcessed and numBytesProcessed that indicate the number
69         *                                              of files and bytes that have been processed in this file (so the first parameter should
70         *                                              only be 1 when the file is finished)
71         * @return                              List structure. Examples:
72         *
73         *   [ success: true, filename: 'abc.fasta', type: 'fasta', numSequences: 200 ]
74         *   [ success: false, filename: 'def.fasta', type: 'fasta', message: 'File is not a valid FASTA file' ]
75         */
76        protected Map countLines( File file, Closure onProgress ) {
77
78                long startTime = System.nanoTime();
79                log.trace "Start counting lines for " + file.getName()
80
81                // Count the number of lines. The real processing is done later on, because
82                // the GROUPS and TAXONOMY files must be handled together
83                long numLines = 0;
84                long bytesProcessed = 0;
85
86                file.eachLine { line ->
87                        numLines++;
88
89                        // Update progress every time 100kB is processed
90                        bytesProcessed += line.size();
91                        if( bytesProcessed > 100000 ) {
92                                onProgress( bytesProcessed, 0 );
93                                bytesProcessed = 0;
94                        }
95                }
96
97                // Update progress and say we're finished
98                onProgress( bytesProcessed, 0 );
99
100                log.trace "Finished counting lines for " + file.getName() + ": " + ( System.nanoTime() - startTime ) / 1000000L
101
102                return [ numLines: numLines, fileSize: file.size() ];
103        }
104
105        /**
106         * Removes all classification sequence objects for a given SequenceData object
107         * @param sequenceData 
108         * @return      True if the deletion has succeeded
109         */
110        public boolean removeClassificationForSequenceData( SequenceData sequenceData ) {
111                if( sequenceData ) {
112                        Sequence.executeUpdate( "DELETE FROM Sequence s WHERE s.sequenceData = ?", [sequenceData] );
113                        return true;
114                } else {
115                        return false;
116                }
117        }
118
119        /**
120         * Store the classification from the given files
121         * @param groupsFile    Filename
122         * @param taxonomyFile
123         * @return
124         */
125        public Map storeClassification( String taxonomyFilename, SequenceData sequenceData, Closure onProgress ) {
126                def hibernateSession = sessionFactory.getCurrentSession();
127                hibernateSession.flush();
128                hibernateSession.clear();
129               
130                File taxonomyFile = fileService.get( taxonomyFilename );
131
132                if( !taxonomyFile ) {
133                        return [ success: false, type: 'taxonomy', filename: groupsFilename, message: 'Taxonomy file doesn\'t exist.' ]
134                }
135               
136                // Open files for reading and create temporary variables
137                def taxonomyReader = taxonomyFile.newReader();
138
139                def parts, sequenceName, classification, sampleName, taxa, taxon, line, groupsLine
140
141                def taxonRegex = /"?([^"(]+)"?(\(\d+\))?/
142                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
143                def unclassified = "unclassified"
144
145                def i = 0;
146                def start = System.currentTimeMillis();
147                def lapTime = start;
148
149                // Create a stateless session in order to speed up inserts
150                StatelessSession statelessSession = sessionFactory.openStatelessSession();
151
152                // For each line in the taxonomy file, read the sequence and save it
153                def discardedSequences = 0;
154                def importedSequences = 0;
155                long bytesProcessed = 0;
156               
157                try {
158                        while( ( line = taxonomyReader.readLine() ) ) {
159                                bytesProcessed += line.size();
160       
161                                // Find the taxon for this line
162                                parts = line.tokenize( "\t" );
163       
164                                if( parts.size() != 2 ) {
165                                        // Skip this line because it is incorrect
166                                        continue;
167                                }
168       
169                                sequenceName = parts[ 0 ];
170                                classification = parts[ 1 ];
171       
172                                // Split classification to check whether the taxon already exists
173                                taxa = classification.tokenize( ';' ).collect { it.replaceAll( taxonRegex, '$1' ) }.findAll { it && it != unclassified }
174                               
175                                // Find taxon or create one if needed
176                                taxon = Taxon.findOrCreateTaxonByPath( taxa, startLevel );
177       
178                                if( taxon ) {
179                                        // Create a new sequence record
180                                        def s = new Sequence()
181                                        s.name = sequenceName
182                                        s.classification = taxon
183                                        s.sequenceData = sequenceData
184       
185                                        statelessSession.insert(s);
186                                } else {
187                                        log.trace( "" + i + " Sequence " + sequenceName + " not imported because it is unclassified." )
188                                }
189       
190                                if( i % 80 == 0 ) {
191                                        onProgress( bytesProcessed );
192                                        bytesProcessed = 0;
193                                }
194       
195                                importedSequences++;
196       
197                                i++;
198                        }
199                       
200                } catch( Exception e ) {
201                        throw e;
202                } finally {
203                        statelessSession.close();
204                        taxonomyReader.close();
205                }
206
207                onProgress( bytesProcessed );
208
209                return [ filename: taxonomyFilename, type: 'taxonomy', success: true, importedSequences: importedSequences, discardedSequences: discardedSequences ]
210        }
211
212        /**
213         * Updates the Classification table with data computed from the raw sequences for the given assay sample
214         * @param assaySample
215         */
216        public void updateClassificationForAssaySample( AssaySample assaySample ) {
217                if( !assaySample )
218                        return;
219                       
220                updateClassificationForAssaySamples( [ assaySample ] );
221        }
222
223        /**
224         * Updates the Classification table with data computed from the raw sequences for the given assay sample
225         * @param assaySample
226         */
227        public void updateClassificationForAssaySamples( List assaySamples ) {
228                if( !assaySamples )
229                        return;
230                       
231                // Clear hibernate session for otherwise the session will be out of sync with the
232                // database
233                def hibernateSession = sessionFactory.getCurrentSession();
234                hibernateSession.flush();
235                hibernateSession.clear();
236
237                // Remove all Classifications for the given assaySample
238                Classification.executeUpdate( "DELETE FROM Classification WHERE assaySample IN (:assaySamples)", [ 'assaySamples': assaySamples ] );
239
240                // Update classification table with data from sequences
241
242                // Unfortunately, Hibernate (through HQL) can't handle bulk insertion if a table has a
243                // composite id (Classification). See http://opensource.atlassian.com/projects/hibernate/browse/HHH-3434
244                // For that reason we execute plain SQL
245                def connection = sessionFactory.currentSession.connection()
246                def statement = connection.createStatement();
247                statement.executeUpdate( connection.nativeSQL( """
248                        INSERT INTO     classification (assay_sample_id, taxon_id, unclassified, total)
249                        SELECT          sa.id, s.classification_id, count( * ),
250                                                (
251                                                        SELECT          count( * )
252                                                        FROM            sequence s2
253                                                        LEFT JOIN       taxonomy t2 ON s2.classification_id = t2.id
254                                                        LEFT JOIN       sequence_data sd2 ON s2.sequence_data_id = sd2.id
255                                                        LEFT JOIN       assay_sample sa2 ON sd2.sample_id = sa2.id
256                                                        WHERE           sa2.id = sa.id AND t2.lft >= t.lft AND t2.rgt <= t.rgt
257                                                )
258                        FROM            sequence s
259                        LEFT JOIN       taxonomy t ON s.classification_id = t.id
260                        LEFT JOIN       sequence_data sd ON s.sequence_data_id = sd.id
261                        LEFT JOIN       assay_sample sa ON sd.sample_id = sa.id
262                        WHERE           sa.id IN (""" + assaySamples.id.join( ", " ) + """ )
263                        GROUP BY        sa.id, s.classification_id
264           """ ) );
265        }
266
267
268        /**
269         * Find a file that matches the parameters and is not already stored
270         * @param inputFiles
271         * @param type
272         * @param numLines
273         * @param alreadyStored
274         * @return
275         */
276        protected Map findMatchingClassificationFile( def inputFiles, String type, long numLines, def alreadyStored ) {
277                return inputFiles.find { it.type == type && it.numLines == numLines && !alreadyStored?.contains( it.filename ) };
278        }
279}
Note: See TracBrowser for help on using the repository browser.