0

I am trying to generate a set of random positive steps that will result in a final location that is close to what I would have gotten from taking a similar number of fixed steps$^1$. I would like the step size to be sampled from gamma distribution$^2$ that has a peak at the fixed step size. If I understand correctly, this means that I want a gamma distribution whose mode is equal to the median.

$$mode = (k − 1) \theta$$ $$median = \theta \gamma^{-1}(k, \frac{\Gamma(k)}{2})$$

Setting these equal to each other and taking $\gamma$ of both sides yields

$$2 \gamma(k, k - 1) = \Gamma(k)$$

Which can be reworked to

$$\gamma(k, k - 1) = \Gamma(k, k - 1)$$

Where $\gamma$ is the lower incomplete gamma, and $\Gamma$ is the upper incomplete gamma function.

There is no closed form solution for $k$ that I am aware of to this problem. Having played with the integrals for a bit, I am unsure how to proceed to find the shape parameter of my distribution that will ensure that the mode and median coincide (getting the scale $\theta$ from the shape is trivial of course).

I am perfectly OK with a numerical solution, but I am not sure how to obtain even that, using commonly defined functions$^3$.

In short, how do I find the shape parameter of a gamma function whose median is equal to its mode?

Related questions, that helped me visualize up to here:

$^1$ I am doing this in Python + numpy, so I want a vector of random values x of size n such that np.cumsum(x) represents a sequence of motonically increasing x-values randomly disturbed from approximately np.linspace(x[0], x[-1], n).

$^2$ Implemented in np.random.gamma.

$^3$ In particular, scipy defines the regularized and non-regularized, complete and incomplete gammas and their inverses.

Mad Physicist
  • 138
  • At the start, it seems you want a member of the gamma family of distributions for which population median equals population mode. At the end, you seem to be dealing with samples. Please clarify. – BruceET Sep 16 '18 at 21:28
  • Not sure this is possible for distributions. Start by setting scale or rate parameter to 1. Then Wikipedia on gamma median suggests an approx. expression for median which you can set equal to expression for mode. Seems not to have possible solution for shape parameter. – BruceET Sep 16 '18 at 21:42
  • @BruceT. I'm never dealing with samples. The problem is reduced to finding the shape parameter of the distribution that will make the upper and lower half equal when the second argument is alpha-1. The scale doesn't seem to figure in the equation at all. – Mad Physicist Sep 16 '18 at 22:15
  • Right about the scale parameter. In my experience with gamma distributions, they tend to have mode > median > mean. // Are you exploring on your own or did a problem ask to to find the shape parameter that makes mode and median equal? – BruceET Sep 16 '18 at 22:31
  • @BruceET. I'm exploring on my own. I've set up the integrals, but my rudimentary attempts at solving them have been fruitless so far. I've played with a bunch of these distributions, and making alpha small enough seems to shift the mode very much to the left. When shape = 0, the distribution looks a lot like a decaying exponential, so there must be some shape for which mode is not greater than median. – Mad Physicist Sep 16 '18 at 23:47
  • @BruceT. I realize that I was using unconventional notation from another question I saw, so I fixed that. – Mad Physicist Sep 17 '18 at 00:02
  • @BruceT. Looking more carefully, I see your point though. It may be that the solution is a value of k less than 1. I hadn't considered that. Will investigate. – Mad Physicist Sep 17 '18 at 00:08
  • Technically, if $k < 1$ there is no mode because there the PDF becomes infinite at 0. If $k = 1,$ dist'n is exponential. For exponential with mean 1, mode is 0 and median is 0.6931472. // What is the motivation for this question? – BruceET Sep 17 '18 at 07:12
  • @BruceT Pretty much exactly what I wrote in the first paragraph and footnote. I'm trying to take a bunch of random steps that appear to be about the same length as the equivalent constant step, and come out to the same total distance in the long run. – Mad Physicist Sep 17 '18 at 07:34

1 Answers1

1

This is in response to your original modeling problem for steps with gamma distributed incremental distances. (I don't understand why you need the mean and mode of the gamma distributed steps to be the same.)

For example, we can use incremental distances $D_i$ identically and independently distributed as $\mathsf{Gamma}(shape=3, rate=3).$ Then $E(D_i) = 1, Var(D_i) = 1/3.$

Look at 50 such steps in time, starting at position $0$ at stage $0.$ At the $i$th stage we move to the right by $D_i.$ Then total distances after $n$ stages are $X_n = \sum_{i=1}^n D_i.$ Below we use R statistical software to simulate 50 stages.

set.seed(2018)       # use this statement to repeat EXACTLY this same simulation
d = rgamma(50, 3, 3) # incremental distances at each stage
x = c(0, cumsum(d))  # successive positions 
head(x)              # first six positions
## 0.0000000 0.6253108 1.0705457 2.0527608 4.0516030 4.7512606
tail(x)              # last six positions
## 43.76088 44.59100 45.67443 46.93781 47.70384 48.97505

The incremental distances average about 1 unit each. After 50 stages the total distance is $48.98 \approx 50$ units. Also, $E(X_{50}) = 50,\, Var(X_{50}) = 50/3,\, SD(X_{50}) = 4.0825,$ and by the Central Limit Theorem $X_{50} \stackrel{aprx}{\sim} \mathsf{Norm}(50, 4.0825),$ so that roughly 95% of the total distances after 50 stages will be within $50 \pm 8.0.$

Here is a plot with step numbers on the horizontal axis and total distance on the vertical axis. Some of the incremental distances (vertical jumps) are less than 1 and some are greater than 1, but their average is about 1.

plot(x, type="s", xlab="Step")

enter image description here

BruceET
  • 51,500
  • Not the mean. I want the median and mode to be the same. The mode is the "apparent" value, while the median ensures that the overall sum ends up the same as the constant steps. – Mad Physicist Sep 17 '18 at 09:00
  • That being said, I like the result you're getting and I'll look into it further. – Mad Physicist Sep 17 '18 at 09:03
  • In my opinion using modes does not make sense, it's means that add. Besides it seems as if you can't have gamma mode = median. What you're doing seems vaguely like queuing processes, random walks, and brownian motion. All of those use means of times or distances. – BruceET Sep 17 '18 at 09:04
  • @BruceT. If this is indeed an XY problem, I'll be more than happy to select your answer. I'll fiddle with it a bit more, because clearly my understanding is severely lacking and you've clarified a few things for me. – Mad Physicist Sep 17 '18 at 14:56
  • Glad for some progress. Not sure what process you're trying to model. Pretty sure mode is not going to be fruitful; no idea what you mean analytically by saying it's the "apparent" value. // Symmetrical dist'ns can have mean = median = mode, e.g BETA(a,a), both shape parameters equal and > 1.But haven't seen it used for this kind of application.Sums of betas are untidy. Potential problem with my $S_n$ is $Var(S_n)$ grows with $n.$ // Maybe Accept to close this Q out. Then tidy up some details, and try more focused version on 'stats.stackexchange'. State problem; say what tried, keep soln open. – BruceET Sep 17 '18 at 17:20
  • I've flagged to close as unclear what you're asking. You are 100% correct. I'm looking for the mean, not the median, and it will never be the mode. – Mad Physicist Sep 17 '18 at 20:43
  • 2
    So yeah, my equations are solving for shape when both the scale and median are 1. Not even what I thought I was asking, but that's what I get for not being careful. I'm going to try to close this but not delete it. You deserve the credit for figuring out what I was actually asking and solving it for me. – Mad Physicist Sep 17 '18 at 21:46
  • To me the important thing is that you found the exchange useful. By voting the 'Accept' you have taken the question off the list of ones awaiting answers. (Now voting to 'close' the question to additional activity, or taking steps to delete it altogether, are would be the hands of others.) – BruceET Sep 17 '18 at 22:03