signalTimeFrequencyFeatureExtractor
Description
Use a signalTimeFrequencyFeatureExtractor object to extract time-frequency
features from a signal. You can use the extracted features to train a machine
learning model or a deep learning network.
Creation
Syntax
Description
sFE = signalTimeFrequencyFeatureExtractorsignalTimeFrequencyFeatureExtractor object with default property values.
sFE = signalTimeFrequencyFeatureExtractor(PropertyName=Value)signalTimeFrequencyFeatureExtractor object. For
example,
sFE = signalTimeFrequencyFeatureExtractor( ...
FrameSize=30,SampleRate=10,SpectralKurtosis=true)signalTimeFrequencyFeatureExtractor object that partitions a signal into frames of 30 samples
each and a rate of 10 Hz to extract the spectral kurtosis feature from the spectrogram of
each frame.Properties
Object Functions
extract | Extract time-domain, frequency-domain, or time-frequency-domain features |
generateMATLABFunction | Create MATLAB function compatible with C/C++ code generation |
getExtractorParameters | Get current parameter values of feature extractor object |
getScalarizationMethods | Get scalarization methods for domain-specific signal features |
setExtractorParameters | Set nondefault parameter values for feature extractor object |
setScalarizationMethods | Set scalarization methods for domain-specific signal features |
Examples
Extract the spectral kurtosis and instantaneous frequency from the spectrogram of a quadratic-swept chirp signal.
Generate a chirp with quadratic instantaneous frequency deviation. The chirp is sampled at 2 Hz for 60 seconds. The instantaneous frequency is 0.5 Hz at t = 0 and crosses 2 Hz at t = 60 seconds.
fs = 2*pi;
t = 0:1/fs:60;
x = chirp(t,0.5,60,2,"quadratic");Create a signalTimeFrequencyFeatureExtractor object and enable the SpectralKurtosis and InstantaneousFrequency time-frequency-domain features.
tfFE = signalTimeFrequencyFeatureExtractor( ...
SpectralKurtosis=true,InstantaneousFrequency=true)tfFE =
signalTimeFrequencyFeatureExtractor with properties:
Properties
FrameSize: []
FrameRate: []
SampleRate: []
IncompleteFrameRule: "drop"
FeatureFormat: "matrix"
Transform: "Spectrogram"
ScalarizationMethod: [1×1 timeFrequencyScalarFeatureOptions]
Enabled Features
SpectralKurtosis, InstantaneousFrequency
Disabled Features
SpectralSkewness, SpectralCrest, SpectralFlatness, SpectralEntropy, TFRidges, InstantaneousBandwidth
MeanEnvelopeEnergy, InstantaneousEnergy, WaveletEntropy, TimeSpectrum, ScaleSpectrum
Extract the features from the chirp signal. Reshape the result into the number of features extracted. Plot the features.
features = extract(tfFE,x); featuresRows = reshape(features,[],2); stackedplot(featuresRows,"*",... DisplayLabels=["Spectral Kurtosis" "Instantaneous Frequency"]) grid on
Set up a three-signal signalDatastore object and a signalTimeFrequencyFeatureExtractor object. Extract the spectral entropy and time-frequency ridges to the signals.
Create a signalDatastore object with three oscillating signals sampled at 5000 Hz for three seconds.
Fs = 5000;
t = 0:1/Fs:3;
members = {chirp(t,0,t(end),Fs/4); ...
2*chirp(t,Fs/20,t(end),Fs/2,"quadratic",[],"convex"); ...
vco(sin(2*pi*t),[0.1 0.4]*Fs,Fs)};
sds = signalDatastore(members,SampleRate=Fs,OutputDataType="single");Create a signalTimeFrequencyFeatureExtractor object. Specify a sample rate. Enable the spectral entropy and time-frequency ridges as features to extract.
tfFE = signalTimeFrequencyFeatureExtractor(SampleRate=Fs, ...
SpectralEntropy=true,TFRidges=true);Set up the extractor parameters for a signalTimeFrequencyFeatureExtractor object. Specify the range between Fs/5 and Fs/2.5 to extract the spectral entropy. Set the penalty parameter to zero to extract the time-frequency ridges.
setExtractorParameters(tfFE,"SpectralEntropy",Range=[Fs/5 Fs/2.5]); setExtractorParameters(tfFE,"TFRidges",Penalty=0);
Extract and plot the time-frequency features for all the signals from a signalDatastore object.
[features,infofeatures] = extract(tfFE,sds);
featureMatrix = cell2mat(features);
spectralEntropy = featureMatrix(:,infofeatures{1}.SpectralEntropy);
TFRidges = featureMatrix(:,infofeatures{1}.TFRidges);
featureNames = ["Spectral Entropy" "Time-Frequency Ridges"];
tiledlayout(numel(members),2)
for id=1:numel(members)
featureValues = {spectralEntropy(id,:) TFRidges(id,:)};
for fid = 1:numel(featureNames)
nexttile
plot(featureValues{fid})
title(featureNames(fid))
subtitle("Signal " + id)
axis tight
end
end
Specify scalarization methods to extract the scalar impulse factor and peak value of the instantaneous energy feature for an oscillating signal.
Set the impulse factor and the peak value as scalarization methods for the instantaneous energy time-frequency domain feature.
opts = scalarFeatureOptions("timefrequency", ... InstantaneousEnergy={'ImpulseFactor';'PeakValue'})
opts =
timeFrequencyScalarFeatureOptions with properties:
SpectralKurtosis: [0×0 string]
SpectralSkewness: [0×0 string]
SpectralCrest: [0×0 string]
SpectralFlatness: [0×0 string]
SpectralEntropy: [0×0 string]
TFRidges: [0×0 string]
InstantaneousBandwidth: [0×0 string]
InstantaneousFrequency: [0×0 string]
InstantaneousEnergy: [2×1 string]
MeanEnvelopeEnergy: [0×0 string]
TimeSpectrum: [0×0 string]
WaveletEntropy: [0×0 string]
ScaleSpectrum: [0×0 string]
All: [0×0 string]
Create a signalTimeFrequencyFeatureExtractor object that returns scalar values for the instantaneous energy feature from the variational mode decomposition of a signal.
tfFE = signalTimeFrequencyFeatureExtractor(Transform="vmd", ... InstantaneousEnergy=true,ScalarizationMethod=opts);
Extract the vector and scalar features for a voltage-controlled oscillating signal.
fs = 10000; t = 0:1/fs:2; x = vco(sawtooth(2*pi*t,0.75),[0.1 0.4]*fs,fs); [features,info] = extract(tfFE,x);
Display the scalar feature values.
% Impulse Factor of the Instantaneous Energy Vector
features(info.InstantaneousEnergyImpulseFactor)ans = 8.7551
% Peak Value of the Instantaneous Energy Vector
features(info.InstantaneousEnergyPeakValue)ans = 1.3385
Extract the time-frequency ridges from a spectrogram using a signalTimeFrequencyFeatureExtractor object.
Load a data file containing an echolocation pulse batsignal, emitted by a big brown bat (Eptesicus fuscus) and measured with a sample rate DT of 7 microseconds. For more information, see the example Find and Track Ridges Using Reassigned Spectrogram. Thanks to Curtis Condon, Ken White, and Al Feng of the Beckman Center at the University of Illinois for the bat data and permission to use it in this example.
load batsignalCreate a MATLAB® timetable using the signal and the time information.
t = (0:length(batsignal)-1)*DT; sg = timetable(seconds(t)',batsignal);
Create an object to extract time-frequency ridges from the spectrogram of a signal.
fs = 1/DT; sTFFE = signalTimeFrequencyFeatureExtractor( ... Transform="spectrogram",SampleRate=fs,TFRidges=true);
Set the leakage, the time resolution, and the overlap percent parameters for the spectrogram analysis method.
lk = 0.9; tRes = 0.00028; % seconds oPercent = 85; setExtractorParameters(sTFFE,"spectrogram", ... Leakage=lk,TimeResolution=tRes,OverlapPercent=oPercent);
Set the number of ridges, number of frequency bins, and penalty parameters for the time-frequency ridges feature.
nRidges = 3; setExtractorParameters(sTFFE,"TFRidges", ... NumRidges=nRidges,NumFrequencyBins=128,Penalty=0.01);
Extract and plot the time-frequency ridges.
tfRidges = reshape(extract(sTFFE,sg),[],nRidges); tStamps = seconds(tRes/2+(1-oPercent/100)*tRes*(0:length(tfRidges)-1)); plot(tStamps,tfRidges/1000) xlabel("Time") ylabel("Frequency (kHz)")
Extract the instantaneous frequency of a signal using the wavelet synchrosqueezed transform (WSST), the time-averaged wavelet spectrum using the continuous wavelet transform (CWT), and the wavelet entropy using the maximal overlap discrete wavelet transform (MODWT).
Load an ECG signal corresponding to record 200 of the MIT-BIH Arrhythmia Database [1]. The workspace variable ecgsig contains the signal, and the variable tm contains the sample times. The sample rate is approximately 400 Hz.
load mit200
Fs = 400;The signal has 10,000 samples. Reshape the signal into a 1000-by-10 matrix. Each column represents one channel in a multichannel signal. Create a signalDatastore object from the matrix.
ecgsig = reshape(ecgsig,[],10); ecgsig = num2cell(ecgsig,1); sds = signalDatastore(ecgsig,SampleRate=Fs);
Set the mean and standard deviation as scalarization methods for the wavelet entropy feature.
opts = scalarFeatureOptions("timefrequency", ... WaveletEntropy=["Mean" "StandardDeviation"]);
Create a signalTimeFrequencyFeatureExtractor object that returns:
Wavelet entropy and the scalar values specified in
optsInstantaneous frequency from a WSST
Time-averaged spectrum from a CWT
You can extract instantaneous frequency from time-frequency signal representations created using methods such as the short-time Fourier transform ("spectrogram"), empirical mode decomposition ("emd"), and the maximal overlap discrete wavelet packet transform ("waveletpacket"). Set the time-frequency analysis method to "synchrosqueezedscalogram". The CWT is the only analysis method you can use to extract the time-averaged wavelet spectrum. The signalTimeFrequencyFeatureExtractor object uses the MODWT, the default analysis method of the wentropy function, to extract the wavelet entropy.
tfFE = signalTimeFrequencyFeatureExtractor( ... SampleRate=Fs, ... Transform="synchrosqueezedscalogram", ... InstantaneousFrequency=true, ... WaveletEntropy=true, ... TimeSpectrum=true, ... ScalarizationMethod=opts);
Extract the features from all the signals.
data = read(sds);
[features,info] = extract(tfFE,sds);
featureMatrix = cell2mat(features);
instFrq = featureMatrix(:,info{1}.InstantaneousFrequency);
timeSpect = featureMatrix(:,info{1}.TimeSpectrum);
waveEnt = featureMatrix(:,info{1}.WaveletEntropy);
waveEntMean = featureMatrix(:,info{1}.WaveletEntropyMean);
waveEntStd = featureMatrix(:,info{1}.WaveletEntropyStandardDeviation);For each signal, compute the mean and standard deviation directly from the extracted wavelet entropy. Confirm they are equal to the extracted scalar values.
[mean(waveEnt,2) waveEntMean]
ans = 10×2
0.7675 0.7675
0.7941 0.7941
0.7533 0.7533
0.7465 0.7465
0.8055 0.8055
0.8141 0.8141
0.8048 0.8048
0.7833 0.7833
0.7999 0.7999
0.7910 0.7910
[std(waveEnt,0,2) waveEntStd]
ans = 10×2
0.1714 0.1714
0.1258 0.1258
0.2063 0.2063
0.1945 0.1945
0.1367 0.1367
0.1123 0.1123
0.1235 0.1235
0.1228 0.1228
0.1123 0.1123
0.1168 0.1168
Display the time-averaged wavelet spectrum of the signals as a waterfall plot.
waterfall(timeSpect) title("Time-Averaged Wavelet Spectrum") xlabel("Frequency Index") ylabel("Channel Index") zlabel("Power")
Display the instantaneous frequency of the signals as a waterfall plot.
waterfall(instFrq) view(30,40) title("Instantaneous Frequency") xlabel("Sample") ylabel("Channel Index") zlabel("Frequency")
Since R2025a
Extract time-domain, frequency-domain, and time-frequency features from healthy bearing vibration signals and faulty bearing vibration signals. While a healthy bearing vibration signal does not have outstanding defects, a faulty bearing vibration signal results from wear-and-tear defects, such as spalls on the gear teeth, eccentricity or gear misalignment, and cracks at the races.
For more information on bearing signal generation and analysis, see Vibration Analysis of Rotating Machinery. To learn more about the feature extraction and model training workflow to identify faulty bearing signals in mechanical systems, see Machine Learning and Deep Learning Classification Using Signal Feature Extraction Objects.
Generate Healthy Bearing Signal
Generate a healthy bearing vibration signal as a sum of three cosine pulses with amplitudes of 0.4 V, 0.2 V, and 1 V, respectively, and frequencies of 22.5 Hz, 8.36 Hz, and 292.5 Hz, respectively, for three seconds and with a sample rate of 20 kHz. Generate Gaussian noise and add it to the signal.
rng("default")
Fs = 20e3;
t = (0:1/Fs:3-1/Fs)';
a = [0.4 0.2 1];
f = [22.5 8.36 292.5];
sClean = cos(2*pi*f.*t)*a';
sHealthy = sClean + 0.2*randn(size(t));Generate Faulty Bearing Signal
Generate a faulty bearing vibration signal by adding a bearing impact signal to the healthy bearing signal. Model each impact as a 3 kHz sinusoid windowed by a Kaiser window. The defect causes a series of 10-millisecond impacts on the bearing.
tImpact = t(t<10e-3)'; xImpact = sin(2*pi*3000*tImpact).*kaiser(length(tImpact),40)'; xImpactBper = 0.33*pulstran(t,0:1/104.5:t(end),xImpact,Fs);
Generate a faulty bearing vibration signal as a sum of the healthy bearing signal, the bearing impact signal, and additive Gaussian noise.
sFaulty = sHealthy + xImpactBper;
Consolidate and Visualize Signals
Bundle the healthy bearing and faulty bearing signals in a signalDatastore object in single precision.
sds = signalDatastore({sHealthy,sFaulty},OutputDataType="single");Plot the power spectrum of the healthy and faulty vibration signals. Observe the peaks that correspond to the bearing impact.
[P,F] = pspectrum([sHealthy sFaulty],Fs); p = plot(F/1000,pow2db(P)); p(1).Marker = "."; xlabel("Frequency (kHz)") ylabel("Power Spectrum (dB)") legend(["Healthy" "Faulty"])
Set Up Feature Extraction Pipeline
Create a signalTimeFeatureExtractor object for time-domain feature extraction.
timeFE = signalTimeFeatureExtractor(SampleRate=Fs,...
RMS=true,ImpulseFactor=true,StandardDeviation=true);Create a signalFrequencyFeatureExtractor object for frequency-domain feature extraction.
freqFE = signalFrequencyFeatureExtractor(SampleRate=Fs, ...
MedianFrequency=true,BandPower=true,PeakAmplitude=true);Create a signalTimeFrequencyFeatureExtractor object to extract time-frequency features from a spectrogram. Set the leakage parameter for the spectrogram to 90%.
timeFreqFE = signalTimeFrequencyFeatureExtractor(SampleRate=Fs, ... SpectralKurtosis=true,SpectralSkewness=true,TFRidges=true); setExtractorParameters(timeFreqFE,"spectrogram",Leakage=0.9);
Extract Multidomain Features
Extract signal features using all three feature extractors for the signals in the signalDatastore object. Concatenate and display the features in a multidomain feature table.
features = cellfun(@(a,b,c) [real(a) real(b) real(c)], ... extract(timeFE,sds),extract(freqFE,sds),extract(timeFreqFE,sds), ... UniformOutput=false); featureMatrix = cell2mat(features); featureTable = array2table(featureMatrix); rows2vars(featureTable)
ans=579×3 table
OriginalVariableNames Var1 Var2
_____________________ _______ _______
{'featureMatrix1' } 0.80115 0.80538
{'featureMatrix2' } 0.80116 0.80539
{'featureMatrix3' } 3.2635 3.1501
{'featureMatrix4' } 292.39 292.41
{'featureMatrix5' } 0.64086 0.64764
{'featureMatrix6' } 0.20977 0.20977
{'featureMatrix7' } 27.474 25.423
{'featureMatrix8' } 35.088 32.666
{'featureMatrix9' } 25.867 24.521
{'featureMatrix10'} 29.091 26.485
{'featureMatrix11'} 36.085 32.242
{'featureMatrix12'} 32.92 31.423
{'featureMatrix13'} 24.421 22.288
{'featureMatrix14'} 26.056 24.772
{'featureMatrix15'} 30.36 28.084
{'featureMatrix16'} 26.464 25.432
⋮
More About
Algorithms
Assume an input signal x sampled at a rate Fs, from
which to extract time-frequency features. When you specify signal framing properties
(FrameSize, FrameRate or FrameOverlapLength, and IncompleteFrameRule), the feature extractor sets up the signal partitioning
operation for x to extract features for each frame. This table shows the
equivalent syntaxes that signalTimeFrequencyFeatureExtractor uses to partition the signal
x into frames of size fl, frame rate
fr or frame overlap length ol, and incomplete frame
rule ifr.
| Frame Specifications | Feature Extractor Object Specification | Signal Framing Operation Equivalency |
|---|---|---|
sFE = signalTimeFrequencyFeatureExtractor( ... FrameSize=fl,FrameRate=fr, ... IncompleteFrameRule=ifr); |
xFrames = framesig(x,fl, ... OverlapLength=fl-fr, ... IncompleteFrameRule=ifr); | |
sFE = signalTimeFrequencyFeatureExtractor( ... FrameSize=fl,FrameOverlapLength=ol, ... IncompleteFrameRule=ifr); |
xFrames = framesig(x,fl, ... OverlapLength=ol, ... IncompleteFrameRule=ifr); |
If you do not specify signal framing properties, signalTimeFrequencyFeatureExtractor
considers x as a single-framed signal.
Given the single-framed input signal x and sample rate
Fs, this table lists the equivalent syntaxes for extracting features
using the signalTimeFrequencyFeatureExtractor object and the individual feature extractor
functions.
| Features | Feature Extractor Object | Individual Feature Extractors |
|---|---|---|
|
sFE = signalTimeFrequencyFeatureExtractor( ... SampleRate=Fs, ... Transform="spectrogram", ... SpectralKurtosis=true, ... SpectralSkewness=true, ... SpectralCrest=true, ... SpectralFlatness=true, ... SpectralEntropy=true, ... TFRidges=true, ... InstantaneousBandwidth=true, ... InstantaneousFrequency=true); features = extract(sFE,x); |
[S,F,T] = pspectrum( ... x,Fs,"spectrogram"); features = [ ... spectralKurtosis(S,F)' ... spectralSkewness(S,F)' ... spectralCrest(S,F)' ... spectralFlatness(S,F)' ... spectralEntropy(S,F)' ... tfridge(S,F)' ... instbw(S,F,T)' ... instfreq(S,F,T)' ... ]; |
|
sFE = signalTimeFrequencyFeatureExtractor( ... SampleRate=Fs, ... Transform="synchrosqueezedspectrogram", ... SpectralKurtosis=true, ... SpectralSkewness=true, ... SpectralCrest=true, ... SpectralFlatness=true, ... SpectralEntropy=true, ... TFRidges=true, ... InstantaneousBandwidth=true, ... InstantaneousFrequency=true); features = extract(sFE,x); |
[S,F,T] = fsst(x,Fs); S = abs(S).^2; F = F/(2*pi); features = [ ... spectralKurtosis(S,F)' ... spectralSkewness(S,F)' ... spectralCrest(S,F)' ... spectralFlatness(S,F)' ... spectralEntropy(S,F)' ... tfridge(S,F)' ... instbw(S,F,T)' ... instfreq(S,F,T)' ... ]; |
|
sFE = signalTimeFrequencyFeatureExtractor( ... SampleRate=Fs, ... Transform="synchrosqueezedscalogram", ... SpectralKurtosis=true, ... SpectralSkewness=true, ... SpectralCrest=true, ... SpectralFlatness=true, ... SpectralEntropy=true, ... TFRidges=true, ... InstantaneousBandwidth=true, ... InstantaneousFrequency=true); features = extract(sFE,x); |
[S,F] = wsst(x,Fs); S = abs(S).^2; T = (0:1:size(S,2)-1)/Fs; features = [ ... spectralKurtosis(S,F)' ... spectralSkewness(S,F)' ... spectralCrest(S,F)' ... spectralFlatness(S,F)' ... spectralEntropy(S,F)' ... tfridge(S,F)' ... instbw(S,F,T)' ... instfreq(S,F,T)' ... ]; |
|
sFE = signalTimeFrequencyFeatureExtractor( ... SampleRate=Fs, ... Transform="emd", ... InstantaneousFrequency=true, ... MeanEnvelopeEnergy=true, ... InstantaneousEnergy=true); features = extract(sFE,x); |
[imf,~,d] = emd(x,MaxNumIMF=10); [~,~,~,F,E] = hht(imf,Fs); m = d.MeanEnvelopeEnergy'; features = [ ... F(:)' ... m zeros(1,10-numel(m)) ... E(:)' ... ]; |
VMD features: |
sFE = signalTimeFrequencyFeatureExtractor( ... SampleRate=Fs,Transform="vmd", ... InstantaneousFrequency=true, ... InstantaneousEnergy=true); features = extract(sFE,x); |
[~,~,~,F,E] = hht(vmd(x),Fs); features = [ ... F(:)' ... E(:)' ... ]; |
|
sFE = signalTimeFrequencyFeatureExtractor( ... SampleRate=Fs, ... Transform="wavelet", ... InstantaneousFrequency=true, ... InstantaneousEnergy=true, ... WaveletEntropy=true); features = extract(sFE,x); |
mra = modwtmra(modwt(x)); imf = permute(mra,[2 1 3]); [~,~,~,F,E] = hht(imf,Fs); features = [ ... F(:)' ... E(:)' ... wentropy(modwt(x))' ... ]; |
|
sFE = signalTimeFrequencyFeatureExtractor( ... SampleRate=Fs, ... Transform="waveletpacket", ... InstantaneousFrequency=true, ... InstantaneousEnergy=true, ... WaveletEntropy=true); features = extract(sFE,x); |
d = modwptdetails(x); imf = permute(d,[2 1 3]); [~,~,~,F,E] = hht(imf,Fs); features = [ ... F(:)' ... E(:)' ... wentropy(modwpt(x))' ... ]; |
|
sFE = signalTimeFrequencyFeatureExtractor( ... SampleRate=Fs, ... Transform="scalogram", ... TimeSpectrum=true, ... ScaleSpectrum=true); features = extract(sFE,x); |
fb = cwtfilterbank( ... SignalLength=numel(x)); features = [ ... timeSpectrum(fb,wt(fb,x))' ... scaleSpectrum(fb,wt(fb,x)) ... ]; |
Note
Because
signalTimeFrequencyFeatureExtractorreuses intermediary representations, the features thatsignalTimeFrequencyFeatureExtractoroutputs might not correspond with the default configuration of features output by the corresponding individual feature extraction functions.To obtain the equivalent syntax for the feature extraction setup based on the properties specified when you create the
signalTimeFrequencyFeatureExtractorobject, usegenerateMATLABFunction.
References
[1] Moody, G.B., and R.G. Mark. “The Impact of the MIT-BIH Arrhythmia Database.” IEEE Engineering in Medicine and Biology Magazine 20, no. 3 (June 2001): 45–50. https://doi.org/10.1109/51.932724.
Extended Capabilities
Version History
Introduced in R2024aSee Also
Functions
Objects
Select a Web Site
Choose a web site to get translated content where available and see local events and offers. Based on your location, we recommend that you select: United States.
You can also select a web site from the following list
Americas
- América Latina (Español)
- Canada (English)
- United States (English)
Europe
- Belgium (English)
- Denmark (English)
- Deutschland (Deutsch)
- España (Español)
- Finland (English)
- France (Français)
- Ireland (English)
- Italia (Italiano)
- Luxembourg (English)
- Netherlands (English)
- Norway (English)
- Österreich (Deutsch)
- Portugal (English)
- Sweden (English)
- Switzerland
- United Kingdom (English)