Generate Random Draw from Gamma Distribution Re: Monte Carlo Simulation in R...
If the number of groups can be set "endogenously" my previous email about the smallest `n' would apply. You can view this as a "waiting time" problem. Here is one approach: x <- round(pmax(1, pmin(rgamma(500, shape=0.067, rate=0.008), 85))) csx <- cumsum(x) ind <- which(csx > 2000)[1] xg <- c(x[1:(ind-1)], 2000 - csx[ind-1]) # group sizes sum(xg) length(xg) # number of groups However, note that the elements of `xg' are not gamma variables (due to min-max constraint), and they are all not independent (due to sum constraint). Hope this helps, Ravi. ---------------------------------------------------------------------------- ------- Ravi Varadhan, Ph.D. Assistant Professor, The Center on Aging and Health Division of Geriatric Medicine and Gerontology Johns Hopkins University Ph: (410) 502-2619 Fax: (410) 614-9625 Email: rvaradhan at jhmi.edu Webpage: http://www.jhsph.edu/agingandhealth/People/Faculty_personal_pages/Varadhan.h tml ---------------------------------------------------------------------------- -------- -----Original Message----- From: r-help-bounces at r-project.org [mailto:r-help-bounces at r-project.org] On Behalf Of Hongwei Dong Sent: Tuesday, November 10, 2009 4:44 PM To: R-help Forum Subject: Re: [R] Generate Random Draw from Gamma Distribution Re: Monte Carlo Simulation in R... By the way, maybe the number of groups can be determined endogenously. It will be better if I do not have to set the total number of groups exogenously. Thanks Garry
On Tue, Nov 10, 2009 at 1:29 PM, Hongwei Dong <pdxdong at gmail.com> wrote:
Sorry for the confusion. Let me put it in this way. Here we have 2000 people and we want to put
them
into 150 groups. The distribution of the group size follows the Gamma distribution with shape parameter 0.067 and scale parameter 0.008. At the same time, the minimum group size is 1, and the largest one should not be bigger than 85. My questions: Can I generate a set of groups that follow the above rules by generating random draws? By the way, I also confused by the rate and scale parameters in R. I did the distribution test in SPSS and got those shape and scale parameters. In SPSS Q-Q plot, the scale parameter is 0.008. I noticed that someone mentioned that this is actually the rate. I'm confused. Thanks a lot. Garry On Tue, Nov 10, 2009 at 12:15 PM, Ravi Varadhan <RVaradhan at jhmi.edu>wrote:
I think he means "rate = 0.008", so he is looking for: rgamma(n, shape=0.067, rate=0.008) Even then his problem is not well-posed. You cannot have both "independent" gamma rv's and have them sum to 2000. Ravi.
----------------------------------------------------------------------------
------- Ravi Varadhan, Ph.D. Assistant Professor, The Center on Aging and Health Division of Geriatric Medicine and Gerontology Johns Hopkins University Ph: (410) 502-2619 Fax: (410) 614-9625 Email: rvaradhan at jhmi.edu Webpage:
http://www.jhsph.edu/agingandhealth/People/Faculty_personal_pages/Varadhan.h
tml<http://www.jhsph.edu/agingandhealth/People/Faculty_personal_pages/Varadh an.html>
----------------------------------------------------------------------------
-------- -----Original Message----- From: r-help-bounces at r-project.org [mailto:r-help-bounces at r-project.org] On Behalf Of David Winsemius Sent: Tuesday, November 10, 2009 2:47 PM To: Hongwei Dong Cc: R-help Forum; Duncan Murdoch Subject: Re: [R] Generate Random Draw from Gamma Distribution Re: Monte Carlo Simulation in R... On Nov 10, 2009, at 2:26 PM, Hongwei Dong wrote:
Exactly! Thanks, Duncan. Let me re-phrase me question like this: 1) X_i values are independent Gammas, with the shape 0.067 and scale 0.008 2) Min(X)=1 and Max(X)=85
You might want to check that your parameterization in in agreement with that used by the rgamma function. Simply using those numbers yields a distribution that does not look as though it would get many qualifying samples. Here are 20 draws without any exclusions outside a range:
> rgamma(20, shape=0.067, scale = 0.008)
[1] 2.213459e-03 2.815705e-05 2.381306e-04 2.264602e-07 1.293713e-07 7.680773e-38 6.441082e-15 6.168961e-13 [9] 5.089033e-06 1.571858e-16 9.869878e-12 1.813121e-13 1.253287e-11 1.852885e-04 4.212802e-07 1.774495e-25 [17] 1.892984e-07 5.927422e-17 1.322638e-12 4.327472e-05 http://finzi.psych.upenn.edu/R/Rhelp02/archive/31459.html
3) SUM(X)=2000 4) Do I also have to define the number of draws? if yes, it could be 250. Based on these restrictions, I want to generate random draw. I'm wondering how I can do this in R. Thanks. Garry On Tue, Nov 10, 2009 at 11:17 AM, Duncan Murdoch <murdoch at stats.uwo.ca>wrote:
On 11/10/2009 1:25 PM, Hongwei Dong wrote:
Hi, Dear R users, I'm wondering if I can do Monte Carlo Simulation in R. My problem is like this: I know variable X follows Gamma distribution with shape parameter 0.067 and scale parameter 0.008. The sum of the X is 2000. I need R help me to simulate a vector of X that satisfies both the probability distribution and the sum. Anyone has a clue to this? Much appreciated.
Your requirements are slightly contradictory or incomplete. Here's one way to fully specify the problem: The X_i values are independent Gammas, with the given shape and scale. You want to simulate from the joint distribution conditional on the event sum(X) == 2000. Is that your problem? I don't know how to do the simulation, but maybe someone else does. Duncan Murdoch
[[alternative HTML version deleted]]
______________________________________________ R-help at r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide
http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code. David Winsemius, MD Heritage Laboratories West Hartford, CT ______________________________________________ R-help at r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
______________________________________________ R-help at r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.