```Date: Mon, 14 Feb 2011 15:15:05 -0800 Reply-To: Dale McLerran Sender: "SAS(r) Discussion" From: Dale McLerran Subject: Re: Linear Mixed Model Simulation: Random Intercepts and Slopes In-Reply-To: Content-Type: text/plain; charset=iso-8859-1 Ryan, In order to generate data in which the random effects are correlated, you can take advantage of the Cholesky decomposition of a symmetric positive definite matrix. The Cholesky decomposition has the property that: _ _ _ _ | a11 0 0 ... 0 | | a11 a21 a31 ... ak1 | | a21 a22 0 ... 0 | | 0 a22 a32 ... ak2 | V = | a31 a32 a33 ... 0 | X | 0 0 a33 ... ak3 | | ... ... ... ... ... | | ... ... ... ... ... | | ak1 ak2 ak3 ... akk | | 0 0 0 ... akk | - - - - for V a symmetrix p.d. matrix. For a 2x2 covariance matrix, we have _ _ V = | a11*a11 0 | | a21*a11 a21*a21 + a22*a22 | - - So, if you want the random intercept variance to be V11, the variance of the slope to be V22, and the covariance to be rho*sqrt(V11*V22), then you can construct V11 = ; V22 = ; rho = ; V21 = rho*sqrt(V11*V22); a11 = sqrt(V11); a21 = V21/a11; a22 = sqrt(V22 - a21*a21); If you have two independent r.v.'s x0j and x1j and you construct u0j = a11*x0j; u1j = a21*x0j + a22*x1j; then Cov(u0j,u1j)=V. I don't have time right now to edit code. But you should be able to apply the above to get what you need. Dale --------------------------------------- Dale McLerran Fred Hutchinson Cancer Research Center mailto: dmclerra@nFoHsCpRaCm.org Ph: (206) 667-2926 Fax: (206) 667-5977 --------------------------------------- --- On Sun, 2/13/11, R B wrote: > From: R B > Subject: Re: Linear Mixed Model Simulation: Random Intercepts and Slopes > To: SAS-L@LISTSERV.UGA.EDU > Date: Sunday, February 13, 2011, 5:09 PM > A minor correction to notation: The > random slope, "Gamma01", > technically should be "Gamma10". The code below should now > be > consistent with standard notation. If anyone happens to > find other > mistakes or perhaps a way to allow the random effects to > be > correlated, I'd welcome the feedback. -Ryan > > /* Generate Data */ > data response; > > seed = 987624; /*seed for random # generator*/ > > do subject = 1 to 1000; /*N=1000 subjects*/ > > Gamma00 = 0.50; /*fixed intercept*/ > Gamma10 = 0.30; /*fixed slope*/ > > u0j = sqrt(0.3)*rannor(seed); /*random intercept error > term*/ > u1j = sqrt(0.3)*rannor(seed); /*random slope error term*/ > > B0J = Gamma00 + u0j; /*random intercept*/ > B1J = Gamma10 + u1j; /*random slope*/ > > do time=1 to 5; /*5 time points*/ > > eij = rannor(seed); /*error term*/ > > y = B0J + B1J*time + eij; /*full equation*/ > > output; > end; > end; > run; > > /* Fit Model */ > proc mixed data=response; > model y = time / solution; > random intercept time / subject = subject; > run; > > On Sun, Feb 13, 2011 at 12:34 PM, R B > wrote: > > Hi, > > > > Below I generate data for a linear mixed model with > random intercept > > and slope terms which are considered a function of > fixed intercept, > > fixed slope, and error terms. As I understand it, the > random slope > > term suggests that the linear effect of time varies > across subjects. I > > believe the correlation between the random intercept > and slope terms > > are assumed to be zero in this simulation. That is, we > are assuming > > that as subject-specific intercepts change, the > subject-specific > > slopes do not change in a consistent way. Often, in > real-world data, > > one does observe such a correlation. I'm not sure how > to adjust the > > simulation code to allow for a non-zero correlation > between these > > random terms. > > > > As always, any thoughts on the simulation and MIXED > code would be > > appreciated. I'm particularly interested in finding > out if there are > > any mistakes in the simulation code. > > > > Ryan > > -- > > > > /* Generate Data */ > > data response; > > > > seed = 987624; /*seed for random # generator*/ > > > > do subject = 1 to 1000; /*N=1000 subjects*/ > > > > Gamma00 = 0.50; /*fixed intercept*/ > > Gamma01 = 0.30; /*fixed slope*/ > > > > u0j = sqrt(0.3)*rannor(seed); /*random intercept > error term*/ > > u1j = sqrt(0.3)*rannor(seed); /*random slope > error term*/ > > > > B0J = Gamma00 + u0j; /*random intercept*/ > > B1J = Gamma01 + u1j; /*random slope*/ > > > > do time=1 to 5; /*5 time points*/ > > > > eij = rannor(seed); /*error term*/ > > > > y = B0J + B1J*time + eij; /*full > equation*/ > > > > output; > > end; > > end; > > run; > > > > /* Fit Model */ > > proc mixed data=response; > > model y = time / solution; > > random intercept time / subject = subject; > > run; > > > ```

Back to: Top of message | Previous page | Main SAS-L page