1 # -*- coding: utf-8 -*-
3 # Copyright (C) 2008-2010 Alberto Gomez-Casado
5 # Massimo Sandal <devicerandom@gmail.com>
6 # W. Trevor King <wking@drexel.edu>
8 # This file is part of Hooke.
10 # Hooke is free software: you can redistribute it and/or
11 # modify it under the terms of the GNU Lesser General Public
12 # License as published by the Free Software Foundation, either
13 # version 3 of the License, or (at your option) any later version.
15 # Hooke is distributed in the hope that it will be useful,
16 # but WITHOUT ANY WARRANTY; without even the implied warranty of
17 # MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
18 # GNU Lesser General Public License for more details.
20 # You should have received a copy of the GNU Lesser General Public
21 # License along with Hooke. If not, see
22 # <http://www.gnu.org/licenses/>.
24 """The ``flatfilt`` module provides :class:`~FlatFiltPlugin` and
25 several associated :class:`~hooke.command.Command`\s for removing flat
26 (featureless) :mod:`~hooke.curve.Curve`\s from
27 :class:`~hooke.playlist.Playlist`\s.
31 :mod:`~hooke.plugin.convfilt` for a convolution-based filter for
32 :class:`~hooke.experiment.VelocityClamp` experiments.
36 from multiprocessing import Queue
38 from numpy import diff
39 from scipy.signal.signaltools import medfilt
41 from ..command import Command, Argument, Success, Failure
42 from ..config import Setting
43 from ..experiment import VelocityClamp
44 from ..plugin import Plugin, argument_to_setting
45 from ..plugin.curve import CurveArgument
46 from ..plugin.playlist import FilterCommand
47 from ..plugin.vclamp import scale
48 from ..util.peak import find_peaks, find_peaks_arguments, Peak, _kwargs
51 class FlatFiltPlugin (Plugin):
52 """Standard-devitiation-based peak recognition and filtering.
55 super(FlatFiltPlugin, self).__init__(name='flatfilt')
56 self._arguments = [ # For Command initialization
57 Argument('median window', type='int', default=7, help="""
58 Median window filter size (in points).
60 Argument('blind window', type='float', default=20e-9, help="""
61 Meters after the contact point where we do not count peaks to avoid
62 non-specific surface interaction.
64 Argument('min peaks', type='int', default=4, help="""
65 Minimum number of peaks for curve acceptance.
67 ] + copy.deepcopy(find_peaks_arguments)
68 # Set flat-filter-specific defaults for the fit_peak_arguments.
69 for key,value in [('cut side', 'both'),
72 ('min deviations', 9.0),
74 ('see double', 10e-9),
76 argument = [a for a in self._arguments if a.name == key][0]
77 argument.default = value
79 Setting(section=self.setting_section, help=self.__doc__)]
80 for argument in self._arguments:
81 self._settings.append(argument_to_setting(
82 self.setting_section, argument))
83 argument.default = None # if argument isn't given, use the config.
84 self._commands = [FlatPeaksCommand(self), FlatFilterCommand(self)]
86 def dependencies(self):
89 def default_settings(self):
93 class FlatPeaksCommand (Command):
94 """Detect peaks in velocity clamp data using noise statistics.
99 Noise analysis on the retraction curve:
101 1) A median window filter (using
102 :func:`scipy.signal.signaltools.medfilt`) smooths the
104 2) The deflection derivative is calculated (using
105 :func:`numpy.diff` which uses forward differencing).
106 3) Peaks in the derivative curve are extracted with
107 :func:`~hooke.plugins.peak.find_peaks`.
109 The algorithm was originally Francesco Musiani's idea.
111 def __init__(self, plugin):
112 config_arguments = [a for a in plugin._arguments
113 if a.name != 'min peaks']
114 # Drop min peaks, since we're not filtering with this
115 # function, just detecting peaks.
116 super(FlatPeaksCommand, self).__init__(
117 name='flat filter peaks',
120 ] + config_arguments,
121 help=self.__doc__, plugin=plugin)
123 def _run(self, hooke, inqueue, outqueue, params):
124 z_data,d_data,params = self._setup(params)
126 while (start_index < len(z_data)
127 and z_data[start_index] < params['blind window']):
129 median = medfilt(d_data[start_index:], params['median window'])
131 peaks = find_peaks(deriv, **_kwargs(params, find_peaks_arguments,
132 argument_input_keys=True))
134 peak.name = 'flat filter of %s' % (params['deflection column name'])
135 peak.index += start_index
138 def _setup(self, params):
139 """Setup `params` from config and return the z piezo and
142 curve = params['curve']
143 if curve.info['experiment'] != VelocityClamp:
144 raise Failure('%s operates on VelocityClamp experiments, not %s'
145 % (self.name, curve.info['experiment']))
146 for col in ['surface z piezo (m)', 'deflection (N)']:
147 if col not in curve.data[0].info['columns']:
150 for block in curve.data:
151 if block.info['name'].startswith('retract'):
155 raise Failure('No retraction blocks in %s.' % curve)
156 z_data = data[:,data.info['columns'].index('surface z piezo (m)')]
157 if 'flattened deflection (N)' in data.info['columns']:
158 params['deflection column name'] = 'flattened deflection (N)'
160 params['deflection column name'] = 'deflection (N)'
161 d_data = data[:,data.info['columns'].index(
162 params['deflection column name'])]
163 for key,value in params.items():
164 if value == None: # Use configured default value.
165 params[key] = self.plugin.config[key]
166 # TODO: better option parser to do this automatically by Argument.type
167 for key in ['blind window', 'median window', 'max cut', 'min deviations', 'min points', 'see double', 'stable']:
168 params[key] = float(params[key])
169 # TODO: convert 'see double' from nm to points
170 return z_data,d_data,params
172 class FlatFilterCommand (FilterCommand):
173 u"""Create a subset playlist of curves with enough flat peaks.
177 This type of filter is of course very raw, and requires relatively
178 conservative settings to safely avoid false negatives (that is, to
179 avoid discarding interesting curves). Using it on the protein
180 unfolding experiments described by Sandal [#sandal2008] it has
181 been found to reduce the data set to analyze by hand by 60-80%.
183 .. [#sandal2008] M. Sandal, F. Valle, I. Tessari, S. Mammi, E. Bergantino,
184 F. Musiani, M. Brucale, L. Bubacco, B. Samorì.
185 "Conformational equilibria in monomeric α-Synuclein at the
186 single molecule level."
188 doi: `10.1371/journal.pbio.0060006 <http://dx.doi.org/10.1371/journal.pbio.0060006>`_
192 FlatCommand : Underlying flat-based peak detection.
194 def __init__(self, plugin):
195 super(FlatFilterCommand, self).__init__(
196 plugin, name='flat filter playlist')
197 self.arguments.extend(plugin._arguments)
199 def filter(self, curve, hooke, inqueue, outqueue, params):
202 filt_command = [c for c in hooke.commands
203 if c.name=='flat filter peaks'][0]
204 filt_command.run(hooke, inq, outq, **params)
206 if not (isinstance(peaks, list) and (len(peaks) == 0
207 or isinstance(peaks[0], Peak))):
208 raise Failure('Expected a list of Peaks, not %s' % peaks)
210 if not isinstance(ret, Success):
212 if params['min peaks'] == None: # Use configured default value.
213 params['min peaks'] = self.plugin.config['min peaks']
214 return len(peaks) >= int(params['min peaks'])