使用深度学习对多光谱图像进行语义分割
此示例说明如何使用 U-Net 对包含七个通道的多光谱图像执行语义分割。
语义分割涉及用类标注图像中的每个像素。语义分割的一个应用是跟踪森林采伐,即森林植被随时间的变化。环保机构需要跟踪森林采伐,以评估和量化某个地区的环境和生态健康。
基于深度学习的语义分割可以通过高分辨率航拍照片精确测量植被覆盖度。一个挑战是区分具有相似视觉特性的类,例如尝试将绿色像素分类为草、灌木或树。为了提高分类准确度,一些数据集包含多光谱图像,这些图像提供关于每个像素的附加信息。例如,Hamlin Beach 国家公园数据集用三个近红外通道为彩色图像补充信息,以提供更清晰的分类。
此示例首先说明如何使用预训练的 U-Net 执行语义分割,然后使用分割结果来计算植被覆盖率。然后,您可以选择使用基于补片的训练方法基于 Hamlin Beach State Parck 数据集训练 U-Net 网络。
下载数据集
此示例使用高分辨率多光谱数据集来训练网络 [1]。图像集是用无人机在纽约 Hamlin Beach 州立公园拍摄的。数据包含已标注的训练集、验证集和测试集,有 18 个对象类标签。数据文件的大小为 3.0 GB。
使用 downloadHamlinBeachMSIData
辅助函数下载数据集的 MAT 文件版本。此函数作为支持文件包含在本示例中。指定 dataDir
作为数据的存放位置。
dataDir = fullfile(tempdir,"rit18_data");
downloadHamlinBeachMSIData(dataDir);
加载该数据集。
load(fullfile(dataDir,"rit18_data.mat")); whos train_data val_data test_data
Name Size Bytes Class Attributes test_data 7x12446x7654 1333663576 uint16 train_data 7x9393x5642 741934284 uint16 val_data 7x8833x6918 855493716 uint16
多光谱图像数据排列为通道数×宽度×高度数组。但是,在 MATLAB® 中,多通道图像排列为宽度×高度×通道数数组。要重构数据以使通道处于第三个维度中,请使用 switchChannelsToThirdPlane
辅助函数。此函数作为支持文件包含在本示例中。
train_data = switchChannelsToThirdPlane(train_data); val_data = switchChannelsToThirdPlane(val_data); test_data = switchChannelsToThirdPlane(test_data);
确认数据结构正确。
whos train_data val_data test_data
Name Size Bytes Class Attributes test_data 12446x7654x7 1333663576 uint16 train_data 9393x5642x7 741934284 uint16 val_data 8833x6918x7 855493716 uint16
将训练数据保存为一个 MAT 文件,将训练标签保存为一个 PNG 文件。这有助于在训练期间使用 imageDatastore
和 pixelLabelDatastore
加载训练数据。
save("train_data.mat","train_data"); imwrite(train_labels,"train_labels.png");
可视化多光谱数据
在此数据集中,RGB 颜色通道是第 3 个、第 2 个和第 1 个图像通道。将训练、验证和测试图像的颜色分量显示为蒙太奇形式。要使图像在屏幕上看起来更亮,请使用 histeq
(Image Processing Toolbox) 函数执行直方图均衡化处理。
figure montage(... {histeq(train_data(:,:,[3 2 1])), ... histeq(val_data(:,:,[3 2 1])), ... histeq(test_data(:,:,[3 2 1]))}, ... BorderSize=10,BackgroundColor="white") title("RGB Component of Training, Validation, and Test Image (Left to Right)")
将训练数据的后三个经过直方图均衡化处理的通道显示为蒙太奇形式。这些通道对应于近红外波段,并基于其热特征突出显示图像的不同分量。例如,第二个通道图像中心附近的树比其他两个通道中的树显示的细节要多。
figure montage(... {histeq(train_data(:,:,4)),histeq(train_data(:,:,5)),histeq(train_data(:,:,6))}, ... BorderSize=10,BackgroundColor="white") title("Training Image IR Channels 1, 2, and 3 (Left to Right)")
通道 7 是指示有效分割区域的掩膜。显示训练、验证和测试图像的掩膜。
figure montage(... {train_data(:,:,7),val_data(:,:,7),test_data(:,:,7)}, ... BorderSize=10,BackgroundColor="white") title("Mask of Training, Validation, and Test Image (Left to Right)")
可视化真实值标签
已标注的图像包含用于分割的真实值数据,每个像素分配给 18 个类之一。获取具有对应 ID 的类的列表。
disp(classes)
0. Other Class/Image Border 1. Road Markings 2. Tree 3. Building 4. Vehicle (Car, Truck, or Bus) 5. Person 6. Lifeguard Chair 7. Picnic Table 8. Black Wood Panel 9. White Wood Panel 10. Orange Landing Pad 11. Water Buoy 12. Rocks 13. Other Vegetation 14. Grass 15. Sand 16. Water (Lake) 17. Water (Pond) 18. Asphalt (Parking Lot/Walkway)
创建一个类名向量。
classNames = [ "RoadMarkings","Tree","Building","Vehicle","Person", ... "LifeguardChair","PicnicTable","BlackWoodPanel",... "WhiteWoodPanel","OrangeLandingPad","Buoy","Rocks",... "LowLevelVegetation","Grass_Lawn","Sand_Beach",... "Water_Lake","Water_Pond","Asphalt"];
在经过直方图均衡化处理的 RGB 训练图像上叠加标签。向图像添加颜色栏。
cmap = jet(numel(classNames)); B = labeloverlay(histeq(train_data(:,:,4:6)),train_labels,Transparency=0.8,Colormap=cmap); figure imshow(B) title("Training Labels") N = numel(classNames); ticks = 1/(N*2):1/N:1; colorbar(TickLabels=cellstr(classNames),Ticks=ticks,TickLength=0,TickLabelInterpreter="none"); colormap(cmap)
执行语义分割
下载预训练的 U-Net 网络。
trainedUnet_url = "https://www.mathworks.com/supportfiles/vision/data/multispectralUnet.mat"; downloadTrainedNetwork(trainedUnet_url,dataDir); load(fullfile(dataDir,"multispectralUnet.mat"));
为了在经过训练的网络上执行语义分割,请对验证数据使用 segmentMultispectralImage
辅助函数。此函数作为支持文件包含在本示例中。segmentMultispectralImage
函数使用 semanticseg
(Computer Vision Toolbox) 函数对图像补片执行分割。需要处理补片,因为图像的大小导致无法一次性处理整个图像。
predictPatchSize = [1024 1024]; segmentedImage = segmentMultispectralImage(val_data,net,predictPatchSize);
为了只提取分割的有效部分,需要将分割的图像乘以验证数据的掩膜通道。
segmentedImage = uint8(val_data(:,:,7)~=0) .* segmentedImage;
figure
imshow(segmentedImage,[])
title("Segmented Image")
语义分割的输出含有噪声。执行图像后处理以去除噪声和杂散像素。使用 medfilt2
(Image Processing Toolbox) 函数从分割中去除椒盐噪声。去除噪声后,可视化分割后的图像。
segmentedImage = medfilt2(segmentedImage,[7,7]);
imshow(segmentedImage,[]);
title("Segmented Image with Noise Removed")
在经过直方图均衡化处理的 RGB 验证图像上叠加分割后的图像。
B = labeloverlay(histeq(val_data(:,:,[3 2 1])),segmentedImage,Transparency=0.8,Colormap=cmap); figure imshow(B) title("Labeled Segmented Image") colorbar(TickLabels=cellstr(classNames),Ticks=ticks,TickLength=0,TickLabelInterpreter="none"); colormap(cmap)
计算植被覆盖率
语义分割结果可用于解决相关生态学问题。例如,植被覆盖的陆地面积百分比是多少?要回答此问题,请找出标注为植被的像素数。标签 ID 2 ("Trees")、13 ("LowLevelVegetation") 和 14 ("Grass_Lawn") 都属于植被类。通过对掩膜图像的关注区域中的像素求和,还可以求得有效像素的总数。
vegetationClassIds = uint8([2,13,14]); vegetationPixels = ismember(segmentedImage(:),vegetationClassIds); validPixels = (segmentedImage~=0); numVegetationPixels = sum(vegetationPixels(:)); numValidPixels = sum(validPixels(:));
通过用植被像素数除以有效像素数来计算植被覆盖率。
percentVegetationCover = (numVegetationPixels/numValidPixels)*100;
fprintf("The percentage of vegetation cover is %3.2f%%.",percentVegetationCover);
The percentage of vegetation cover is 51.72%.
该示例的其余部分说明如何基于 Hamlin Beach 数据集训练 U-Net。
为训练创建随机补片提取数据存储
使用随机补片提取数据存储将训练数据馈送给网络。此数据存储从一个包含真实值图像的图像数据存储和包含像素标签的像素标签数据存储中提取多个对应的随机补片。补片是防止大图像耗尽内存和有效增加可用训练数据量的常用方法。
首先从 imageDatastore
中的 "train_data.mat"
加载训练图像。由于 MAT 文件格式是非标准图像格式,您必须使用 MAT 文件读取器来读取图像数据。您可以使用 MAT 文件读取器辅助函数 matRead6Channels
,它从训练数据中提取前六个通道,省略包含掩膜的最后一个通道。此函数作为支持文件包含在本示例中。
imds = imageDatastore("train_data.mat",FileExtensions=".mat",ReadFcn=@matRead6Channels);
创建一个 pixelLabelDatastore
(Computer Vision Toolbox) 以存储包含 18 个已标注区域的标签补片。
pixelLabelIds = 1:18;
pxds = pixelLabelDatastore("train_labels.png",classNames,pixelLabelIds);
从图像数据存储和像素标签数据存储创建一个 randomPatchExtractionDatastore
(Image Processing Toolbox)。每个小批量包含 16 个大小为 256×256 像素的补片。在轮的每次迭代中提取一千个小批量。
dsTrain = randomPatchExtractionDatastore(imds,pxds,[256,256],PatchesPerImage=16000);
随机补片提取数据存储 dsTrain
在一轮训练的每次迭代中向网络提供小批量数据。预览数据存储以探查数据。
inputBatch = preview(dsTrain); disp(inputBatch)
InputImage ResponsePixelLabelImage __________________ _______________________ {256×256×6 uint16} {256×256 categorical} {256×256×6 uint16} {256×256 categorical} {256×256×6 uint16} {256×256 categorical} {256×256×6 uint16} {256×256 categorical} {256×256×6 uint16} {256×256 categorical} {256×256×6 uint16} {256×256 categorical} {256×256×6 uint16} {256×256 categorical} {256×256×6 uint16} {256×256 categorical}
创建 U-Net 网络层
此示例使用 U-Net 网络的一种变体。U-Net 中的初始卷积层序列与最大池化层交叠,从而逐步降低输入图像的分辨率。这些层后跟一系列使用上采样算子散布的卷积层,从而会连续增加输入图像的分辨率 [2]。U-Net 的名称源于网络可以绘制成形似字母 U 的对称形状。
此示例对 U-Net 进行了修改以在卷积中使用零填充,从而使卷积的输入和输出具有相同的大小。使用辅助函数 createUnet
创建一个具有几个预选超参数的 U-Net。此函数作为支持文件包含在本示例中。
inputTileSize = [256,256,6]; lgraph = createUnet(inputTileSize); disp(lgraph.Layers)
58×1 Layer array with layers: 1 'ImageInputLayer' Image Input 256×256×6 images with 'zerocenter' normalization 2 'Encoder-Section-1-Conv-1' 2-D Convolution 64 3×3×6 convolutions with stride [1 1] and padding [1 1 1 1] 3 'Encoder-Section-1-ReLU-1' ReLU ReLU 4 'Encoder-Section-1-Conv-2' 2-D Convolution 64 3×3×64 convolutions with stride [1 1] and padding [1 1 1 1] 5 'Encoder-Section-1-ReLU-2' ReLU ReLU 6 'Encoder-Section-1-MaxPool' 2-D Max Pooling 2×2 max pooling with stride [2 2] and padding [0 0 0 0] 7 'Encoder-Section-2-Conv-1' 2-D Convolution 128 3×3×64 convolutions with stride [1 1] and padding [1 1 1 1] 8 'Encoder-Section-2-ReLU-1' ReLU ReLU 9 'Encoder-Section-2-Conv-2' 2-D Convolution 128 3×3×128 convolutions with stride [1 1] and padding [1 1 1 1] 10 'Encoder-Section-2-ReLU-2' ReLU ReLU 11 'Encoder-Section-2-MaxPool' 2-D Max Pooling 2×2 max pooling with stride [2 2] and padding [0 0 0 0] 12 'Encoder-Section-3-Conv-1' 2-D Convolution 256 3×3×128 convolutions with stride [1 1] and padding [1 1 1 1] 13 'Encoder-Section-3-ReLU-1' ReLU ReLU 14 'Encoder-Section-3-Conv-2' 2-D Convolution 256 3×3×256 convolutions with stride [1 1] and padding [1 1 1 1] 15 'Encoder-Section-3-ReLU-2' ReLU ReLU 16 'Encoder-Section-3-MaxPool' 2-D Max Pooling 2×2 max pooling with stride [2 2] and padding [0 0 0 0] 17 'Encoder-Section-4-Conv-1' 2-D Convolution 512 3×3×256 convolutions with stride [1 1] and padding [1 1 1 1] 18 'Encoder-Section-4-ReLU-1' ReLU ReLU 19 'Encoder-Section-4-Conv-2' 2-D Convolution 512 3×3×512 convolutions with stride [1 1] and padding [1 1 1 1] 20 'Encoder-Section-4-ReLU-2' ReLU ReLU 21 'Encoder-Section-4-DropOut' Dropout 50% dropout 22 'Encoder-Section-4-MaxPool' 2-D Max Pooling 2×2 max pooling with stride [2 2] and padding [0 0 0 0] 23 'Mid-Conv-1' 2-D Convolution 1024 3×3×512 convolutions with stride [1 1] and padding [1 1 1 1] 24 'Mid-ReLU-1' ReLU ReLU 25 'Mid-Conv-2' 2-D Convolution 1024 3×3×1024 convolutions with stride [1 1] and padding [1 1 1 1] 26 'Mid-ReLU-2' ReLU ReLU 27 'Mid-DropOut' Dropout 50% dropout 28 'Decoder-Section-1-UpConv' 2-D Transposed Convolution 512 2×2×1024 transposed convolutions with stride [2 2] and cropping [0 0 0 0] 29 'Decoder-Section-1-UpReLU' ReLU ReLU 30 'Decoder-Section-1-DepthConcatenation' Depth concatenation Depth concatenation of 2 inputs 31 'Decoder-Section-1-Conv-1' 2-D Convolution 512 3×3×1024 convolutions with stride [1 1] and padding [1 1 1 1] 32 'Decoder-Section-1-ReLU-1' ReLU ReLU 33 'Decoder-Section-1-Conv-2' 2-D Convolution 512 3×3×512 convolutions with stride [1 1] and padding [1 1 1 1] 34 'Decoder-Section-1-ReLU-2' ReLU ReLU 35 'Decoder-Section-2-UpConv' 2-D Transposed Convolution 256 2×2×512 transposed convolutions with stride [2 2] and cropping [0 0 0 0] 36 'Decoder-Section-2-UpReLU' ReLU ReLU 37 'Decoder-Section-2-DepthConcatenation' Depth concatenation Depth concatenation of 2 inputs 38 'Decoder-Section-2-Conv-1' 2-D Convolution 256 3×3×512 convolutions with stride [1 1] and padding [1 1 1 1] 39 'Decoder-Section-2-ReLU-1' ReLU ReLU 40 'Decoder-Section-2-Conv-2' 2-D Convolution 256 3×3×256 convolutions with stride [1 1] and padding [1 1 1 1] 41 'Decoder-Section-2-ReLU-2' ReLU ReLU 42 'Decoder-Section-3-UpConv' 2-D Transposed Convolution 128 2×2×256 transposed convolutions with stride [2 2] and cropping [0 0 0 0] 43 'Decoder-Section-3-UpReLU' ReLU ReLU 44 'Decoder-Section-3-DepthConcatenation' Depth concatenation Depth concatenation of 2 inputs 45 'Decoder-Section-3-Conv-1' 2-D Convolution 128 3×3×256 convolutions with stride [1 1] and padding [1 1 1 1] 46 'Decoder-Section-3-ReLU-1' ReLU ReLU 47 'Decoder-Section-3-Conv-2' 2-D Convolution 128 3×3×128 convolutions with stride [1 1] and padding [1 1 1 1] 48 'Decoder-Section-3-ReLU-2' ReLU ReLU 49 'Decoder-Section-4-UpConv' 2-D Transposed Convolution 64 2×2×128 transposed convolutions with stride [2 2] and cropping [0 0 0 0] 50 'Decoder-Section-4-UpReLU' ReLU ReLU 51 'Decoder-Section-4-DepthConcatenation' Depth concatenation Depth concatenation of 2 inputs 52 'Decoder-Section-4-Conv-1' 2-D Convolution 64 3×3×128 convolutions with stride [1 1] and padding [1 1 1 1] 53 'Decoder-Section-4-ReLU-1' ReLU ReLU 54 'Decoder-Section-4-Conv-2' 2-D Convolution 64 3×3×64 convolutions with stride [1 1] and padding [1 1 1 1] 55 'Decoder-Section-4-ReLU-2' ReLU ReLU 56 'Final-ConvolutionLayer' 2-D Convolution 18 1×1×64 convolutions with stride [1 1] and padding [0 0 0 0] 57 'Softmax-Layer' Softmax softmax 58 'Segmentation-Layer' Pixel Classification Layer Cross-entropy loss
选择训练选项
使用具有动量的随机梯度下降 (SGDM) 优化来训练网络。使用 trainingOptions
函数指定 SGDM 的超参数设置。
训练深度网络是很费时间的。通过指定高学习率可加快训练速度。然而,这可能会导致网络的梯度爆炸或不受控制地增长,阻碍网络训练成功。要将梯度保持在有意义的范围内,请通过将 "GradientThreshold"
指定为 0.05
来启用梯度裁剪,并指定 "GradientThresholdMethod"
以使用梯度的 L2-范数。
initialLearningRate = 0.05; maxEpochs = 150; minibatchSize = 16; l2reg = 0.0001; options = trainingOptions("sgdm",... InitialLearnRate=initialLearningRate, ... Momentum=0.9,... L2Regularization=l2reg,... MaxEpochs=maxEpochs,... MiniBatchSize=minibatchSize,... LearnRateSchedule="piecewise",... Shuffle="every-epoch",... GradientThresholdMethod="l2norm",... GradientThreshold=0.05, ... Plots="training-progress", ... VerboseFrequency=20);
训练网络或下载预训练网络
要训练网络,请将以下代码中的 doTraining
变量设置为 true
。使用 trainNetwork
函数训练模型。
在 GPU 上(如果有)进行训练。使用 GPU 需要 Parallel Computing Toolbox™ 和支持 CUDA® 的 NVIDIA® GPU。有关详细信息,请参阅GPU Computing Requirements (Parallel Computing Toolbox)。在 NVIDIA Titan X 上训练大约需要 20 个小时。
doTraining = false; if doTraining net = trainNetwork(dsTrain,lgraph,options); modelDateTime = string(datetime("now",Format="yyyy-MM-dd-HH-mm-ss")); save(fullfile(dataDir,"multispectralUnet-"+modelDateTime+".mat"),"net"); end
评估分割准确度
分割验证数据。
segmentedImage = segmentMultispectralImage(val_data,net,predictPatchSize);
将分割后的图像和真实值标签另存为 PNG 文件。该示例使用这些文件来计算准确度度量。
imwrite(segmentedImage,"results.png"); imwrite(val_labels,"gtruth.png");
使用 pixelLabelDatastore
(Computer Vision Toolbox) 加载分割结果和真实值。
pxdsResults = pixelLabelDatastore("results.png",classNames,pixelLabelIds); pxdsTruth = pixelLabelDatastore("gtruth.png",classNames,pixelLabelIds);
使用 evaluateSemanticSegmentation
(Computer Vision Toolbox) 函数衡量语义分割的全局准确度。
ssm = evaluateSemanticSegmentation(pxdsResults,pxdsTruth,Metrics="global-accuracy");
Evaluating semantic segmentation results ---------------------------------------- * Selected metrics: global accuracy. * Processed 1 images. * Finalizing... Done. * Data set metrics: GlobalAccuracy ______________ 0.90411
全局准确度分数表明,正确分类的像素稍大于 90%。
参考资料
[1] Kemker, R., C. Salvaggio, and C. Kanan."High-Resolution Multispectral Dataset for Semantic Segmentation."CoRR, abs/1703.01918. 2017.
[2] Ronneberger, O., P. Fischer, and T. Brox."U-Net:Convolutional Networks for Biomedical Image Segmentation."CoRR, abs/1505.04597. 2015.
[3] Kemker, Ronald, Carl Salvaggio, and Christopher Kanan."Algorithms for Semantic Segmentation of Multispectral Remote Sensing Imagery Using Deep Learning."ISPRS Journal of Photogrammetry and Remote Sensing, Deep Learning RS Data, 145 (November 1, 2018):60-77. https://doi.org/10.1016/j.isprsjprs.2018.04.014.
另请参阅
trainingOptions
| trainNetwork
| randomPatchExtractionDatastore
(Image Processing Toolbox) | pixelLabelDatastore
(Computer Vision Toolbox) | semanticseg
(Computer Vision Toolbox) | evaluateSemanticSegmentation
(Computer Vision Toolbox) | imageDatastore
| histeq
(Image Processing Toolbox) | unetLayers
(Computer Vision Toolbox)
相关主题
- Getting Started with Semantic Segmentation Using Deep Learning (Computer Vision Toolbox)
- 使用深度学习进行语义分割
- 使用扩张卷积进行语义分割
- Datastores for Deep Learning