</pre><img vspace="5" hspace="5" src="lungbayesdemo_09.png"> <h2>References<a name="27"></a></h2>         <p>[1] Pearl J., "Probabilistic Reasoning in Intelligent Systems", Morgan Kaufmann, San Mateo, California, 1988.</p>         <p>[2] Neapolitan R., "Learning Bayesian Networks", Pearson Prentice Hall, Upper Saddle River, New Jersey, 2004.</p>         <p class="footer"><br>            Published with MATLAB&reg; 7.6<br></p>      </div>      <!--##### SOURCE BEGIN #####%% Modeling Lung Cancer Diagnosis Using Bayesian Network Inference
% This demo illustrates a simple Bayesian Network example for exact
% probabilistic inference using Pearl's message-passing algorithm. 

%% Introduction
% Bayesian networks (or belief networks) are probabilistic graphical models
% representing a set of variables and their dependencies.
% The graphical nature of Bayesian networks and the ability of describing
% uncertainty of complex relationships in a compact manner provide a method
% for modelling almost any type of data.

%%
% Consider the following example, representing a simplified model to help
% diagnose the patients arriving at a respiratory clinic. A history of smoking
% has a direct influence on both whether or not a patient has bronchitis
% and whether or not a patient has lung cancer. In turn, the presence or
% absence of lung cancer has direct influence on the results of a chest
% x-ray test. We are interested in doing probabilistic inference involving
% features that are not directly related, and for which the conditional
% probability cannot be readily computed using a simple application of the
% Bayes' theorem.

%% Creating the Bayesian Network 
% A Bayesian network consists of a direct-acyclic graph (DAG) in which
% every node represents a variable and every edge represents a dependency
% between variables. We construct this graph by specifying an adjacency
% matrix where the element on row _i_ and column _j_ contains the number of
% edges directed from node _i_ to node _j_. The variables of the models are
% specified by the graph's nodes: |S| (smoking history), |B| (bronchitis),
% |L| (lung cancer) and |X| (chest x-ray). The variables are discrete and
% can take only two values: true (|t|) or false (|f|). 

%=== setup
adj = [0 1 1 0; 0 0 0 0; 0 0 0 1; 0 0 0 0]; % adjacency matrix
nodeNames = {'S', 'B', 'L', 'X'};           % nodes 
S = 1; B = 2; L = 3; X = 4;                 % node identifiers
n = numel(nodeNames);                       % number of nodes
t = 1; f  = 2;                              % true and false
values = cell(1,n);                         % values assumed by variables
for i = 1:numel(nodeNames)
    values{i} = [t f];                      
end

%%
% In addition to the graph structure, we need to specify the
% parameters of the model, namely the conditional probability
% distribution. For discrete variables, this distribution can be
% represented as a table (Conditional Probability Table, |CPT|), which lists
% the probability that a node takes on each of its value, given the value
% combinations of its parents.

%=== Conditional Probability Table
CPT{S} = [.2 .8];
CPT{B}(:,t) = [.25 .05] ; CPT{B}(:,f) = 1 - CPT{B}(:,t);
% CPT{L}(:,t) = [.03 .0005]; CPT{L}(:,f) = 1 - CPT{L}(:,t);
CPT{L}(:,t) = [.3 .005]; CPT{L}(:,f) = 1 - CPT{L}(:,t);
CPT{X}(:,t) = [.6 .02]; CPT{X}(:,f) = 1 - CPT{X}(:,t);

%% Visualizing the Bayesian Network as a Graph
% We can visualize the network structure using the |biograph|
% object. The properties of nodes and edges can be changed as
% desidered.

%=== draw the network
nodeLabels = {'Smoking', 'Bronchitis', 'Lung Cancer', 'Abnormal Xrays'};
bg = biograph(adj, nodeLabels, 'arrowsize', 4);
set(bg.Nodes, 'shape', 'ellipse');
bgInViewer = view(bg);

%=== save as figure
bgFig = figure;
copyobj(bgInViewer.hgAxes,bgFig)

%=== annotate using the CPT
[xp, xn] = find(adj);     % xp = parent id, xn = node id
pa(xn) = xp;              % parents
pa(1) = 1;                % root is parent of itself

s1 = cell(1,n); s2 = cell(1,n); pos = zeros(n,2);
for i = 2:n
    pos(i,:) = bgInViewer.Nodes(i).Position;
    s1{i} = sprintf('P(%s|%s=t) = %f', nodeNames{i}, nodeNames{pa(i)}, CPT{i}(1,t));
    s2{i} = sprintf('P(%s|%s=f) = %f', nodeNames{i}, nodeNames{pa(i)}, CPT{i}(2,t));
end

pos(1,:) = bgInViewer.Nodes(1).Position; % root
s1{1} = sprintf('P(%s=t) = %f', nodeNames{1}, CPT{1}(1));
s2{1} = ' ';

text(pos(:,1)+2, pos(:,2)-10, s1)
text(pos(:,1)+2, pos(:,2)-15, s2)


%% Initializing the Bayesian Network
% The process of computing the probability distribution of variables given
% specific evidence is called probabilistic inference. By exploiting local
% independencies among nodes, Pearls [1] developed a message-passing algorithm
% for exact inference in singly-connected networks. The algorithm can
% compute the conditional probability of any variable given any set of
% evidence by propagation of beliefs between neighboring nodes. For more
% information about the message-passing algorithm see [2].
% We can create and initiate a Bayesian network for the example under
% consideration as follows:

root = find(sum(adj,1)==0); % root is any node with no parent
[nodes, edges] = bnMsgPassCreate(adj, values, CPT);
[nodes, edges] = bnMsgPassInitiate(nodes, edges, root)

%%
% The algorithm parameters, including the conditional probability of each
% node given the evidence, are stored in the fields of the MATLAB structures
% |nodes| and |edges|. Using the function |customnodedraw|, we can
% visualize the distribution of the conditional probability given an empty
% set of evidence in a series of pie charts, as shown below. 

%=== conditional probability given the empty set []
for i = 1:n
    disp(['P(' nodeNames{i}, '|[]) = ' num2str(nodes(i).P(1))]);
end

%=== assign relevant info to each node handle
nodeHandles = bgInViewer.Nodes;
for i = 1:n
    nodeHandles(i).UserData.Distribution = [nodes(i).P];
end

%=== draw customized nodes
bgInViewer.ShowTextInNodes = 'none';

set(nodeHandles, 'shape','circle')
bgInViewer.CustomNodeDrawFcn = @(node) customnodedraw(node);
%bgInViewer.Scale = .7
bgInViewer.dolayout


%%
% Suppose we are interested in evaluating the likelihood that a patient
% with bronchitis has lung cancer. We instantiate |B=t| (true) and we
% update the network as follows:

%=== inference with B = t
evNode = B; 
evValue = t;
[n1, e1, A1, a1] = bnMsgPassUpdate(nodes, edges, [], [], evNode, evValue);

for i = 1:n
    disp(['P(' nodeNames{i}, '|B=t) = ' num2str(n1(i).P(1))]);
end

%== plot and compare
figure(); subplot(2,1,1); 
x = cat(1,nodes.P);
bar(x, 'stacked'); set(gca, 'xticklabel', nodeNames);
ylabel('Probability');
title('Initialized network with empty evidence set')
legend({'true', 'false'}, 'location', 'SouthEastOutside')
hold on; subplot(2,1,2);
x1 = cat(1,n1.P); 
bar(x1, 'stacked'); set(gca, 'xticklabel', nodeNames);
ylabel('Probability');
title('Updated network with evidence B=true')
legend({'true', 'false'}, 'location', 'SouthEastOutside')

%%
% With the observation that the patient has bronchitis (|B = t|), the probability of
% a true condition for all other nodes has increased. In particular, the
% probability of smoking history increases because smoking is one leading
% cause of chronic bronchitis. In turn, because smoking is also associated
% with lung cancer, the probability of lung cancer increases and so does the
% probability of an abnormal chest x-ray test.

%% 
% Suppose the patient has not been evaluated for bronchitis but the chest
% x-ray shows some abnormalities. We instantiate |X = t| and we
% intialize again the network with the new evidence.

evNode = X; 
evValue = t;
[n2, e2, A2, a2] = bnMsgPassUpdate(nodes, edges, [], [], evNode, evValue);

for i = 1:n
    disp(['P(' nodeNames{i}, '|X=t) = ' num2str(n2(i).P(1))]);
end

%% 
% Given the observed abnormal x-ray results, the probability of lung cancer
% has increased significantly because of the direct dependency of node |X| (x-rays) on
% node |L| (lung cancer).

%%
% Finally, suppose the patient has both been diagnosed with bronchitis and
% received positive results for his/her chest x-ray. We update the previous
% state of the network (|X| = |t|) with the new evidence(|B| = |t|), as shown below: 

evNode = B; 
evValue = t;
[n3, e3, A3, a3] = bnMsgPassUpdate(n2, e2, A2, a2, evNode, evValue);

for i = 1:n
    disp(['P(' nodeNames{i}, '|B=t,X=t) = ' num2str(n3(i).P(1))]);
end

%%
% We can compare the three situations by plotting the probabilities as bar
% charts. Evidence of bronchitis and evidence of abnormal x-rays increase
% the probability of lung cancer and smoking history, one indirectly and
% the other directly.

figure(); subplot(3,1,1); 
bar(x1,  'stacked'); set(gca, 'xticklabel', nodeNames);
ylabel('Probability'); title('Bronchitis diagnosis');
legend({'true', 'false'}, 'location', 'SouthEastOutside')
hold on; subplot(3,1,2);
x2 = cat(1,n2.P);
bar(x2,'stacked'); set(gca, 'xticklabel', nodeNames);
ylabel('Probability'); title('Abnormal x-rays');
legend({'true', 'false'}, 'location', 'SouthEastOutside')
hold on; subplot(3,1,3);
x3 = cat(1,n3.P);
bar(x3, 'stacked'); set(gca, 'xticklabel', nodeNames);
ylabel('Probability'); title('Bronchitis and abnormal x-ray');
legend({'true', 'false'}, 'location', 'SouthEastOutside')
%%
% We can now compare the effect of a positive versus negative bronchitis
% diagnosis in presence of abnormal x-ray results. We instantiate |B = f|
% (false) and compare with previous estimates |B = t| (true).

evNode = B; 
evValue = f;
[n4, e4, A4, a4] = bnMsgPassUpdate(n2, e2, [], [], evNode, evValue);

for i = 1:n
    disp(['P(' nodeNames{i}, '|B=f,X=t) = ' num2str(n4(i).P(1))]);
end

figure();
bar3([n3(S).P(:,t) n4(S).P(:,t); n3(L).P(:,t) n4(L).P(:,t)]);
colormap(summer); zlabel('Probability'); 
set(gca,'xticklabel',{'Smoking','Lung Cancer'},'yticklabel', {'With Bronchitis', 'Without Bronchitis'});