深層学習を使用したハイパースペクトルイメージの分類

この例では次を使用します。

この例では、分類用のカスタムスペクトル畳み込みニューラルネットワーク (CSCNN) を使用してハイパースペクトルイメージを分類する方法を説明します。

この例では、Image Processing Toolbox™ Hyperspectral Imaging Library が必要となります。Image Processing Toolbox Hyperspectral Imaging Library はアドオンエクスプローラーからインストールできます。アドオンのインストールの詳細については、アドオンの取得と管理を参照してください。Image Processing Toolbox Hyperspectral Imaging Library は MATLAB® Online™ および MATLAB® Mobile™ ではサポートされないため、デスクトップの MATLAB® が必要となります。

ハイパースペクトルイメージングでは、可視スペクトルを含む紫外線から遠赤外線までの範囲のさまざまな波長でオブジェクトの空間特徴およびスペクトル特徴を測定します。可視スペクトルの赤、緑、青の部分を感知する 3 種類のセンサーのみを使用するカラーイメージングと異なり、ハイパースペクトルイメージは多数のチャネルを含むことができます。このため、ハイパースペクトルイメージでは、RGB イメージでは同一に見えるオブジェクトの違いを有効にすることができます。

この例では、材質ごとに異なるスペクトルシグネチャに基づいて 16 種類の植生と地形を分類する学習を行う CSCNN を使用します。この例では、CSCNN を学習させる方法を説明しますが、分類の実行に使用できる事前学習済みのネットワークも提供します。

ハイパースペクトルデータセットの読み込み

この例では、Image Processing Toolbox™ Hyperspectral Imaging Library に含まれているインディアンパインデータセットを使用します。このデータセットは、サイズ 145×145 ピクセル、カラーチャネル数 220 のハイパースペクトルイメージ 1 つで構成されます。また、このデータセットには、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;

学習データの前処理

関数hyperpcaを使用してスペクトルバンドの数を 30 まで減らします。この関数は、主成分分析 (PCA) を実行し、最もユニークなシグネチャを持つスペクトルバンドを選択します。

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

イメージデータを正規化します。

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

補助関数 createImagePatchesFromHypercube を使用してハイパースペクトルイメージを、サイズ 25×25 ピクセル、チャネル数 30 のパッチに分割します。この関数は、この例にサポートファイルとして添付されています。また、この関数は、パッチごとにラベルを 1 つ返します。これは、中心ピクセルのラベルです。

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);

数値ラベルを categorical に変換します。

patchLabels = categorical(patchLabels);

パッチをランダムに学習用とテスト用のデータセットに分割します。

[trainingIdx,valIdx,testIdx] = dividerand(numCubes,0.3,0,0.7);
dataInputTrain = patchesLabeled(trainingIdx,:,:,:);
dataLabelTrain = patchLabels(trainingIdx,1);
dataInputTest = patchesLabeled(testIdx,:,:,:);
dataLabelTest = patchLabels(testIdx,1);

入力データを転置します。

dataInputTransposeTrain = permute(dataInputTrain,[2 3 4 1]); 
dataInputTransposeTest = permute(dataInputTest,[2 3 4 1]);

学習データとテストデータのバッチを読み取るデータストアを作成します。

dsTrain = augmentedImageDatastore(inputSize,dataInputTransposeTrain,dataLabelTrain);
dsTest = augmentedImageDatastore(inputSize,dataInputTransposeTest,dataLabelTest);

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);

ディープネットワークデザイナーを使用してネットワークを可視化します。

deepNetworkDesigner(lgraph)

学習オプションの指定

必要なネットワークパラメーターを指定します。この例では、初期学習率を 0.001、バッチサイズを 256 とし、Adam 最適化を使用してネットワークに 100 エポック学習させます。

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=dsTest, ...
    ValidationFrequency=100);

ネットワークの学習

この例では既定で、インディアンパインデータセット用の事前学習済み分類器がダウンロードされます。この事前学習済みのネットワークを使用することで、学習の完了を待たずにインディアンパインデータセットを分類できます。

ネットワークに学習させるには、次のコードで変数 doTraining を true に設定します。自分でネットワークに学習させる場合は、CUDA 対応 NVIDIA™ GPU の使用を強くお勧めします。GPU を使用するには Parallel Computing Toolbox™ が必要です。サポートされる GPU デバイスの詳細については、GPU 計算の要件 (Parallel Computing Toolbox)を参照してください。

doTraining = false;
if doTraining
    net = trainNetwork(dsTrain,lgraph,options);
    modelDateTime = string(datetime("now",Format="yyyy-MM-dd-HH-mm-ss"));
    save("trainedIndianPinesCSCNN-"+modelDateTime+".mat","net");
else
    dataDir = pwd;
    trainedNetwork_url = "https://ssd.mathworks.com/supportfiles/image/data/trainedIndianPinesCSCNN.mat";
    downloadTrainedNetwork(trainedNetwork_url,pwd);
    load(fullfile(dataDir,"trainedIndianPinesCSCNN.mat"));
end

Downloading pretrained network.
This can take several minutes to download...
Done.

学習済み CSCNN を使用したハイパースペクトルイメージの分類

テストデータセットに対する分類の精度を計算します。ここでの精度は、すべてのクラスにおいて正しいピクセル分類が行われる割合です。

predictionTest = classify(net,dsTest);
accuracy = sum(predictionTest == dataLabelTest)/numel(dataLabelTest);
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)

参考

深層学習を使用したハイパースペクトル イメージの分類

ハイパースペクトル データセットの読み込み