В этом примере показано, как классифицировать отказы на акустические записи воздушных компрессоров с помощью сети рассеивания вейвлета и машины опорных векторов (SVM). Пример обеспечивает возможность использовать графический процессор, чтобы ускориться, расчет рассеивания вейвлета преобразовывают. Если вы хотите использовать графический процессор, у вас должны быть Parallel Computing Toolbox™ и поддерживаемый графический процессор. Смотрите Поддержку графического процессора Релизом (Parallel Computing Toolbox) для деталей.
Набор данных состоит из акустических записей, собранных на одноступенчатом воздушном компрессоре типа оплаты [1]. Данные производятся на уровне 16 кГц. Технические требования воздушного компрессора следующие:
Область значений Давления воздуха: 0-500 lb/m2, 0-35 Kg/cm2
Асинхронный двигатель: 5 л. с., 415 В, 5:00, 50 Гц, 1440 об/мин
Датчик давления: введите PR-15, область значений 100-213 фунтов на квадратный дюйм
Каждая запись представляет одно из 8 состояний, которое включает здоровое состояние и 7 дефектных состояний. 7 дефектных состояний:
Отказ утечки вставила клапан (LIV)
Отказ клапана выхода утечки (LOV)
Отказ обратного клапана (NRV)
Поршневой кольцевой отказ
Отказ маховика
Отказ пояса наездника
Подшипник отказа
Загрузите набор данных и разархивируйте файл данных в папке, где у вас есть разрешение записи. Этот пример принимает, что вы загружаете данные во временной директории, определяемой как tempdir
в MATLAB®. Если вы принимаете решение использовать другую папку, замените той папкой tempdir
в следующем. Записи хранятся как .wav файлы в папках, названных по имени их соответствующего состояния.
url = 'https://www.mathworks.com/supportfiles/audio/AirCompressorDataset/AirCompressorDataset.zip'; downloadFolder = fullfile(tempdir,'AirCompressorDataSet'); if ~exist(fullfile(tempdir,'AirCompressorDataSet'),'dir') loc = websave(downloadFolder,url); unzip(loc,fullfile(tempdir,'AirCompressorDataSet')) end
Используйте audioDatastore
управлять доступом к данным. Каждая подпапка содержит только записи обозначенного класса. Используйте имена папок в качестве меток класса.
datasetLocation = fullfile(tempdir,'AirCompressorDataSet','AirCompressorDataset'); ads = audioDatastore(datasetLocation,'IncludeSubfolders',true,... 'LabelSource','foldernames');
Исследуйте количество примеров в каждом классе. Существует 225 записей в каждом классе для в общей сложности 1 800 записей.
countcats(ads.Labels)
ans = 8×1
225
225
225
225
225
225
225
225
Разделите данные в наборы обучающих данных и наборы тестов. Используйте 80% данных для обучения и протяните остающиеся 20% для тестирования. Переставьте данные, однажды разделяющие.
rng default
ads = shuffle(ads);
[adsTrain,adsTest] = splitEachLabel(ads,0.8,0.2);
Проверьте, что количество примеров в каждом классе является ожидаемым номером.
uniqueLabels = unique(adsTrain.Labels); tblTrain = countEachLabel(adsTrain); tblTest = countEachLabel(adsTest); H = bar(uniqueLabels,[tblTrain.Count, tblTest.Count],'stacked'); legend(H,["Training Set","Test Set"],'Location','NorthEastOutside')
Выберите некоторые случайные примеры из набора обучающих данных для графического вывода. Каждая запись имеет 50 000 выборок, произведенных на уровне 16 кГц.
idx = randperm(numel(adsTrain.Files),8); Fs = 16e3; for n = 1:numel(idx) x = audioread(adsTrain.Files{idx(n)}); t = (0:size(x,1)-1)/Fs; subplot(4,2,n); plot(t,x); if n == 7 || n == 8 xlabel('Seconds'); end title(string(adsTrain.Labels(idx(n)))); end
Создайте сеть рассеивания вейвлета на основе характеристик данных. Установите шкалу инвариантности составлять 0,5 секунды.
N = 5e4; Fs = 16e3; IS = 0.5; sn = waveletScattering('SignalLength',N,'SamplingFrequency',Fs,... 'InvarianceScale',0.5);
С этими сетевыми настройками существует 330 рассеивающихся путей и 25 раз окна на пример. Вы видите это со следующим кодом.
[~,numpaths] = paths(sn); Ncfs = numCoefficients(sn); sum(numpaths)
ans = 330
Ncfs
Ncfs = 25
Обратите внимание, что это уже представляет 6-кратное сокращение размера данных для каждой записи. Мы уменьшали размер данных с 50 000 выборок до 8 250. В этом случае мы получаем дальнейшее сокращение подвыборкой окна времени на коэффициент 6, уменьшая размер каждого примера на коэффициент 25.
Получите функции рассеивания вейвлета наборов обучающих данных и наборов тестов. Если у вас есть подходящий графический процессор и Parallel Computing Toolbox, можно установить useGPU
к true
ускорять рассеивающееся преобразование. Функциональный helperBatchScatFeatures
получает рассеивающееся преобразование и подпроизводит окна времени на коэффициент 6.
batchsize = 64; useGPU = false; scTrain = []; while hasdata(adsTrain) sc = helperBatchScatFeatures(adsTrain,sn,N,batchsize,useGPU); scTrain = cat(3,scTrain,sc); end
Повторите процесс для протянутого набора тестов.
scTest = []; while hasdata(adsTest) sc = helperBatchScatFeatures(adsTest,sn,N,batchsize,useGPU); scTest = cat(3,scTest,sc); end scTest = gather(scTest);
Определите количество рассеивающихся путей и количество окон времени. Существует 330 рассеивающихся путей и пять раз окна.
[Npaths,numTimeWindows,~] = size(scTrain);
Измените рассеивающиеся преобразования для наборов обучающих данных и наборов тестов так, чтобы каждая строка была окном времени через все 330 рассеивающихся путей. Поскольку существует пять раз, когда окна в этом подпроизведенном рассеивании вейвлета преобразовывают, размер TrainFeatures
5*1440-by-330. Точно так же размер TestFeatures
5*360-by-330.
TrainFeatures = permute(scTrain,[2 3 1]); TrainFeatures = reshape(TrainFeatures,[],Npaths,1); TestFeatures = permute(scTest,[2 3 1]); TestFeatures = reshape(TestFeatures,[],Npaths,1);
Для того, чтобы подбирать нашу модель SVM, реплицируйте метки так, чтобы была метка для каждой строки TrainFeatures
.
trainLabels = adsTrain.Labels; numTrainSignals = numel(trainLabels); trainLabels = repmat(trainLabels,1,numTimeWindows); trainLabels = reshape(trainLabels',numTrainSignals*numTimeWindows,1);
Используйте классификатор машины опорных векторов (SVM) мультикласса с ядром кубического полинома. Соответствуйте SVM к обучающим данным.
template = templateSVM(... 'KernelFunction', 'polynomial', ... 'PolynomialOrder', 3, ... 'KernelScale', 'auto', ... 'BoxConstraint', 1, ... 'Standardize', true); classificationSVM = fitcecoc(... TrainFeatures, ... trainLabels, ... 'Learners', template, ... 'Coding', 'onevsall','ClassNames',uniqueLabels);
Определите точность перекрестной проверки модели SVM с помощью 5-кратной перекрестной проверки.
partitionedModel = crossval(classificationSVM, 'KFold', 5);
validationAccuracy = (1 - kfoldLoss(partitionedModel))*100
validationAccuracy = single
99.5278
Точность перекрестной проверки составляет более чем 99%, когда каждый раз окно отдельно классифицируется. Эта эффективность превосходна, но фактическая точность перекрестной проверки на самом деле выше. Поскольку существует пять окон на пример, мы должны использовать все пять раз окна, чтобы присвоить метку класса. Существуют многочисленные способы сделать это, но здесь мы используем голосование простого большинства по этим пяти окнам. Если нет никакого большинства, мы присваиваем класс "NoUniqueMode" и полагаем что ошибка классификации.
Получите предсказания класса из модели перекрестной проверки и исследуйте точность с помощью голосования простого большинства.
validationPredictions = kfoldPredict(partitionedModel); TrainVotes = helperMajorityVote(validationPredictions,adsTrain.Labels,uniqueLabels); crossvalAccuracy = sum(TrainVotes == adsTrain.Labels)/numel(adsTrain.Labels)*100; crossvalAccuracy
crossvalAccuracy = 100
Точность перекрестной проверки составляет близкие 100%. Примените модель к протянутому набору тестов. Используйте голосование простого большинства, чтобы присвоить метку класса.
predLabels = predict(classificationSVM,TestFeatures); [TestVotes,TestCounts] = helperMajorityVote(predLabels,adsTest.Labels,uniqueLabels); testAccuracy = sum(TestVotes == adsTest.Labels)/numel(adsTest.Labels)*100; testAccuracy
testAccuracy = 100
Точность на протянутом наборе тестов составляет близкие 100%, указывающих, что наша модель делает вывод хорошо к невидимым данным. Постройте график беспорядка. Обратите внимание на то, что каждый тестовый пример произвел класс большинства (уникальный режим).
figure confusionchart(adsTest.Labels,TestVotes)
В этом примере мы использовали рассеивание вейвлета, преобразовывают с SVM, чтобы классифицировать отказы на воздушный компрессор. Рассеивающееся преобразование обеспечило устойчивый набор функций, которыми SVM смог достигнуть превосходной перекрестной проверки и проведения испытаний.
[1] Verma, Нищел К., Рахул Кумар Севэкула, Сонэл Диксит и Аль Сэлур. “Интеллектуальный основанный на условии Контроль Используя Акустические Сигналы для Воздушных Компрессоров”. Транзакции IEEE на Надежности 65, № 1 (март 2016): 291–309. https://doi.org/10.1109/TR.2015.2459684.
helperBatchScatFeatures - Эта функция возвращает время вейвлета, рассеивая матрицу функции для данного входного сигнала. Рассеивающиеся функции подпроизводятся на коэффициент 6. Если useGPU
установлен в true
, рассеивающееся преобразование вычисляется на графическом процессоре.
function sc = helperBatchScatFeatures(ds,sn,N,batchsize,useGPU) % This function is only intended to support examples in the Wavelet % Toolbox. It may be changed or removed in a future release. % Read batch of data from audio datastore batch = helperReadBatch(ds,N,batchsize); if useGPU batch = gpuArray(batch); end % Obtain scattering features S = sn.featureMatrix(batch,'transform','log'); gather(batch); S = gather(S); % Subsample the features sc = S(:,1:6:end,:); end
helperReadBatch - Эта функция читает пакеты заданного размера от datastore и возвращает выходной параметр в одинарной точности. Каждый столбец выхода является отдельным сигналом от datastore. Выход может иметь меньше столбцов, чем batchsize, если datastore не имеет достаточного количества записей.
function batchout = helperReadBatch(ds,N,batchsize) % This function is only in support of Wavelet Toolbox examples. It may % change or be removed in a future release. % % batchout = readReadBatch(ds,N,batchsize) where ds is the Datastore and % ds is the Datastore % batchsize is the batchsize kk = 1; while(hasdata(ds)) && kk <= batchsize tmpRead = read(ds); batchout(:,kk) = cast(tmpRead(1:N),'single'); %#ok<AGROW> kk = kk+1; end end
helperMajorityVote - Эта функция получает решение большинством голосов для меток класса. Если нет никакого большинства, "NoUniqueMode" возвращен и обработан как ошибка.
function [ClassVotes,ClassCounts] = helperMajorityVote(predLabels,origLabels,classes) % This function is in support of wavelet scattering examples only. It may % change or be removed in a future release. % Make categorical arrays if the labels are not already categorical predLabels = categorical(predLabels); origLabels = categorical(origLabels); % Expects both predLabels and origLabels to be categorical vectors Npred = numel(predLabels); Norig = numel(origLabels); Nwin = Npred/Norig; predLabels = reshape(predLabels,Nwin,Norig); ClassCounts = countcats(predLabels); [mxcount,idx] = max(ClassCounts); ClassVotes = classes(idx); tmpsum = sum(ClassCounts == mxcount); ClassVotes(tmpsum > 1) = categorical({'NoUniqueMode'}); end