Calculate confidence interval for gamma function

Hi!

I am supposed to compute the confidence intervals for a gamma distributed sample (samplesize n=5, mean value mu=1, standard deviation=2 and alfa=1, beta=1) and repeat the simulation 10 000 times with a for-loop. The goal is to compute how many of the confidence intervals that contains the true mean value when the sample size is 5 and compare it to a samplesize of 100.

When using the following code I get that all of the confidence intervals contains the mean value when n=100 and that 97,45% of the intervals contain the mean value when n=5 which is too large. The answer is supposed to be less than 95% when n=5 and about 95% when n=100. When calculating the CI I should assume normal distribution?

Can anyone see what I am doing wrong?

Thank you!

clc

clear

a=1;

b=1;

n = 5; % samplesize

nu = n - 1; % degrees of freedom

conf = 0.95; %konfidensgraden är 0.95

alpha = 1 - conf;

pLo = alpha/2; %lower limit

pUp = 1 - alpha/2; %upper limit

number=0; % number of samples containing the true mean value before the loop

for i=1:10000

g = gamrnd(a,b,n,1);

gmean = mean(g); % calculated mean value

se = std(g); % calculated standard deviation

crit = tinv([pLo pUp], nu); %

ci = gmean + crit*se % calculating the confidence interval

if ci(1)<1<ci(2) % checking if the true mean value (1) is withing the CI

number=number+1 %number of times the mean value is within the CI

end

end

proportion=number/10000

Answer 1

[Star Strider]

I am not certain what you are permitted to use. I would use fitdist then paramci. With respect to assuming a distribution of the means, consider the Central Limit Theorem.

[Elin Wildlock]

Thanks! We are supposed to use mean, std, tinv and calculate it as an "usual confidence interval".

[Star Strider]

My pleasure!

The t-distribution approximates the normal distrbution with greater than about 30 degrees-of-freedom. The Central Limit Theorem would seem to apply here.

If my Answer helped you solve your problem, please Accept it!

.

Answer 2

[Paul]

A few comments:

With a = 1, b = 1, the variance of the distribution is 1, not 4 (as implied by your statement that the standard deviation = 2)

>> var(gamrnd(a,b,10000,1))

ans =

    0.9924

I think your missing a sqrt(n) term in your caculation. Change to:

ci = gmean + crit*se/sqrt(n); % calculating the confidence interval

Finally, this line, in general, doesn't do what you think

if ci(1)<1<ci(2) % checking if the true mean value (1) is withing the CI

For example, this expression evaluates to false

>> -2 < -1.5 < -1

ans =

  logical

   0

So you should change that to

if ci(1) <= 1 && 1 <= ci(2)  % replace 1 with a*b to be more general?

After I made these two changes, the result came in at ~0.95 for n=100 and ~0.88 for n = 5. I don't know for sure what's going on with that second result but suspect it has something to do with the small sample size. Is there a different CI formula for small sample sizes?

[Elin Wildlock]

Hi, Paul and thank you so much for taking your time. No we were supposed to use the same formula for both to see the difference when we increase the sample size.