Chapter 9

Chapter 9

9.7 Write a function that takes a sample from a Weibull distribution and estimates its parameters

log(log(1/1-y)) = (k * log(x)) - (k * log(lambda))

Given this, we should be able to take the log(log(1/1-y)) and the log(x) to get a line with slope equal to the shape parameter k and intercept k * log(lambda).

I first did this in R to iterate on the transformation. Here is the code:

# shape = k
shape = 2
# scale = lambda
scale = 1

weibull.rsample = sort(rweibull(1e5, shape, scale))
weibull.ecdf = (1:length(weibull.rsample))/length(weibull.rsample)
df = data.frame(cbind(x=weibull.rsample, y=weibull.ecdf))
df$xp = log(df$x)
df$yp = log(log(1/(1-df$y)))

# exclude y = 1
model = lm(df$yp ~ df$xp, subset=seq(1,length(df$yp)-1))
summary(model)
# -shape * log(scale) = ~0
intercept = coef(model)[1]
# shape ~ 2
slope = coef(model)[2]

And then in clojure, I wrote a test:

(facts :least-squares
  ; use model fitting of a Weibull distribution to test the correct resolution
  ; of the shape and scale parameters from a liner transformation.
  (let [k 0.5
        lambda 1.0
        w (repeatedly 100000 #(random/weibullvariate k lambda))
        [x y] (random/weibull-line w)
        model (cor/least-squares x y)]
    ; the slope should be k
    (:slope model) => (roughly k 0.1)
    ; the intercept should be 0
    (:intercept model) => (roughly (* k (Math/log lambda)) 0.006)))