Simplified Monte-Carlo simulations (STB-20: ssi6.1) ---------------------------------- ^simul^ progname^, r^eps^(^#^)^ [ ^a^rgs^(^whatever^) d^ots ] Description ----------- ^simul^ eases the programming task of performing Monte-Carlo type simulations. progname is the name of a program that performs a single simulation. Typing "^simul^ progname^, reps(^#^)^" iterates progname for # replications and collects the results. ^simul^ calls progname two ways. At the outset, ^simul^ issues "progname ^?^" and expects progname to set the global macro $S_1 to contain a list of variable names under which results are to be stored. Thereafter, ^simul^ issues straight "progname" calls and expects progname to perform a single simulation and to store the results using ^post^. Details of ^post^ can be found in ^help post^, but enough information is provided below to use ^post^ successfully. Options ------- ^reps(^#^)^ is not optional -- it specifies the number of replications to be performed. ^args(^whatever^)^ specifies any arguments to be passed to progname on invocation. The query call is then of the form "progname ^?^ whatever" and subsequent calls of the form "progname whatever". ^dots^ requests a dot be placed on the screen at the beginning of every call to progname, thus providing entertainment during a long simulation. Remarks ------- progname must have the following outline: ^program define^ progname ^if "`1'"=="?" {^ ^global S_1 "^variable names^"^ ^exit^ ^}^ perform single simulation ^post^ results ^end^ There must be the same number of results following the ^post^ command as variable names following the ^global S_1^ command. Example ------- Make a data set containing means and variances of 100-observation samples from a log-normal distribution. Perform the experiment 10,000 times: ^program define lnsim^ ^version 3.1^ ^if "`1'"=="?" {^ ^global S_1 "mean var"^ ^exit^ ^}^ ^drop _all^ ^set obs 100^ ^gen z = exp(invnorm(uniform()))^ ^summarize z^ ^post _result(3) _result(4)^ ^end^ It is instructive to compare this program to the same example in ^help post^. In any case, ^lnsim^ can then be executed 10,000 times by typing: . ^simul lnsim, reps(10000)^ Technical note: debugging a simulation --------------------------------------- Before executing our ^lnsim^ simulator, we can verify it works by typing: . ^lnsim ?^ . ^display "$S_1"^ This verifies that the program sets the global macro $S_1 correctly on a query call. We can then try executing ^lnsim^: . ^postfile $S_1 using myfile^ . ^lnsim^ The ^postfile^ command opens a file for the posted results (see ^help post^). Invoking ^lnsim^ should perform a single simulation. We could then close the result file and examine it: . ^postclos^ . ^use myfile, clear^ Passing arguments to the simulation ----------------------------------- Consider a more complicated problem: Let's experiment with estimating Y = a + b*X + u when the true model has a=1, b=2, and u=2*(z+c*x) and z is N(0,1). We will keep the parameter estimates and standard errors and experiment with varying c. X will be fixed across the experiments but will originally be generated as N(0,1). We begin by interactively making the true data: . ^drop _all^ . ^set obs 100^ . ^gen x = invnorm(uniform())^ . ^gen true_y = 1 + 2*x^ . ^save truth^ Our program is: Passing arguments to the simulation, continued ---------------------------------------------- ^program define hetero^ ^version 3.1^ ^if "`1'"=="?" {^ ^global S_1 "a se_a b se_b"^ ^exit^ ^}^ ^use truth, clear ^ ^gen y = true_y + 2*(invnorm(uniform()) + `1'*x)^ ^regress y x^ ^post _b[_cons] _se[_cons] _b[x] _se[x]^ ^end^ Note the use of ^`1'^ in our statement for generating y. ^`1'^ is an argument. If the argument's value is 2, then the last part of the statement is equivalent to "^2*x^". We can run 10,000 simulations, setting the argument to 2 by typing: . ^simul hetero, args(2) reps(10000)^ Passing arguments to the simulation, continued ---------------------------------------------- We can set the argument to 2 by typing: . ^simul hetero, args(2) reps(10000)^ We might then analyze that data and try a different experiment, this time setting the argument to 1.5: . ^simul hetero, args(1.5) reps(10000)^ A more efficient implementation ------------------------------- Our simulation program is: ^program define hetero^ ^version 3.1^ ^if "`1'"=="?" {^ ^global S_1 "a se_a b se_b"^ ^exit^ ^}^ ^use truth, clear ^ ^gen y = true_y + 2*(invnorm(uniform()) + `1'*x)^ ^regress y x^ ^post _b[_cons] _se[_cons] _b[x] _se[x]^ ^end^ Actually, we can do better because the query call can also be used to initial- ize ourselves. In particular we could load the file truth once at the outset rather than reloading it every replication. The faster version reads: A more efficient implementation, continued ------------------------------------------ ^program define hetero^ ^version 3.1^ ^if "`1'"=="?" {^ ^use truth, clear^ (load the data just once) ^global S_1 "a se_a b se_b"^ ^exit^ ^}^ ^gen y = true_y + 2*(invnorm(uniform()) + `1'*x)^ ^regress y x^ ^post _b[_cons] _se[_cons] _b[x] _se[x]^ ^drop y^ (because we will recreate it next time) ^end^ Assume we plan on replicating the experiment 10,000 times. In our original draft, we used the truth data once per replication, meaning we performed 10,000 uses. The new version does this only once, saving 9,999 unnecessary uses. Also see -------- STB: ssi6.1 (STB-20) On-line: ^help^ for ^post^