Question

Use the Monte Carlo method to estimate the mean, variance, and 2.5th and 97.5th percentiles of a gamma distribution G(2, 4). Compare the estimates with the true values and comment how the sample size affects the estimates.

Answer 1

Given information: shape = 2, scale = 4 for the gamma distribution under study (rate = 1/scale = 1/4)

Let us take different number of samples to study effect of sample size on estimates

R-code and Solution :

ss = c(50,100,200,500,1000,10000,100000)

sh = 2

sc = 4

i = 0

means = list()

vars = list()

quant_025 = list()

quant_975 = list()

for(s in ss)

{

i = i+1

sample_series = rgamma(n = s, shape = sh, scale = sc)

means[i] = mean(sample_series)

vars[i] = var(sample_series)

quant_025[i] = qgamma(p = c(0.025), shape = sh, scale = sc)

quant_975[i] = qgamma(p = c(0.975), shape = sh, scale = sc)

}

i = 0

true_mean = sh*sc

# 8

true_variance = sh*sc^2

# 32

# Means:

## True mean = 8

## as.matrix(means)

# [1] 9.82338

# [2] 8.970444

# [3] 7.859794

# [4] 8.204468

# [5] 7.613365

# [6] 8.055619

# [7] 7.990096 --> As sample size increases, the estimated mean tends to the theoretical mean

# Variances:

## True variance = 32

## as.matrix(vars)

# [1] 35.64035

# [2] 54.33689

# [3] 29.01417

# [4] 32.25097

# [5] 27.00697

# [6] 32.92603

# [7] 32.10989 --> As sample size increases, the estimated variance tends to the theoretical variance

# Quantiles:

## as.matrix(quant_025) and as.matrix(quant_975)

## 2.5% 97.5%

# [1] 0.9688371 22.28657

# [2] 0.9688371 22.28657

# [3] 0.9688371 22.28657

# [4] 0.9688371 22.28657

# [5] 0.9688371 22.28657

# [6] 0.9688371 22.28657

# [7] 0.9688371 22.28657

R-code and Solution

Explanation: