Mercurial > hg > camir-aes2014

annotate toolboxes/FullBNT-1.0.7/bnt/CPDs/@gaussian_CPD/Old/gaussian_CPD.m @ 0:e9a9cd732c1e tip

Find changesets by keywords (author, files, the commit message), revision number or hash, or revset expression.
first hg version after svn
author wolffd
date Tue, 10 Feb 2015 15:05:51 +0000
parents
children
rev   line source
wolffd@0 1 function CPD = gaussian_CPD(varargin)
wolffd@0 2 % GAUSSIAN_CPD Make a conditional linear Gaussian distrib.
wolffd@0 3 %
wolffd@0 4 % To define this CPD precisely, call the continuous (cts) parents (if any) X,
wolffd@0 5 % the discrete parents (if any) Q, and this node Y. Then the distribution on Y is:
wolffd@0 6 % - no parents: Y ~ N(mu, Sigma)
wolffd@0 7 % - cts parents : Y|X=x ~ N(mu + W x, Sigma)
wolffd@0 8 % - discrete parents: Y|Q=i ~ N(mu(i), Sigma(i))
wolffd@0 9 % - cts and discrete parents: Y|X=x,Q=i ~ N(mu(i) + W(i) x, Sigma(i))
wolffd@0 10 %
wolffd@0 11 % CPD = gaussian_CPD(bnet, node, ...) will create a CPD with random parameters,
wolffd@0 12 % where node is the number of a node in this equivalence class.
wolffd@0 13 %
wolffd@0 14 % The list below gives optional arguments [default value in brackets].
wolffd@0 15 % (Let ns(i) be the size of node i, X = ns(X), Y = ns(Y) and Q = prod(ns(Q)).)
wolffd@0 16 %
wolffd@0 17 % mean - mu(:,i) is the mean given Q=i [ randn(Y,Q) ]
wolffd@0 18 % cov - Sigma(:,:,i) is the covariance given Q=i [ repmat(eye(Y,Y), [1 1 Q]) ]
wolffd@0 19 % weights - W(:,:,i) is the regression matrix given Q=i [ randn(Y,X,Q) ]
wolffd@0 20 % cov_type - if 'diag', Sigma(:,:,i) is diagonal [ 'full' ]
wolffd@0 21 % tied_cov - if 1, we constrain Sigma(:,:,i) to be the same for all i [0]
wolffd@0 22 % clamp_mean - if 1, we do not adjust mu(:,i) during learning [0]
wolffd@0 23 % clamp_cov - if 1, we do not adjust Sigma(:,:,i) during learning [0]
wolffd@0 24 % clamp_weights - if 1, we do not adjust W(:,:,i) during learning [0]
wolffd@0 25 % cov_prior_weight - weight given to I prior for estimating Sigma [0.01]
wolffd@0 26 %
wolffd@0 27 % e.g., CPD = gaussian_CPD(bnet, i, 'mean', [0; 0], 'clamp_mean', 'yes')
wolffd@0 28 %
wolffd@0 29 % For backwards compatibility with BNT2, you can also specify the parameters in the following order
wolffd@0 30 % CPD = gaussian_CPD(bnet, self, mu, Sigma, W, cov_type, tied_cov, clamp_mean, clamp_cov, clamp_weight)
wolffd@0 31 %
wolffd@0 32 % Sometimes it is useful to create an "isolated" CPD, without needing to pass in a bnet.
wolffd@0 33 % In this case, you must specify the discrete and cts parents (dps, cps) and the family sizes, followed
wolffd@0 34 % by the optional arguments above:
wolffd@0 35 % CPD = gaussian_CPD('self', i, 'dps', dps, 'cps', cps, 'sz', fam_size, ...)
wolffd@0 36
wolffd@0 37
wolffd@0 38 if nargin==0
wolffd@0 39 % This occurs if we are trying to load an object from a file.
wolffd@0 40 CPD = init_fields;
wolffd@0 41 clamp = 0;
wolffd@0 42 CPD = class(CPD, 'gaussian_CPD', generic_CPD(clamp));
wolffd@0 43 return;
wolffd@0 44 elseif isa(varargin{1}, 'gaussian_CPD')
wolffd@0 45 % This might occur if we are copying an object.
wolffd@0 46 CPD = varargin{1};
wolffd@0 47 return;
wolffd@0 48 end
wolffd@0 49 CPD = init_fields;
wolffd@0 50
wolffd@0 51 CPD = class(CPD, 'gaussian_CPD', generic_CPD(0));
wolffd@0 52
wolffd@0 53
wolffd@0 54 % parse mandatory arguments
wolffd@0 55 if ~isstr(varargin{1}) % pass in bnet
wolffd@0 56 bnet = varargin{1};
wolffd@0 57 self = varargin{2};
wolffd@0 58 args = varargin(3:end);
wolffd@0 59 ns = bnet.node_sizes;
wolffd@0 60 ps = parents(bnet.dag, self);
wolffd@0 61 dps = myintersect(ps, bnet.dnodes);
wolffd@0 62 cps = myintersect(ps, bnet.cnodes);
wolffd@0 63 fam_sz = ns([ps self]);
wolffd@0 64 else
wolffd@0 65 disp('parsing new style')
wolffd@0 66 for i=1:2:length(varargin)
wolffd@0 67 switch varargin{i},
wolffd@0 68 case 'self', self = varargin{i+1};
wolffd@0 69 case 'dps', dps = varargin{i+1};
wolffd@0 70 case 'cps', cps = varargin{i+1};
wolffd@0 71 case 'sz', fam_sz = varargin{i+1};
wolffd@0 72 end
wolffd@0 73 end
wolffd@0 74 ps = myunion(dps, cps);
wolffd@0 75 args = varargin;
wolffd@0 76 end
wolffd@0 77
wolffd@0 78 CPD.self = self;
wolffd@0 79 CPD.sizes = fam_sz;
wolffd@0 80
wolffd@0 81 % Figure out which (if any) of the parents are discrete, and which cts, and how big they are
wolffd@0 82 % dps = discrete parents, cps = cts parents
wolffd@0 83 CPD.cps = find_equiv_posns(cps, ps); % cts parent index
wolffd@0 84 CPD.dps = find_equiv_posns(dps, ps);
wolffd@0 85 ss = fam_sz(end);
wolffd@0 86 psz = fam_sz(1:end-1);
wolffd@0 87 dpsz = prod(psz(CPD.dps));
wolffd@0 88 cpsz = sum(psz(CPD.cps));
wolffd@0 89
wolffd@0 90 % set default params
wolffd@0 91 CPD.mean = randn(ss, dpsz);
wolffd@0 92 CPD.cov = 100*repmat(eye(ss), [1 1 dpsz]);
wolffd@0 93 CPD.weights = randn(ss, cpsz, dpsz);
wolffd@0 94 CPD.cov_type = 'full';
wolffd@0 95 CPD.tied_cov = 0;
wolffd@0 96 CPD.clamped_mean = 0;
wolffd@0 97 CPD.clamped_cov = 0;
wolffd@0 98 CPD.clamped_weights = 0;
wolffd@0 99 CPD.cov_prior_weight = 0.01;
wolffd@0 100
wolffd@0 101 nargs = length(args);
wolffd@0 102 if nargs > 0
wolffd@0 103 if ~isstr(args{1})
wolffd@0 104 % gaussian_CPD(bnet, self, mu, Sigma, W, cov_type, tied_cov, clamp_mean, clamp_cov, clamp_weights)
wolffd@0 105 if nargs >= 1 & ~isempty(args{1}), CPD.mean = args{1}; end
wolffd@0 106 if nargs >= 2 & ~isempty(args{2}), CPD.cov = args{2}; end
wolffd@0 107 if nargs >= 3 & ~isempty(args{3}), CPD.weights = args{3}; end
wolffd@0 108 if nargs >= 4 & ~isempty(args{4}), CPD.cov_type = args{4}; end
wolffd@0 109 if nargs >= 5 & ~isempty(args{5}) & strcmp(args{5}, 'tied'), CPD.tied_cov = 1; end
wolffd@0 110 if nargs >= 6 & ~isempty(args{6}), CPD.clamped_mean = 1; end
wolffd@0 111 if nargs >= 7 & ~isempty(args{7}), CPD.clamped_cov = 1; end
wolffd@0 112 if nargs >= 8 & ~isempty(args{8}), CPD.clamped_weights = 1; end
wolffd@0 113 else
wolffd@0 114 CPD = set_fields(CPD, args{:});
wolffd@0 115 end
wolffd@0 116 end
wolffd@0 117
wolffd@0 118 % Make sure the matrices have 1 dimension per discrete parent.
wolffd@0 119 % Bug fix due to Xuejing Sun 3/6/01
wolffd@0 120 CPD.mean = myreshape(CPD.mean, [ss ns(dps)]);
wolffd@0 121 CPD.cov = myreshape(CPD.cov, [ss ss ns(dps)]);
wolffd@0 122 CPD.weights = myreshape(CPD.weights, [ss cpsz ns(dps)]);
wolffd@0 123
wolffd@0 124 CPD.init_cov = CPD.cov; % we reset to this if things go wrong during learning
wolffd@0 125
wolffd@0 126 % expected sufficient statistics
wolffd@0 127 CPD.Wsum = zeros(dpsz,1);
wolffd@0 128 CPD.WYsum = zeros(ss, dpsz);
wolffd@0 129 CPD.WXsum = zeros(cpsz, dpsz);
wolffd@0 130 CPD.WYYsum = zeros(ss, ss, dpsz);
wolffd@0 131 CPD.WXXsum = zeros(cpsz, cpsz, dpsz);
wolffd@0 132 CPD.WXYsum = zeros(cpsz, ss, dpsz);
wolffd@0 133
wolffd@0 134 % For BIC
wolffd@0 135 CPD.nsamples = 0;
wolffd@0 136 switch CPD.cov_type
wolffd@0 137 case 'full',
wolffd@0 138 ncov_params = ss*(ss-1)/2; % since symmetric (and positive definite)
wolffd@0 139 case 'diag',
wolffd@0 140 ncov_params = ss;
wolffd@0 141 otherwise
wolffd@0 142 error(['unrecognized cov_type ' cov_type]);
wolffd@0 143 end
wolffd@0 144 % params = weights + mean + cov
wolffd@0 145 if CPD.tied_cov
wolffd@0 146 CPD.nparams = ss*cpsz*dpsz + ss*dpsz + ncov_params;
wolffd@0 147 else
wolffd@0 148 CPD.nparams = ss*cpsz*dpsz + ss*dpsz + dpsz*ncov_params;
wolffd@0 149 end
wolffd@0 150
wolffd@0 151
wolffd@0 152
wolffd@0 153 clamped = CPD.clamped_mean & CPD.clamped_cov & CPD.clamped_weights;
wolffd@0 154 CPD = set_clamped(CPD, clamped);
wolffd@0 155
wolffd@0 156 %%%%%%%%%%%
wolffd@0 157
wolffd@0 158 function CPD = init_fields()
wolffd@0 159 % This ensures we define the fields in the same order
wolffd@0 160 % no matter whether we load an object from a file,
wolffd@0 161 % or create it from scratch. (Matlab requires this.)
wolffd@0 162
wolffd@0 163 CPD.self = [];
wolffd@0 164 CPD.sizes = [];
wolffd@0 165 CPD.cps = [];
wolffd@0 166 CPD.dps = [];
wolffd@0 167 CPD.mean = [];
wolffd@0 168 CPD.cov = [];
wolffd@0 169 CPD.weights = [];
wolffd@0 170 CPD.clamped_mean = [];
wolffd@0 171 CPD.clamped_cov = [];
wolffd@0 172 CPD.clamped_weights = [];
wolffd@0 173 CPD.init_cov = [];
wolffd@0 174 CPD.cov_type = [];
wolffd@0 175 CPD.tied_cov = [];
wolffd@0 176 CPD.Wsum = [];
wolffd@0 177 CPD.WYsum = [];
wolffd@0 178 CPD.WXsum = [];
wolffd@0 179 CPD.WYYsum = [];
wolffd@0 180 CPD.WXXsum = [];
wolffd@0 181 CPD.WXYsum = [];
wolffd@0 182 CPD.nsamples = [];
wolffd@0 183 CPD.nparams = [];
wolffd@0 184 CPD.cov_prior_weight = [];