wolffd@0: function marginal = marginal_nodes(engine, query)
wolffd@0: % MARGINAL_NODES Compute the marginal on the specified query nodes (gaussian)
wolffd@0: % marginal = marginal_nodes(engine, query)
wolffd@0: 
wolffd@0: % Compute sum_{Hsum} Pr(Hkeep, Hsum | o)
wolffd@0: H = engine.hnodes;
wolffd@0: bnet = bnet_from_engine(engine);
wolffd@0: ns = bnet.node_sizes;
wolffd@0: Hkeep = myintersect(H, query);
wolffd@0: Hsum = mysetdiff(H, Hkeep);
wolffd@0: 
wolffd@0: [marginal.mu, marginal.Sigma] = marginalize_gaussian(engine.Hmu, engine.HSigma, Hkeep, Hsum, ns);
wolffd@0: marginal.domain = query;
wolffd@0: marginal.T = 1;
wolffd@0: