% Solutions to Exercises in Chapter 14
% An Engineer’s Guide to MATLAB, 2nd ed.
% ©Prentice Hall 2005 Saddle River NJ
% Exercise 14.1
%(a)
disp(['probability of 8 calls = ' num2str(poisspdf(8, 5))])
disp(['probability of 2 calls = ' num2str(poisspdf(2, 5))])
%(b)
disp(['probability of a busy signal = ' num2str(1-poisscdf(10, 5))])
Answers
probability of 8 calls = 0.065278
probability of 2 calls = 0.084224
probability of a busy signal = 0.013695
% Exercise 14.2
disp(['Probability that at least 12 withstand load = ' num2str(1-binocdf(11, 15, 0.7))])
Answers
Probability that at least 12 withstand load = 0.29687
% Exercise 14.3
dat=[88.4 93.2 87.4 94.3 93.0 94.3 89.0 90.5 90.8 93.1 92.8 91.9;...
92.6 93.2 89.2 94.8 93.3 94.0 93.2 91.7 91.5 92.0 90.7 93.8];
for k=1:2
el(k)=var(dat(k,:))+(mean(dat(k,:))-92)^2;
end
disp(['L1 = ' num2str(el(1)) ' L2 = ' num2str(el(2))])
Answers
L1 = 5.5904 L2 = 2.6936
% Exercise 14.4
1
,p=0.2; n=24; m=0:n;
disp(['Expected value = ' num2str(n*p) ' variance = ' num2str(n*p*(1-p))])
disp(['Probability that no more than 2 are underweight = ' num2str(binocdf(2, n, p))])
disp(['Probability that none is underweight = ' num2str(binopdf(0, n, p)) ])
plot([m; m], [zeros(1, n+1); binopdf(m, n, p)], 'k', m, binocdf(m, n, p), 'k*')
Answers
expected value = 4.8 variance = 3.84
probability that no more than 2 are underweight = 0.11452
probability that none is underweight = 0.0047224
1
0.9
0.8
0.7
0.6
0.5
0.4
0.3
0.2
0.1
0
0 5 10 15 20 25
% Exercise 14.5
function Exercise14_5
%(a)
phat=[0.6 0.25 0.15]; x=[32 14 9];
n=sum(x); e=n*phat;
x2test=sum((x-e).^2./e);
disp(['pvalue = ' num2str(1-chi2cdf(x2test, length(x)-1))])
%(b)
xo=[3 14 20 25 14 6 2 0 1];
n=sum(xo); lambda =3; cutoff=5;
eo=[n*poisspdf(0:length(xo)-1, lambda) n*(1-poisscdf(length(xo)-1, lambda))];
[vee, vxx]=chicombine(eo, [xo 0], cutoff);
disp([' xi ei'])
disp([num2str(vee') repmat(' ', length(vee), 1) num2str(vxx')])
2
,x2test=sum((vxx-vee).^2./vee);
disp(['pvalue = ' num2str(1-chi2cdf(x2test, length(vee)-1))])
function [ve, vx]=chicombine(eo, xo, cutoff)
eon=zeros(1, length(eo)); xon=zeros(1, length(eo));
n=1; nend=length(eo); cnt=1;
while n<=nend
z=cumsum(eo(n:nend));
ind=find(z<cutoff);
nn=n+max(ind);
if max(ind)>0
if nn<=nend
eon(cnt)=sum(eo(n:nn));
xon(cnt)=sum(xo(n:nn));
n=n+ind+1;
cnt=cnt+1;
else
eon(cnt-1)=eon(cnt-1)+sum(eo(n:nend));
xon(cnt-1)=xon(cnt-1)+sum(xo(n:nend));
n=n+ind+1;
cnt=cnt+1;
end
else
eon(cnt)=eo(n);
xon(cnt)=xo(n);
cnt=cnt+1;
n=n+1;
end
end
ve=eon(1:max(find(eon)));
vx=xon(1:max(find(eon)));
Answers
(a)
pvalue = 0.94979
(b)
xi ei
16.9276 17
19.0436 20
19.0436 25
14.2827 14
3
, 8.5696 6
7.13303 3
pvalue = 0.40592
% Exercise 14.6
figure(1)
t=linspace(0, 1.5);
r=1-expcdf(t);
plot(t, exppdf(t)./r, 'k-', t, r, 'k--')
title('Hazard rate and reliability for exponential distribution')
legend('Hazard rate','Reliability')
figure(2)
b=[0.5 1 2 4];
for k=1:4
subplot(2, 2, k)
r=1-weibcdf(t, 1, b(k));
plot(t, weibpdf(t, 1, b(k))./r, 'k-', t, r, 'k--')
title(['Weibull: beta = ' num2str(b(k))])
legend('Hazard rate','Reliability')
end
Answers
Weibull: beta = 0.5 Weibull: beta = 1
5 1.2
Hazard rate and reliability for exponential distribution Hazard rate Hazard rate
1.1 Reliability Reliability
4 1
Hazard rate
Reliability
1 3 0.8
2 0.6
0.9
1 0.4
0.8
0 0.2
0 0.5 1 1.5 0 0.5 1 1.5
0.7
Weibull: beta = 2 Weibull: beta = 4
0.6 3 14
Hazard rate Hazard rate
2.5 12
Reliability Reliability
0.5 10
2
8
1.5
0.4 6
1
4
0.3
0.5 2
0 0
0.2 0 0.5 1 1.5 0 0.5 1 1.5
0 0.5 1 1.5
% Exercise 14.7
dat=[1.55 3.05 3.65 5.20 7.75 10.45 10.85 10.90 12.65 15.25 ...
15.70 16.35 17.70 17.95 19.45 19.80 20.05 32.75 35.45 49.35];
4
% An Engineer’s Guide to MATLAB, 2nd ed.
% ©Prentice Hall 2005 Saddle River NJ
% Exercise 14.1
%(a)
disp(['probability of 8 calls = ' num2str(poisspdf(8, 5))])
disp(['probability of 2 calls = ' num2str(poisspdf(2, 5))])
%(b)
disp(['probability of a busy signal = ' num2str(1-poisscdf(10, 5))])
Answers
probability of 8 calls = 0.065278
probability of 2 calls = 0.084224
probability of a busy signal = 0.013695
% Exercise 14.2
disp(['Probability that at least 12 withstand load = ' num2str(1-binocdf(11, 15, 0.7))])
Answers
Probability that at least 12 withstand load = 0.29687
% Exercise 14.3
dat=[88.4 93.2 87.4 94.3 93.0 94.3 89.0 90.5 90.8 93.1 92.8 91.9;...
92.6 93.2 89.2 94.8 93.3 94.0 93.2 91.7 91.5 92.0 90.7 93.8];
for k=1:2
el(k)=var(dat(k,:))+(mean(dat(k,:))-92)^2;
end
disp(['L1 = ' num2str(el(1)) ' L2 = ' num2str(el(2))])
Answers
L1 = 5.5904 L2 = 2.6936
% Exercise 14.4
1
,p=0.2; n=24; m=0:n;
disp(['Expected value = ' num2str(n*p) ' variance = ' num2str(n*p*(1-p))])
disp(['Probability that no more than 2 are underweight = ' num2str(binocdf(2, n, p))])
disp(['Probability that none is underweight = ' num2str(binopdf(0, n, p)) ])
plot([m; m], [zeros(1, n+1); binopdf(m, n, p)], 'k', m, binocdf(m, n, p), 'k*')
Answers
expected value = 4.8 variance = 3.84
probability that no more than 2 are underweight = 0.11452
probability that none is underweight = 0.0047224
1
0.9
0.8
0.7
0.6
0.5
0.4
0.3
0.2
0.1
0
0 5 10 15 20 25
% Exercise 14.5
function Exercise14_5
%(a)
phat=[0.6 0.25 0.15]; x=[32 14 9];
n=sum(x); e=n*phat;
x2test=sum((x-e).^2./e);
disp(['pvalue = ' num2str(1-chi2cdf(x2test, length(x)-1))])
%(b)
xo=[3 14 20 25 14 6 2 0 1];
n=sum(xo); lambda =3; cutoff=5;
eo=[n*poisspdf(0:length(xo)-1, lambda) n*(1-poisscdf(length(xo)-1, lambda))];
[vee, vxx]=chicombine(eo, [xo 0], cutoff);
disp([' xi ei'])
disp([num2str(vee') repmat(' ', length(vee), 1) num2str(vxx')])
2
,x2test=sum((vxx-vee).^2./vee);
disp(['pvalue = ' num2str(1-chi2cdf(x2test, length(vee)-1))])
function [ve, vx]=chicombine(eo, xo, cutoff)
eon=zeros(1, length(eo)); xon=zeros(1, length(eo));
n=1; nend=length(eo); cnt=1;
while n<=nend
z=cumsum(eo(n:nend));
ind=find(z<cutoff);
nn=n+max(ind);
if max(ind)>0
if nn<=nend
eon(cnt)=sum(eo(n:nn));
xon(cnt)=sum(xo(n:nn));
n=n+ind+1;
cnt=cnt+1;
else
eon(cnt-1)=eon(cnt-1)+sum(eo(n:nend));
xon(cnt-1)=xon(cnt-1)+sum(xo(n:nend));
n=n+ind+1;
cnt=cnt+1;
end
else
eon(cnt)=eo(n);
xon(cnt)=xo(n);
cnt=cnt+1;
n=n+1;
end
end
ve=eon(1:max(find(eon)));
vx=xon(1:max(find(eon)));
Answers
(a)
pvalue = 0.94979
(b)
xi ei
16.9276 17
19.0436 20
19.0436 25
14.2827 14
3
, 8.5696 6
7.13303 3
pvalue = 0.40592
% Exercise 14.6
figure(1)
t=linspace(0, 1.5);
r=1-expcdf(t);
plot(t, exppdf(t)./r, 'k-', t, r, 'k--')
title('Hazard rate and reliability for exponential distribution')
legend('Hazard rate','Reliability')
figure(2)
b=[0.5 1 2 4];
for k=1:4
subplot(2, 2, k)
r=1-weibcdf(t, 1, b(k));
plot(t, weibpdf(t, 1, b(k))./r, 'k-', t, r, 'k--')
title(['Weibull: beta = ' num2str(b(k))])
legend('Hazard rate','Reliability')
end
Answers
Weibull: beta = 0.5 Weibull: beta = 1
5 1.2
Hazard rate and reliability for exponential distribution Hazard rate Hazard rate
1.1 Reliability Reliability
4 1
Hazard rate
Reliability
1 3 0.8
2 0.6
0.9
1 0.4
0.8
0 0.2
0 0.5 1 1.5 0 0.5 1 1.5
0.7
Weibull: beta = 2 Weibull: beta = 4
0.6 3 14
Hazard rate Hazard rate
2.5 12
Reliability Reliability
0.5 10
2
8
1.5
0.4 6
1
4
0.3
0.5 2
0 0
0.2 0 0.5 1 1.5 0 0.5 1 1.5
0 0.5 1 1.5
% Exercise 14.7
dat=[1.55 3.05 3.65 5.20 7.75 10.45 10.85 10.90 12.65 15.25 ...
15.70 16.35 17.70 17.95 19.45 19.80 20.05 32.75 35.45 49.35];
4
