Вычисление скалограммы с помощью графического процессора

Наиболее эффективным способом вычисления скалограмм на GPU является использование cwtfilterbank. Для вычисления скалограмм на GPU выполняются следующие шаги:

При первом использовании WT способ, cwtfilterbank кэширует вейвлет-фильтры на GPU. В результате достигается существенная экономия времени при вычислениях при получении скалограмм множества сигналов с использованием одного и того же банка фильтров и одного и того же типа данных. Ниже приведен пример рекомендуемого рабочего процесса. В качестве примера можно использовать образец гитарной музыки, содержащий 661 500 образцов.

[y,fs] = audioread('guitartune.wav');
plot(y)
grid on

Поскольку большинство графических процессоров NVIDIA значительно эффективнее при использовании данных с одной, а не двойной точностью, приведение сигнала к одной точности.

y = single(y);

Создайте набор фильтров в соответствии с длиной сигнала и частотой дискретизации. Для глубокого обучения частота дискретизации обычно не требуется и поэтому может быть исключена.

fb = cwtfilterbank('SignalLength',length(y),'SamplingFrequency',fs);

Наконец, переместите сигнал в графический процессор с помощью gpuArray и вычисляют CWT данных. Постройте график результирующей скалограммы.

[cfs,f] = fb.wt(gpuArray(y));
t = 0:1/fs:(length(y)*1/fs)-1/fs;
imagesc(t,f,abs(cfs))
axis xy
ylabel('Hz')
xlabel('Seconds')

Использовать gather для возврата коэффициентов CWT и любых других выходных сигналов в CPU.

cfs = gather(cfs);
f = gather(f);

Чтобы продемонстрировать эффективность, полученную при использовании GPU, время вычисления CWT на GPU и CPU для одного и того же сигнала. Время вычислений GPU, сообщенное здесь, получают с использованием NVIDIA Titan V с вычислительной способностью 7,0.

ygpu = gpuArray(y);
fgpu = @()fb.wt(ygpu);
Tgpu = gputimeit(fgpu)

Tgpu = 0.2658

Повторите то же самое измерение на CPU и проверьте отношение GPU к времени CPU, чтобы увидеть сокращение времени вычислений.

fcpu = @()fb.wt(y);
Tcpu = timeit(fcpu)

Tcpu = 3.7088

Tcpu/Tgpu

ans = 13.9533

Скалограммы в глубоком обучении

Обычным применением CWT в глубоком обучении является использование скалограммы сигнала в качестве входного «изображения» для глубокого CNN. Это обязательно требует вычисления нескольких скалограмм, по одной для каждого сигнала в обучающих, валидационных и тестовых наборах. Хотя GPU часто используются для ускорения обучения и вывода в глубокой сети, также полезно использовать GPU для ускорения извлечения любых функций или предварительной обработки данных, необходимых для повышения надежности глубокой сети.

Для иллюстрации этого следующий раздел применяет вейвлет-скалограммы к классификации электрокардиограмм человека (ЭКГ). Скалограммы используются с теми же данными, обработанными во временных рядах классификации с использованием вейвлет-анализа и глубокого обучения. В этом примере обучение передаче с помощью GoogLeNet и SqueeEcNet использовалось для классификации форм сигналов ЭКГ на одну из трех категорий. Описание данных и способ их получения повторяются здесь для удобства.

Описание и загрузка данных ЭКГ

Данные ЭКГ получены от трех групп людей: лиц с сердечной аритмией (ARR), лиц с застойной сердечной недостаточностью (CHF) и лиц с нормальными синусовыми ритмами (СМП). Всего существует 162 записей ЭКГ из трех баз данных PhysioNet: база данных аритмии MIT-BIH [2] [3], база данных нормального синусового ритма MIT-BIH [3] и база данных застойной сердечной недостаточности BIDMC [1] [3]. Более конкретно, 96 записей от людей с аритмией, 30 записей от людей с застойной сердечной недостаточностью и 36 записей от людей с нормальными синусовыми ритмами. Цель состоит в том, чтобы обучить модель различать ARR, CHF и NSR.

Эти данные можно получить из репозитория MathWorks GitHub. Чтобы загрузить данные с веб-сайта, нажмите Code и выбрать Download ZIP. Сохранить файл physionet_ECG_data-main.zip в папке, в которой имеется разрешение на запись. В инструкциях для этого примера предполагается, что файл загружен во временный каталог. tempdir, в MATLAB. Измените последующие инструкции по распаковке и загрузке данных, если вы решили загрузить данные в папку, отличную от tempdir.

После загрузки данных из GitHub распакуйте файл во временном каталоге.

unzip(fullfile(tempdir,'physionet_ECG_data-main.zip'),tempdir)

Распаковка создает папку physionet-ECG_data-main во временном каталоге. Эта папка содержит текстовый файл README.md и ECGData.zip. ECGData.zip файл содержит

ECGData.mat содержит данные, используемые в этом примере. Текстовый файл, Modified_physionet_data.txt, требуется политикой копирования PhysioNet и предоставляет исходные атрибуты для данных, а также описание этапов предварительной обработки, применяемых к каждой записи ЭКГ.

Расстегнуть молнию ECGData.zip в physionet-ECG_data-main. Загрузите файл данных в рабочую область MATLAB.

unzip(fullfile(tempdir,'physionet_ECG_data-main','ECGData.zip'),...
    fullfile(tempdir,'physionet_ECG_data-main'))
load(fullfile(tempdir,'physionet_ECG_data-main','ECGData.mat'))

ECGData - структурный массив с двумя полями: Data и Labels. Data поле представляет собой матрицу 162 на 65536, где каждая строка представляет собой запись ЭКГ, дискретизированную при 128 герц. Labels представляет собой массив диагностических меток типа 162 на 1, по одной для каждой строки Data. Тремя диагностическими категориями являются: 'ARR', 'CHF', и 'NSR'. Используйте функцию помощника, helperRandomSplit, чтобы разделить данные на наборы обучения и проверки с 80% данных, выделенных для обучения, и 20% для проверки. Преобразуйте диагностические метки ЭКГ в категории.

[trainData, validationData, trainLabels, validationLabels] = helperRandomSplit(80,ECGData);
trainLabels = categorical(trainLabels);
validationLabels = categorical(validationLabels);

Есть 130 записей в trainData набор и 32 записи при проверкеData. По конструкции данные обучения содержат 80,25% (130/162) данных. Напомним, что класс ARR представляет 59,26% данных (96/162), класс CHF - 18,52% (30/162), а класс NSR - 22,22% (36/162). Проверьте процент каждого класса в обучающих и тестовых наборах. Проценты в каждом из них согласуются с общими процентами классов в наборе данных.

Ctrain = countcats(trainLabels)./numel(trainLabels).*100

Ctrain = 3×1

   59.2308
   18.4615
   22.3077

Cvalid = countcats(validationLabels)./numel(validationLabels).*100

Cvalid = 3×1

   59.3750
   18.7500
   21.8750

Скалограммы с глубоким CNN - данные ЭКГ

Вычисление скалограммы на GPU

Вычислите скалограммы как для обучающих, так и для проверочных наборов. Набор useGPU кому true для использования графического процессора и false для вычисления скалограмм на CPU. Чтобы смягчить влияние больших входных матриц на CNN и создать больше примеров обучения и проверки, helperECGScalograms разбивает каждый сигнал ЭКГ на четыре неперекрывающихся сегмента по 16384 выборки каждый и вычисляет скалограммы для всех четырех сегментов. Реплицируйте метки в соответствии с расширенным набором данных. В этом случае получают оценку затраченного времени вычисления.

frameLength = 16384;
useGPU = true;
tic;
Xtrain = helperECGScalograms(trainData,frameLength,useGPU);

Computing scalograms...
Processed 50 files out of 130
Processed 100 files out of 130
...done

T = toc;
sprintf('Elapsed time is %1.2f seconds',T)

ans = 
'Elapsed time is 4.22 seconds'

trainLabels = repelem(trainLabels,4);

С помощью графического процессора Titan V было вычислено 502 скалограммы приблизительно за 4,2 секунды. Настройка useGPU кому false и повторение вышеприведенного вычисления демонстрирует скорость, полученную при использовании графического процессора. В этом случае использование ЦП потребовало 33,3 секунды для вычисления скалограмм. Вычисление GPU выполнялось более чем в 7 раз быстрее.

Повторите тот же процесс для данных проверки.

useGPU = true;
Xvalid = helperECGScalograms(validationData,frameLength,useGPU);

Computing scalograms...
...done

validationLabels = repelem(validationLabels,4);

Затем настройте глубокий CNN для обработки наборов обучения и проверки. Используемая здесь простая сеть не оптимизирована. Этот CNN используется только для иллюстрации сквозного рабочего процесса для случаев, когда скалограммы помещаются в память.

sz = size(Xtrain);
specSize = sz(1:2);
imageSize = [specSize 1];
dropoutProb = 0.3;

layers = [
    imageInputLayer(imageSize)
    
    convolution2dLayer(3,12,'Padding','same')
    batchNormalizationLayer
    reluLayer
    
    maxPooling2dLayer(2,'Stride',2)
    
   
    convolution2dLayer(3,20,'Padding','same')
    batchNormalizationLayer
    reluLayer
    
    maxPooling2dLayer(2,'Stride',2)
    
   
    
    convolution2dLayer(3,32,'Padding','same')
    batchNormalizationLayer
    reluLayer
    dropoutLayer(dropoutProb)
    fullyConnectedLayer(3)
    softmaxLayer
    classificationLayer];

Используйте следующие варианты обучения.

options = trainingOptions('sgdm',...
    'InitialLearnRate', 1e-4,...
    'LearnRateDropPeriod',18,...
    'MiniBatchSize', 20,...
    'MaxEpochs',25,...
    'L2Regularization',1e-1,...
    'Plots', 'training-progress',...
    'Verbose',false,...
    'Shuffle','every-epoch',...
    'ExecutionEnvironment','auto',...
    'ValidationData',{Xvalid,validationLabels});

Обучите сеть и измерьте ошибку проверки.

trainNetwork(Xtrain,trainLabels,layers,options);

Несмотря на то, что используемый здесь простой CNN не оптимизирован, точность проверки постоянно находится в диапазоне от 80 до 90 процентов. Это сравнимо с точностью проверки, достигнутой с более мощным и оптимизированным SqueeEcNet, показанным в примере классификации временных рядов с использованием вейвлет-анализа и глубокого обучения. Кроме того, это гораздо более эффективное использование скалограммы, потому что в этом примере скалограммы должны были быть масштабированы как изображения RGB, совместимые с SqueeEnet, сохранены на диск в соответствующем формате изображения, а затем переданы в глубокую сеть с помощью imageDatastore.

Распознавание речевых цифр - вычисления графического процессора с использованием хранилища данных Transform

В этом разделе показано, как ускорить вычисление скалограммы с помощью графического процессора в преобразованном рабочем процессе хранилища данных.

Данные

Клонируйте или загрузите набор данных Free Spoken Digit Dataset (FSDD), доступный по адресу https://github.com/Jakobovski/free-spoken-digit-dataset. FSDD - это открытый набор данных, который со временем может расти. В этом примере используется версия, зафиксированная на 01/29/2019, которая состоит из 2000 записей английских цифр от 0 до 9, полученных из четырех динамиков. Два из говорящих в этой версии являются носителями американского английского языка и два говорящих - неосновными носителями английского языка с бельгийским французским и немецким акцентом. Данные дискретизируются на частоте 8000 Гц.

Другие подходы к этому набору данных, включая вейвлет-рассеяние, см. в разделе Распознавание разговорных цифр с вейвлет-рассеянием и глубоким обучением.

Использовать audioDatastore управление доступом к данным и обеспечение произвольного разделения записей на обучающие и тестовые наборы. Установите location свойство для расположения папки записей FSDD на компьютере, например:

pathToRecordingsFolder = '/home/user/free-spoken-digit-dataset/recordings';
location = pathToRecordingsFolder;

Пункт audioDatastore в это место.

ads = audioDatastore(location);

Вспомогательная функция helpergenLabels создает категориальный массив меток из файлов FSDD. Перечислите классы и количество примеров в каждом классе.

ads.Labels = helpergenLabels(ads);
summary(ads.Labels)

     0      200 
     1      200 
     2      200 
     3      200 
     4      200 
     5      200 
     6      200 
     7      200 
     8      200 
     9      200 

Преобразованное хранилище данных

Сначала разбейте FSDD на учебные и тестовые наборы. Выделите 80% данных обучающему набору и сохраните 20% для тестового набора.

rng default;
ads = shuffle(ads);
[adsTrain,adsTest] = splitEachLabel(ads,0.8);
countEachLabel(adsTrain)

ans=10×2 table
    Label    Count
    _____    _____

      0       160 
      1       160 
      2       160 
      3       160 
      4       160 
      5       160 
      6       160 
      7       160 
      8       160 
      9       160 

Затем создайте банк фильтров CWT и преобразованные хранилища данных как для обучения, так и для тестирования с помощью функции помощника. helperDigitScalogram. Преобразованное хранилище данных преобразует каждую запись в сигнал длиной 8192, вычисляет скалограмму на GPU и собирает данные обратно в CPU.

reset(gpuDevice(1))
fb = cwtfilterbank('SignalLength',8192);
adsSCTrain = transform(adsTrain,@(audio,info)helperDigitScalogram(audio,info,fb),'IncludeInfo',true);
adsSCTest = transform(adsTest,@(audio,info)helperDigitScalogram(audio,info,fb),'IncludeInfo',true);

Глубокий CNN

Создайте глубокий CNN для обучения с преобразованным хранилищем данных, adscTrain. Как и в первом примере, сеть не оптимизирована. Смысл состоит в том, чтобы показать рабочий процесс с использованием скалограмм, вычисленных на GPU для данных, не имеющих памяти.

numClasses = 10;

dropoutProb = 0.2;
numF = 12;
layers = [
    imageInputLayer([101 8192 1])

    convolution2dLayer(5,numF,'Padding','same')
    batchNormalizationLayer
    reluLayer

    maxPooling2dLayer(3,'Stride',2,'Padding','same')

    convolution2dLayer(3,2*numF,'Padding','same')
    batchNormalizationLayer
    reluLayer

    maxPooling2dLayer(3,'Stride',2,'Padding','same')

    convolution2dLayer(3,4*numF,'Padding','same')
    batchNormalizationLayer
    reluLayer

    maxPooling2dLayer(3,'Stride',2,'Padding','same')

    convolution2dLayer(3,4*numF,'Padding','same')
    batchNormalizationLayer
    reluLayer
    convolution2dLayer(3,4*numF,'Padding','same')
    batchNormalizationLayer
    reluLayer

    maxPooling2dLayer(2)

    dropoutLayer(dropoutProb)
    fullyConnectedLayer(numClasses)
    softmaxLayer
    classificationLayer('Classes',categories(ads.Labels));
   ]; 

Настройка параметров обучения для сети.

miniBatchSize = 25;
options = trainingOptions('adam', ...
    'InitialLearnRate',1e-4, ...
    'MaxEpochs',30, ...
    'MiniBatchSize',miniBatchSize, ...
    'Shuffle','every-epoch', ...
    'Plots', 'training-progress',...
    'Verbose',false,...
    'ExecutionEnvironment','gpu');

Обучение сети.

trainedNet = trainNetwork(adsSCTrain,layers,options);

В этом случае тренировка была завершена за 25 минут и 10 секунд. Если вы прокомментируете вызов gpuArray в helperDigitScalogram и использовать CPU для получения скалограмм, время обучения значительно увеличивается. В этом случае наблюдалось увеличение от 25 минут и 10 секунд до 45 минут и 38 секунд.

Используйте обученную сеть для прогнозирования цифровых меток для тестового набора.

Ypredicted = classify(trainedNet,adsSCTest,'ExecutionEnvironment','CPU');
cnnAccuracy = sum(Ypredicted == adsTest.Labels)/numel(Ypredicted)*100

cnnAccuracy = 96.2500

Время вывода с использованием GPU составляло приблизительно 22 секунды. С помощью CPU время вывода удвоилось до 45 секунд.

Производительность обученной сети на тестовых данных близка к 96%. Это сравнимо с показателями распознавания речевых цифр с вейвлет-рассеянием и глубоким обучением.

Резюме

В этом примере показано, как использовать графический процессор для ускорения вычисления скалограммы. В примере представлен оптимальный рабочий процесс для эффективного вычисления скалограмм как для данных в памяти, так и для данных из памяти, считанных с диска с помощью преобразованных хранилищ данных.

Ссылки

function Labels = helpergenLabels(ads)
% This function is only for use in Wavelet Toolbox examples. It may be
% changed or removed in a future release.
tmp = cell(numel(ads.Files),1);
expression = "[0-9]+_";
for nf = 1:numel(ads.Files)
    idx = regexp(ads.Files{nf},expression);
    tmp{nf} = ads.Files{nf}(idx);
end
Labels = categorical(tmp);

end

function X = helperECGScalograms(data,window,useGPU)

disp("Computing scalograms...");
Nsig = size(data,1);
Nsamp = size(data,2);
Nsegment = Nsamp/window;

fb = cwtfilterbank('SignalLength',window,'Voices',10);
Ns = length(fb.Scales);
X = zeros([Ns,window,1,Nsig*Nsegment],'single');
start = 0;
if useGPU
    data = gpuArray(single(data'));
else
    data = single(data');
end
for ii = 1:Nsig
    ts = data(:,ii);
    ts = reshape(ts,window,Nsegment);
    ts = (ts-mean(ts))./max(abs(ts));
    
    for kk = 1:size(ts,2)
        cfs = fb.wt(ts(:,kk));
        X(:,:,1,kk+start) = gather(abs(cfs));
        
    end
    start = start+Nsegment;
    
    if mod(ii,50) == 0
        disp("Processed " + ii + " files out of " + Nsig)
    end
    
end

disp("...done");
data = gather(data);

end

function [x,info] = helperReadSPData(x,info)
% This function is only for use Wavelet Toolbox examples. It may change or
% be removed in a future release.

N = numel(x);
if N > 8192
    x = x(1:8192);
elseif N < 8192
    pad = 8192-N;
    prepad = floor(pad/2);
    postpad = ceil(pad/2);
    x = [zeros(prepad,1) ; x ; zeros(postpad,1)];
end
x = x./max(abs(x));

end

function [dataout,info] = helperDigitScalogram(audioin,info,fb)
audioin = single(audioin);
audioin = gpuArray(audioin);
audioin = helperReadSPData(audioin);
cfs = gather(abs(fb.wt(audioin)));
audioin = gather(audioin);
dataout = {cfs,info.Label};
end