STA 506 Assignment 4: Method of Moment Estimation

1 Problem 2

For this problem, we are given a data set that consists of recovery times (in days) after a specific knee surgery. This data comes from a log-logistic distribution, also known as a Fisk distribution. We will estimate alpha and beta two ways: first is via the method of moments, and the second is via bootstrapping.

time = c(8.23, 12.74, 14.83, 16.61, 18.16, 19.55, 20.80, 21.94, 23.00, 23.98, 24.89, 25.75, 26.56, 27.34, 28.08, 28.79, 29.48, 30.15, 30.81, 31.45, 32.08, 32.70, 33.31, 33.92, 34.53, 35.13, 35.73, 36.33, 36.93, 37.53, 38.14, 38.75, 39.37, 40.00, 40.64, 41.29, 41.95, 42.63, 43.33, 
44.05, 44.79, 45.56, 46.36, 47.20, 48.08, 49.02, 50.03, 51.12, 52.32, 53.65)

1.1 Part A

First we will estimate alpha and beta via method of moments. In order to find our moments, we will turn to the internet.

The first moment of the distribution is: \[ \mu_1 = \frac{\alpha\pi/\beta}{\sin(\pi/\beta)} \]

and the second moment is: \[ \mu_2 = \frac{2\pi\alpha^2/\beta}{\sin(2\pi/\beta)} \]

Next, let’s caclulate the sample moments:

m1.samp <- mean(time)
m2.samp <- mean(time^2)
m2.sampvar <- var(time)

m1.samp

[1] 34.1922

m2.samp

[1] 1288.845

m2.sampvar

[1] 122.1819

Now we can create a function to estimate beta and alpha. Starting with beta:

## Estimating Beta

beta.function <- function(beta) {
 m1.act <- (pi/beta)/sin(pi/beta)
 m2.act <- (2*pi/beta)/sin(2*pi/beta)
 
alpha <- m1.samp/m1.act
 
var.theory <- alpha^2*(m2.act-m1.act^2)
var.theory - m2.sampvar
 
 
}

beta.est <- uniroot(beta.function, interval = c(2.01, 20)) 
 # Beta needs to be greater than 2
beta.est$root

[1] 5.95252

We get an estimate for alpha to be 5.9525.

Next, let’s use the estimate for beta to estimate alpha:

alpha.est <- m1.samp/((pi/beta.est$root)/sin(pi/beta.est$root))
alpha.est

[1] 32.62681

We get an alpha of 32.62681.

To summarize so far: Through the method of moments, we estimated the parameters of the distribution from which our sample of times may have came from. Our parameter estimates are: \[ \alpha = 32.62681, \beta = 5.95252 \]

1.2 Part B

Next, we estimate alpha and beta via bootstrapping.