Page 236 - Computational Statistics Handbook with MATLAB

Chapter 6: Monte Carlo Methods for Inferential Statistics 223

which is the sample second central moment. We write our own simple func-
tion called mom (included in the Computational Statistics Toolbox) to estimate
this.
% This function will calculate the sample 2nd
% central moment for a given sample vector x.
function mr = mom(x)
n = length(x);
mu = mean(x);
mr = (1/n)*sum((x-mu).^2);
We use this function as an input argument to bootstrp to get the bootstrap-
t confidence interval. The MATLAB code given below also shows how to get
the bootstrap estimate of standard error for each bootstrap sample. First we
load the data and get the observed value of the statistic.

load forearm
n = length(forearm);
alpha = 0.1;
B = 1000;
thetahat = mom(forearm);
Now we get the bootstrap replicates using the function bootstrp. One of
the optional output arguments from this function is a matrix of indices for the
resamples. As shown below, each column of the output bootsam contains
the indices to a bootstrap sample. We loop through all of the bootstrap sam-
ples to estimate the standard error of the bootstrap replicate using that resa-
mple.

% Get the bootstrap replicates and samples.
[bootreps, bootsam] = bootstrp(B,'mom',forearm);
% Set up some storage space for the SE’s.
sehats = zeros(size(bootreps));
% Each column of bootsam contains indices
% to a bootstrap sample.
for i = 1:B
% Extract the sample from the data.
xstar = forearm(bootsam(:,i));
bvals(i) = mom(xstar);
% Do bootstrap using that sample to estimate SE.
sehats(i) = std(bootstrp(25,'mom',xstar));
end
zvals = (bootreps - thetahat)./sehats;
Then we get the estimate of the standard error that we need for the endpoints
of the interval.
% Estimate the SE using the bootstrap.
SE = std(bootreps);

© 2002 by Chapman & Hall/CRC