Statalist


[Date Prev][Date Next][Thread Prev][Thread Next][Date index][Thread index]

Re: st: Drawing 3 normally distributed variable with a given correlationstructure to an existing variable


From   Jeph Herrin <[email protected]>
To   [email protected]
Subject   Re: st: Drawing 3 normally distributed variable with a given correlationstructure to an existing variable
Date   Mon, 24 Sep 2007 13:07:26 -0400

Perhaps I'm missing something, but isn't this a
job for -drawnorm- ?

. matrix C= <your corr matrix>
. drawnorm x1 x2 x3 x4, n(100000) corr(C)

hth,
Jeph



Rachel wrote:
I posted a simpler version of this question a few months ago and
received several useful answers.  Basically, I have an existing
variable x1 distributed as N(0,1), and need to generate three
additional variables x2, x3, x4, such that x2,x3, and x4 are
distributed N(0,1) and such that there is a given 4x4 correlation
matrix for the four variables.

 -corr2data- doesn't work for this purpose, so it has to be coded
manually.  One excellent answer was from Paulo, who suggested a way of
doing this for 2 newly drawn variables.  I am trying to expand this to
3 variables but have been unsuccessful.  Could Paulo or someone else
suggest a way to do this?

Here is Paulo's original code:

clear
set obs 100000
local p12=0.5
local p13=0.5
local p23=0.5
local vx1=1
gen double x1=sqrt(`vx1')*invnorm(uniform())
local vu1=`vx1'*(1-`p12'^2)/`p12'^2
gen double x2=x1+sqrt(`vu1')*invnorm(uniform())
local vx2=`vx1'+`vu1'
local covx1x2=sqrt(`vx1')*sqrt(`vx2')*`p12'
local alpha=(sqrt(`vx1')*`p13')/(sqrt(`vx2')*`p23')
local beta=(`covx1x2'-`alpha'*`vx2')/(`alpha'*`covx1x2'-`vx1')
local part1=`beta'^2*`vx1'+`vx2'+2*`covx1x2'
local part2=(`beta'*`vx1'+`covx1x2')^2
local part3=`vx1'*`p13'^2
local vu2=`part2'/`part3'-`part1'
gen double x3=`beta'*x1+x2+sqrt(`vu2')*invnorm(uniform())
corr x1 x2 x3
*
*   For searches and help try:
*   http://www.stata.com/support/faqs/res/findit.html
*   http://www.stata.com/support/statalist/faq
*   http://www.ats.ucla.edu/stat/stata/


*
*   For searches and help try:
*   http://www.stata.com/support/faqs/res/findit.html
*   http://www.stata.com/support/statalist/faq
*   http://www.ats.ucla.edu/stat/stata/



© Copyright 1996–2024 StataCorp LLC   |   Terms of use   |   Privacy   |   Contact us   |   What's new   |   Site index