%% ================== Matlab code fragment 2.1 ==================
p = 0.8;
s = 0;
for k = 5:2:13  % loop from 5 to 13 by steps of 2
   s = s + p^k;
end
fprintf('The answer is %g\n',(1-p)*s)
%% ================== Matlab code fragment 2.2 ==================
p = 0.8;
pvec = (1-p)*p.^[5:2:13];
fprintf('The answer is %g\n',sum(pvec))
%% ================== Matlab code fragment 2.3 ==================
lambda = 10;
k = [0:15];  % k = [ 0 1 2 ... 15 ]
pvec = exp(-lambda)*lambda.^k./gamma(k+1);
fprintf('The answer is %g\n',1-sum(pvec))
%% ================== Matlab code fragment 2.4 ==================
p = ones(1,101)/101;	% p(i) = P(X=i) = 1/101, i = 0,...,100
k=[0:100];
i = find(cos(2*pi*k/10) <= 1/2);
fprintf('The answer is %g\n',sum(p(i)))
%% ================== Matlab code fragment 2.5 ==================
format rat % print numbers as ratios of small integers
pY = sum(P)
%% ================== Matlab code fragment 2.6 ==================
for i = 1:3
   for j = 1:5
      gmat(i,j) = g(i,j);
   end
end
prob = sum(P(find(gmat<6)))
%% ================== Matlab code fragment 2.7 ==================
s = 0;
for i = 1:3
   for j = 1:5
      s = s + i*j*P(i,j);
   end
end
[n,d] = rat(s); % to express answer fraction
fprintf('E[XY] = %i/%i = %g\n',n,d,s)
%% ================== Matlab code fragment 2.8 ==================
format rat
pY = sum(P)
y = [ 1 2 3 4 5 ]
mY = y*pY'
varY = ((y-mY).^2)*pY'
pX = sum(P')
x = [ 1 2 3 ]
mX = x*pX'
varX = ((x-mX).^2)*pX'
rhoxy = (s-mX*mY)/sqrt(varX*varY)
%% ================== Matlab code fragment 2.9 ==================
% function [x1,y1] = allpairs(x,y)
% lx = length(x);
% ly = length(y);
% x1 = kron(ones(ly,1),x);
% y1 = kron(y,ones(lx,1));
%% ================== Matlab code fragment 2.10 ==================
x = [ 1 2 3 ]'
y = [ 10 20 30 ]'
[x1,y1] = allpairs(x,y);
pairs = [x1 y1]
allsums = x1+y1
%% ================== Matlab code fragment 2.11 ==================
i = [1:3]';
j = [1:5]';
[x,y]=allpairs(i,j);
prob = sum(P(find(x.*y<6)))
