Простой пример/прецедент для BNT gaussian_CPD?

В конце концов я понял это, поэкспериментировав с BNT в командной строке MATLAB. Вот как я определил свою сеть классификатора, используя узлы gaussian_CPD:

DAG = false(5); DAG(1, 2:5) = true
bnet = mk_bnet(DAG, [2 1 1 2 1], 'discrete', 1);
bnet.CPD{1} = tabular_CPD(bnet, 1, 'prior_type', 'dirichlet');
for node = 2:5
   bnet.CPD{node} = gaussian_CPD(bnet, node);
end
bnet

DAG =

     0     1     1     1     1
     0     0     0     0     0
     0     0     0     0     0
     0     0     0     0     0
     0     0     0     0     0

bnet = 

               equiv_class: [1 2 3 4 5]
                    dnodes: 1
                  observed: []
                     names: {}
                    hidden: [1 2 3 4 5]
               hidden_bitv: [1 1 1 1 1]
                       dag: [5x5 logical]
                node_sizes: [2 1 1 2 1]
                    cnodes: [2 3 4 5]
                   parents: {[1x0 double]  [1]  [1]  [1]  [1]}
    members_of_equiv_class: {[1]  [2]  [3]  [4]  [5]}
                       CPD: {[1x1 tabular_CPD]  [1x1 gaussian_CPD]  [1x1 gaussian_CPD]  [1x1 gaussian_CPD]  [1x1 gaussian_CPD]}
             rep_of_eclass: [1 2 3 4 5]
                     order: [1 5 4 3 2]

Чтобы обучить его, я использовал свой оригинальный классификатор, который помог мне пометить набор из 300 образцов, а затем просто прогнал 2/3 из них через обучающий алгоритм.

bnet = learn_params(bnet, lsamples);
CPD = struct(bnet.CPD{1}); % Peek inside CPD{1}
dispcpt(CPD.CPT);

1 : 0.6045 
2 : 0.3955

Выходные данные из dispcpt дают приблизительное представление о разбивке между назначениями классов в помеченных образцах в обучающем наборе.

Чтобы протестировать новый классификатор, я провел последнюю 1/3 результатов как по исходной, так и по новой байесовской сети. Вот код, который я использовал для новой сети:

engine = jtree_inf_engine(bnet);
evidence = cell(1, 5);
tresults = cell(3, length(tsamples));
tresults(3, :) = tsamples(1, :);
for i = 1:length(tsamples)
    evidence(2:5) = tsamples(2:5, i);
    marginal = marginal_nodes(enter_evidence(engine, evidence), 1);
    tresults{1, i} = find(marginal.T == max(marginal.T)); % Generic decision point
    tresults{2, i} = marginal.T(1);
end
tresults(:, 1:8)

ans = 

    [         2]    [     1]    [         2]    [         2]    [         2]    [     1]    [     1]    [     1]
    [1.8437e-10]    [0.9982]    [3.3710e-05]    [3.8349e-04]    [2.2995e-11]    [0.9997]    [0.9987]    [0.5116]
    [         2]    [     1]    [         2]    [         2]    [         2]    [     1]    [     1]    [     2]

Затем, чтобы выяснить, есть ли какое-либо улучшение, я построил наложенные ROC-диаграммы. Как оказалось, моя первоначальная сеть работала достаточно хорошо, поэтому было трудно с уверенностью сказать, работала ли обученная сеть, использующая гауссовские CPD, лучше. Печать областей под ROC-кривыми показала, что новая сеть действительно работает немного лучше. (base area — исходная сеть, а area — новая.)

conf = cell2mat(tresults(2,:));
hit = cell2mat(tresults(3,:)) == 1;
[~, ~, basearea] = plotROC(baseconf, basehit, 'r')
hold all;
[~, ~, area] = plotROC(conf, hit, 'b')
hold off;

basearea =

    0.9371

area =

    0.9555

Я публикую это здесь, чтобы в следующий раз, когда мне нужно будет это сделать, я смог найти ответ... Надеюсь, кто-то еще найдет это полезным.

Ниже я представляю полный пример, который иллюстрирует, как построить наивную байесовскую сеть с помощью BNT Toolbox. Я использую часть набора данных cars. Он содержит как дискретные, так и непрерывные атрибуты.

Просто для удобства я использую пару функций, для которых требуется набор инструментов Statistical.

Начнем с подготовки набора данных:

%# load dataset
D = load('carsmall');

%# keep only features of interest
D = rmfield(D, {'Mfg','Horsepower','Displacement','Model'});

%# filter the rows to keep only two classes
idx = ismember(D.Origin, {'USA' 'Japan'});
D = structfun(@(x)x(idx,:), D, 'UniformOutput',false);
numInst = sum(idx);

%# replace missing values with mean
D.MPG(isnan(D.MPG)) = nanmean(D.MPG);

%# convert discrete attributes to numeric indices 1:mx
[D.Origin,~,gnOrigin] = grp2idx( cellstr(D.Origin) );
[D.Cylinders,~,gnCylinders] = grp2idx( D.Cylinders );
[D.Model_Year,~,gnModel_Year] = grp2idx( D.Model_Year );

Далее мы строим нашу графическую модель:

%# info about the nodes
nodeNames = fieldnames(D);
numNodes = numel(nodeNames);
node = [nodeNames num2cell((1:numNodes)')]';
node = struct(node{:});
dNodes = [node.Origin node.Cylinders node.Model_Year];
cNodes = [node.MPG node.Weight node.Acceleration];
depNodes = [node.MPG node.Cylinders node.Weight ...
            node.Acceleration node.Model_Year];

vals = cell(1,numNodes);
vals(dNodes) = cellfun(@(f) unique(D.(f)), nodeNames(dNodes), 'Uniform',false);
nodeSize = ones(1,numNodes);
nodeSize(dNodes) = cellfun(@numel, vals(dNodes));

%# DAG
dag = false(numNodes);
dag(node.Origin, depNodes) = true;

%# create naive bayes net
bnet = mk_bnet(dag, nodeSize, 'discrete',dNodes, 'names',nodeNames, ...
    'observed',depNodes);
for i=1:numel(dNodes)
    name = nodeNames{dNodes(i)};
    bnet.CPD{dNodes(i)} = tabular_CPD(bnet, node.(name), ...
        'prior_type','dirichlet');
end
for i=1:numel(cNodes)
    name = nodeNames{cNodes(i)};
    bnet.CPD{cNodes(i)} = gaussian_CPD(bnet, node.(name));
end

%# visualize the graph
[~,~,h] = draw_graph(bnet.dag, nodeNames);
hTxt = h(:,1); hNodes = h(:,2);
set(hTxt(node.Origin), 'FontWeight','bold', 'Interpreter','none')
set(hNodes(node.Origin), 'FaceColor','g')
set(hTxt(depNodes), 'Color','k', 'Interpreter','none')
set(hNodes(depNodes), 'FaceColor','y')

Теперь разбиваем данные на обучение/тестирование:

%# build samples as cellarray
data = num2cell(cell2mat(struct2cell(D)')');

%# split train/test: 1/3 for testing, 2/3 for training
cv = cvpartition(D.Origin, 'HoldOut',1/3);
trainData = data(:,cv.training);
testData = data(:,cv.test);
testData(1,:) = {[]};    %# remove class

Наконец, мы изучаем параметры из обучающего набора и предсказываем класс тестовых данных:

%# training
bnet = learn_params(bnet, trainData);

%# testing
prob = zeros(nodeSize(node.Origin), sum(cv.test));
engine = jtree_inf_engine(bnet);         %# Inference engine
for i=1:size(testData,2)
    [engine,loglik] = enter_evidence(engine, testData(:,i));
    marg = marginal_nodes(engine, node.Origin);
    prob(:,i) = marg.T;

end
[~,pred] = max(prob);
actual = D.Origin(cv.test)';

%# confusion matrix
predInd = full(sparse(1:numel(pred),pred,1));
actualInd = full(sparse(1:numel(actual),actual,1));
conffig(predInd, actualInd);             %# confmat

%# ROC plot and AUC
figure
[~,~,auc] = plotROC(max(prob), pred==actual, 'b')
title(sprintf('Area Under the Curve = %g',auc))
set(findobj(gca, 'type','line'), 'LineWidth',2)

Результаты, достижения:

и мы можем извлечь CPT и среднее/сигму в каждом узле:

cellfun(@(x)dispcpt(struct(x).CPT), bnet.CPD(dNodes), 'Uniform',false)
celldisp(cellfun(@(x)struct(x).mean, bnet.CPD(cNodes), 'Uniform',false))
celldisp(cellfun(@(x)struct(x).cov, bnet.CPD(cNodes), 'Uniform',false))