基于深度学习的多光谱图像语义分割

这个例子使用了:

打开实时脚本

这个例子展示了如何使用U-Net对具有七个通道的多光谱图像进行语义分割。

语义分割涉及到用一个类标记图像中的每个像素。语义分割的一个应用是跟踪森林砍伐，即森林覆盖随时间的变化。环境机构跟踪森林砍伐，以评估和量化一个地区的环境和生态健康。

基于深度学习的语义分割可以从高分辨率航拍照片中精确测量植被覆盖。一个挑战是区分具有相似视觉特征的类，例如试图将一个绿色像素分类为草、灌木或树。为了提高分类的准确性，一些数据集包含多光谱图像，提供关于每个像素的附加信息。例如，Hamlin海滩州立公园的数据集用三个近红外通道补充了彩色图像，提供了更清晰的分类。

这个例子展示了如何使用基于深度学习的语义分割技术从一组多光谱图像中计算一个地区的植被覆盖百分比。

下载数据

本例使用高分辨率多光谱数据集来训练网络[1］．这组照片是用无人机在纽约哈姆林海滩州立公园上空拍摄的。数据包含有标记的训练、验证和测试集，有18个对象类标签。数据文件的大小为~3.0 GB。

下载数据集的mat文件版本downloadHamlinBeachMSIDatahelper函数。该函数作为支持文件附加到示例中。万博1manbetx

imageDir = tempdir;url =“http://www.cis.rit.edu/ ~ rmk6217 / rit18_data.mat”；downloadHamlinBeachMSIData (url, imageDir);

检查培训数据

将数据集加载到工作区中。

负载(fullfile (imageDir“rit18_data”，“rit18_data.mat”));

检查数据的结构。

谁train_dataval_datatest_data

名称大小字节类属性test_data 7x12446x7654 1333663576 uint16 train_data 7x9393x5642 741934284 uint16 val_data 7x8833x6918 855493716 uint16

所述多光谱图像数据排列为numChannels——- - - - - -宽度——- - - - - -高度数组。然而，在MATLAB®中，多通道图像的排列方式为宽度——- - - - - -高度——- - - - - -numChannels数组。要重塑数据，使通道处于第三维度，可以使用helper函数，switchChannelsToThirdPlane．该函数作为支持文件附加到示例中。万博1manbetx

train_data = switchChannelsToThirdPlane(train_data);val_data = switchChannelsToThirdPlane(val_data);test_data = switchChannelsToThirdPlane(test_data);

确认数据具有正确的结构。

谁train_dataval_datatest_data

名称大小字节类属性test_data 12446x7654x7 1333663576 uint16 train_data 9393x5642x7 741934284 uint16 val_data 8833x6918x7 855493716 uint16

RGB颜色通道是第三、第二和第一个图像通道。以蒙太奇的形式显示训练、验证和测试图像的颜色组件。要使图像在屏幕上显示得更亮，可以使用histeq函数。

图蒙太奇(.．.{histeq(train_data(:，:，[3 2 1]))，.．.Histeq (val_data(:，:，[3 2 1]))，.．.Histeq (test_data(:，:，[3 2 1]))}，.．.“BorderSize”10“写成BackgroundColor”，“白色”)标题(“训练图像(左)、验证图像(中)和测试图像(右)的RGB组件”）

以蒙太奇的形式显示训练数据的最后三个直方图均衡通道。这些通道对应于近红外波段，并根据其热特征突出显示图像的不同组成部分。例如，靠近第二个通道图像中心的树比其他两个通道中的树显示出更多的细节。

图蒙太奇(.．.{histeq (train_data (:: 4)),.．.histeq (train_data (:: 5)),.．.histeq (train_data (:: 6))},.．.“BorderSize”10“写成BackgroundColor”，“白色”)标题(“训练图像的红外通道1(左)、2(中)和3(右)”）

通道7是表示有效分割区域的掩码。显示训练、验证和测试图像的掩码。

图蒙太奇(.．.{train_data (:: 7).．.val_data (:: 7).．.test_data (:: 7)},.．.“BorderSize”10“写成BackgroundColor”，“白色”)标题(“训练图像(左)、验证图像(中)和测试图像(右)的掩码”）

标记的图像包含分割的地面真相数据，每个像素分配给18个类中的一个。获取类及其对应id的列表。

disp(类)

0.其他类/图像边界道路标记树3。建设4。交通工具(小汽车、卡车或公共汽车)6人。救生员座椅7。野餐桌黑色木质面板 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 = [“路标”，“树”，“建筑”，“汽车”，“人”，.．.“LifeguardChair”，“PicnicTable”，“BlackWoodPanel”，.．.“WhiteWoodPanel”，“OrangeLandingPad”，“浮”，“石头”，.．.“LowLevelVegetation”，“Grass_Lawn”，“Sand_Beach”，.．.“Water_Lake”，“Water_Pond”，“沥青”];

将标签覆盖在直方图均衡化的RGB训练图像上。为图像添加一个颜色条。

cmap = jet(numel(classNames));B = labeloverlay(histq (train_data(:，:，4:6))，train_labels，“透明”, 0.8,“Colormap”,提出);图imshow(B)标题(“培训”标签N = numel(classNames);ticks = 1/(N*2):1/N:1;colorbar (“TickLabels”cellstr(类名),“滴答”蜱虫,“TickLength”0,“TickLabelInterpreter”，“没有”）;colormap城市规划机构(cmap)

将培训数据保存为MAT文件，培训标签保存为PNG文件。

保存(“train_data.mat”，“train_data”）;imwrite (train_labels“train_labels.png”）;

创建随机补丁提取训练数据存储

使用随机补丁提取数据存储将训练数据提供给网络。该数据存储从包含地面真实图像和像素标签数据的图像数据存储和像素标签数据存储中提取多个相应的随机补丁。打补丁是一种常见的技术，可以防止大图像的内存耗尽，并有效增加可用的训练数据量。

从存储训练图像开始“train_data.mat”在一个imageDatastore．因为MAT文件格式是非标准图像格式，所以必须使用MAT文件读取器才能读取图像数据。你可以使用辅助MAT文件阅读器，matReader，该方法从训练数据中提取前六个通道，并删除包含掩码的最后一个通道。该函数作为支持文件附加到示例中。万博1manbetx

imds = imageDatastore(“train_data.mat”，“FileExtensions”，“.mat”，“ReadFcn”, @matReader);

创建一个pixelLabelDatastore存储包含18个标记区域的标签补丁。

pixelLabelIds = 1:18;pxds = pixelLabelDatastore(“train_labels.png”一会,pixelLabelIds);

创建一个randomPatchExtractionDatastore从图像数据存储和像素标签数据存储。每个小批包含16个大小为256 × 256像素的补丁。每次迭代提取1000个小批次。

dsTrain = randomPatchExtractionDatastore(imds,pxds，[256,256]，“PatchesPerImage”, 16000);

随机补丁提取数据存储dsTrain在历的每个迭代中向网络提供小批量数据。预览数据存储以查看数据。

inputBatch =预览(dsTrain);disp (inputBatch)

InputImage ResponsePixelLabelImage  __________________ _______________________ { 256×256×6 uint16}{256×256分类}{256×256×6 uint16}{256×256分类}{256×256×6 uint16}{256×256分类}{256×256×6 uint16}{256×256分类}{256×256×6 uint16}{256×256分类}{256×256×6 uint16}{256×256分类}{256×256×6 uint16}{256×256分类}{256×256×6 uint16}{256×256分类}

创建U-Net网络层

本例使用的是U-Net网络的变体。在U-Net中，初始的卷积层序列中穿插着最大池化层，依次降低了输入图像的分辨率。这些层之后是一系列的卷积层，其中穿插着上采样算子，依次增加输入图像的分辨率[2］．U- net的名字来源于这样一个事实:这个网络可以画成一个像字母U一样的对称形状。

这个示例修改U-Net，使其在卷积中使用零填充，以便卷积的输入和输出具有相同的大小。使用助手函数，createUnet，以创建一个带有一些预选超参数的U-Net。该函数作为支持文件附加到示例中。万博1manbetx

inputTileSize = [256,256,6];lgraph = createUnet(inputTileSize);disp (lgraph.Layers)

58×1带有图层的图层数组:1“ImageInputLayer”图像输入图像256×256×6”zerocenter正常化2 Encoder-Section-1-Conv-1卷积64 3×3×6旋转步[1]和填充[1 1 1 1]3‘Encoder-Section-1-ReLU-1 ReLU ReLU 4 Encoder-Section-1-Conv-2卷积64 3×3×64旋转步[1]和填充[1 1 1 1]5‘Encoder-Section-1-ReLU-2 ReLU ReLU 6“Encoder-Section-1-MaxPool”马克斯池2×2马克斯池步(2 - 2)和填充[0 0 0 0]7“Encoder-Section-2-Conv-1”卷积128 3×3×64旋转步[1]和填充(1 1 1)8“Encoder-Section-2-ReLU-1”ReLU ReLU 9 Encoder-Section-2-Conv-2卷积128 3×3×128旋转步[1]和填充[1 1 1 1]10 ' Encoder-Section-2-ReLU-2 ReLU ReLU 11“Encoder-Section-2-MaxPool”马克斯池2×2马克斯池步(2 - 2)和填充[0 0 0 0]12 Encoder-Section-3-Conv-1卷积256 3×3×128旋转步[1]和填充[1 1 1 1]13“Encoder-Section-3-ReLU-1”ReLU ReLU 14 Encoder-Section-3-Conv-2卷积256 3×3×256旋转步[1]和填充[1 1 1 1]15 ' Encoder-Section-3-ReLU-2 ReLU ReLU 16“Encoder-Section-3-MaxPool”马克斯池2×2马克斯池步(2 - 2)和填充[0 0 0 0]17 Encoder-Section-4-Conv-1卷积512 3×3×256旋转步[1]和填充[1 1 1 1]18 ' Encoder-Section-4-ReLU-1 ReLU ReLU 19 Encoder-Section-4-Conv-2卷积512 3×3×512旋转步[1]和填充[1 1 1 1] 20 ' Encoder-Section-4-ReLU-2 ReLU ReLU 21“Encoder-Section-4-DropOut”辍学50%辍学22“Encoder-Section-4-MaxPool”马克斯池2×2马克斯池步(2 - 2)和填充[0 0 0 0]23 Mid-Conv-1卷积1024 3×3×512旋转步[1]和填充[1 1 1 1]24的Mid-ReLU-1 ReLU ReLU 25 Mid-Conv-2卷积1024 3×3×1024旋转步[1]和填充[1 1 1 1]26 ' Mid-ReLU-2 ReLU ReLU 27“Mid-DropOut”辍学辍学28 50%Decoder-Section-1-UpConv转置卷积512 2×2×1024转置运算跨2[2]和裁剪[0 0 0 0]29的Decoder-Section-1-UpReLU ReLU ReLU 30的Decoder-Section-1-DepthConcatenation深度连接深度连接2输入31 Decoder-Section-1-Conv-1卷积512 3×3×1024旋转步[1]和填充[1 1 1 1]32的Decoder-Section-1-ReLU-1 ReLU ReLU 33 Decoder-Section-1-Conv-2卷积512 3×3×512旋转步[1]和填充[1 1 1 1]34 ' Decoder-Section-1-ReLU-2 ReLU ReLU 35 Decoder-Section-2-UpConv转置卷积256 2×2×512转置运算跨2[2]和裁剪[0 0 0 0]36的Decoder-Section-2-UpReLU ReLU ReLU 37 Decoder-Section-2-DepthConcatenation的深度连接深度连接2输入38 Decoder-Section-2-Conv-1卷积256 3×3×512旋转步[1]和填充[1 1 1 1]39的Decoder-Section-2-ReLU-1 ReLU ReLU 40“Decoder-Section-2-Conv-2”卷积256 3×3×256旋转步[1]和填充(1 1 1)41 ' Decoder-Section-2-ReLU-2 ReLU ReLU 42 Decoder-Section-3-UpConv转置卷积128 2×2×256转置运算跨2[2]和裁剪[0 0 0 0]43的Decoder-Section-3-UpReLU ReLU ReLU 44的Decoder-Section-3-DepthConcatenation深度连接深度连接2输入45‘Decoder-Section-3-Conv-1卷积128 3×3×256旋转步[1]和填充[1 1 1 1]46‘Decoder-Section-3-ReLU-1 ReLU ReLU 47 Decoder-Section-3-Conv-2卷积128 3×3×128旋转步[1]和填充[1 1 1 1]48 ' Decoder-Section-3-ReLU-2 ReLU ReLU 49 Decoder-Section-4-UpConv转置卷积64 2×2×128转置运算跨2[2]和种植50[0 0 0 0]的Decoder-Section-4-UpReLU ReLU ReLU 51的Decoder-Section-4-DepthConcatenation深度连接深度连接2输入52 Decoder-Section-4-Conv-1卷积64 3×3×128旋转步[1]和填充[1 1 1 1]53的Decoder-Section-4-ReLU-1 ReLU ReLU 54 Decoder-Section-4-Conv-2卷积64 3×3×64旋转步[1]和填充[1 1 1 1]55的Decoder-Section-4-ReLU-2 ReLU ReLU 56 Final-ConvolutionLayer的卷积18 1×1×64旋转步[1]和填充[0 0 0 0]57 ' Softmax-Layer Softmax Softmax 58 Segmentation-Layer的像素分类层叉损失

选择培训项目

利用随机梯度下降动量(SGDM)优化训练网络。方法指定SGDM的超参数设置trainingOptions(深度学习工具箱)函数。

训练一个深度网络是很耗时的。通过指定较高的学习率来加速培训。然而，这可能导致网络的梯度爆发或不受控制地增长，阻止网络训练成功。要将渐变保持在有意义的范围内，请通过指定来启用渐变剪切“GradientThreshold”作为0.05，并指定“GradientThresholdMethod”用梯度的l2范数。

initialLearningRate = 0.05;maxEpochs = 150;minibatchSize = 16;L2reg = 0.0001;options = trainingOptions(“个”，.．.“InitialLearnRate”initialLearningRate,.．.“动量”, 0.9,.．.“L2Regularization”l2reg,.．.“MaxEpochs”maxEpochs,.．.“MiniBatchSize”minibatchSize,.．.“LearnRateSchedule”，“分段”，.．.“洗牌”，“every-epoch”，.．.“GradientThresholdMethod”，“l2norm”，.．.“GradientThreshold”, 0.05,.．.“阴谋”，“训练进步”，.．.“VerboseFrequency”, 20);

训练网络

默认情况下，该示例为该数据集下载预先训练的U-Net版本downloadTrainedUnethelper函数。该函数作为支持文件附加到示例中。万博1manbetx预先训练的网络使您可以运行整个示例，而不必等待训练完成。

要训练网络，请设置doTraining变量的真正的．的方法训练模型trainNetwork(深度学习工具箱)函数。

如果有GPU，请使用GPU进行训练。使用GPU需要并行计算工具箱™和CUDA®支持的NVIDIA®GPU。有关更多信息，请参见GPU支万博1manbetx持按版本划分(并行计算工具箱)．在NVIDIA Titan X上训练大约需要20个小时。

doTraining = false;如果doTraining [net,info] = trainNetwork(dsTrain,lgraph,options);modelDateTime = string(datetime)“现在”，“格式”，“yyyy-MM-dd-HH-mm-ss”));保存(strcat (“multispectralUnet——”modelDateTime,“时代——”num2str (maxEpochs),“.mat”)，“净”）;其他的trainedUnet_url =“//www.tianjin-qmedu.com/万博1manbetxsupportfiles/vision/data/multispectralUnet.mat”；downloadTrainedUnet (trainedUnet_url imageDir);负载(fullfile (imageDir“trainedUnet”，“multispectralUnet.mat”));结束

您现在可以使用U-Net对多光谱图像进行语义分割。

根据测试数据预测结果

要在经过训练的网络上执行转发，可以使用helper函数，segmentImage，使用验证数据集。该函数作为支持文件附加到示例中。万博1manbetxsegmentImage方法对图像补丁执行分割semanticseg函数。

predictPatchSize = [1024 1024];segmentedImage = segmentmage (val_data,net,predictPatchSize);

为了只提取分割的有效部分，将分割后的图像乘以验证数据的掩码通道。

segmentedImage = uint8(val_data(:，:，7)~=0) .* segmentedImage;图imshow(segmentedImage，[]) title(“分割图像”）

语义分割的输出是有噪声的。执行后图像处理，以去除噪声和杂散像素。使用medfilt2函数从分割中去除盐和胡椒噪声。可视化去噪后的分割图像。

segmentedImage = medfilt2(segmentedImage，[7,7]);imshow (segmentedImage []);标题(“去噪分割图像”）

将分割后的图像叠加在直方图均衡化RGB验证图像上。

B = labeloverlay(histq (val_data(:，:，[3 2 1]))，segmentedImage，“透明”, 0.8,“Colormap”,提出);图imshow(B)标题(“标记验证图像”) colorbar (“TickLabels”cellstr(类名),“滴答”蜱虫,“TickLength”0,“TickLabelInterpreter”，“没有”）;colormap城市规划机构(cmap)

将分割后的图像和ground truth标签保存为PNG文件。这些将用于计算精度指标。

imwrite (segmentedImage“results.png”）;imwrite (val_labels“gtruth.png”）;

量化分割精度

创建一个pixelLabelDatastore用于分割结果和地面真值标签。

pxdsResults = pixelLabelDatastore(“results.png”一会,pixelLabelIds);pxdsTruth = pixelLabelDatastore(“gtruth.png”一会,pixelLabelIds);

使用的方法来测量语义分割的全局精度evaluateSemanticSegmentation函数。

ssm = evaluateSemanticSegmentation(pxdsResults,pxdsTruth，“指标”，“global-accuracy”）;

评估语义分割结果  ---------------------------------------- * 所选指标:全球精度。*处理1张图像。*完成……完成了。*数据集指标:GlobalAccuracy ______________ 0.90698

全局精度评分表明，只有超过90%的像素被正确分类。

计算植被覆盖范围

本例的最终目标是计算多光谱图像中的植被覆盖范围。

找出标记植被的像素数。标签id 2(“Trees”)、13(“LowLevelVegetation”)和14(“Grass_Lawn”)是植被类。还可以通过对掩码图像的ROI内的像素相加来找到有效像素的总数。

vegetationClassIds = uint8([2,13,14]);vegetationPixels = ismember(segmentedImage(:)，vegetationClassIds);validPixels = (segmentedImage~=0);numVegetationPixels = sum(vegetationPixels(:));numValidPixels = sum(validPixels(:));

用植被像素数除以有效像素数计算植被覆盖的百分比。

percentVegetationCover = (numVegetationPixels/numValidPixels)*100;流(“植被覆盖率为%3.2f%%。”, percentVegetationCover);

植被覆盖率为51.72%。

参考文献

R. Kemker, C. Salvaggio和C. Kanan。用于语义分割的高分辨率多光谱数据集。, abs / 1703.01918。2017.

Ronneberger, O. P. Fischer, T. Brox。U-Net:生物医学图像分割的卷积网络。, abs / 1505.04597。2015.

另请参阅

外部网站

https://github.com/rmkemker/RIT-18