n = 1000
mygamma = 1 # log link ==> 0 is the answer
delta = 0
y = delta + mygamma / rnorm(n)^2 # This is Levy(mygamma, delta)
# Cf. Table 1.1 of Nolan for Levy(1,0)
sum(y > 1) / length(y) # Should be 0.6827
sum(y > 2) / length(y) # Should be 0.5205
fit = vglm(y ~ 1, levy(delta=delta), trace=TRUE) # 1 parameter
fit = vglm(y ~ 1, levy(idelta=delta, igamma=mygamma),
trace=TRUE) # 2 parameters
coef(fit, matrix=TRUE)
Coef(fit)
summary(fit)
weights(fit, type="w")[1:4,]
Run the code above in your browser using DataLab