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

本例使用：

此示例示出了如何培养的U净卷积神经网络与七个通道进行多光谱图像的语义分割：三个彩色通道，三近红外信道，和一个掩模。

该示例展示了如何训练一个U-Net网络，并提供了一个预训练的U-Net网络。如果您选择训练U-Net网络，强烈建议使用支持cuda的NVIDIA™GPU，计算能力3.0或更高(需要并行计算工具箱™)。

介绍

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

深学习基于语义分割可以产生从高分辨率的航空照片植被覆盖的精确测量。一个挑战是相似的视觉特性区分类别，如尝试将绿色像素草，灌木，或树木分类。为了提高分级精度，某些数据集包含提供关于每个像素的附加信息多光谱图像。例如，哈姆林海滩州立公园数据与提供的类的更清楚的分离近红外信道设置补充剂的彩色图像。

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

下载数据

此示例使用高分辨率多光谱数据组对网络进行训练[1]．这组图像是由一架无人机在纽约哈姆林海滩州立公园上空拍摄的。该数据包含带有18个对象类标签的标记训练、验证和测试集。数据文件大小为~3.0 GB。

下载数据集的MAT-file版本downloadHamlinBeachMSIData助手功能。这个函数作为支持文件附加到示例中。金宝app

imageDir = tempdir;url =“http://www.cis.rit.edu/ ~ rmk6217 / rit18_data.mat”下载Hamlinbeachmsidata（url，imageDir）；

此外，使用。下载该数据集的U-Net预训练版本downloadTrainedUnet助手功能。这个函数作为支持文件附加到示例中。金宝app预训练的模型使您无需等待培训，完全运行整个例子。

trainedUnet_url ='//www.tatmou.com/金宝appsupportfiles/vision/data/multispectralUnet.mat'；downloadTrainedUnet (trainedUnet_url imageDir);

检查训练数据

加载数据设置到工作区。

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

检查数据的结构。

谁train_dataval_data测试数据

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

多光谱图像数据被布置为numChannels-借-宽度-借-高度数组。而在MATLAB®中，多通道图像排列为宽度-借-高度-借-numChannels数组。要重塑数据以使通道处于三维，请使用helper函数，switchChannelsToThirdPlane．这个函数作为支持文件附加到示例中。金宝app

列车数据=道岔通道至第三车道（列车数据）；车辆数据=道岔通道至第三车道（车辆数据）；测试数据=道岔通道至第三车道（测试数据）；

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

谁train_dataval_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)))},......“边界大小”10，“背景色”，“白色”） 标题（“训练图像(左)、验证图像(中)和测试图像(右)的RGB组件”）

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

数字蒙太奇（......{histeq (train_data (:: 4)),......histeq (train_data (:: 5)),......histeq (train_data (:: 6))},......“边界大小”10，“背景色”，“白色”） 标题（“训练图像的IR通道1(左)、2(中)和3(右)”）

通道7的掩模，指示有效分割区域。显示为培训，验证和测试图像掩码。

数字蒙太奇（......{train_data（：，：，7），......val_data（：，：，7），......test_data (:: 7)},......“边界大小”10，“背景色”，“白色”） 标题（“训练图像遮罩(左)、验证图像(中)和测试图像(右)”）

将标记的图像包含地面实况数据的分割，与分配给18类中的一个的每个像素。获取类的列表以及它们相应的ID。

disp(类)

0，其它类/图像边界1.路标2.树3.建筑4.车辆（汽车，卡车，或总线）5.人6.救生员椅7.野餐表8.黑色人造板9.白人造板10橙色着落垫11水浮标12.岩石13.其他植被14.草15.砂16.将水（湖）17.将水（池）18沥青（停车场/走道）

创建类名的载体。

一会= [“路标”，“树”，“建筑”，“汽车”，“人”，......“救生员发”，“PicnicTable”，“BlackWoodPanel”，......“WhiteWoodPanel”，“OrangeLandingPad”，“浮标”，“岩石”，......“LowLevelVegetation”，“Grass_Lawn”，“Sand_Beach”，......“Water_Lake”，“Water_Pond”，“沥青”];

在直方图均衡化的RGB训练图像上叠加标签。给图像添加一个颜色栏。

cmap=jet（numel（类名））；B=labeloverlay（histeq（列车数据（：，：，4:6）），列车标签，“透明”, 0.8,“颜色表”，CMAP）;图题（“培训”标签) imshow(B) N = numel(classNames);蜱虫= 1 / (N * 2): 1 / N: 1;colorbar ('TickLabels'，cellstr（类名），“蜱”蜱虫,'TickLength'，0，“TickLabelInterpreter”，“没有”）;颜色表（CMAP）

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

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

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

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

从存储训练图像开始“train_data.mat”在一个图像数据存储．因为MAT文件格式是非标准的图像格式，所以必须使用MAT文件阅读器来读取图像数据。你可以使用辅助MAT文件读取器，matReader，即提取物从训练数据和遗漏含有该掩模的最后一个信道的前六个通道。这个函数作为支持文件附加到示例中。金宝app

imd = imageDatastore (“train_data.mat”，“文件扩展名”，“.mat”，“ReadFcn”，@matReader）；

创建一个像素标签数据库存储包含18个标记区域的标签修补程序。

pixelLabelIds = 1:18;pxds = pixelLabelDatastore (“train_labels.png”，类名，pixelLabelIds）;

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

pxds dsTrain = randomPatchExtractionDatastore (imd, [256256],'PatchesPerImage'，16000）;

随机补丁抽取数据存储dsTrain在epoch的每次迭代中向网络提供小批量的数据。预览数据存储以探索数据。

inputBatch=预览（dsTrain）；显示（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 categorical} {256×256×6 uint16} {256×256 categorical}

创建掌中宽带网络层

本例使用U-Net网络的一种变体。在U-Net中，最初的卷积层序列被最大池化层点缀，依次降低了输入图像的分辨率。这些层之后是一系列交错上采样算子的卷积层，依次提高了输入图像的分辨率[2]网络的名称来源于这样一个事实，即网络可以画成一个对称的形状，如字母U。

本例修改了U-Net，使其在卷积中使用零填充，以便卷积的输入和输出具有相同的大小。使用辅助函数，createUnet，以创建带有几个预选超参数的U形网络。此函数作为支持文件附加到示例中。金宝app

inputTileSize = [256,256,6];lgraph = createUnet（inputTileSize）;DISP（lgraph.Layers）

58x1 Layer array with layers:256 x256x6 ImageInputLayer的图像输入图像的zerocenter正常化64 3 x3x6 Encoder-Section-1-Conv-1的卷积运算与步幅[1]和填充[1 1 1 1]3‘Encoder-Section-1-ReLU-1 ReLU ReLU 4 Encoder-Section-1-Conv-2卷积64 3 x3x64旋转步[1]和填充(1 1 1)5‘Encoder-Section-1-ReLU-2 ReLU ReLU 6“Encoder-Section-1-MaxPool”马克斯池2 x2马克斯池步(2 - 2)和填充[0 0 0 0]7 Encoder-Section-2-Conv-1卷积128 3 x3x64旋转步[1]和填充(1 1 1)8“Encoder-Section-2-ReLU-1”ReLU ReLU 9 Encoder-Section-2-Conv-2卷积128 3 x3x128曲线玲珑步[1]和填充[1 1 1 1]10“Encoder-Section-2-ReLU-2”ReLU ReLU 11“Encoder-Section-2-MaxPool”马克斯池2 x2马克斯池步[2 2]和填充[0 0 0 0]12 Encoder-Section-3-Conv-1卷积256 3 x3x128旋转步[1]和填充[1 1 1 1]13的Encoder-Section-3-ReLU-1 ReLU ReLU 14“Encoder-Section-3-Conv-2”Convolution 256 3x3x256 convolutions with stride [1 1] and padding [1 1 1 1] 15 'Encoder-Section-3-ReLU-2' ReLU ReLU 16 'Encoder-Section-3-MaxPool' Max Pooling 2x2 max pooling with stride [2 2] and padding [0 0 0 0] 17 'Encoder-Section-4-Conv-1' Convolution 512 3x3x256 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' Convolution 512 3x3x512 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' Max Pooling 2x2 max pooling with stride [2 2] and padding [0 0 0 0] 23 'Mid-Conv-1' Convolution 1024 3x3x512 convolutions with stride [1 1] and padding [1 1 1 1] 24 'Mid-ReLU-1' ReLU ReLU 25 'Mid-Conv-2' Convolution 1024 3x3x1024 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' Transposed Convolution 512 2x2x1024 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' Convolution 512 3x3x1024 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' Convolution 512 3x3x512 convolutions with stride [1 1] and padding [1 1 1 1] 34 'Decoder-Section-1-ReLU-2' ReLU ReLU 35 'Decoder-Section-2-UpConv' Transposed Convolution 256 2x2x512 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' Convolution 256 3x3x512 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' Convolution 256 3x3x256 convolutions with stride [1 1] and padding [1 1 1 1] 41 'Decoder-Section-2-ReLU-2' ReLU ReLU 42 'Decoder-Section-3-UpConv' Transposed Convolution 128 2x2x256 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' Convolution 128 3x3x256 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' Convolution 128 3x3x128 convolutions with stride [1 1] and padding [1 1 1 1] 48 'Decoder-Section-3-ReLU-2' ReLU ReLU 49 'Decoder-Section-4-UpConv' Transposed Convolution 64 2x2x128 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' Convolution 64 3x3x128 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' Convolution 64 3x3x64 convolutions with stride [1 1] and padding [1 1 1 1] 55 'Decoder-Section-4-ReLU-2' ReLU ReLU 56 'Final-ConvolutionLayer' Convolution 18 1x1x64 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)优化。方法指定SGDM的超参数设置trainingOptions（深学习工具箱）函数。

训练深层网络非常耗时。通过指定较高的学习率来加速训练。但是，这可能会导致网络的梯度爆炸或无法控制地增长，从而阻止网络成功训练。若要将梯度保持在有意义的范围内，请通过指定来启用梯度剪裁“GradientThreshold”作为0．05,并指定'GradientThresholdMethod'使用梯度的l2范数。

initialLearningRate = 0.05;maxEpochs = 150;minibatchSize = 16;l2reg = 0.0001;选择= trainingOptions ('SGDM'，......“InitialLearnRate”initialLearningRate,......“动量”, 0.9,......'L2Regularization'l2reg,......“MaxEpochs”maxEpochs,......'MiniBatchSize'，小批量，......'LearnRateSchedule'，“分段”，......“洗牌”，“every-epoch”，......'GradientThresholdMethod'，“l2norm”，......“GradientThreshold”，0.05，......“阴谋”，“训练进度”，......“详细频率”, 20);

培训网络

在配置训练选项和随机补丁提取数据存储后，使用trainNetwork（深学习工具箱）函数。要训练网络，设置doTraining的参数真正的．强烈建议使用支持cuda的NVIDIA™GPU，其计算能力为3.0或更高。

如果你保留doTraining参数，如下所示假，则示例返回一个预先训练的U-Net网络。

注：在NVIDIA上培训大约需要20小时™ 泰坦X和可能需要更长的时间取决于您的GPU硬件。

doTraining = false;如果doTraining modelDateTime = datestr（现在，“dd-mmm-yyyy-HH-MM-SS”）;(网络,信息)= trainNetwork (dsTrain、lgraph选项);保存([“multispectralUnet -”modelDateTime的时代,num2str（maxEpochs）“.mat”]，“净”，“选项”）;其他的负载(fullfile (imageDir“trainedUnet”，“multispectralUnet.mat”））;结束

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

预测测试数据的结果

要在经过训练的网络上执行向前传球，请使用辅助功能，segmentImage与验证数据集。这个函数作为支持文件附加到示例中。金宝appsegmentImage的方法对图像块进行分割semanticseg函数。

predictPatchSize = [1024 1024];segmentedImage = segmentImage (val_data,净,predictPatchSize);

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

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

语义分割的输出有噪声。执行图像后处理以去除噪声和杂散像素。使用医疗过滤器2函数从分割除去椒盐噪声。可视化移除了噪声的分割图像。

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

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

B = labeloverlay（histeq（val_data（：，：，[3 2 1]）），segmentedImage，“透明”, 0.8,“颜色表”，CMAP）;图imshow（B）标题（“标记的验证图像”）彩条（'TickLabels'，cellstr（类名），“蜱”蜱虫,'TickLength'，0，“TickLabelInterpreter”，“没有”）;颜色表（CMAP）

将分割图像和地面真相标签保存为PNG文件。这些将用于计算精度指标。

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

Quantify的分割精度

创建一个像素标签数据库对于分割结果和地面实况标签。

pxdsResults = pixelLabelDatastore (“results.png”，类名，pixelLabelIds）;pxdsTruth = pixelLabelDatastore（“gtruth.png”，类名，pixelLabelIds）;

用该方法对语义切分的全局精度进行度量evaluateSemanticSegmentation函数。

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

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

全局精度分数表明，90%以上的像素分类正确。

植被覆盖的程度计算

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

找到标记植被的像素数。标签id 2(“Trees”)、13(“LowLevelVegetation”)和14(“Grass_Lawn”)是植被类。并通过对掩模图像感兴趣区域内像素的求和得到有效像素的总数。

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

将植被像素数除以有效像素数，计算植被覆盖率百分比。

percentVegetationCover =（numVegetationPixels / numValidPixels）* 100;fprintf中（“植被覆盖率为%3.2f%%”, percentVegetationCover);

植被覆盖率为51.72%。

工具书类

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

[2] Ronneberger，O.，P.费，和T.是BrOx。“U-网：卷积网络的医学图像分割。”CORR，ABS / 1505.04597。2015年。

另请参阅

外部网站

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