Main Content

PointNet による深層学習を使用した点群の分類

この例では、点群分類用の PointNet ネットワークに学習させる方法を説明します。

点群データは、LiDAR、レーダー、深度カメラといったさまざまなセンサーで取得されます。これらのセンサーは、シーン内にあるオブジェクトの 3 次元位置情報を取得します。この情報は、自動運転や拡張現実におけるさまざまなアプリケーションで役立ちます。たとえば、車両と歩行者を区別することは、自律型車両の経路を計画する際に非常に重要となります。しかし、オブジェクト当たりのデータのスパース性、オブジェクトのオクルージョン、センサーのノイズの問題があるため、点群データを使用してロバストな分類器に学習させることは困難です。これらのさまざまな課題に対処するため、点群データからロバストな特徴表現を直接学習する深層学習手法が提案されています。点群分類のための独創性に富んだ深層学習手法のひとつとして、PointNet があります [1]。

この例では、University of Sydney によって作成された Sydney Urban Objects データセットで PointNet 分類器に学習させます [2]。このデータ セットは、LiDAR センサーを使用して都市環境から取得された点群データの集合を提供します。データセットにはラベル付けされた 100 個のオブジェクトが含まれており、各オブジェクトは、自動車、歩行者、バスなどの 14 のカテゴリに属しています。

データ セットの読み込み

Sydney Urban Objects データセットを、一時ディレクトリにダウンロードして解凍します。

downloadDirectory = tempdir;
datapath = downloadSydneyUrbanObjects(downloadDirectory);

この例の最後にリストされている補助関数 loadSydneyUrbanObjectsData を使用して、ダウンロードした学習データ セットと検証データ セットを読み込みます。最初の 3 つのデータ フォルダーを学習に使用し、4 番目のデータ フォルダーを検証に使用します。

foldsTrain = 1:3;
foldsVal = 4;
dsTrain = loadSydneyUrbanObjectsData(datapath,foldsTrain);
dsVal = loadSydneyUrbanObjectsData(datapath,foldsVal);

いずれかの学習サンプルを読み取り、pcshow を使用してそのサンプルを表示します。

data = read(dsTrain);
ptCloud = data{1,1};
label = data{1,2};

figure
pcshow(ptCloud.Location,[0 1 0],"MarkerSize",40,"VerticalAxisDir","down")
xlabel("X")
ylabel("Y")
zlabel("Z")
title(label)

データセット内のラベルの分布をさらに理解するため、ラベルを読み取って各ラベルに割り当てられている点の数をカウントします。

dsLabelCounts = transform(dsTrain,@(data){data{2} data{1}.Count});
labelCounts = readall(dsLabelCounts);
labels = vertcat(labelCounts{:,1});
counts = vertcat(labelCounts{:,2});

次に、ヒストグラムを使用して、クラスの分布を表示します。

figure
histogram(labels)

ラベルのヒストグラムを見ると、データセットのバランスがとれておらず、自動車と歩行者にバイアスがかかっているのがわかります。これが、ロバストな分類器の学習を妨げています。クラスの不均衡の問題に対処するには、頻度の少ないクラスをオーバーサンプリングします。Sydney Urban Objects データセットの場合、頻度の少ないクラスに対応するファイルを複製すれば、クラスの不均衡の問題に簡単に対処できます。

ラベルごとにファイルをグループ化し、クラス当たりの観測値の数をカウントします。そして、クラス当たりの観測値の数が目的の値となるように、この例の最後にリストされている補助関数 randReplicateFiles を使用してファイルをランダムにオーバーサンプリングします。

rng(0)
[G,classes] = findgroups(labels);
numObservations = splitapply(@numel,labels,G);
desiredNumObservationsPerClass = max(numObservations);
files = splitapply(@(x){randReplicateFiles(x,desiredNumObservationsPerClass)},dsTrain.Files,G);
files = vertcat(files{:});
dsTrain.Files = files;

データ拡張

クラスの不均衡に対処するためにファイルを複製すると、学習データの大部分が同一となるため、ネットワークが過適合になる可能性が高くなります。この影響を打ち消すには、補助関数 transformaugmentPointCloud を使用して学習データにデータ拡張を適用します。このデータ拡張では、点群がランダムに回転され、点がランダムに削除され、ガウス ノイズを使用してランダムなジッターが点に与えられます。

dsTrain = transform(dsTrain,@augmentPointCloud);

拡張された学習サンプルのいずれかをプレビューします。

data = preview(dsTrain);
ptCloud = data{1,1};
label = data{1,2};

figure
pcshow(ptCloud.Location,[0 1 0],"MarkerSize",40,"VerticalAxisDir","down")
xlabel("X")
ylabel("Y")
zlabel("Z")
title(label)

学習済みネットワークのパフォーマンスを測定するためのデータは、元のデータセットを表すものでなければならないため、データ拡張は検証データまたはテスト データには適用されないことに注意してください。

データ前処理

学習と予測に使用する点群を準備するには、2 つの前処理手順が必要です。

まず、学習中にバッチ処理できるように、各点群から使用する固定の点数を選択します。最適な点の数は、データセットや、オブジェクトの形状を正確に取得するのに必要な点の数によって異なります。適切な点数を選択できるように、クラス当たりの最小点数、最大点数、平均点数を計算します。

minPointCount = splitapply(@min,counts,G);
maxPointCount = splitapply(@max,counts,G);
meanPointCount = splitapply(@(x)round(mean(x)),counts,G);

stats = table(classes,numObservations,minPointCount,maxPointCount,meanPointCount)
stats=14×5 table
       classes        numObservations    minPointCount    maxPointCount    meanPointCount
    ______________    _______________    _____________    _____________    ______________

    4wd                      15               140              1955              751     
    building                 15               193              8455             2708     
    bus                      11               126             11767             2190     
    car                      64                52              2377              528     
    pedestrian              107                20               297              110     
    pillar                   15                80               751              357     
    pole                     15                13               253               90     
    traffic lights           36                38               352              161     
    traffic sign             40                18               736              126     
    tree                     24                53              2953              470     
    truck                     9               445              3013             1376     
    trunk                    42                32               766              241     
    ute                      12                90              1380              580     
    van                      28                91              5809             1125     

クラス当たりの点数は、クラス内およびクラス間で大きく異なるため、すべてのクラスに当てはまる値を選択するのは困難です。経験則は、多くの点を処理することによる計算コストの増加を抑えつつ、十分な数の点を選択してオブジェクトの形状を適切に取得できるようにすることです。この 2 つの側面のトレードオフとして 1024 という値が適しています。経験的解析に基づいて最適な点数を選択することもできます。ただし、これについてはこの例では扱いません。関数 transform を使用して、学習セットと検証セットに含まれる 1024 個の点を選択します。

numPoints = 1024;
dsTrain = transform(dsTrain,@(data)selectPoints(data,numPoints));
dsVal = transform(dsVal,@(data)selectPoints(data,numPoints));

最後の前処理手順では、データ値の範囲が大きく異なることを考慮して、点群データを 0 から 1 の範囲で正規化します。たとえば、LiDAR センサーの近くにあるオブジェクトは、遠くにあるオブジェクトと比べて値が小さくなります。このような違いが存在すると、学習中にネットワークが収束しなくなる可能性があります。transform を使用して、学習セットと検証セットに含まれる点群データを正規化します。

dsTrain = transform(dsTrain,@preprocessPointCloud);
dsVal = transform(dsVal,@preprocessPointCloud);

拡張して前処理した学習データをプレビューします。

data = preview(dsTrain);
figure
pcshow(data{1,1},[0 1 0],"MarkerSize",40,"VerticalAxisDir","down");
xlabel("X")
ylabel("Y")
zlabel("Z")
title(data{1,2})

PointNet モデルの定義

PointNet 分類モデルは、2 つのコンポーネントで構成されています。1 つ目のコンポーネントは、学習によって疎な点群データを密な特徴ベクトルに符号化する点群符号化器です。2 つ目のコンポーネントは、符号化された各点群のカテゴリカル クラスを予測する分類器です。

PointNet 符号化器モデルは、さらに 4 つのモデルで構成され、その後ろに max 演算が続きます。

  1. 入力変換モデル

  2. 共有 MLP モデル

  3. 特徴変換モデル

  4. 共有 MLP モデル

共有 MLP モデルは、一連の畳み込み演算、バッチ正規化演算、および ReLU 演算を使用して実装されます。畳み込み演算は、重みが入力点群全体で共有されるように構成されます。変換モデルは、共有 MLP、および各点群に適用される学習可能な変換行列で構成されています。共有 MLP と max 演算によって、PointNet 符号化器が点の処理順序に左右されなくなり、変換モデルが向きの変化に左右されなくなります。

PointNet 符号化器のモデル パラメーターの定義

共有 MLP と変換モデルは、入力チャネルの数と隠れチャネルのサイズによってパラメーター化されます。この例で使用されている値は、Sydney Urban Objects データセットに関するこれらのハイパーパラメーターの調整によって選択されたものです。PointNet を別のデータセットに適用する場合は、ハイパーパラメーターの調整を別途行わなければならないことに注意してください。

入力変換モデルについて、入力チャネルのサイズを 3 に設定し、隠れチャネルのサイズを 64、128、256 に設定し、この例の最後にリストされている補助関数 initializeTransform を使用してモデル パラメーターを初期化します。

inputChannelSize = 3;
hiddenChannelSize1 = [64,128];
hiddenChannelSize2 = 256;
[parameters.InputTransform, state.InputTransform] = initializeTransform(inputChannelSize,hiddenChannelSize1,hiddenChannelSize2);

最初の共有 MLP モデルについて、入力チャネルのサイズを 3 に設定し、隠れチャネルのサイズを 64 に設定し、この例の最後にリストされている補助関数 initializeSharedMLP を使用してモデル パラメーターを初期化します。

inputChannelSize = 3;
hiddenChannelSize = [64 64];
[parameters.SharedMLP1,state.SharedMLP1] = initializeSharedMLP(inputChannelSize,hiddenChannelSize);

特徴変換モデルについて、入力チャネルのサイズを 64 に設定し、隠れチャネルのサイズを 64、128、256 に設定し、この例の最後にリストされている補助関数 initializeTransform を使用してモデル パラメーターを初期化します。

inputChannelSize = 64;
hiddenChannelSize1 = [64,128];
hiddenChannelSize2 = 256;
[parameters.FeatureTransform, state.FeatureTransform] = initializeTransform(inputChannelSize,hiddenChannelSize,hiddenChannelSize2);

2 番目の共有 MLP モデルについて、入力チャネルのサイズを 64 に設定し、隠れチャネルのサイズを 64 に設定し、この例の最後にリストされている関数 initializeSharedMLP を使用してモデル パラメーターを初期化します。

inputChannelSize = 64;
hiddenChannelSize = 64;
[parameters.SharedMLP2,state.SharedMLP2] = initializeSharedMLP(inputChannelSize,hiddenChannelSize);

PointNet 分類器のモデル パラメーターの定義

PointNet 分類器モデルは、共有 MLP、全結合演算、およびソフトマックス活性化で構成されています。分類器モデルについて、入力のサイズを 64 に設定し、隠れチャネルのサイズを 512 および 256 に設定し、この例の最後にリストされている補助関数 initalizeClassifier を使用してモデル パラメーターを初期化します。

inputChannelSize = 64;
hiddenChannelSize = [512,256];
numClasses = numel(classes);
[parameters.ClassificationMLP, state.ClassificationMLP] = initializeClassificationMLP(inputChannelSize,hiddenChannelSize,numClasses);

PointNet 関数の定義

この例の最後のモデル関数の節にリストされている関数 pointnetClassifier を作成し、PointNet モデルの出力を計算します。この関数モデルは、点群データ、学習可能なモデル パラメーター、モデルの状態、およびモデルが学習と予測どちらの出力を返すかを指定するフラグを入力に取ります。このネットワークは、入力点群を分類するための予測を返します。

モデル勾配関数の定義

この例のモデル勾配関数の節にリストされている関数 modelGradients を作成します。この関数は、モデル パラメーター、モデルの状態、および入力データのミニバッチを入力に取り、モデルの学習可能なパラメーターと対応する損失の損失勾配を返します。

学習オプションの指定

学習を 10 エポック行い、バッチ サイズ 128 でデータを読み込みます。初期学習率を 0.002 に設定し、L2 正則化係数を 0.01 に設定します。

numEpochs = 10;
learnRate = 0.002;
miniBatchSize = 128;
l2Regularization = 0.01;
learnRateDropPeriod = 15;
learnRateDropFactor = 0.5;

Adam 最適化のオプションを初期化します。

gradientDecayFactor = 0.9;
squaredGradientDecayFactor = 0.999;

PointNet の学習

カスタム学習ループを使用してモデルに学習させます。

学習の開始時に、データをシャッフルします。

それぞれの反復で次を行います。

  • データのバッチを読み取ります。

  • モデル勾配を評価します。

  • L2 重み正則化を適用します。

  • adamupdate を使用してモデル パラメーターを更新します。

  • 学習の進行状況プロットを更新します。

各エポックの最後で、検証データセットと照らし合わせてモデルを評価し、混同メトリクスを収集して学習の進行に伴う分類精度を測定します。

learnRateDropPeriod で指定されたエポックが完了した後に、learnRateDropFactor で指定された係数だけ学習率を減少させます。

Adam オプティマイザーで使用されるパラメーターの勾配と要素単位の勾配の 2 乗の移動平均を初期化します。

avgGradients = [];
avgSquaredGradients = [];

doTraining が true の場合、モデルに学習させます。そうでない場合は、事前学習済みのネットワークを読み込みます。

この学習は、12 GB の GPU メモリ搭載の NVIDIA Titan X で検証済みであることに注意してください。GPU のメモリがこれより少ない場合、学習時にメモリ不足が発生する可能性があります。これが発生した場合は、miniBatchSize を小さくします。このネットワークに学習させるには、約 5 分かかります。GPU ハードウェアによっては、時間がかかる場合があります。

doTraining = false;

if doTraining
    
    % Create a minibatchqueue to batch data from training and validation
    % datastores. Use the batchData function, listed at the end of the
    % example, to batch the point cloud data and one-hot encode the label 
    % data.
    numOutputsFromDSRead = 2;
    mbqTrain = minibatchqueue(dsTrain,numOutputsFromDSRead,...
        "MiniBatchSize", miniBatchSize,...
        "MiniBatchFcn",@batchData,...
        "MiniBatchFormat",["SCSB" "BC"]);
    
    mbqVal = minibatchqueue(dsVal,numOutputsFromDSRead,...
        "MiniBatchSize", miniBatchSize,... 
        "MiniBatchFcn",@batchData,...
        "MiniBatchFormat",["SCSB" "BC"]);
 
    % Use the configureTrainingProgressPlot function, listed at the end of the
    % example, to initialize the training progress plot to display the training
    % loss, training accuracy, and validation accuracy.
    [lossPlotter, trainAccPlotter,valAccPlotter] = initializeTrainingProgressPlot;
    
    numClasses = numel(classes);
    iteration = 0;
    start = tic;
    for epoch = 1:numEpochs
        
        % Shuffle data every epoch.
        shuffle(mbqTrain);
      
        % Iterate through data set.
        while hasdata(mbqTrain)
            iteration = iteration + 1;
            
            % Read next batch of training data.
            [XTrain, YTrain] = next(mbqTrain);            
            
            % Evaluate the model gradients and loss using dlfeval and the
            % modelGradients function.
            [gradients, loss, state, acc] = dlfeval(@modelGradients,XTrain,YTrain,parameters,state);
            
            % L2 regularization.
            gradients = dlupdate(@(g,p) g + l2Regularization*p,gradients,parameters);
            
            % Update the network parameters using the Adam optimizer.
            [parameters, avgGradients, avgSquaredGradients] = adamupdate(parameters, gradients, ...
                avgGradients, avgSquaredGradients, iteration,...
                learnRate,gradientDecayFactor, squaredGradientDecayFactor);
            
            % Update the training progress.
            D = duration(0,0,toc(start),"Format","hh:mm:ss");
            title(lossPlotter.Parent,"Epoch: " + epoch + ", Elapsed: " + string(D))
            addpoints(lossPlotter,iteration,double(gather(extractdata(loss))))
            addpoints(trainAccPlotter,iteration,acc);
            drawnow
        end
        
        % Evaluate the model on validation data.
        cmat = sparse(numClasses,numClasses);
        while hasdata(mbqVal)
            
            % Read next batch of validation data.
            [XVal, YVal] = next(mbqVal);

            % Compute label predictions.
            isTraining = false;
            YPred = pointnetClassifier(XVal,parameters,state,isTraining);
            
            % Choose prediction with highest score as the class label for
            % XTest.
            [~,YValLabel] = max(YVal,[],1);
            [~,YPredLabel] = max(YPred,[],1);
            
            % Collect confusion metrics.
            cmat = aggreateConfusionMetric(cmat,YValLabel,YPredLabel);
        end
        
        % Update training progress plot with average classification accuracy.
        acc = sum(diag(cmat))./sum(cmat,"all");
        addpoints(valAccPlotter,iteration,acc);
        
        % Upate the learning rate.
        if mod(epoch,learnRateDropPeriod) == 0
            learnRate = learnRate * learnRateDropFactor;
        end
        
        % Reset training and validation data queues.
        reset(mbqTrain);
        reset(mbqVal);
    end

else
    % Download pretrained model parameters, model state, and validation
    % results.
    pretrainedURL = 'https://ssd.mathworks.com/supportfiles/vision/data/pointnetSydneyUrbanObjects.zip'; 
    pretrainedResults = downloadPretrainedPointNet(pretrainedURL);
   
    parameters = pretrainedResults.parameters;
    state = pretrainedResults.state;
    cmat = pretrainedResults.cmat;
    
    % Move model parameters to the GPU if possible and convert to a dlarray.
    parameters = prepareForPrediction(parameters,@(x)dlarray(toDevice(x,canUseGPU)));
    
    % Move model state to the GPU if possible.
    state = prepareForPrediction(state,@(x)toDevice(x,canUseGPU));
end

検証用の混同行列を表示します。

figure
chart = confusionchart(cmat,classes);

学習精度と検証精度の平均値を計算します。

acc = sum(diag(cmat))./sum(cmat,"all")
acc = 0.5742

Sydney Urban Objects データセットは学習サンプルの数が少ないため、検証精度を 60% より高くすることは困難です。このモデルは、補助関数 augmentPointCloudData で定義された拡張を行わない場合、学習データに簡単に過適合してしまいます。PointNet 分類器のロバスト性を向上させるには、さらに学習を行う必要があります。

PointNet を使用した点群データの分類

pcread を使用して点群データを読み込み、学習で使用するのと同じ関数を使用して点群を前処理し、その結果を dlarray に変換します。

ptCloud = pcread("car.pcd");
X = preprocessPointCloud(ptCloud);
dlX = dlarray(X{1},"SCSB");

モデル関数 pointnetClassifier を使用して、点群のラベルを予測します。

YPred = pointnetClassifier(dlX,parameters,state,false);
[~,classIdx] = max(YPred,[],1);

点群、および予測されたラベルのうちスコアが最も高いラベルを表示します。

figure
pcshow(ptCloud.Location,[0 1 0],"MarkerSize",40,"VerticalAxisDir","down")
title(classes(classIdx))

モデル勾配関数

関数 modelGradients は、データ dlX のミニバッチ、対応するターゲット dlY、および学習可能なパラメーターを入力に取り、学習可能なパラメーターと対応する損失の損失勾配を返します。この損失には、PointNet 符号化器によって予測された特徴変換行列をほぼ直交にするために設計された正則化項が含まれています。勾配を計算するため、学習ループの中で関数 dlfeval を使用して、関数 modelGradients を評価します。

function [gradients, loss, state, acc] = modelGradients(X,Y,parameters,state)

% Execute the model function.
isTraining = true;
[YPred,state,dlT] = pointnetClassifier(X,parameters,state,isTraining);

% Add regularization term to ensure feature transform matrix is
% approximately orthogonal.
K = size(dlT,1);
B = size(dlT, 4);
I = repelem(eye(K),1,1,1,B);
dlI = dlarray(I,"SSCB");
treg = mse(dlI,pagemtimes(dlT,permute(dlT,[2 1 3 4])));
factor = 0.001;

% Compute the loss.
loss = crossentropy(YPred,Y) + factor*treg;

% Compute the parameter gradients with respect to the loss. 
gradients = dlgradient(loss, parameters);

% Compute training accuracy metric.
[~,YTest] = max(Y,[],1);
[~,YPred] = max(YPred,[],1);
acc = gather(extractdata(sum(YTest == YPred)./numel(YTest)));

end

PointNet 分類器関数

関数 pointnetClassifier は、点群データ dlX、学習可能なモデル パラメーター、モデルの状態、およびモデルが学習と予測どちらの出力を返すかを指定するフラグ isTraining を入力に取ります。その後、この関数は、PointNet の符号化器と多層パーセプトロンを呼び出し、分類特徴量を抽出します。学習中は、各パーセプトロン演算が完了した後に、ドロップアウトが適用されます。最後のパーセプトロンが完了した後、fullyconnect 演算によって分類特徴量とクラス数がマッピングされ、ソフトマックス活性化によって出力がラベルの確率分布に正規化されます。確率分布、更新されたモデル状態、および PointNet 符号化器によって予測された特徴変換行列が出力として返されます。

function [dlY,state,dlT] = pointnetClassifier(dlX,parameters,state,isTraining)

% Invoke the PointNet encoder.
[dlY,state,dlT] = pointnetEncoder(dlX,parameters,state,isTraining);

% Invoke the classifier.
p = parameters.ClassificationMLP.Perceptron;
s = state.ClassificationMLP.Perceptron;
for k = 1:numel(p) 
     
    [dlY, s(k)] = perceptron(dlY,p(k),s(k),isTraining);
      
    % If training, apply inverted dropout with a probability of 0.3.
    if isTraining
        probability = 0.3; 
        dropoutScaleFactor = 1 - probability;
        dropoutMask = ( rand(size(dlY), "like", dlY) > probability ) / dropoutScaleFactor;
        dlY = dlY.*dropoutMask;
    end
    
end
state.ClassificationMLP.Perceptron = s;

% Apply final fully connected and softmax operations.
weights = parameters.ClassificationMLP.FC.Weights;
bias = parameters.ClassificationMLP.FC.Bias;
dlY = fullyconnect(dlY,weights,bias);
dlY = softmax(dlY);
end

PointNet 符号化器関数

関数 pointnetEncoder は、入力変換、共有 MLP、特徴変換、2 番目の共有 MLP、および max 演算を使用して入力 dlX を処理し、max 演算の結果を返します。

function [dlY,state,T] = pointnetEncoder(dlX,parameters,state,isTraining)
% Input transform.
[dlY,state.InputTransform] = dataTransform(dlX,parameters.InputTransform,state.InputTransform,isTraining);

% Shared MLP.
[dlY,state.SharedMLP1.Perceptron] = sharedMLP(dlY,parameters.SharedMLP1.Perceptron,state.SharedMLP1.Perceptron,isTraining);

% Feature transform.
[dlY,state.FeatureTransform,T] = dataTransform(dlY,parameters.FeatureTransform,state.FeatureTransform,isTraining);

% Shared MLP.
[dlY,state.SharedMLP2.Perceptron] = sharedMLP(dlY,parameters.SharedMLP2.Perceptron,state.SharedMLP2.Perceptron,isTraining);

% Max operation.
dlY = max(dlY,[],1);
end

共有多層パーセプトロン関数

共有多層パーセプトロン関数 perceptron は、一連のパーセプトロン演算を使用して入力 dlX を処理し、最後のパーセプトロンの結果を返します。

function [dlY,state] = sharedMLP(dlX,parameters,state,isTraining)
dlY = dlX;
for k = 1:numel(parameters) 
    [dlY, state(k)] = perceptron(dlY,parameters(k),state(k),isTraining);
end
end

パーセプトロン関数

パーセプトロン関数は、畳み込み、バッチ正規化、ReLU 演算を使用して入力 dlX を処理し、ReLU 演算の出力を返します。

function [dlY,state] = perceptron(dlX,parameters,state,isTraining)
% Convolution.
W = parameters.Conv.Weights;
B = parameters.Conv.Bias;
dlY = dlconv(dlX,W,B);

% Batch normalization. Update batch normalization state when training.
offset = parameters.BatchNorm.Offset;
scale = parameters.BatchNorm.Scale;
trainedMean = state.BatchNorm.TrainedMean;
trainedVariance = state.BatchNorm.TrainedVariance;
if isTraining
    [dlY,trainedMean,trainedVariance] = batchnorm(dlY,offset,scale,trainedMean,trainedVariance);
    
    % Update state.
    state.BatchNorm.TrainedMean = trainedMean;
    state.BatchNorm.TrainedVariance = trainedVariance;
else
    dlY = batchnorm(dlY,offset,scale,trainedMean,trainedVariance);
end

% ReLU.
dlY = relu(dlY);
end

データ変換関数

関数 dataTransform は、共有 MLP、max 演算、およびもうひとつの共有 MLP を使用して入力 dlX を処理し、変換行列 T を予測します。この変換行列は、バッチ行列乗算演算を使用して入力 dlX に適用されます。この関数は、バッチ行列乗算の結果と変換行列を返します。

function [dlY,state,T] = dataTransform(dlX,parameters,state,isTraining)

% Shared MLP.
[dlY,state.Block1.Perceptron] = sharedMLP(dlX,parameters.Block1.Perceptron,state.Block1.Perceptron,isTraining);

% Max operation.
dlY = max(dlY,[],1);

% Shared MLP.
[dlY,state.Block2.Perceptron] = sharedMLP(dlY,parameters.Block2.Perceptron,state.Block2.Perceptron,isTraining);

% Transform net (T-Net). Apply last fully connected operation as W*X to
% predict tranformation matrix T.
dlY = squeeze(dlY); % N-by-B
T = parameters.Transform * stripdims(dlY); % K^2-by-B

% Reshape T into a square matrix.
K = sqrt(size(T,1));
T = reshape(T,K,K,1,[]); % [K K 1 B]
T = T + eye(K);

% Apply to input dlX using batch matrix multiply. 
[C,B] = size(dlX,[3 4]);  % [M 1 K B]
dlX = reshape(dlX,[],C,1,B); % [M K 1 B]
Y = pagemtimes(dlX,T);
dlY = dlarray(Y,"SCSB");
end

モデル パラメーター初期化関数

関数 initializeTransform

関数 initializeTransform は、2 つの共有 MLP のチャネル サイズと隠れチャネルの数を入力として受け取り、初期化されたパラメーターを struct で返します。このパラメーターは、He の重み初期化を使用して初期化されます [3]。

function [params,state] = initializeTransform(inputChannelSize,block1,block2)
[params.Block1,state.Block1] = initializeSharedMLP(inputChannelSize,block1);
[params.Block2,state.Block2] = initializeSharedMLP(block1(end),block2);

% Parameters for the transform matrix.
params.Transform = dlarray(zeros(inputChannelSize^2,block2(end)));
end

関数 initializeSharedMLP

関数 initializeSharedMLP は、チャネル サイズと隠れチャネル サイズを入力として受け取り、初期化されたパラメーターを struct で返します。このパラメーターは、He の重み初期化を使用して初期化されます。

function [params,state] = initializeSharedMLP(inputChannelSize,hiddenChannelSize)
weights = initializeWeightsHe([1 1 inputChannelSize hiddenChannelSize(1)]);
bias = zeros(hiddenChannelSize(1),1,"single");
p.Conv.Weights = dlarray(weights);
p.Conv.Bias = dlarray(bias);

p.BatchNorm.Offset = dlarray(zeros(hiddenChannelSize(1),1,"single"));
p.BatchNorm.Scale = dlarray(ones(hiddenChannelSize(1),1,"single"));

s.BatchNorm.TrainedMean = zeros(hiddenChannelSize(1),1,"single");
s.BatchNorm.TrainedVariance = ones(hiddenChannelSize(1),1,"single");

params.Perceptron(1) = p;
state.Perceptron(1) = s;

for k = 2:numel(hiddenChannelSize)
    weights = initializeWeightsHe([1 1 hiddenChannelSize(k-1) hiddenChannelSize(k)]);
    bias = zeros(hiddenChannelSize(k),1,"single");
    p.Conv.Weights = dlarray(weights);
    p.Conv.Bias = dlarray(bias);
    
    p.BatchNorm.Offset = dlarray(zeros(hiddenChannelSize(k),1,"single"));
    p.BatchNorm.Scale = dlarray(ones(hiddenChannelSize(k),1,"single"));
    
    s.BatchNorm.TrainedMean = zeros(hiddenChannelSize(k),1,"single");
    s.BatchNorm.TrainedVariance = ones(hiddenChannelSize(k),1,"single");

    params.Perceptron(k) = p;
    state.Perceptron(k) = s;
end
end

関数 initializeClassificationMLP

関数 initializeClassificationMLP は、チャネル サイズ、隠れチャネル サイズ、クラスの数を入力として受け取り、初期化されたパラメーターを struct で返します。共有 MLP は、He の重み初期化を使用して初期化され、最終的な全結合演算はランダムなガウス値を使用して初期化されます。

function [params,state] = initializeClassificationMLP(inputChannelSize,hiddenChannelSize,numClasses)
[params,state] = initializeSharedMLP(inputChannelSize,hiddenChannelSize);

weights = initializeWeightsGaussian([numClasses hiddenChannelSize(end)]);
bias = zeros(numClasses,1,"single");
params.FC.Weights = dlarray(weights);
params.FC.Bias = dlarray(bias);
end

関数 initializeWeightsHe

関数 initializeWeightsHe は、He の初期化を使用してパラメーターを初期化します。

function x = initializeWeightsHe(sz)
fanIn = prod(sz(1:3));
stddev = sqrt(2/fanIn);
x = stddev .* randn(sz);
end

関数 initializeWeightsGaussian

関数 initializeWeightsGaussian は、ガウスの初期化を使用して標準偏差 0.01 でパラメーターを初期化します。

function x = initializeWeightsGaussian(sz)
x = randn(sz,"single") .* 0.01;
end

データ前処理関数

関数 preprocessPointCloudData

関数 preprocessPointCloudData は、入力データから X、Y、Z の点データを抽出し、そのデータを 0 ~ 1 の範囲で正規化します。この関数は、正規化された X、Y、Z データを返します。

function data = preprocessPointCloud(data)

if ~iscell(data)
    data = {data};
end

numObservations = size(data,1);
for i = 1:numObservations
    % Scale points between 0 and 1.
    xlim = data{i,1}.XLimits;
    ylim = data{i,1}.YLimits;
    zlim = data{i,1}.ZLimits;
    
    xyzMin = [xlim(1) ylim(1) zlim(1)];
    xyzDiff = [diff(xlim) diff(ylim) diff(zlim)];
    
    data{i,1} = (data{i,1}.Location - xyzMin) ./ xyzDiff;
end
end

関数 selectPoints

関数 selectPoints は、目的の数だけ点をサンプリングします。目的の数より多い点が点群に含まれている場合、この関数は pcdownsample を使用してランダムに点を選択します。そうでない場合、この関数はデータを複製して目的の数の点を生成します。

function data = selectPoints(data,numPoints) 
% Select the desired number of points by downsampling or replicating
% point cloud data.
numObservations = size(data,1);
for i = 1:numObservations    
    ptCloud = data{i,1};
    if ptCloud.Count > numPoints
        percentage = numPoints/ptCloud.Count;
        data{i,1} = pcdownsample(ptCloud,"random",percentage);   
    else    
        replicationFactor = ceil(numPoints/ptCloud.Count);
        ind = repmat(1:ptCloud.Count,1,replicationFactor);
        data{i,1} = select(ptCloud,ind(1:numPoints));
    end 
end
end

データ拡張関数

関数 augmentPointCloudData は、z 軸を中心とした点群の回転、点の 30% のドロップ、ガウス ノイズによる点の位置へのジッター付加をランダムに実行します。

function data = augmentPointCloud(data)
   
numObservations = size(data,1);
for i = 1:numObservations
    
    ptCloud = data{i,1};
    
    % Rotate the point cloud about "up axis", which is Z for this data set.
    tform = randomAffine3d(...
        "XReflection", true,...
        "YReflection", true,...
        "Rotation",@randomRotationAboutZ);
    
    ptCloud = pctransform(ptCloud,tform);
    
    % Randomly drop out 30% of the points.
    if rand > 0.5
        ptCloud = pcdownsample(ptCloud,'random',0.3);
    end
    
    if rand > 0.5
        % Jitter the point locations with Gaussian noise with a mean of 0 and 
        % a standard deviation of 0.02 by creating a random displacement field.
        D = 0.02 * randn(size(ptCloud.Location));
        ptCloud = pctransform(ptCloud,D);   
    end
    
    data{i,1} = ptCloud;
end
end

function [rotationAxis,theta] = randomRotationAboutZ()
rotationAxis = [0 0 1];
theta = 360*rand;
end

サポート関数

関数 aggregateConfusionMetric

関数 aggregateConfusionMetric は、予測される結果 YPred と期待される結果 YTest をインクリメンタルに混同行列に追加します。

function cmat = aggreateConfusionMetric(cmat,YTest,YPred)
YTest = gather(extractdata(YTest));
YPred = gather(extractdata(YPred));
[m,n] = size(cmat);
cmat = cmat + full(sparse(YTest,YPred,1,m,n));
end

関数 initializeTrainingProgressPlot

関数 initializeTrainingProgressPlot は、学習損失、学習精度、および検証精度を表示する 2 つのプロットを構成します。

function [plotter,trainAccPlotter,valAccPlotter] = initializeTrainingProgressPlot()
% Plot the loss, training accuracy, and validation accuracy.
figure

% Loss plot
subplot(2,1,1)
plotter = animatedline;
xlabel("Iteration")
ylabel("Loss")

% Accuracy plot
subplot(2,1,2)
trainAccPlotter = animatedline('Color','b');
valAccPlotter = animatedline('Color','g');
legend('Training Accuracy','Validation Accuracy','Location','northwest');
xlabel("Iteration")
ylabel("Accuracy")
end

関数 replicateFiles

関数 replicateFiles は、一連のファイルをランダムにオーバーサンプリングし、numDesired 個の要素をもつ一連のファイルを返します。

function files = randReplicateFiles(files,numDesired)
n = numel(files);
ind = randi(n,numDesired,1);
files = files(ind);
end

関数 downloadSydneyUrbanObjects

関数 downloadSydneyUrbanObjects は、データセットをダウンロードして一時ディレクトリに保存します。

function datapath = downloadSydneyUrbanObjects(dataLoc)

if nargin == 0
    dataLoc = pwd;
end

dataLoc = string(dataLoc);

url = "http://www.acfr.usyd.edu.au/papers/data/";
name = "sydney-urban-objects-dataset.tar.gz";

datapath = fullfile(dataLoc,'sydney-urban-objects-dataset');
if ~exist(datapath,'dir')
    disp('Downloading Sydney Urban Objects data set...');
    untar(url+name,dataLoc);
end

end

関数 loadSydneyUrbanObjectsData

関数 The loadSydneyUrbanObjectsData は、Sydney Urban Objects データ セットから点群とラベル データを読み込むためのデータストアを作成します。

function ds = loadSydneyUrbanObjectsData(datapath,folds)

if nargin == 0
    return;
end

if nargin < 2
    folds = 1:4;
end

datapath = string(datapath);
path = fullfile(datapath,'objects',filesep);

% Add folds to datastore.
foldNames{1} = importdata(fullfile(datapath,'folds','fold0.txt'));
foldNames{2} = importdata(fullfile(datapath,'folds','fold1.txt'));
foldNames{3} = importdata(fullfile(datapath,'folds','fold2.txt'));
foldNames{4} = importdata(fullfile(datapath,'folds','fold3.txt'));
names = foldNames(folds);
names = vertcat(names{:});

fullFilenames = append(path,names);
ds = fileDatastore(fullFilenames,'ReadFcn',@extractTrainingData,'FileExtensions','.bin');
end

関数 batchData

関数 batchData は、データを順番に並べてバッチを作成し、GPU で処理できるようにそのデータを GPU に移動します。

function [X,Y] = batchData(ptCloud,labels)
X = cat(4,ptCloud{:});
labels = cat(1,labels{:});
Y = onehotencode(labels,2);
end

関数 extractTrainingData

関数 extractTrainingData は、Sydney Urban Objects データ セットから点群とラベル データを抽出します。

function dataOut = extractTrainingData(fname)

[pointData,intensity] = readbin(fname);

[~,name] = fileparts(fname);
name = string(name);
name = extractBefore(name,'.');
name = replace(name,'_',' ');

labelNames = ["4wd","building","bus","car","pedestrian","pillar",...
    "pole","traffic lights","traffic sign","tree","truck","trunk","ute","van"];

label = categorical(name,labelNames);

dataOut = {pointCloud(pointData,'Intensity',intensity),label};

end

関数 readbin

関数 readbin は、Sydney Urban Object バイナリ ファイルから点群データを読み取ります。

function [pointData,intensity] = readbin(fname)
% readbin Read point and intensity data from Sydney Urban Object binary
% files.

% names = ['t','intensity','id',...
%          'x','y','z',...
%          'azimuth','range','pid']
% 
% formats = ['int64', 'uint8', 'uint8',...
%            'float32', 'float32', 'float32',...
%            'float32', 'float32', 'int32']

fid = fopen(fname, 'r');
c = onCleanup(@() fclose(fid));

fseek(fid,10,-1); % Move to the first X point location 10 bytes from beginning
X = fread(fid,inf,'single',30);
fseek(fid,14,-1);
Y = fread(fid,inf,'single',30);
fseek(fid,18,-1);
Z = fread(fid,inf,'single',30);

fseek(fid,8,-1);
intensity = fread(fid,inf,'uint8',33);

pointData = [X,Y,Z];
end

関数 downloadPretrainedPointNet

関数 downloadPretrainedPointNet は、事前学習済みの PointNet モデルをダウンロードします。

function data = downloadPretrainedPointNet(pretrainedURL)
% Download and load a pretrained pointnet model.
if ~exist('pointnetSydneyUrbanObjects.mat', 'file')
    if ~exist('pointnetSydneyUrbanObjects.zip', 'file')
        disp('Downloading pretrained detector (5 MB)...');
        websave('pointnetSydneyUrbanObjects.zip', pretrainedURL);
    end
    unzip('pointnetSydneyUrbanObjects.zip');
end
data = load("pointnetSydneyUrbanObjects.mat");
end

関数 prepareForPrediction

関数 prepareForPrediction は、入れ子にされた構造体データにユーザー定義関数を適用するのに使用されます。これは、モデル パラメーターと状態データを GPU に移動するのに使用されます。

function p = prepareForPrediction(p,fcn)

for i = 1:numel(p)
    p(i) = structfun(@(x)invoke(fcn,x),p(i),'UniformOutput',0);
end

    function data = invoke(fcn,data)
        if isstruct(data)
            data = prepareForPrediction(data,fcn);
        else
            data = fcn(data);
        end
    end
end

% Move data to the GPU.
function x = toDevice(x,useGPU)
if useGPU
    x = gpuArray(x);
end
end

参考文献

[1] Charles, R. Qi, Hao Su, Mo Kaichun, and Leonidas J. Guibas. "PointNet: Deep Learning on Point Sets for 3D Classification and Segmentation." In 2017 IEEE Conference on Computer Vision and Pattern Recognition (CVPR), 77–85. Honolulu, HI: IEEE, 2017. https://doi.org/10.1109/CVPR.2017.16.

[2] de Deuge, Mark, Alastair Quadras, Calvin Hung, and Bertrand Douillard. "Unsupervised Feature Learning for Classification of Outdoor 3D Scans." In Australasian Conference on Robotics and Automation 2013 (ACRA 13). Sydney, Australia: ACRA, 2013.

[3] He, Kaiming, Xiangyu Zhang, Shaoqing Ren, and Jian Sun. "Delving Deep into Rectifiers: Surpassing Human-Level Performance on ImageNet Classification." In 2015 IEEE International Conference on Computer Vision (ICCV), 1026–34. Santiago, Chile: IEEE, 2015. https://doi.org/10.1109/ICCV.2015.123.

関連するトピック