|
1 | 1 | %% Windsor House Price data. |
2 | | -% This file creates Figures 8.1-8.7. |
| 2 | +% This file creates Figures 8.1-8.7 |
| 3 | +% and Table 8.1 |
3 | 4 |
|
4 | 5 | %% Beginning of code |
5 | 6 | prin=0; |
6 | 7 |
|
7 | 8 | load hprice.txt; |
8 | 9 |
|
9 | 10 | % setup parameters |
10 | | -n=size(hprice,1); |
| 11 | + |
11 | 12 | y=hprice(:,1); |
12 | 13 | X=hprice(:,2:5); |
13 | 14 |
|
| 15 | +n=size(hprice,1); |
| 16 | +p=size(X,2)+1; |
| 17 | + |
14 | 18 | bayes=struct; |
15 | | -n0=5; |
| 19 | +n0=6; |
16 | 20 | bayes.n0=n0; |
17 | 21 |
|
18 | | -% set \beta components |
| 22 | +% beta prior setting |
19 | 23 | beta0=zeros(5,1); |
20 | 24 | beta0(2,1)=10; |
21 | 25 | beta0(3,1)=5000; |
22 | 26 | beta0(4,1)=10000; |
23 | 27 | beta0(5,1)=10000; |
24 | 28 | bayes.beta0=beta0; |
25 | 29 |
|
26 | | -% \tau |
| 30 | +% tau0 prior setting |
27 | 31 | s02=1/4.0e-8; |
28 | 32 | tau0=1/s02; |
29 | 33 | bayes.tau0=tau0; |
|
34 | 38 | R(3,3)=.15; |
35 | 39 | R(4,4)=.6; |
36 | 40 | R(5,5)=.6; |
| 41 | +diagR=diag(R); |
37 | 42 | R=inv(R); |
38 | 43 | bayes.R=R; |
39 | 44 |
|
40 | 45 | %% Create initial yXplot |
41 | 46 | % Not given in the book |
42 | 47 | % yXplot(y,X) |
43 | 48 |
|
44 | | -%% Create Figures 8.1 and 8.4 |
| 49 | +%% Create Table 8.1 |
| 50 | +% Windsor House Price data: prior and least squares estimates of parameters |
45 | 51 |
|
46 | | -outBA=FSRB(y,X,'bayes',bayes', 'plots',1,'xlim',[280 n]); |
47 | | -dout=n-length(outBA.ListOut); |
| 52 | +betaLS=regstats(y,X,"linear",{'beta', 'mse'}); |
| 53 | +betaOLS=betaLS.beta; |
| 54 | +tauOLS=1/(betaLS.mse*(n-p)/n); |
48 | 55 |
|
49 | | -fig=findobj(0,'tag','fsr_yXplot'); |
50 | | -figure(fig(1)) |
51 | | -set(gcf,'Tag','fsr_yXplot81') |
| 56 | +T1=[[beta0; tau0], [diagR; n0], [betaOLS; tauOLS]]; |
| 57 | +nam=["beta"+(0:4) "tau"]; |
| 58 | +T1table=array2table(T1,'RowNames',nam,'VariableNames',["mean" "V0jj" "Least squares"]); |
| 59 | +disp('Table 8.1') |
| 60 | +disp(T1table) |
52 | 61 |
|
53 | | -if prin ==1 |
54 | | - print -depsc h4.eps; |
55 | | -else |
56 | | - set(gcf,'Name', 'Figure 8.4'); |
57 | | - sgtitle('Figure 8.4') |
58 | | -end |
| 62 | +%% Create Figure 8.1 |
| 63 | + |
| 64 | +outFSRB=FSRB(y,X,'bayes',bayes', 'plots',1,'xlim',[280 n]); |
| 65 | +dout=n-length(outFSRB.ListOut); |
| 66 | + |
| 67 | +fig=findobj(0,'tag','fsr_yXplot'); |
| 68 | +delete(fig(1)) |
59 | 69 |
|
60 | 70 | fig=findobj(0,'tag','pl_fsr'); |
61 | 71 | figure(fig(1)) |
|
234 | 244 |
|
235 | 245 | drawnow |
236 | 246 |
|
| 247 | +%% Figure 8.4 will be created at the end of the file |
| 248 | +% Note that Figure 8.4 is created later because we need both the output of FSR |
| 249 | +% (frequentist forward search) and FSRB (Bayesian forward search) |
| 250 | + |
237 | 251 | %% Create Figure 8.5 (frequentist analysis) |
238 | 252 | % Monitoring of 95 per cent and 99 per cent confidence intervals |
239 | 253 | % of beta and sigma2 |
|
317 | 331 |
|
318 | 332 | drawnow |
319 | 333 |
|
| 334 | + |
320 | 335 | %% Create Figure 8.6 |
321 | 336 | figure |
322 | 337 | n0=250; |
|
345 | 360 | set(gcf,'Name', 'Figure 8.7'); |
346 | 361 | sgtitle('Figure 8.7') |
347 | 362 | end |
| 363 | + |
| 364 | + |
| 365 | +%% Create Figure 8.4 |
| 366 | +% Note that plot is created later because we need both the output of FSR |
| 367 | +% (frequentist forward search) and FSRB (Bayesian forward search) |
| 368 | +group=ones(n,1); |
| 369 | +commonOutliers=intersect(outFSR.outliers,outFSRB.outliers); |
| 370 | +onlyFSRB=setdiff(outFSRB.outliers,outFSR.outliers); |
| 371 | +group(commonOutliers)=2; |
| 372 | +group(onlyFSRB)=3; |
| 373 | +yXplot(y,X,'group',group) |
| 374 | +legend(["Normal observations" "Outliers found both by FSR and FSRB" "Outlier found just by FSRB"]) |
| 375 | + |
| 376 | +if prin ==1 |
| 377 | + print -depsc h4.eps; |
| 378 | +else |
| 379 | + set(gcf,'Name', 'Figure 8.4'); |
| 380 | + sgtitle('Figure 8.4') |
| 381 | +end |
348 | 382 | %InsideREADME |
349 | 383 |
|
350 | 384 |
|
|
0 commit comments