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

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