5 Bootstrapping

The bootstrap function makes use of both the simulate and fitmodel functions and performs parametric bootstrapping from a given point in the parameter space. It takes four arguments; the name of the migration model, a parameter vector specifying at which point in the parameter space to simulate from, the observed gene frequency matrix, and, optionally, the number of bootstrap replicates to produce (which by default is 100).

The distribution of the data, conditioned on the observed gene frequency mean at each loci, is simulated. This is why the observed gene frequency matrix is needed by bootstrap (to compute the means). 50 bootstrap replicates from the maximum likelihood estimates obtained in section 4 is produced by the following call: ¹

   > boots <- bootstrap(steppingstone,fit1$par,p,nboots=50)

If you get the feeling that this call has caused your computer to hang, you can set the optional argument dots=T. Then, at least, some small dots start to appear every now and then.

The bootstrap function returns a matrix with rows corresponding to each bootstrap replicate of the parameter vector:

   > boots
                [,1]       [,2]
    [1,] 0.002350949 0.12999907
    [2,] 0.002001756 0.09736528
    [3,] 0.001058648 0.14776179
    [4,] 0.009843646 0.11854703
     .        .           .
     .        .           .

The variance-covariance matrix of the replicates is thus easy to compute using S-Plus' function

   > var(boots)
                 [,1]          [,2]
   [1,]  5.317043e-05 -0.0000885599
   [2,] -8.855990e-05  0.0007198054

Taking the square root of the diagonal elements gives us an estimate of the standard errors

   > sqrt(diag(var(boots)))
   [1] 0.007291805 0.026829190

which aren't to large -- the coefficient of variation for $m_0$ is only 0.2. The uncertainty in $u$ is quite big though.