GMR.m 1.21 KB
Newer Older
Sylvain Calinon's avatar
Sylvain Calinon committed
1 2 3 4 5 6 7 8 9 10 11
function [expData, expSigma, H] = GMR(model, DataIn, in, out)
%Gaussian mixture regression (GMR)
%Sylvain Calinon, Danilo Bruno, 2014

nbData = size(DataIn,2);
nbVarOut = length(out);

MuTmp = zeros(nbVarOut,model.nbStates);
expData = zeros(nbVarOut,nbData);
expSigma = zeros(nbVarOut,nbVarOut,nbData);
for t=1:nbData
Milad Malekzadeh's avatar
Milad Malekzadeh committed
12 13 14
	%Compute activation weight
	%See Eq. (3.0.5) in doc/TechnicalReport.pdf
	for i=1:model.nbStates
Milad Malekzadeh's avatar
Milad Malekzadeh committed
15
		H(i,t) = model.Priors(i) * gaussPDF(DataIn(:,t), model.Mu(in,i), model.Sigma(in,in,i));
Milad Malekzadeh's avatar
Milad Malekzadeh committed
16 17 18 19
	end
	H(:,t) = H(:,t)/sum(H(:,t));
	%Compute expected conditional means
	%See Eq. (3.0.3) in doc/TechnicalReport.pdf
Sylvain Calinon's avatar
Sylvain Calinon committed
20 21 22 23 24
	for i=1:model.nbStates
		MuTmp(:,i) = model.Mu(out,i) + model.Sigma(out,in,i)/model.Sigma(in,in,i) * (DataIn(:,t)-model.Mu(in,i));
		expData(:,t) = expData(:,t) + H(i,t) * MuTmp(:,i);
	end
	%Compute expected conditional covariances
Milad Malekzadeh's avatar
Milad Malekzadeh committed
25
	%See Eq. (3.0.4) in doc/TechnicalReport.pdf
Sylvain Calinon's avatar
Sylvain Calinon committed
26 27
	for i=1:model.nbStates
		SigmaTmp = model.Sigma(out,out,i) - model.Sigma(out,in,i)/model.Sigma(in,in,i) * model.Sigma(in,out,i);
Milad Malekzadeh's avatar
Milad Malekzadeh committed
28
		expSigma(:,:,t) = expSigma(:,:,t) + H(i,t) * (SigmaTmp + MuTmp(:,i)*MuTmp(:,i)');
Milad Malekzadeh's avatar
Milad Malekzadeh committed
29
		for j=1:model.nbStates
Milad Malekzadeh's avatar
Milad Malekzadeh committed
30
			expSigma(:,:,t) = expSigma(:,:,t) - H(i,t)*H(j,t) * (MuTmp(:,i)*MuTmp(:,j)');
Milad Malekzadeh's avatar
Milad Malekzadeh committed
31
		end
Sylvain Calinon's avatar
Sylvain Calinon committed
32 33 34
	end
end