diff options
Diffstat (limited to 'sourcecodes/bnt-master/BNT/CPDs/@mlp_CPD')
8 files changed, 404 insertions, 0 deletions
diff --git a/sourcecodes/bnt-master/BNT/CPDs/@mlp_CPD/CVS/Entries b/sourcecodes/bnt-master/BNT/CPDs/@mlp_CPD/CVS/Entries new file mode 100644 index 00000000..1cef220a --- /dev/null +++ b/sourcecodes/bnt-master/BNT/CPDs/@mlp_CPD/CVS/Entries @@ -0,0 +1,6 @@ +/convert_to_table.m/1.1.1.1/Wed May 29 15:59:54 2002// +/maximize_params.m/1.1.1.1/Wed May 29 15:59:54 2002// +/mlp_CPD.m/1.1.1.1/Wed May 29 15:59:54 2002// +/reset_ess.m/1.1.1.1/Wed May 29 15:59:54 2002// +/update_ess.m/1.1.1.1/Wed May 29 15:59:54 2002// +D diff --git a/sourcecodes/bnt-master/BNT/CPDs/@mlp_CPD/CVS/Repository b/sourcecodes/bnt-master/BNT/CPDs/@mlp_CPD/CVS/Repository new file mode 100644 index 00000000..9fadd7fe --- /dev/null +++ b/sourcecodes/bnt-master/BNT/CPDs/@mlp_CPD/CVS/Repository @@ -0,0 +1 @@ +FullBNT/BNT/CPDs/@mlp_CPD diff --git a/sourcecodes/bnt-master/BNT/CPDs/@mlp_CPD/CVS/Root b/sourcecodes/bnt-master/BNT/CPDs/@mlp_CPD/CVS/Root new file mode 100644 index 00000000..f3bd14a6 --- /dev/null +++ b/sourcecodes/bnt-master/BNT/CPDs/@mlp_CPD/CVS/Root @@ -0,0 +1 @@ +:ext:nsaunier@bnt.cvs.sourceforge.net:/cvsroot/bnt diff --git a/sourcecodes/bnt-master/BNT/CPDs/@mlp_CPD/convert_to_table.m b/sourcecodes/bnt-master/BNT/CPDs/@mlp_CPD/convert_to_table.m new file mode 100644 index 00000000..7e25d072 --- /dev/null +++ b/sourcecodes/bnt-master/BNT/CPDs/@mlp_CPD/convert_to_table.m @@ -0,0 +1,80 @@ +function T = convert_to_table(CPD, domain, evidence) +% CONVERT_TO_TABLE Convert a mlp CPD to a table, incorporating any evidence +% T = convert_to_table(CPD, domain, evidence) + +self = domain(end); +ps = domain(1:end-1); % self' parents +%cps = myintersect(ps, cnodes); % self' continous parents +cnodes = domain(CPD.cpndx); +cps = myintersect(ps, cnodes); +odom = domain(~isemptycell(evidence(domain))); % obs nodes in the net +assert(myismember(cps, odom)); % !ALL the CTS parents must be observed! +ns(cps)=1; +dps = mysetdiff(ps, cps); % self' discrete parents +dobs = myintersect(dps, odom); % discrete obs parents + +% Extract the params compatible with the observations (if any) on the discrete parents (if any) + +if ~isempty(dobs), + dvals = cat(1, evidence{dobs}); + ns_eff= CPD.sizes; % effective node sizes + ens=ns_eff; + ens(dobs) = 1; + S=prod(ens(dps)); + subs = ind2subv(ens(dps), 1:S); + mask = find_equiv_posns(dobs, dps); + for i=1:length(mask), + subs(:,mask(i)) = dvals(i); + end + support = subv2ind(ns_eff(dps), subs)'; +else + ns_eff= CPD.sizes; + support=[1:prod(ns_eff(dps))]; +end + +W1=[]; b1=[]; W2=[]; b2=[]; + +W1 = CPD.W1(:,:,support); +b1= CPD.b1(support,:); +W2 = CPD.W2(:,:,support); +b2= CPD.b2(support,:); +ns(odom) = 1; +dpsize = prod(ns(dps)); % overall size of the self' discrete parents + +x = cat(1, evidence{cps}); +ndata=size(x,2); + +if ~isempty(evidence{self}) % + app=struct(CPD); % + ns(self)=app.mlp{1}.nout; % pump up self to the original dimension if observed + clear app; % +end % + +T =zeros(dpsize, ns(self)); % +for i=1:dpsize % + W1app = W1(:,:,i); % + b1app = b1(i,:); % + W2app = W2(:,:,i); % + b2app = b2(i,:); % for each of the dpsize combinations of self'parents values + z = tanh(x(:)'*W1app + ones(ndata, 1)*b1app); % we tabulate the corrisponding glm model + a = z*W2app + ones(ndata, 1)*b2app; % (element of the cell array CPD.glim) + appoggio = normalise(exp(a)); % + T(i,:)=appoggio; % + W1app=[]; W2app=[]; b1app=[]; b2app=[]; % + z=[]; a=[]; appoggio=[]; % +end % + +if ~isempty(evidence{self}) + appoggio=[]; % + appoggio=zeros(1,ns(self)); % + r = evidence{self}; %...if self is observed => in output there's only the probability of the 'true' class + for i=1:dpsize % + appoggio(i)=T(i,r); % + end + T=zeros(dpsize,1); + for i=1:dpsize + T(i,1)=appoggio(i); + end + clear appoggio; + ns(self) = 1; +end diff --git a/sourcecodes/bnt-master/BNT/CPDs/@mlp_CPD/maximize_params.m b/sourcecodes/bnt-master/BNT/CPDs/@mlp_CPD/maximize_params.m new file mode 100644 index 00000000..19d0a1be --- /dev/null +++ b/sourcecodes/bnt-master/BNT/CPDs/@mlp_CPD/maximize_params.m @@ -0,0 +1,34 @@ +function CPD = maximize_params(CPD, temp) +% MAXIMIZE_PARAMS Find ML params of an MLP using Scaled Conjugated Gradient (SCG) +% CPD = maximize_params(CPD, temperature) +% temperature parameter is ignored + +if ~adjustable_CPD(CPD), return; end +options = foptions; + +% options(1) >= 0 means print an annoying message when the max. num. iter. is reached +if CPD.verbose + options(1) = 1; +else + options(1) = -1; +end +%options(1) = CPD.verbose; + +options(2) = CPD.wthresh; +options(3) = CPD.llthresh; +options(14) = CPD.max_iter; + +dpsz=length(CPD.mlp); + +for i=1:dpsz + mask=[]; + mask=find(CPD.eso_weights(:,:,i)>0); % for adapting the parameters we use only positive weighted example + if ~isempty(mask), + CPD.mlp{i} = netopt_weighted(CPD.mlp{i}, options, CPD.parent_vals(mask',:), CPD.self_vals(mask',:,i), CPD.eso_weights(mask',:,i), 'scg'); + + CPD.W1(:,:,i)=CPD.mlp{i}.w1; % update the parameters matrix + CPD.b1(i,:)=CPD.mlp{i}.b1; % + CPD.W2(:,:,i)=CPD.mlp{i}.w2; % update the parameters matrix + CPD.b2(i,:)=CPD.mlp{i}.b2; % + end +end diff --git a/sourcecodes/bnt-master/BNT/CPDs/@mlp_CPD/mlp_CPD.m b/sourcecodes/bnt-master/BNT/CPDs/@mlp_CPD/mlp_CPD.m new file mode 100644 index 00000000..7e9d3f6b --- /dev/null +++ b/sourcecodes/bnt-master/BNT/CPDs/@mlp_CPD/mlp_CPD.m @@ -0,0 +1,139 @@ +function CPD = mlp_CPD(bnet, self, nhidden, w1, b1, w2, b2, clamped, max_iter, verbose, wthresh, llthresh) +% MLP_CPD Make a CPD from a Multi Layer Perceptron (i.e., feedforward neural network) +% +% We use a different MLP for each discrete parent combination (if there are any discrete parents). +% We currently assume this node (the child) is discrete. +% +% CPD = mlp_CPD(bnet, self, nhidden) +% will create a CPD with random parameters, where self is the number of this node and nhidden the number of the hidden nodes. +% The params are drawn from N(0, s*I), where s = 1/sqrt(n+1), n = length(X). +% +% CPD = mlp_CPD(bnet, self, nhidden, w1, b1, w2, b2) allows you to specify the params, where +% w1 = first-layer weight matrix +% b1 = first-layer bias vector +% w2 = second-layer weight matrix +% b2 = second-layer bias vector +% These are assumed to be the same for each discrete parent combination. +% If any of these are [], random values will be created. +% +% CPD = mlp_CPD(bnet, self, nhidden, w1, b1, w2, b2, clamped) allows you to prevent the params from being +% updated during learning (if clamped = 1). Default: clamped = 0. +% +% CPD = mlp_CPD(bnet, self, nhidden, w1, b1, w2, b2, clamped, max_iter, verbose, wthresh, llthresh) +% alllows you to specify params that control the M step: +% max_iter - the maximum number of steps to take (default: 10) +% verbose - controls whether to print (default: 0 means silent). +% wthresh - a measure of the precision required for the value of +% the weights W at the solution. Default: 1e-2. +% llthresh - a measure of the precision required of the objective +% function (log-likelihood) at the solution. Both this and the previous condition must +% be satisfied for termination. Default: 1e-2. +% +% For learning, we use a weighted version of scaled conjugated gradient in the M step. + +if nargin==0 + % This occurs if we are trying to load an object from a file. + CPD = init_fields; + CPD = class(CPD, 'mlp_CPD', discrete_CPD(0,[])); + return; +elseif isa(bnet, 'mlp_CPD') + % This might occur if we are copying an object. + CPD = bnet; + return; +end +CPD = init_fields; + +assert(myismember(self, bnet.dnodes)); +ns = bnet.node_sizes; + +ps = parents(bnet.dag, self); +dnodes = mysetdiff(1:length(bnet.dag), bnet.cnodes); +dps = myintersect(ps, dnodes); +cps = myintersect(ps, bnet.cnodes); +dpsz = prod(ns(dps)); +cpsz = sum(ns(cps)); +self_size = ns(self); + +% discrete/cts parent index - which ones of my parents are discrete/cts? +CPD.dpndx = find_equiv_posns(dps, ps); +CPD.cpndx = find_equiv_posns(cps, ps); + +CPD.mlp = cell(1,dpsz); +for i=1:dpsz + CPD.mlp{i} = mlp(cpsz, nhidden, self_size, 'softmax'); + if nargin >=4 & ~isempty(w1) + CPD.mlp{i}.w1 = w1; + end + if nargin >=5 & ~isempty(b1) + CPD.mlp{i}.b1 = b1; + end + if nargin >=6 & ~isempty(w2) + CPD.mlp{i}.w2 = w2; + end + if nargin >=7 & ~isempty(b2) + CPD.mlp{i}.b2 = b2; + end + W1app(:,:,i)=CPD.mlp{i}.w1; + W2app(:,:,i)=CPD.mlp{i}.w2; + b1app(i,:)=CPD.mlp{i}.b1; + b2app(i,:)=CPD.mlp{i}.b2; +end +if nargin < 8, clamped = 0; end +if nargin < 9, max_iter = 10; end +if nargin < 10, verbose = 0; end +if nargin < 11, wthresh = 1e-2; end +if nargin < 12, llthresh = 1e-2; end + +CPD.self = self; +CPD.max_iter = max_iter; +CPD.verbose = verbose; +CPD.wthresh = wthresh; +CPD.llthresh = llthresh; + +% sufficient statistics +% Since MLP is not in the exponential family, we must store all the raw data. +% +CPD.W1=W1app; % Extract all the parameters of the node for handling discrete obs parents +CPD.W2=W2app; % +nparaW=[size(W1app) size(W2app)]; % +CPD.b1=b1app; % +CPD.b2=b2app; % +nparab=[size(b1app) size(b2app)]; % + +CPD.sizes=bnet.node_sizes(:); % used in CPD_to_table to pump up the node sizes + +CPD.parent_vals = []; % X(l,:) = value of cts parents in l'th example + +CPD.eso_weights=[]; % weights used by the SCG algorithm + +CPD.self_vals = []; % Y(l,:) = value of self in l'th example + +% For BIC +CPD.nsamples = 0; +CPD.nparams=prod(nparaW)+prod(nparab); +CPD = class(CPD, 'mlp_CPD', discrete_CPD(clamped, ns([ps self]))); + +%%%%%%%%%%% + +function CPD = init_fields() +% This ensures we define the fields in the same order +% no matter whether we load an object from a file, +% or create it from scratch. (Matlab requires this.) + +CPD.mlp = {}; +CPD.self = []; +CPD.max_iter = []; +CPD.verbose = []; +CPD.wthresh = []; +CPD.llthresh = []; +CPD.approx_hess = []; +CPD.W1 = []; +CPD.W2 = []; +CPD.b1 = []; +CPD.b2 = []; +CPD.sizes = []; +CPD.parent_vals = []; +CPD.eso_weights=[]; +CPD.self_vals = []; +CPD.nsamples = []; +CPD.nparams = []; diff --git a/sourcecodes/bnt-master/BNT/CPDs/@mlp_CPD/reset_ess.m b/sourcecodes/bnt-master/BNT/CPDs/@mlp_CPD/reset_ess.m new file mode 100644 index 00000000..ba7a7101 --- /dev/null +++ b/sourcecodes/bnt-master/BNT/CPDs/@mlp_CPD/reset_ess.m @@ -0,0 +1,12 @@ +function CPD = reset_ess(CPD) +% RESET_ESS Reset the Expected Sufficient Statistics for a CPD (mlp) +% CPD = reset_ess(CPD) + +CPD.W1 = []; +CPD.W2 = []; +CPD.b1 = []; +CPD.b2 = []; +CPD.parent_vals = []; +CPD.eso_weights=[]; +CPD.self_vals = []; +CPD.nsamples = 0; \ No newline at end of file diff --git a/sourcecodes/bnt-master/BNT/CPDs/@mlp_CPD/update_ess.m b/sourcecodes/bnt-master/BNT/CPDs/@mlp_CPD/update_ess.m new file mode 100644 index 00000000..353a0b5c --- /dev/null +++ b/sourcecodes/bnt-master/BNT/CPDs/@mlp_CPD/update_ess.m @@ -0,0 +1,131 @@ +function CPD = update_ess(CPD, fmarginal, evidence, ns, cnodes, hidden_bitv) +% UPDATE_ESS Update the Expected Sufficient Statistics of a CPD (MLP) +% CPD = update_ess(CPD, family_marginal, evidence, node_sizes, cnodes, hidden_bitv) +% +% fmarginal = overall posterior distribution of self and its parents +% fmarginal(i1,i2...,ik,s)=prob(Pa1=i1,...,Pak=ik, self=s| X) +% +% => 1) prob(self|Pa1,...,Pak)=fmarginal/prob(Pa1,...,Pak) with prob(Pa1,...,Pak)=sum{s,fmarginal} +% [self estimation -> CPD.self_vals] +% 2) prob(Pa1,...,Pak) [SCG weights -> CPD.eso_weights] +% +% Hidden_bitv is ignored + +% Written by Pierpaolo Brutti + +if ~adjustable_CPD(CPD), return; end + +dom = fmarginal.domain; +cdom = myintersect(dom, cnodes); +assert(~any(isemptycell(evidence(cdom)))); +ns(cdom)=1; + +self = dom(end); +ps=dom(1:end-1); +dpdom=mysetdiff(ps,cdom); + +dnodes = mysetdiff(1:length(ns), cnodes); + +ddom = myintersect(ps, dnodes); % +if isempty(evidence{self}), % if self is hidden in what follow we must + ddom = myintersect(dom, dnodes); % consider its dimension +end % + +odom = dom(~isemptycell(evidence(dom))); +hdom = dom(isemptycell(evidence(dom))); % hidden parents in domain + +dobs = myintersect(ddom, odom); +dvals = cat(1, evidence{dobs}); +ens = ns; % effective node sizes +ens(dobs) = 1; + +dpsz=prod(ns(dpdom)); +S=prod(ens(ddom)); +subs = ind2subv(ens(ddom), 1:S); +mask = find_equiv_posns(dobs, ddom); +for i=1:length(mask), + subs(:,mask(i)) = dvals(i); +end +supportedQs = subv2ind(ns(ddom), subs); + +Qarity = prod(ns(ddom)); +if isempty(ddom), + Qarity = 1; +end +fullm.T = zeros(Qarity, 1); +fullm.T(supportedQs) = fmarginal.T(:); + +% For dynamic (recurrent) net------------------------------------------------------------- +% ---------------------------------------------------------------------------------------- +high=size(evidence,1); % slice height +ss_ns=ns(1:high); % single slice nodes sizes +pos=self; % +slice_num=0; % +while pos>high, % + slice_num=slice_num+1; % find active slice + pos=pos-high; % pos=self posistion into a single slice +end % + +last_dim=pos-1; % +if isempty(evidence{self}), % + last_dim=pos; % +end % last_dim=last reshaping dimension +reg=dom-slice_num*high; +dex=myintersect(reg(find(reg>=0)), [1:last_dim]); % +rs_dim=ss_ns(dex); % reshaping dimensions + +if slice_num>0, + act_slice=[]; past_ancest=[]; % + act_slice=slice_num*high+[1:high]; % recover the active slice nodes + % past_ancest=mysetdiff(ddom, act_slice); + past_ancest=mysetdiff(ps, act_slice); % recover ancestors contained into past slices + app=ns(past_ancest); + rs_dim=[app(:)' rs_dim(:)']; % +end % +if length(rs_dim)==1, rs_dim=[1 rs_dim]; end % +if size(rs_dim,1)~=1, rs_dim=rs_dim'; end % + +fullm.T=reshape(fullm.T, rs_dim); % reshaping the marginal + +% ---------------------------------------------------------------------------------------- +% ---------------------------------------------------------------------------------------- + +% X = cts parent, R = discrete self + +% 1) observations vector -> CPD.parents_vals ------------------------------------------------- +x = cat(1, evidence{cdom}); + +% 2) weights vector -> CPD.eso_weights ------------------------------------------------------- +if isempty(evidence{self}) % R is hidden + sum_over=length(rs_dim); + app=sum(fullm.T, sum_over); + pesi=reshape(app,[dpsz,1]); + clear app; +else + pesi=reshape(fullm.T,[dpsz,1]); +end + +assert(approxeq(sum(pesi),1)); + +% 3) estimate (if R is hidden) or recover (if R is obs) self'value---------------------------- +if isempty(evidence{self}) % R is hidden + app=mk_stochastic(fullm.T); % P(self|Pa1,...,Pak)=fmarginal/prob(Pa1,...,Pak) + app=reshape(app,[dpsz ns(self)]); % matrix size: prod{j,ns(Paj)} x ns(self) + r=app; + clear app; +else + r = zeros(dpsz,ns(self)); + for i=1:dpsz + if pesi(i)~=0, r(i,evidence{self}) = 1; end + end +end +for i=1:dpsz + if pesi(i) ~=0, assert(approxeq(sum(r(i,:)),1)); end +end + +CPD.nsamples = CPD.nsamples + 1; +CPD.parent_vals(CPD.nsamples,:) = x(:)'; +for i=1:dpsz + CPD.eso_weights(CPD.nsamples,:,i)=pesi(i); + CPD.self_vals(CPD.nsamples,:,i) = r(i,:); +end |