How can I calculate the median in a gamma distribution?

Question

I have the below distribution, which I want to calculate the median of:

x=0:0.01:10;
x=[x' x' x' x' x'];
a=ones(1,1001)';
a=[a*2 a*4 a*6 a*8 a*10];
b=2;
f = gampdf(x,a,b);
plot(x,f)
grid on

Visualization of Distribution

You need to find the value `m` for which the integral from 0 to `m` give you 0.5. You can do that with `c = cumsum(f)*dx` where `dx = 0.01` in your case. After it's just a matter of using `find(c>0.5, 1, 'first')`. — TwistedSim, May 01 '18 at 14:46

score 1 · Accepted Answer · answered May 01 '18 at 16:00

1

User TwistedSim answered my question.

You need to find the value m for which the integral from 0 to m give you 0.5. You can do that with c = cumsum(f)*dx where dx = 0.01 in your case. After it's just a matter of using find(c>0.5, 1, 'first').

answered May 01 '18 at 16:00

Valerii Bosiak

53
7

Mh... apart from the fact that this solution only provides an answer for the first distribution (there are `5` here, one in each column of `f`)... `72` doesn't seem like a good median value to me. I suppose it's an index... but yet an index to something that doesn't sound right. A correct median for `a=2` and `b=2` should be approximately `3.35`. – Tommaso Belluzzo May 01 '18 at 17:57

score 0 · Answer 2 · answered May 01 '18 at 21:47

dx = 0.04;
b=2;

x=kron([0:dx:10]',ones(1,5));
a=kron(ones(size(x,1),1),[2:2:10]);
f = gampdf(x,a,b);
cf = cumsum(f)*dx;
[i,j] = find(cf(1:end-1,:)<0.5 & cf(2:end,:)>=0.5);

cflow  = cf(sub2ind(size(x),i  ,j));
cfhigh = cf(sub2ind(size(x),i+1,j));

xm = x(i)+dx/2 + (0.5-cflow)./(cfhigh-cflow) * dx
fm = gampdf(xm',a(1,:),b)';

plot(x,f)
hold on; plot([xm xm]',[fm fm]','*'); hold off
grid on

How can I calculate the median in a gamma distribution?

2 Answers2