Классификация гиперспектральных изображений с помощью глубокого обучения

В этом примере используются:

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

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

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

Загрузить набор гиперспектральных данных

В этом примере используется набор данных Indian Pines, включенный в библиотеку гиперспектральных изображений Toolbox™ обработки изображений. Набор данных состоит из одного гиперспектрального изображения размером 145 на 145 пикселей с 220 цветовыми каналами. Набор данных также содержит изображение этикетки с 16 классами, такими как «Альфальфа», «Кукуруза», «Травяные пастбища», «Травяные деревья» и «Каменные стальные башни».

Чтение гиперспектрального изображения с помощью hypercube функция.

hcube = hypercube('indian_pines.dat');

Визуализация ложной цветовой версии изображения с помощью colorize функция.

rgbImg = colorize(hcube,'method','rgb');
imshow(rgbImg)

Загрузите метки истинности земли и укажите количество классов.

gtLabel = load('indian_pines_gt.mat');
gtLabel = gtLabel.indian_pines_gt;
numClasses = 16;

Предварительные данные обучения

Уменьшите число спектральных диапазонов до 30 с помощью hyperpca функция. Эта функция выполняет анализ основных компонентов (PCA) и выбирает спектральные полосы с наиболее уникальными сигнатурами.

dimReduction = 30;
imageData = hyperpca(hcube,dimReduction);

Нормализуйте данные изображения.

sd = std(imageData,[],3);
imageData = imageData./sd;

Разбейте гиперспектральное изображение на фрагменты размером 25 на 25 пикселей с 30 каналами с помощью createImagePatchesFromHypercube функция помощника. Эта функция присоединена к примеру как вспомогательный файл. Функция также возвращает одну метку для каждого фрагмента, который является меткой центрального пикселя.

windowSize = 25;
inputSize = [windowSize windowSize dimReduction];
[allPatches,allLabels] = createImagePatchesFromHypercube(imageData,gtLabel,windowSize);

indianPineDataTransposed = permute(allPatches,[2 3 4 1]);
dsAllPatches = augmentedImageDatastore(inputSize,indianPineDataTransposed,allLabels);

Не все кубы в этом наборе данных имеют метки. Однако для обучения сети требуются помеченные данные. Выберите только помеченные кубы для обучения. Подсчитайте количество доступных помеченных исправлений.

patchesLabeled = allPatches(allLabels>0,:,:,:);
patchLabels = allLabels(allLabels>0);
numCubes = size(patchesLabeled,1);

Преобразование числовых меток в категориальные.

patchLabels = categorical(patchLabels);

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

[trainingIndex,validationIndex,testInd] = dividerand(numCubes,0.3,0.7,0);
dataInputTrain = patchesLabeled(trainingIndex,:,:,:);
dataLabelTrain = patchLabels(trainingIndex,1);

Транспонируйте входные данные.

dataInputTransposeTrain = permute(dataInputTrain,[2 3 4 1]); 
dataInputVal = patchesLabeled(validationIndex,:,:,:);
dataLabelVal = patchLabels(validationIndex,1);
dataInputTransposeVal = permute(dataInputVal,[2 3 4 1]);

Создание хранилищ данных, считывающих пакеты обучающих и тестовых данных.

imdsTrain = augmentedImageDatastore(inputSize,dataInputTransposeTrain,dataLabelTrain);
imdsTest = augmentedImageDatastore(inputSize,dataInputTransposeVal,dataLabelVal);

Создание сети классификации CSCNN

Определите архитектуру CSCNN.

layers = [
    image3dInputLayer(inputSize,'Name','Input','Normalization','None')
    convolution3dLayer([3 3 7],8,'Name','conv3d_1')
    reluLayer('Name','Relu_1')
    convolution3dLayer([3 3 5],16,'Name','conv3d_2')
    reluLayer('Name','Relu_2')
    convolution3dLayer([3 3 3],32,'Name','conv3d_3')
    reluLayer('Name','Relu_3')
    convolution3dLayer([3 3 1],8,'Name','conv3d_4')
    reluLayer('Name','Relu_4')
    fullyConnectedLayer(256,'Name','fc1')
    reluLayer('Name','Relu_5')
    dropoutLayer(0.4,'Name','drop_1')
    fullyConnectedLayer(128,'Name','fc2')
    dropoutLayer(0.4,'Name','drop_2')
    fullyConnectedLayer(numClasses,'Name','fc3')
    softmaxLayer('Name','softmax')
    classificationLayer('Name','output')];
lgraph = layerGraph(layers);

Визуализация сети с помощью Deep Network Designer.

deepNetworkDesigner(lgraph)

Укажите параметры обучения

Укажите требуемые параметры сети. В этом примере обучайте сеть в течение 100 периодов с начальной скоростью обучения 0,001, размером партии 256 и оптимизацией Адама.

numEpochs = 100;
miniBatchSize = 256;
initLearningRate = 0.001;
momentum = 0.9;
learningRateFactor = 0.01;

options = trainingOptions('adam', ...
    'InitialLearnRate',initLearningRate, ...
    'LearnRateSchedule','piecewise', ...
    'LearnRateDropPeriod',30, ...
    'LearnRateDropFactor',learningRateFactor, ...
    'MaxEpochs',numEpochs, ...
    'MiniBatchSize',miniBatchSize, ...
    'GradientThresholdMethod','l2norm', ...
    'GradientThreshold',0.01, ...
    'VerboseFrequency',100, ...
    'ValidationData',imdsTest, ...
    'ValidationFrequency',100);

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

По умолчанию в примере загружается предварительно подготовленный классификатор для набора данных Indian Pines с помощью downloadTrainedIndianPinesCSCNN функция помощника. Эта функция присоединена к примеру как вспомогательный файл. Предварительно обученная сеть позволяет классифицировать набор данных Indian Pines, не дожидаясь завершения обучения.

Для обучения сети установите doTraining переменная в следующем коде true. Если вы хотите обучить сеть, настоятельно рекомендуется использовать графический процессор CUDA с поддержкой NVIDIA™. Для использования графического процессора требуется Toolbox™ параллельных вычислений. Дополнительные сведения о поддерживаемых GPU-устройствах см. в разделе Поддержка GPU по выпуску (Parallel Computing Toolbox).

doTraining = false;
if doTraining
    modelDateTime = string(datetime('now','Format',"yyyy-MM-dd-HH-mm-ss"));
    net = trainNetwork(imdsTrain,lgraph,options);
    save(strcat("IndianPinesCSCNN-",modelDateTime,"-Epoch-",num2str(numEpochs),".mat"),'net');
else
    dataFolder = fullfile(tempdir,"indianPines");
    trainedHyperspectralCSCNN_url = 'https://ssd.mathworks.com/supportfiles/image/data/trainedIndianPinesCSCNN.mat';
    downloadTrainedIndianPinesCSCNN(trainedHyperspectralCSCNN_url,dataFolder);
    load(fullfile(dataFolder,'trainedIndianPinesCSCNN.mat'));
end

Pretrained CSCNN network already exists.

Классифицировать гиперспектральное изображение с помощью обученного CSCNN

Вычислите точность классификации для набора тестовых данных. Здесь точность - это доля правильной классификации пикселей по всем классам.

predictionTest = classify(net,imdsTest);
accuracy = sum(predictionTest == dataLabelVal)/numel(dataLabelVal);
disp(['Accuracy of the test data = ', num2str(accuracy)])

Accuracy of the test data = 0.99359

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

prediction = classify(net,dsAllPatches);
prediction = double(prediction);

Сеть обучается только на помеченных патчах. Поэтому предсказанная классификация немеченых пикселей бессмысленна. Найдите немаркированные исправления и установите для метки значение 0.

patchesUnlabeled = find(allLabels==0);
prediction(patchesUnlabeled) = 0;

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

[m,n,d] = size(imageData);
indianPinesPrediction = reshape(prediction,[n m]);
indianPinesPrediction = indianPinesPrediction';

Отображение истинной и предсказанной классификации.

cmap = parula(numClasses);

figure
tiledlayout(1,2,"TileSpacing","Tight")
nexttile
imshow(gtLabel,cmap)
title("Ground Truth Classification")

nexttile
imshow(indianPinesPrediction,cmap)
colorbar
title("Predicted Classification")

Чтобы выделить неправильно классифицированные пикселы, отобразите составное изображение истинности земли и прогнозируемых меток. Серые пикселы обозначают идентичные метки, а цветные - разные метки.

figure
imshowpair(gtLabel,indianPinesPrediction)

См. также

Связанные темы

Семантическая сегментация мультиспектральных изображений с использованием глубокого обучения

Документация