How to calculate average Dice score across the set of test volumes (Deep Learning, U Net Model)

Dear all,

Based on Segment3DLungTumorUsingVnet as attached, below is the command to Calculate the average Dice score across the set of test volumes. (line 122-127)

meanDiceBackground = mean(diceResult(:,1));
disp(['Average Dice score of background across ',num2str(j), ...
    ' test volumes = ',num2str(meanDiceBackground)])
meanDiceTumor = mean(diceResult(:,2));
disp(['Average Dice score of tumor across ',num2str(j), ...
    ' test volumes = ',num2str(meanDiceTumor)])

Here is my coding for Unet Model Deep Learning.

Anyone can help me to Calculate the average Dice score across the set of test volume?

%% first, read the image data and labelled images
clc
clear all; close all;
dataSetDir = fullfile('C:\Users\Akmal\Desktop\I-131 256 28.02.2020\I-131 SPECT NEMA VALIDATION 01112019 256X256 26.09.2021 petang');
imageDir = fullfile(dataSetDir,'Image');
labelDir = fullfile(dataSetDir,'PixelLabelData');
imds = imageDatastore(imageDir);
% view data set images origional
% figure
% for i = 1:23
%     subplot(5,5,i)
%     I = readimage(imds,i);
%     imshow(I)
%     title('training labels')
% end
%% train the data. if network already, then just drag it into command window
classNames = ["foreground" "background"];
labelIDs = [1 2];
pxds = pixelLabelDatastore(labelDir, classNames, labelIDs);
imds1 = imageDatastore(labelDir);
% figure
% for i = 1:5
%     subplot(3,3,i)
%     I = readimage(imds1,i);
%     imshow(I)
%     title('training labels')
% end
ds = pixelLabelImageDatastore(imds,pxds);
tbl = countEachLabel(pxds)
totalNumberOfPixels = sum(tbl.PixelCount);
frequency = tbl.PixelCount / totalNumberOfPixels;
inverseFrequency = 1./frequency
%    layerf = pixelClassificationLayer(...
%       'Classes',tbl.Name,'ClassWeights',inverseFrequency)
%
layerf=pixelClassificationLayer("Name","Segmentation-Layer")
lgraph = layerGraph();
tempLayers = [
    
imageInputLayer([512/2 512/2 1],"Name","ImageInputLayer")
convolution2dLayer([4 4],64,"Name","Encoder-Stage-1-Conv-1","Padding","same","WeightsInitializer","he")
reluLayer("Name","Encoder-Stage-1-ReLU-1")
convolution2dLayer([4 4],64,"Name","Encoder-Stage-1-Conv-2","Padding","same","WeightsInitializer","he")
reluLayer("Name","Encoder-Stage-1-ReLU-2")];
lgraph = addLayers(lgraph,tempLayers);
tempLayers = [
    maxPooling2dLayer([2 2],"Name","Encoder-Stage-1-MaxPool","Stride",[4 4])
    convolution2dLayer([4 4],128,"Name","Encoder-Stage-2-Conv-1","Padding","same","WeightsInitializer","he")
    reluLayer("Name","Encoder-Stage-2-ReLU-1")
    convolution2dLayer([4 4],128,"Name","Encoder-Stage-2-Conv-2","Padding","same","WeightsInitializer","he")
    
    reluLayer("Name","Encoder-Stage-2-ReLU-2")];
lgraph = addLayers(lgraph,tempLayers);
tempLayers = [
    maxPooling2dLayer([2 2],"Name","Encoder-Stage-2-MaxPool","Stride",[4 4])
    convolution2dLayer([4 4],256,"Name","Encoder-Stage-3-Conv-1","Padding","same","WeightsInitializer","he")
    reluLayer("Name","Encoder-Stage-3-ReLU-1")
    convolution2dLayer([4 4],256,"Name","Encoder-Stage-3-Conv-2","Padding","same","WeightsInitializer","he")
    reluLayer("Name","Encoder-Stage-3-ReLU-2")];
lgraph = addLayers(lgraph,tempLayers);
tempLayers = [
    dropoutLayer(0.5,"Name","Encoder-Stage-3-DropOut")
    maxPooling2dLayer([2 2],"Name","Encoder-Stage-3-MaxPool","Stride",[4 4])
    convolution2dLayer([4 4],512,"Name","Bridge-Conv-1","Padding","same","WeightsInitializer","he")
    reluLayer("Name","Bridge-ReLU-1")
    convolution2dLayer([4 4],512,"Name","Bridge-Conv-2","Padding","same","WeightsInitializer","he")
    reluLayer("Name","Bridge-ReLU-2")
    dropoutLayer(0.5,"Name","Bridge-DropOut")
    transposedConv2dLayer([2 2],256,"Name","Decoder-Stage-1-UpConv","BiasLearnRateFactor",2,"Stride",[4 4],"WeightsInitializer","he")
    reluLayer("Name","Decoder-Stage-1-UpReLU")];
lgraph = addLayers(lgraph,tempLayers);
tempLayers = [
    depthConcatenationLayer(2,"Name","Decoder-Stage-1-DepthConcatenation")
    convolution2dLayer([4 4],256,"Name","Decoder-Stage-1-Conv-1","Padding","same","WeightsInitializer","he")
    reluLayer("Name","Decoder-Stage-1-ReLU-1")
    convolution2dLayer([4 4],256,"Name","Decoder-Stage-1-Conv-2","Padding","same","WeightsInitializer","he")
    reluLayer("Name","Decoder-Stage-1-ReLU-2")
    transposedConv2dLayer([2 2],128,"Name","Decoder-Stage-2-UpConv","BiasLearnRateFactor",2,"Stride",[4 4],"WeightsInitializer","he")
    reluLayer("Name","Decoder-Stage-2-UpReLU")];
lgraph = addLayers(lgraph,tempLayers);
tempLayers = [
    depthConcatenationLayer(2,"Name","Decoder-Stage-2-DepthConcatenation")
    convolution2dLayer([4 4],128,"Name","Decoder-Stage-2-Conv-1","Padding","same","WeightsInitializer","he")
    reluLayer("Name","Decoder-Stage-2-ReLU-1")
    convolution2dLayer([4 4],128,"Name","Decoder-Stage-2-Conv-2","Padding","same","WeightsInitializer","he")
    reluLayer("Name","Decoder-Stage-2-ReLU-2")
    transposedConv2dLayer([2 2],64,"Name","Decoder-Stage-3-UpConv","BiasLearnRateFactor",2,"Stride",[4 4],"WeightsInitializer","he")
    reluLayer("Name","Decoder-Stage-3-UpReLU")];
lgraph = addLayers(lgraph,tempLayers);
tempLayers = [
    depthConcatenationLayer(2,"Name","Decoder-Stage-3-DepthConcatenation")
    convolution2dLayer([4 4],64,"Name","Decoder-Stage-3-Conv-1","Padding","same","WeightsInitializer","he")
    reluLayer("Name","Decoder-Stage-3-ReLU-1")
    convolution2dLayer([4 4],64,"Name","Decoder-Stage-3-Conv-2","Padding","same","WeightsInitializer","he")
    reluLayer("Name","Decoder-Stage-3-ReLU-2")
    convolution2dLayer([1 1],3,"Name","Final-ConvolutionLayer","Padding","same","WeightsInitializer","he")
    softmaxLayer("Name","Softmax-Layer")
    pixelClassificationLayer("Name","Segmentation-Layer")
    
    ];
lgraph = addLayers(lgraph,tempLayers);
% clean up helper variable
clear tempLayers;
lgraph = connectLayers(lgraph,"Encoder-Stage-1-ReLU-2","Encoder-Stage-1-MaxPool");
lgraph = connectLayers(lgraph,"Encoder-Stage-1-ReLU-2","Decoder-Stage-3-DepthConcatenation/in2");
lgraph = connectLayers(lgraph,"Encoder-Stage-2-ReLU-2","Encoder-Stage-2-MaxPool");
lgraph = connectLayers(lgraph,"Encoder-Stage-2-ReLU-2","Decoder-Stage-2-DepthConcatenation/in2");
lgraph = connectLayers(lgraph,"Encoder-Stage-3-ReLU-2","Encoder-Stage-3-DropOut");
lgraph = connectLayers(lgraph,"Encoder-Stage-3-ReLU-2","Decoder-Stage-1-DepthConcatenation/in2");
lgraph = connectLayers(lgraph,"Decoder-Stage-1-UpReLU","Decoder-Stage-1-DepthConcatenation/in1");
lgraph = connectLayers(lgraph,"Decoder-Stage-2-UpReLU","Decoder-Stage-2-DepthConcatenation/in1");
lgraph = connectLayers(lgraph,"Decoder-Stage-3-UpReLU","Decoder-Stage-3-DepthConcatenation/in1");
% lgraph = connectLayers(lgraph,'relu12','skipConv1');
% lgraph = connectLayers(lgraph,'Encoder-Stage-2-Conv-2','add22/in2');
% lgraph = connectLayers(lgraph,'relu22','');
% Plot Layers
figure,plot(lgraph);
imageSize = [256 256 1];
numClasses = 2;
encoderDepth = 3;
lgraph = unetLayers(imageSize,numClasses,'EncoderDepth',encoderDepth)
% split data
[imdsTrain, imdsVal, pxdsTrain, pxdsVal] = partitionCamVidData2(imds,pxds);
pximds = pixelLabelImageDatastore(imdsTrain,pxdsTrain);
pximdsVal = pixelLabelImageDatastore(imdsVal,pxdsVal);
options1 = trainingOptions('adam', ...
    'InitialLearnRate',1e-3, ...
    'MaxEpochs',100, ...
    'LearnRateDropFactor',5e-1, ...
    'LearnRateDropPeriod',10, ...
    'ValidationData',pximdsVal,...
    'ValidationFrequency',3, ...
    'LearnRateSchedule','piecewise', ...
    'MiniBatchSize',4,'Plots','training-progress');
net1 = trainNetwork(pximds,lgraph,options1);
function [imdsTrain, imdsTest, pxdsTrain, pxdsTest] = partitionCamVidData2(imds,pxds)
% Partition CamVid data by randomly selecting 60% of the data for training. The
% rest is used for testing.
% Set initial random state for example reproducibility.
rng(0); 
numFiles = numel(imds.Files);
shuffledIndices = randperm(numFiles);
% Use 60% of the images for training.
N = round(0.60 * numFiles);
trainingIdx = shuffledIndices(1:N);
% Use the rest for testing.
testIdx = shuffledIndices(N+1:end);
% Create image datastores for training and test.
trainingImages = imds.Files(trainingIdx);
testImages = imds.Files(testIdx);
imdsTrain = imageDatastore(trainingImages);
imdsTest = imageDatastore(testImages);
% Extract class and label IDs info.
classes = pxds.ClassNames;
labelIDs = 1:numel(pxds.ClassNames);
% Create pixel label datastores for training and test.
trainingLabels = pxds.Files(trainingIdx);
testLabels = pxds.Files(testIdx);
pxdsTrain = pixelLabelDatastore(trainingLabels, classes, labelIDs);
pxdsTest = pixelLabelDatastore(testLabels, classes, labelIDs);
end

Answer 1

[Avadhoot]

Hi Mohd Akmal Masud,

As I have inferred from your question, you want to calculate the Dice scores across all the samples in your test set. You have already calculated the dice scores for a few samples as you have illustrated in the beginning. You can calculate the Dice score across the whole test set by creating a simple function as follows:

function calculateAverageDice(net, voldsTest, pxdsTest)
    % Initialize Dice score vector
    numTestImages = numel(voldsTest.Files);
    diceScores = zeros(numTestImages, 2); % Assuming two classes: foreground and background
    
    % Perform inference on test images and compute Dice score
    for i = 1:numTestImages
        % Read test image and ground truth
        testImage = readimage(voldsTest, i);
        groundTruth = readimage(pxdsTest, i);
        
        % Perform segmentation
        predictedLabel = semanticseg(testImage, net);
        
        % Compute Dice score for each class
        diceScores(i, :) = dice(groundTruth, predictedLabel);
    end
    
    % Calculate the average Dice score across all test images
    meanDiceForeground = mean(diceScores(:, 1));
    meanDiceBackground = mean(diceScores(:, 2));
    
    % Display the results
    disp(['Average Dice score of foreground across ', num2str(numTestImages), ...
          ' test volumes = ', num2str(meanDiceForeground)]);
    disp(['Average Dice score of background across ', num2str(numTestImages), ...
          ' test volumes = ', num2str(meanDiceBackground)]);
end

This function takes 3 inputs as follows:

1. The neural network (UNet in your case).
2. The test image datastore (you have defined this as voldsTest)
3. The ground truth image datastore (defined here as pxdsTest)

The predicted label is generated by using the "semanticseg" function on your neural network. The mean values are calculated for the whole test set using the "mean" function. As your code also calculates two different means for background and foreground each, the same is followed here. As the end of the function, both the means are displayed on the screen.

You can call this function by using the following line:

calculateAverageDice(net, voldsTest, pxdsTest);

I hope it helps.