%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% % % Q1.1 % %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% b = 100; n=20; niter = 30; MLE = zeros(niter,1); MAP = zeros(niter,1); l = 0.2 for a=1:niter X=exprnd(1/l,n,50); MLE(a) = mean((1./mean(X)-l).^2); MAP(a) = mean(((n+a-1)./(sum(X)+b)-l).^2); end plot(1:niter,MLE,1:niter,MAP) legend(['MLE', 'MAP']) b = 100; a=30 niter = 1000; MLE = zeros(niter,1); MAP = zeros(niter,1); l = 0.2 for n=1:niter X=exprnd(1/l,n,50); MLE(n) = mean((1./mean(X)'-l).^2); MAP(n) = mean(((n+a-1)./(sum(X)'+b)-l).^2); end semilogy(1:niter,MLE,1:niter,MAP) legend(['MLE', 'MAP']) %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% % % Q1.2 % %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% load('estimate.mat') niters = 100; alpha = 1; alphas = zeros(niters,1); alpha0 = 1; alphas0 = zeros(niters,1); obj = zeros(niters,1); n = size(xs,1); for i=1:niters g0 = n*log(alpha0/mean(xs)) + n- n*digamma(alpha0) + sum(log(xs)) - sum(xs)/mean(xs); gp = n/alpha0 - n*trigamma(alpha0); alpha0 -= g0/gp; %Newton's method alphas0(i) = alpha0; g = n*log(alpha/mean(xs)) + n- n*digamma(alpha) + sum(log(xs)) - sum(xs)/mean(xs); alpha += 0.01*g; %GD alphas(i) = alpha; end plot(1:niters,alphas,1:niters,alphas0) legend("GD","Newton") %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% % % Q2.2 % %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% function [ model ] = NaiveBayes(XTrain, yTrain) yTrain = indicator(yTrain); yCounts = sum(yTrain); model.Py = yCounts/sum(yCounts); % For each p(x|y) compute mu, sigma^2 estimates, indexed by [y][x] model.mean = zeros(size(yTrain,2),size(XTrain,2)); model.sigma = zeros(size(yTrain,2),size(XTrain,2)); for i=1:size(yTrain,2) Xs = XTrain(logical(yTrain(:,i)),:); model.mean(i,:) = mean(Xs); model.sigma(i,:) = std(Xs); end end function [predictions] = NaiveBayesClassify(model, XTest) log_p = zeros(size(XTest,1),size(model.mean,1)); for i=1:size(XTest,1) for k=1:size(model.mean,1) log_p(i,k) = log(model.Py(k)) + sum(log(normpdf(XTest(i,:),model.mean(k,:),model.sigma(k,:)))); end end [max_values predictions] = max(log_p,[],2); end function [ model ] = AvgNaiveBayes(XTrain, yTrain) yTrain = indicator(yTrain); yCounts = sum(yTrain); model.Py = yCounts/sum(yCounts); model.mean = zeros(size(yTrain,2),size(XTrain,2)); model.sigma = zeros(size(yTrain,2),size(XTrain,2)); for i=1:size(yTrain,2) Xs = XTrain(logical(yTrain(:,i)),:); model.mean(i,:) = mean(Xs); model.sigma(i,:) = sqrt(var(Xs)); end model.x_mean = mean(XTrain); model.x_sigma = std(XTrain); log_px = zeros(size(XTrain,2),1); [yMax yIndices] = max(yTrain,[],2); log_cond_px = zeros(size(XTrain,2),1); for j=1:size(XTrain,2) log_px(j) = sum(log(normpdf(XTrain(:,j),model.x_mean(j),model.x_sigma(j)))); log_cond_px(j) = sum(log(normpdf(XTrain(:,j),model.mean(yIndices,j),model.sigma(yIndices,j)))); end model.log_px = log_px; model.log_cond_px = log_cond_px; end function [predictions] = AvgNaiveBayesClassify_T(model, XTest, Beta) log_p = zeros(size(XTest,1),size(model.mean,1)); for i=1:size(XTest,1) for k=1:size(model.mean,1) f = model.log_px + log(normpdf(XTest(i,:),model.x_mean,model.x_sigma))'; g = model.log_cond_px + log(normpdf(XTest(i,:),model.mean(k,:),model.sigma(k,:)))'; m = max(f,g); log_p(i,k) = log(model.Py(k)) + sum(m + log(exp(f-m) + (1/Beta)*exp(g-m))); end end [max_values predictions] = max(log_p,[],2); end