classdef bayesnet < graphicalmodel
properties
ec
end
methods
function obj = bayesnet(factors,varargin)
obj = obj@graphicalmodel(factors,varargin{:});
obj.g = graph(sparse(obj.length(),obj.length()));
for i=1:obj.length()
obj.g([factors{i}.parents],factors{i}.child) = true;
end
if isempty(obj.ec), obj.ec = 1:obj.length(); end
for i=unique(obj.ec)
eclass = find(obj.ec == i);
ess = obj.factors{eclass(1)}.essclone();
for j=eclass
obj.factors{j}.ess = ess;
end
end
end
function write(obj,filename,type,varargin)
switch type
case 'hugin'
bn2hugin(obj,filename,varargin{:});
case 'dot'
bn2dot(obj, filename,varargin{:});
otherwise
error('unrecognized format');
end
end
function data = simulate(obj,nsamples)
sorted = obj.g.topological_sort();
M = numel(sorted);
data = nan(nsamples,M);
domains = cellfun(@(x)(x.domain),obj.factors,'UniformOutput',false);
for i=1:nsamples
bevidence = false(M,1);
for j=sorted
dim = bevidence(domains{j});
val = data(i,domains{j}); val = val(dim);
data(i,j) = obj.factors{j}.sample(val);
bevidence(j) = true;
end
end
end
function l = loglik(obj,data)
l = 0;
for d=1:size(data,1)
for i=1:obj.length()
f = obj.factors{i};
l = l + f.loglik(data(d,i),data(d,f.parents()));
end
end
l = l ./ size(data,1);
end
function d = dim(obj)
d = 0;
for i=1:length(obj)
d = d + dim(obj.factors{i});
end
end
function s = bic(obj,data)
s = obj.loglik(data) - dim(obj)*log(length(obj))/2;
end
function s = aic(obj,data)
s = obj.loglik(data) - dim(obj);
end
end
methods (Access = private)
function bn2hugin(obj,filename,varargin)
fid = fopen(strcat(filename,'.net'), 'wt');
fprintf(fid, 'net{}\n\n');
for i=1:obj.length()
if isa(obj.factors{i},'gaussian_cpd')
fprintf(fid, 'continuous node C%d {\n\tposition = (%d %d);\n',i,round(rand*1000),round(rand*500));
if obj.factors{i}.name
fprintf(fid,'\tlabel = "%s";\n',obj.factors{i}.name);
end
fprintf(fid,'}\n\n');
elseif isa(obj.factors{i},'discrete_cpd')
fprintf(fid, 'node C%d \n{\n\tposition = (%d %d);\n\tstates = ( ',i,round(rand*1000),round(rand*500));
if ~isempty(obj.factors{i}.statenames)
for j=1:obj.size(i), fprintf(fid, '"%s" ',obj.factors{i}.statenames{j}); end
else
for j=1:obj.size(i), fprintf(fid, '"%d" ',j); end
end
fprintf(fid, ');\n');
if obj.factors{i}.name
fprintf(fid,'\tlabel = "%s";\n',obj.factors{i}.name);
end
fprintf(fid,'}\n\n');
else
class(obj.factors{i})
error('unsupported distribution');
end
end
for i=1:obj.length()
f = obj.factors{i};
fprintf(fid, 'potential (C%d',i);
if length(f.domain) > 1, fprintf(fid,' |'); end
dparents = sort(f.dparents,'descend');
for j=1:length(f.dparents)
fprintf(fid, ' C%d', dparents(j));
end
cparents = sort(f.cparents,'descend');
for j=1:length(f.cparents)
fprintf(fid, ' C%d', cparents(j));
end
fprintf(fid,')\n{\n\tdata =\n ');
szc = cumprod(f.dsize);
if isa(f,'continuous_cpd')
mu = f.mu;
beta = f.beta;
sigma2 = f.sigma2;
for j=1:numel(mu)
szc = szc(szc > 1);
for k=1:length(szc)
if mod(j-1,szc(k)) == 0, fprintf(fid,'( '); end
end
fprintf(fid,'normal ( %f',mu(j));
for c = 1:length(beta{j})
if beta{j}(c) >= 0, fprintf(fid,' + '); else fprintf(fid,' -'); end
fprintf(fid,'%f * C%d',abs(beta{j}(c)),f.cparents(c));
end
fprintf(fid,', %f ) ',sigma2(j));
for k=1:length(szc)
if mod(j,szc(k)) == 0, fprintf(fid,') '); end
end
end
else
data = f.p;
for j=1:numel(data)
for k=1:length(szc)
if mod(j-1,szc(k)) == 0, fprintf(fid,'('); end
end
fprintf(fid,'%f ',data(j));
for k=1:length(szc)
if mod(j,szc(k)) == 0, fprintf(fid,')'); end
end
end
end
fprintf(fid,';\n}\n\n');
end
fclose(fid);
end
function bn2dot(obj,filename,varargin)
ext = 'ps';
v = varargin2struct(varargin);
if isfield(v,'extension'), ext = v.extension; end
fid = fopen(strcat(filename,'.dot'), 'wt');
fprintf(fid, 'digraph G{\n');
for i=1:length(obj.factors)
if ~isempty(obj.factors{i}.name)
fprintf(fid,'%d [label="%s"]\n',i,obj.factors{i}.name);
else
fprintf(fid,'%d\n',i);
end
end
for i=1:length(obj.factors)
if isa(obj.factors{i},'gaussian_cpd')
b = obj.factors{i}.betas;
for j=1:length(obj.factors{i}.cparents)
signs = zeros(1,numel(b));
for k=1:numel(b)
signs(k) = sign(b{k}(j));
end
if all(signs == 1)
fprintf(fid,'%d -> %d [color=red];\n',j,i);
elseif all(signs == -1)
fprintf(fid,'%d -> %d [color=blue];\n',j,i);
else
fprintf(fid,'%d -> %d;\n',j,i);
end
end
else
for j=obj.factors{i}.cparents
fprintf(fid,'%d -> %d;\n',j,i);
end
end
for j=obj.factors{i}.dparents
fprintf(fid,'%d -> %d;\n',j,i);
end
end
fprintf(fid, '}\n');
fclose(fid);
if ispc, shell = 'dos'; else shell = 'unix'; end
cmdline = strcat([shell '(''dot -T' ext ' ' strcat(filename,'.dot') ' -o ' strcat(filename,'.',ext) ''')']);
r = eval(cmdline);
end
end
methods (Static)
function obj = random(N,varargin)
sz = 2*ones(1,N);
continuous = [];
g = [];
for i=1:2:length(varargin)
switch varargin{i},
case 'graph', g = varargin{i+1};
case 'size', sz = varargin{i+1};
case 'continuous', continuous = varargin{i+1};
end
end
sz(continuous) = 1;
if isempty(continuous), continuous = find(sz == 1); end
if isempty(g)
g = sparse(N,N);
for i=2:N
parents = randperm(i-1);
parents = parents(1:min(i-1,1+floor(2*rand)));
if ~ismember(i,continuous)
parents = setdiff(parents,continuous);
end
g(parents,i) = 1;
end
end
factors = cell(1,N);
for i=1:N
parents = find(g(:,i))';
if isempty(parents)
parents = [];
end
if ismember(i,continuous)
cparents = intersect(parents,continuous);
if isempty(cparents), cparents = []; end
dparents = setdiff(parents,continuous);
if isempty(dparents)
dparents = [];
szd = [1 1];
else
szd = sz(dparents);
if length(szd) == 1, szd = [szd 1]; end
end
mu = randn(szd);
beta = cell(szd);
if ~isempty(cparents)
for j=1:prod(szd)
beta{j} = randn(numel(cparents),1);
end
end
sigma1 = rand(szd);
factors{i} = gaussian_cpd(i,cparents,dparents,mu,beta,sigma1);
else
if isempty(parents)
szd = [sz(i) 1];
else
szd = sz([i parents]);
end
factors{i} = multinomial_cpd(i,parents,rand(szd));
end
end
obj = bayesnet(factors);
end
end
end
Input argument "factors" is undefined.
Error in ==> bayesnet>bayesnet.bayesnet at 26
obj = obj@graphicalmodel(factors,varargin{:});