st: RE: random number generator for gamma

[jhilbe@aol.com]

I get the Stata Digest, which comes a day later. Unfortunately I have 
been so busy lately that
I rarely get to actually review the queries and responses. I did check 
this morning, however,
so have a chance to explain a bit about the rndgamx.

The pseudo-random number generator you are using is part of a suite of 
pseudo-random
number generators that Walter Linde-Zwirble and I wrote about 15 years 
ago. These were
used by many members of the Stata community until Stata provided its 
new (with version
11) group of pseudo-random number generators last year. Let's just call 
them random number
generators for now. I suggest using Stata's official random number 
generators,
unless you want a short way to create synthetic data sets for members 
of the generalized
linear models family. That's what the random number generators that 
have a final x in the
name are designed to do -- such as rndgamx.

In short, the -rndgam- command is to generate gamma random numbers with 
a specficied mean
and scale.  -rndgamx- is for constructing synthetic GLM gamma models; 
eg canonical inverse link
or log linked gamma models.

Although I suggest using --rgamma-- instead of --rndgam--,  from what I 
have observed the
two generators yield similar results. Know that our random number 
generators are based on
the rejection method; Stata's generators are not. If you don't have 
version 11, however,
I think you need to stick with rndgam.

Let me explain the point of rndgamx, and all of our generators that end 
in x. Again, they were
designed to use for constructing synthetic models. The mean value is 
NOT a constant, as it is for
-rndgam-, but rather a variable having a variety of values, based on 
predictor values.

Let's use randgamx to construct a synthetic log-gamma model with 
specified coefficients as:
intercept=1, coefficient of x1=.5, and coefficient of x2 = -1.75. The 
variables x1 and x2
will be uniformally distributed.

The code:

. set obs 100000
obs was 0, now 100000

. gen x1 = runiform()
. gen x2 = runiform()
. gen xb = 2 + .5*x1 - 1.75*x2
. gen mu=exp(xb)
. rndgamx mu, s(1)    #  GLM gamma models have scale=1
( Generating . )
Variable xg created.

. glm xg x1 x2, fam(gam) link(log) nolog

Generalized linear models                          No. of obs      =    
100000
Optimization     : ML                              Residual df     =    
99997
                                                   Scale parameter =   
.995001
Deviance         =  114525.1513                    (1/df) Deviance =  
1.145286
Pearson          =  99497.11697                    (1/df) Pearson  =   
.995001     <=  point estimate of scale

Variance function: V(u) = u^2                      [Gamma]
Link function    : g(u) = ln(u)                    [Log]

                                                   AIC             =   
4.74273
Log likelihood   = -237133.5085                    BIC             =  
-1036733

-------------------------------------------------------------------------
-----
            |                 OIM
          xg |      Coef.   Std. Err.      z    P>|z|     [95% Conf. 
Interval]
-------------+-----------------------------------------------------------
-----
          x1 |   .5002173   .0109491    45.69   0.000     .4787575    
.5216771
          x2 |  -1.758486   .0109075  -161.22   0.000    -1.779864   
-1.737108
       _cons |   2.000245   .0083384   239.88   0.000     1.983903    
2.016588
-------------------------------------------------------------------------
-----

See how well the log-gamma response, given x1 and x2, can be fit using 
Stata's
glm command with a gamma family and log link.

I have not used or really thought of  -- at least for 10 years -- using 
-rndgamx-
in manner other than for construcing synthetic GLM models. It was not 
designed for
any other purpose.  I suggest using either -rndgam- or Stata's -rgamma- 
commands for
generating gamma data with specified mean and scale parameters. .

Joseph Hilbe

=========================================================================

Date: Wed, 29 Sep 2010 16:26:50 +0100
From: Nick Cox <n.j.cox@durham.ac.uk>
Subject: st: RE: random number generator for gamma

Joe Hilbe will no doubt answer this. Meanwhile official Stata has 
-rgamma()-.

Nick
n.j.cox@durham.ac.uk

Tomas Lind

I´ve been using the -rnd- ado-package from Hilbe/Linde-Zwirble (update 
from STB-28: sg44).

The following syntax produces sensible results (mean and sd) when I 
generate
20000 Gamma distributed random numbers with shape=11.11 and scale=0.18

clear *
rndgam 20000 11.11  0.18
summarize xg    // mean=2.0 sd=0.60

However, with -rndgamx- I don´t get the same values. I have probably 
missed
something fundamental. Someone know what goes on?

clear *
set obs 20000
gen m=2     // mean=2
rndgamx  m , s(11.11)
summarize xg    // mean=247 sd=74






*
*   For searches and help try:
*   http://www.stata.com/help.cgi?search
*   http://www.stata.com/support/statalist/faq
*   http://www.ats.ucla.edu/stat/stata/