I'm learning about gamma distributions and I want to make sure my reasoning is correct for an example I found here.
Suppose you are fishing and you expect to get a fish once every 1/2 hour. Compute the probability that you will have to wait between 2 to 4 hours before you catch 4 fish.
One fish every 1/2 hour means we would expect to get θ = 1 / 0.5 = 2 fish every hour on average. Using θ = 2 and k = 4.
I re-created this in
scipy and got the same answer
import scipy.stats as stats stats.gamma.cdf(4, a=4, scale=2) - stats.gamma.cdf(2, a=4, scale=2) >> 0.12388838262529915
Now I wanted to pay with the results and see if I could get answers which were intuitive.
Say instead of catching 1 fish every half hour, you instead catch 5 fish every half hour so θ = scale = 10 (10 fish on average per hour). And I want to see the probability I have to wait between wait between 4 and 6 hours to catch 50 fish.
My shape parameter here is 50 (I believe), and I use the same formula from above but with tweaking a few things.
stats.gamma.cdf(6, a=50, scale=10) - stats.gamma.cdf(4, a=50, scale=10) >> 1.4758744479975356e-76
I have two questions:
- Does the above answer seem correct? What is the intuition behind it?
- My original reasoning was if I'm catching on average 10 fish an hour, it would be probable (at least higher than my answer) I would catch 50 somewhere between 4 and 6 hours. Is this incorrect?
Any guidance would be appreciated.
I believe the intuition is (scale * a) = average waiting time. So in my example 50*10 = 500 hours for waiting time ... which would be inline with my results almost zero for the second example.
So I was inputting the scale wrong based on the assumptions to my question.
So, if I revise the inputs and choose 1/5 as my scale. The average waiting time is 10 hours to catch 50 fish. Now, I get a probably (below) which is more in line with expectations.
scale = 1/5 (10 / 50 = 1/5) ... and I want to get the probability I have to wait between wait between 12 and 8 hours to catch 50 fish →
stats.gamma.cdf(12, a=50, loc=0, scale=1/5) - stats.gamma.cdf(8, a=50, loc=0, scale=1/5) 0.8452582522469133
Does this seem correct?
I think that answer makes sense. No?
I believe I had my rate/scale parameters mixed up, when estimating the distribution. I was under the impression if I had the scale/rate, I could input it directly into the distribution parameter. My shape parameter (however) was correct.
This is not the case.
As well, choosing scale or rate, will change how the parameters are chosen. If I chose to use alpha and beta as my parameters in my gamma distribution (X ~ GAMMA(alpha, beta) I need to calculate them using these formulas:
mean = alpha/beta and
var = alpha/beta^2
So, using a modified example from above ... if I'm getting 10 fish per hour I get
10 = alpha/beta and if my variance of fish is 4 then
Now I can solve for alpha and beta using
sympy as shown below.
from sympy import symbols, Eq, solve # 10 fish per hour on AVG (1/10 is the rate) # variance of 4 fish per hour … not this is NOT std a, b = symbols('a b' ) eq1 = Eq(a/b - 10) # Use the equation for E[X] eq2 = Eq((a)/b**2 - 4) # Use the equation for V[X] solve((eq1, eq2), (a, b)) # [(25, 5/2)]
Therefore I would use alpha=25 and beta=5/2 for my inputs, and solve the same as above.
Hope this helps someone in the future who is struggling with questions such as this.