% The waste incinerator emissions example from Lauritzen (1992),% "Propogation of Probabilities, Means and Variances in Mixed Graphical Association Models", % JASA 87(420): 1098--1108% node numbersF = 1; W = 2; E = 3; B = 4; C = 5; D = 6; Min = 7; Mout = 8; L = 9;% node sizes - all cts nodes are scalar, all discrete nodes are binaryns = ones(1, 9);dnodes = [F W B];ns(dnodes) = 2;% topology (p 1099, fig 1)dag = zeros(9);dag(F,E)=1;dag(W,[E Min D]) = 1;dag(E,D)=1;dag(B,[C D])=1;dag(D,[L Mout])=1;dag(Min,Mout)=1;% params (p 1102)bnet = mk_bnet(dag, ns, dnodes);bnet.CPD{B} = tabular_CPD(bnet, B, [0.85 0.15]'); % 1=stable, 2=unstablebnet.CPD{F} = tabular_CPD(bnet, F, [0.95 0.05]'); % 1=intact, 2=defectbnet.CPD{W} = tabular_CPD(bnet, W, [2/7 5/7]'); % 1=industrial, 2=householdcov_type = 'diag';tied_cov = 0;mu = reshape([-3.9 -0.4 -3.2 -0.5], [1 2 2]);Sigma = reshape([0.00002 0.0001 0.00002 0.0001], [1 1 2 2]);bnet.CPD{E} = gaussian_CPD(bnet, E, cov_type, tied_cov, mu, Sigma);mu = reshape([6.5 6.0 7.5 7.0], [1 2 2]);Sigma = reshape([0.03 0.04 0.1 0.1], [1 1 2 2]);weights = reshape([1 1 1 1], [1 1 2 2]);bnet.CPD{D} = gaussian_CPD(bnet, D, cov_type, tied_cov, mu, Sigma, weights);mu = reshape([-2 -1], [1 2]);Sigma = reshape([0.1 0.3], [1 1 2]);bnet.CPD{C} = gaussian_CPD(bnet, C, cov_type, tied_cov, mu, Sigma);bnet.CPD{L} = gaussian_CPD(bnet, L, cov_type, tied_cov, 3, 0.25, -0.5);mu = reshape([0.5 -0.5], [1 2]);Sigma = reshape([0.1 0.005], [1 1 2]);bnet.CPD{Min} = gaussian_CPD(bnet, Min, cov_type, tied_cov, mu, Sigma);bnet.CPD{Mout} = gaussian_CPD(bnet, Mout, cov_type, tied_cov, 0, 0.002, [1 1]);%engine = jtree_inf_engine(bnet);engine = cond_gauss_inf_engine(bnet);% Initial valuesevidence = cell(1,9);marginals = infer_all(engine, evidence);% Compare results with line 1 of table on p1107stol = 0.1; % the std only seems to agree to 1dp, presumably due to rounding errormtol = 0.01; % mean toleranceptol = 0.01; % discrete toleranceapproxeq(marginals{F}.T(1), 0.95, ptol)approxeq(marginals{W}.T(2), 0.71, ptol)approxeq(marginals{B}.T(1), 0.85, ptol)approxeq(marginals{Min}.mu, -0.21, mtol)approxeq(sqrt(marginals{Min}.Sigma), 0.46, stol)approxeq(marginals{Mout}.mu, 2.83, mtol)approxeq(sqrt(marginals{Mout}.Sigma), 0.86, stol)sqrt(marginals{Mout}.Sigma)% Add evidence (p 1105, top right)evidence{W} = 1; % industrialevidence{L} = 1.1;evidence{C} = -0.9;marginals = infer_all(engine, evidence);