# NOT RUN {
mu = numeric(8)
dfun = function(x) dgamma(x,10,5)
mu[1] = integrate(function(x) x*dfun(x),0,Inf)$val
for (j in 2:8) mu[j] = integrate(function(x) (x-mu[1])^j*dfun(x),0,Inf)$val
Sigma_fun(mu)
# }
Run the code above in your browser using DataLab