index previous next

iSignal: Interactive smoothing, differentiation, and signal analysis (version 3.4)

[Smoothing]  [Differentiation]  [Peak Sharpening]  [Background subtraction]  [Peak measurement]  [Curve Fitting]  [Frequency Spectrum]  [Examples]  [Matlab/Octave command-line function]  [Have a question? Email me]

iSignal is a Matlab function, written as a single self-contained m-file, for performing smoothing, differentiation, peak sharpening (resolution enhancement), Fourier frequency spectrum, least-squares peak fitting, and other useful functions on time-series data. Using simple keystrokes, you can adjust the signal processing parameters continuously while observing the effect on your signal dynamically. Click here to download the ZIP file "iSignal3.zip" that also includes some sample data for testing. You can also download it from the Matlab File Exchange. Version 3.5, August, 2014. T. C. O'Haver (toh@umd.edu). The syntax is:

pY=isignal(Data);
or
[pY,Spectrum]=isignal(Data,xcenter,xrange,SmoothMode,SmoothWidth,ends,DerivativeMode,Sharpen,Sharp1,Sharp2,SlewRate,MedianWidth,SpectrumMode);

"Data" may be a 2-column matrix with the independent variable (x-values) in the first column and dependent variable (y values) in the second column, or separate x and y vectors, or a single y-vector (in which case the data points are plotted against their index numbers on the x axis). Only the first argument (Data) is required; all the others are optional. Returns the processed independent axis ('pY') vector (and, in the Spectrum Mode, the frequency spectrum matrix, 'Spectrum') as the output arguments. Plots the data in the figure window, the lower half of the window showing the entire signal, and the upper half showing a selected portion controlled by the pan and zoom keystrokes, with the initial pan and zoom settings optionally controlled by input arguments 'xcenter' and 'xrange', respectively. Other keystrokes allow you to control the smooth type, width, and ends treatment, the derivative order (0th through 5th), and peak sharpening. (The initial values of all these parameters can be passed to the function via the optional input arguments SmoothMode, SmoothWidth, ends, DerivativeMode, Sharpen, Sharp1, Sharp2, SlewRate, and MedianWidth. See the examples below). Press K to see all the keyboard commands. Click to enlarge
 
Smoothing   
Reference: http://terpconnect.umd.edu/~toh/spectrum/Smoothing.html
The S key (or input argument "SmoothMode") cycles through five smoothing modes:
    If SmoothMode=0, the signal is not smoothed.
    If SmoothMode=1, rectangular (sliding-average or boxcar)
    If SmoothMode=2, triangular (2 passes of sliding-average)
    If SmoothMode=3, pseudo-Gaussian (3 passes of sliding-average)
    If SmoothMode=4, Savitzky-Golay smooth (thanks to Diederick).

The A and Z keys (or optional input argument SmoothWidth) control the SmoothWidth, w. Example shown in the figure on the right.

The X key toggles "ends" between 0 and 1. This determines how the "ends" of the signal (the first w/2 points and the last w/2 points) are handled when smoothing:
    If ends=0, the ends are zero.
    If ends=1, the ends are smoothed with progressively smaller smooths the closer to the end.
Generally, ends=1 is best, except in some cases using the derivative mode when ends=0 result in better vertical centering of the signal.

Note: when you are smoothing peaks, you can easily measure the effect of smoothing on peak height and width by turning on peak measure mode (press P) and then press S to cycle through the smooth modes.
Click to enlarge
 

There are two special functions for removing or reducing sharp spikes in signals: the M key, which implements a median filter (it asks you to enter the spike width, e.g. 1,2, 3... points) and the ~ key, which limits the maximum rate of change.


Differentiation  
Reference: http://terpconnect.umd.edu/~toh/spectrum/Differentiation.html
The
D / Shift-D keys (or optional input argument DerivativeMode) increase/decrease the derivative order. The default is 0. Careful optimization of smoothing of derivatives is critical for acceptable signal-to-noise ratio. Example shown in the figure on the right. In iSignal version 1.9 and above, in SmoothModes 1 through 3, derivatives are computed with respect to the independent variable (x-values), corrected for non-uniform x axis intervals. In SmoothMode 4 (Savitzky-Golay) the derivatives are computed by the Savitzky-Golay algorithm.
 
Peak sharpening (resolution enhancement)
Reference: http://terpconnect.umd.edu/~toh/spectrum/ResolutionEnhancement.html
The E key (or optional input argument "Sharpen") turns off and on peak sharpening (resolution enhancement). The sharpening strength is controlled by the F and V keys (or optional input argument "Sharp1") and B and G keys (or optional argument "Sharp2"). The optimum values depend on the peak shape and width. For peaks of Gaussian shape, a reasonable value for Sharp1 is PeakWidth2/25 and for Sharp2 is PeakWidth4/800 (or PeakWidth2/6 and PeakWidth4/700 for Lorentzian peaks), where PeakWidth is the full-width at half maximum of the peaks expressed in number of data points.  However, you don't need to do the math yourself; iSignal can calculate sharpening and smoothing settings for Gaussian and for Lorentzian peak shapes using the Y and U keys, respectively. Just isolate a single typical peak in the upper window using the pan and zoom keys, then press Y for Gaussian or U for Lorentzian peaks.  (The optimum settings depends on the width of the peak, so if your signal has peaks of widely different widths, one setting will not be optimum for all the peaks). You can fine-tune the sharpening with the F/V and G/B keys and the smoothing with the A/Z keys.

You can expect a decrease in peak width (and corresponding increase in peak height) of about 20% - 50%, depending on the shape of the peak (the peak area is largely unaffected by sharpening).  Excessive sharpening leads to baseline artifacts and increased noise. iSignal allows you to experimentally determine the values of these parameters that give the best trade-off between sharpening, noise, and baseline artifacts, for your purposes. 
You can easily measure the effect of sharpening quantitatively by turning on peak measure mode (press P) and then press E to toggle the sharpen mode off and on. Note: In iSignal version 1.9 and above, only the Savitzky-Golay smooth mode is used for peak sharpening. Example shown in the figure on the right.

Signal measurement
The cursor keys control the position of the green cursor (left and right arrow keys) and the distance between the dotted red cursors (up and down arrow keys) that mark the selected range displayed in the upper graph window. The label under the top graph window shows the value of the signal (y) at the green cursor, the peak-to-peak (min and max) signal range, the area under the signal, and the standard deviation within the selected range (the dotted cursors). Pressing the Q key prints out a table of the signal information in the command window.


Signal is measured by placing the green
center cursor on top of the peak

Noise is the standard deviation measured
on a flat portion of the baseline

Signal-to-noise ratio (SNR) measurement of a signal with very high SNR.  Left: The peak height of the largest signal peak is measured by placing the green center cursor on the largest peak; peak-to-peak signal=66,000. Right: The noise is measured on a flat portion of the baseline: standard deviation of noise=0.01, therefore the SNR=66,000/.01 = 6,600,000

Frequency Spectrum mode
Version 2.8 and later adds a Frequency Spectrum mode, toggled on and off by the Shift-S key, that computes the Fourier frequency spectrum of the segment of the signal displayed in the upper window and displays it in the lower window (temporarily replacing the full-signal display). Use the pan and zoom keys to adjust the region of the signal to be viewed. Press Shift-A to cycle through four plot modes (linear, semilog X, semilog Y, or log-log) and press Shift-X to toggle between a frequency on the x axis and time on the x-axis. The frequency (or time) and the harmonic number of the largest component of the spectrum is displayed on the left. All signal processing functions remain active in the frequency spectrum mode  (smooth, derivative, etc) so you can observe the effect of these functions on the frequency spectrum of the signal, as in the animated figure on the right. Press Shift-S again to return to the normal mode. To start off in the spectrum mode, set the 13th input argument, SpectrumMode, to 1. To save the spectrum as a new variable, call isignal with the output arguments [pY,Spectrum]:
>> x=0:.1:60; y=sin(x)+sin(10.*x);
>> [pY,Spectrum]=isignal([x;y],30,30,4,3,1,0,0,1,0,0,0,1);
>> plot(Spectrum) or isignal(Spectrum); or ipf(Spectrum); or ipeak(Spectrum); or whatever.
In version 3.4, Shift-Z toggles on and off peak detection and labeling on the frequency/time spectrum. You can adjust the peak detection in lines 1970-1973; see
http://terpconnect.umd.edu/~toh/spectrum/PeakFindingandMeasurement.htm. The Shift-W command displays the 3D waterfall spectrum, useful for signals whose frequency spectrum varies over the duration of the signal. You are asked to choose the number of segments into which to divide the signal (that is, the number of spectra) and the type of 3D display (mesh, contour, surface, etc).
Background subtraction
There are two ways to subtract the background from the signal: automatic and manual. To select an automatic baseline correction mode, press the T key repeatedly; it cycles thorough four modes: No baseline correction, linear baseline subtraction, quadratic baseline subtraction, flat baseline correction, then back to no baseline correction.  When baseline mode is linear, a straight-line baseline connecting the two ends of the signal segment in the upper panel will be automatically subtracted. When baseline mode is quadratic, a parabolic baseline connecting the two ends of the signal segment in the upper panel will be automatically subtracted. The baseline is calculated by computing a linear (or quadratic) least-squares fit to the signal in the first 1/10th of the points and the last 1/10th of the points. Try to adjust the pan and zoom to include some of the baseline at the beginning and end of the segment in the upper window, allowing the automatic baseline subtract gets a good reading of the baseline. The flat baseline mode is used only for peak fitting. The calculation of the signal amplitude, peak-to-peak signal, and peak area are all recalculated based on the baseline-subtracted signal in the upper window. If you are measuring peaks superimposed on a background, the use of the autozero mode will have a big effect on the measured peak height, width, and area, but very little effect on the peak x-axis position, as demonstrated by these two figures.

Click to
                enlarge


In addition to the four autozero baseline subtraction modes for peak measurement, a piecewise linear baseline can be subtracted manually from the entire signal in one operation. The Backspace key starts background correction operation. In the command window, type in the number of background points to click and press the Enter key. The cursor changes to crosshairs; click along the presumed background in the figure window, starting to the left of the x axis and placing the last click to the right of the x axis. When the last point is clicked, the linearly interpolated baseline between those points is subtracted from the signal. To restore the original background (i.e. to start over), press the '\' key (just below the backspace key).

Peak and valley measurement
The P key toggles off and on the "peak parabola" mode, which attempts to measure the one peak (or valley) that is centered in the upper window under the green cursor by superimposing a least-squares best-fit parabola, in red, on the center portion of the signal displayed in the upper window. (Zoom in so that the red overlays just the top of the peak or the bottom of the valley as closely as possible). The peak position, height, and "Gaussian width" are measured by least-squares curve fitting of a parabola to the isolated peak. The "RSquared" value is the coefficient of determination; the closer to 1.000 the better. The peak widths will most accurate if the peaks are Gaussian; other shapes, and very noisy peaks of any shape, will give only approximate results. However, the position, height, and area values are pretty good for any peak shape as long as the "RSquared" value is at least 0.99. The "total area" is measured by the trapezoidal method over the entire selected segment displayed in the upper window. The "SNR" is the signal-to-noise-ratio of the peak under the green cursor; it's the ratio of the peak height to the standard deviation of the residuals between the data and the best-fit line in red (version 2.4 and above). Example shown in the figure on the right. If the peaks are superimposed on a non-zero background, subtract the background before measuring peaks, either by using the autozero mode (T key) or the multi-point background subtraction (backspace key). Press the R key to print out the peak parameters in the command window.

In version 2.2 and above, peak area is measured two ways: the "Gaussian area" and the "Total area". The "Gaussian area"  is the area under the Gaussian that is best fit to the center portion of the signal displayed in the upper window, marked in red.  The "Total area" is the area by the trapezoidal method over the entire selected segment displayed in the upper window. If the portion of the signal displayed in the upper window is a pure Gaussian with no noise and a zero baseline, then the two measures should agree almost exactly.  If the peak is not Gaussian in shape, then the total area is likely to be more accurate, as long as the peak is well separated from other peaks.  If the peaks are overlapped, but have a known shape, then peak fitting (Shift-F) will give more accurate peak areas.

Peak fitting
Version 2.6 adds an iterative curve fitting method performed by peakfit.m. This is the most accurate method for the measurements of the areas of highly overlapped peaks. Press the Shift-F key, then select the desired peak shape from the menu, enter the number of peaks and the number of repeat trial fits.  A graph of the fit is displayed in Figure window 2 and a table of results is printed out in the command window. Note: for more flexible curve fitting, use ipf.m, which allows manual optimization of peak groupings and start positions. See Peak Fitter for details.

Saving the results
To save the processed signal to the disc as a x.y matrix in a.mat file, press the 'o' key, the type in the desired file name, then press Enter or click Save.

Other keystroke controls
The L key toggles off and on the Overlay mode, which shows the original signal as a dotted line overlayed on the current processed signal, for the purposes of comparison. The tab key restores the original signal and cursor settings. The ";" key sets the selected region to zero (to eliminate artifacts and spikes). The "-" (minus sign) key is used to negate the signal (flip + for -).  Press H to toggle display of semilog y plot in the lower window, which is useful for signals with very wide dynamic range, as in the example in the figures below (zero and negative points are ignored in the log plot).

Linear y-axis mode

Log y mode (H key)

The C key condenses the signal by the specified factor n, replacing each group of n points with their average (n must be an integer, such as 2,3, 4, etc). The I key replaces the signal with a linearly interpolated version containing m data points. This can be used either to increase or decrease the x-axis interval of the signal or to convert unevenly spaced values to evenly spaced values. After pressing C or I, you must type in the value of n or m respectively.

In version 2.96 you can press Shift-C, then click on the graph to print out the x,y coordinates of that point.  This works on both the upper and lower panels, and on the frequency spectrum as well.

In version 3.2, Shift-P to plays the segment of the signal displayed in the upper window as audio via computer's sound output. You are asked for the sampling rate, but you can put any number here - the larger the number the shorted and higher-pitched will be the sound - or you can simply press Enter to set the rate to 8000 Hz.

In version 3.5, Shift-Ctrl-F transfers the current signal to Interactive Peak Fitter (ipf.m) and Shift-Ctrl-P  transfers the current signal to Interactive Peak Detector (iPeak.m), if those functions are installed in your Matlab path.

Press K to see all the keyboard commands.

EXAMPLE 1: Single input argument; data in a two columns of a matrix [x;y] or in a single y vector
              >> isignal(y);
              >> isignal([x;y]);
 
EXAMPLE 2: Two input arguments. Data in separate x and y vectors.
              >> isignal(x,y);
 
EXAMPLE 3: Three or four input arguments. The last two arguments specify the initial values of pan (xcenter) and zoom (xrange) in the last two input arguments. Using data in the ZIP file:          
              >> load data.mat

              >> isignal(DataMatrix,180,40); or
              >> isignal(x,y,180,40);
 
EXAMPLE 4: As above, but additionally specifies initial values of SmoothMode, SmoothWidth, ends, and DerivativeMode in the last four input arguments. 
              >> isignal(DataMatrix,180,40,2,9,0,1);
 
EXAMPLE 5: As above, but additionally specifies initial values of the peak sharpening parameters Sharpen, Sharp1, and Sharp2 in the last three input arguments.  Press the E key to toggle sharpening on and off for comparison.
              >> isignal(DataMatrix,180,40,4,19,0,0,1,51,6000);
 
EXAMPLE 6:  Using the built-in "humps" function:
              >> x=[0:.005:2];y=humps(x);Data=[x;y];

                        4th derivative of the peak at x=0.9:
              >> isignal(Data,0.9,0.5,1,3,1,4); (shown on the right)

                        Peak sharpening applied to the peak at x=0.3:
              >> isignal(Data,0.3,0.5,1,3,1,0,1,220,5400);  (shown on the right)
                        (Press 'E' key to toggle sharpening ON/OFF to compare)

EXAMPLE 7: Measurement of peak area.  This example generates four Gaussian peaks, all with the exact same peak height (1.00) and area (1.77). The first peak (at x=4) is isolated, the second peak (x=9) is slightly overlapped with the third one, and the last two peaks (at x= 13 and 15) are strongly overlapped.  To measure the area under a peak using the perpendicular drop method, position the dotted red marker lines at the minimum between the overlapped peaks.  

>> x=[0:.01:20];
>> y=exp(-(x-4).^2)+exp(-(x-9).^2)+exp(-(x-13).^2)+exp(-(x-15).^2);
>> isignal(x,y); 

Greater accuracy in peak area measurement using iSignal can be obtained by using the peak sharpening function to reduce the overlap between the peaks. This reduces the peak widths, increases the peak heights, but has no effect on the peak areas.  

EXAMPLE 8: Single peak with random spikes (shown in the figure on the right). Compare smoothing vs spike filter (M key) and slew rate limit (~ key) to remove spikes. 

>> x=-5:.01:5;
>> y=exp(-(x).^2);for n=1:1000,if randn()>2,y(n)=rand()+y(n);,end,end;

>> isignal(x,y);
EXAMPLE 9:
The demo function isignaldemo2 creates a test signal containing three peaks with heights 1, 2, and 3, but with equal widths, superimposed on a very strong curved baseline, plus added random white noise. The objective is to extract a measure that is proportional to the peak height but independent of the baseline strength. Suggested approaches: (a) Use automatic or manual baseline subtraction to remove the baseline; or (b) use differentiation (with smoothing) to suppress the baseline.

EXAMPLE 10:   Direct entry into frequency spectrum mode, plotting returned frequency spectrum.
  >> x=0:.1:60; y=sin(x)+sin(10.*x);
  >> [pY,SpectrumOut]=isignal([x;y],30,30,4,3,1,0,0,1,0,0,0,1);
  >> plot(SpectrumOut)


KEYBOARD CONTROLS (Version 3.5):
   Pan left and right..........Coarse pan: < and >
                                Fine pan: left and right cursor arrows
                                Nudge: [ and ]
   Zoom in and out.............Coarse zoom: / and " 

                                Fine zoom: up and down cursor arrows
   Resets pan and zoom.........ESC
  
   Select entire signal........Ctrl-A
   Adjust smooth width.........A,Z  (A=>more, Z=>less)

   Adjust smooth type..........S (None, Rectangular, Triangle, Gaussian, Savitzky-Golay)
   Toggle smooth ends..........X (0=ends zeroed  1=ends smoothed (slower)
   Adjust derivative order.....D/Shift-D Increase/Decrease derivative order
   Toggle peak sharpening......E (0=OFF 1=ON)
   Sharpening for Gaussian.....Y  Set sharpen settings for Gaussian
   Sharpening for Lorentzian...U  Set sharpen settings for Lorentzian
   Adjust sharp1...............F,V  F=>sharper, V=>less sharpening
   Adjust sharp2   ............G,B  G=>sharper, B=>less sharpening
   Slew rate limit (0=OFF).....~ Largest allowed y change between points
   Spike filter width (0=OFF)..M  
median filter eliminates sharp spikes (New)
   Toggle peak parabola........P  fits parabola to center, labels vertex
   Fits peak in upper window...Shift-F (Asks for shape, number of peaks, etc)
   Spectrum mode on/off........Shift-S (Shift-A and Shift-X to change axes)
   Peak labels on spectrum.....Shift-Z in spectrum mode
   Display Waterfall spectrum..Shift-W  Allows choice of mesh, surf, contour, or pcolor
   Click graph to print x,y....Shift-C  Click graph to print coordinates
   Print peak report...........R  prints position, height, width, area
   Toggle overlay mode.........L Overlays original signal as dotted line
   Toggle log y mode...........H  semilog plot in lower window
   Cycle baseline mode.........T  none, linear, quadratic, or flat baseline mode
   Restores original signal....Tab key resets to original signal and modes
   Baseline subtraction........Backspace, then click baseline at multiple points
   Restore background..........\  to cancel previous background subtraction
   Invert signal...............-  Invert (negate) the signal (flip + and -)
   Remove offset...............0  (zero) set minimum signal to zero
   Sets region to zero.........;  sets selected region to zero.
   Absolute value..............+  Computes absolute value of entire signal
   Condense signal.............C  Condense oversampled signal by factor N
   Interpolate signal..........i  Interpolate (re-sample) to N points
   Print keyboard commands.....K  prints this list
   Print signal report.........Q  prints signal info and current settings
   Print isignal arguments.....W  prints isignal (current arguments)
   Save output to disk.........O as .mat file with processed signal matrix and
(in spectrum mode) the frequency spectrum.
   Play signal as sound........Shift-P  Play the signal segment in the upper panel
through computer sound system 
   Set sound sample rate.......Shift-R for the Shift-P command.
   Switch to ipf.m.............
Shift-Ctrl-F  Transfer current signal to Interactive Peak Fitter
   Switch to iPeak.............Shift-Ctrl-P  Transfer current signal to Interactive Peak Detector

Earlier versions
Earlier versions of iSignal as also available:1.5 1.61.7, 1.81.9, 2.0, 2.4, 2.5, 2.6, 2.7, 2.8, 2.9.

ProcessSignal, a Matlab/Octave command-line function that performs smoothing and differentiation on the time-series data set x,y (column or row vectors). Type "help ProcessSignal". Returns the processed signal as a vector that has the same shape as x, regardless of the shape of y. The syntax is Processed=ProcessSignal(x,y,DerivativeMode,w,type,ends,Sharpen,factor1,factor2,SlewRate,MedianWidth)


August, 2014
Tom O'Haver
Professor Emeritus
Department of Chemistry and Biochemistry
The University of Maryland at College Park
toh@umd.edu
http://www.wam.umd.edu/~toh


index previous next

Number of unique visits since May 17, 2008:

counter create hit