February 21, 2012 --- Class 13 --- Review of Biased Estimators,
Application of the Jackknife Method
Activities:
Biased Esimators Using Mathematica
We quickly reviewed the example of biased estimators in the
Mathematica handout. We made use of the Mathematica notebooks in
~sg/Documents: feb16_v1.nb, feb16_v2.nb and October_10_2006.nb
on the Apple cluster in Swain 246.
We found that for small samples sizes, the integrals, which
represent the average value for 1/(mean value) we get numbers
larger than 2, even though for an infinite sample we expect 2.
We plotted the values obtained for different samples sizes vs
1/(sample size) and found that the curve was linearly approaching
2 as 1/(sample size) approaches zero.
Bias reduction
Recognizing that a plot of the result vs. 1/n shows the error is
linear in 1/n, where n is the sample size:
S_n = A +e/n
where A is the limit of n -> infinity and e determines the size of
the error. Considering samples sizes n and n-1, it is easy to
solve for A.
A = n S_n - (n-1) S_{n-1}
The jackknife method allows us to resample our n values to also
get samples of size (n-1). If we have a sample of size n, we use
it to compute the statistic S_n. Instead of getting a new sample
with n-1 new variables, we throw one of the values out of our
sample of size n to create a sample of size n-1 and use that to
compute S_{n-1}. However, we can throw any of the values out, so
we might as well take the average over throwing out each of the n
values. Let
S'^j_{n-1}
represent the statistic on the sample of size n-1 from which the jth
element has been eliminated. Our estimate for S_{n-1} is
(1/n) \sum_{j=1}^n S'^j_{n-1}
Biased Estimators and the Jackknife method
A program ~sg/jackknife/jackknife.c implements the jackknife method
of bias reduction. The jackknife program can take up to two
arguments. The first is the sample size, the second, optional
argument, is the number of samples to create. On each line of
output is the sample estimate of 1/(sum of x_i), and the jackknife
value (which is corrected for bias). We examined the code an how
it implements the formulae displayed on the board.
To analyze this method, we pick a bunch of samples sizes to study
and generate a large number of samples for each size.
To see the difference between the naive estimate and the bias
corrected values, two programs hist and blockerr are useful.
(Another simple technique is to use axis to create a scatterplot
using the m 0 option so the points are not connected. A diagonal line
can easily be added to the plot. We see that the points in the
scatterplot are mostly on one side of the diagonal. The jackknife
values tend to be smaller than the naive value. This is most dramatic
for small sample sizes.)
blockerr is a simple program for blocking numbers together and
reporting the mean value and its error. Being able to block the
data is useful for correlated data. We will deal with this
complication later. For now, just use a blocksize of 1.
The shell script ~sg/jackknife/runem.csh is set up to make 10000
samples for each sample size. The range of samples sizes is from 2
to 25. Using awk and the program blockerr, the results are
averaged and place in a file results.
Here is a shell script to produce results for a variety of
sample sizes
#!/bin/csh
# run through different samples sizes
foreach i (2 4 6 8 10 12 14 18 20 25)
jackknife $i 10000 >samples_$i
echo -n Naive $i " " >>! results
awk '{print $1}' samples_$i | blockerr 1 |awk '{print $8, $10}' >>results
echo -n Jackknife $i " " >>! results
awk '{print $2}' samples_$i | blockerr 1 |awk '{print $8, $10}' >>results
end
The results file can easily be displayed with the command:
awk '{print $2,$3,$4}' results |axis e |xplot
The upper set of values is from the naive estimate, and the lower
is the jackknife corrected values. One immediately sees that the
jackknife values much more quickly approach 2, the asymptotic
value.
The file ~sg/jackknife/jackknife.ax on Swain cluster, was
prepared before class, and includes both the above results
and the analytic results derived in the last
class. The crosses are from doing the integral numerically with
Mathematica. One can see excellent agreement between the
statistical approach here and the analytic results.
The following supplementary material was not covered in class.
Other approaches to removing the sample bias include using the
basic idea behind the jackknife bias reduction of the Mathematica
generated results, and applying a "second order" jacknife technique
on the jackknifed results to remove the 1/n^2 errors.
(Assuming there are error terms e/n + e'/n^2, derive a formula
involving S_n, S_(n-1) and S_(n-2) that removes both error terms to
get a second order bias reduced result.)
Applying the jackknife to the Mathematica results works
quite nicely as you can see by the following command:
cat jackknife.ax j.math.ax |axis |xplot
The lower crosses are the results from Mathematica. The commands
for doing the analysis can be found in ~sg/jackknife/bias2.math