function scg_unstable()% the objective of this script is to test if the stable conditonal gaussian% inference can handle the numerical instability problem described on% page.151 of 'Probabilistic networks and expert system' by Cowell, Dawid, Lauritzen and% Spiegelhalter, 1999.A = 1; Y = 2;n = 2;ns = ones(1, n);dnodes = [A];cnodes = Y;ns = [2 1];dag = zeros(n);dag(A, Y) = 1;bnet = mk_bnet(dag, ns, dnodes);bnet.CPD{A} = tabular_CPD(bnet, A, [0.5 0.5]'); bnet.CPD{Y} = gaussian_CPD(bnet, Y, 'mean', [0 1], 'cov', [1e-5 1e-6]);evidence = cell(1, n);pot_type = 'cg';potYgivenA = convert_to_pot(bnet.CPD{Y}, pot_type, [A Y], evidence);potA = convert_to_pot(bnet.CPD{A}, pot_type, A, evidence);potYandA = multiply_by_pot(potYgivenA, potA);potA2 = marginalize_pot(potYandA, A);thresh = 1; % 0dp[g,h,K] = extract_can(potA);assert(approxeq(g(:)', [-0.693147 -0.693147], thresh))[g,h,K] = extract_can(potYgivenA);assert(approxeq(g(:)', [4.83752 -499994], thresh))assert(approxeq(h(:)', [0 1e6]))assert(approxeq(K(:)', [1e5 1e6]))[g,h,K] = extract_can(potYandA);assert(approxeq(g(:)', [4.14437 -499995], thresh))assert(approxeq(h(:)', [0 1e6]))assert(approxeq(K(:)', [1e5 1e6]))[g,h,K] = extract_can(potA2);%assert(approxeq(g(:)', [-0.69315 -1]))gassert(approxeq(g(:)', [-0.69315 -0.69315]))if 0pot_type = 'scg';spotYgivenA = convert_to_pot(bnet.CPD{Y}, pot_type, [A Y], evidence);spotA = convert_to_pot(bnet.CPD{A}, pot_type, A, evidence);spotYandA = direct_combine_pots(spotYgivenA, spotA); spotA2 = marginalize_pot(spotYandA, A);spotA=struct(spotA);spotA2=struct(spotA2);for i=1:2  assert(approxeq(spotA2.scgpotc{i}.p, spotA.scgpotc{i}.p))  assert(approxeq(spotA2.scgpotc{i}.A, spotA.scgpotc{i}.A))  assert(approxeq(spotA2.scgpotc{i}.B, spotA.scgpotc{i}.B))  assert(approxeq(spotA2.scgpotc{i}.C, spotA.scgpotc{i}.C))endend%%%%%%%%%%%function [g,h,K] = extract_can(pot)pot = struct(pot);D = length(pot.can);g = zeros(1, D);h = zeros(1, D);K = zeros(1, D);for i=1:D  S = struct(pot.can{i});  g(i) = S.g;  if length(S.h) > 0    h(i) = S.h;    K(i) = S.K;  endend