[R-sig-ME] Fw: Fw: simulating a linear random intercept model with exogeneity assumption

mervat mohamed mervat_moh2006 at yahoo.com
Thu Jul 13 13:24:30 CEST 2017


dear ir. Thierry Onkelinx
thank you very much  ir. Thierry Onkelinx i appreciate your help 
mervat



On Wednesday, July 12, 2017, 4:02:03 PM GMT+2, Thierry Onkelinx <thierry.onkelinx at inbo.be> wrote:

Dear Mervat,
Here is how I would write it. No need for the multivariate rnorm since you set all covariances are zero, making them independent. All rnorm() are by default independent so need for a special command to make them independent. 
library(dplyr)grp <- 10nindiv <- 20b0 <- 2b <- 0.2u <- rnorm(grp, mean = 0, sd = 0.05)e <- rnorm(grp * nindiv, mean = 0, sd = 1)expand.grid(  group = seq_len(grp),  individual = seq_len(nindiv)) %>%  mutate(    X = rnorm(n(), mean = 0, sd = 1),    id = interaction(group, individual),    mu = b0 + b * X + u[group],    y = mu + e  )
Best regards,
ir. Thierry Onkelinx
Instituut voor natuur- en bosonderzoek / Research Institute for Nature and Forest 
team Biometrie & Kwaliteitszorg / team Biometrics & Quality Assurance 
Kliniekstraat 25
1070 Anderlecht
Belgium

To call in the statistician after the experiment is done may be no more than asking him to perform a post-mortem examination: he may be able to say what the experiment died of. ~ Sir Ronald Aylmer Fisher
The plural of anecdote is not data. ~ Roger Brinner 
The combination of some data and an aching desire for an answer does not ensure that a reasonable answer can be extracted from a given body of data. ~ John Tukey
2017-07-12 10:51 GMT+02:00 mervat mohamed <mervat_moh2006 at yahoo.com>:

dear ir. Thierry Onkelinx
thank you for your advise , my program code is as follows:
## The code of the program

# yj[i] = b0j+ b1*xj[i]+ej[i]
# b0j = b0 + u0j, u0j ~ N(0,1)              random intercept
# b1                                                    fixed slope
# -- load the used package
library(lme4)

set.seed(1)

    n=200000                                ##  the total size of the simulated population
    grp=100                                   ##  the number of groups of the simulated population
    nindiv=2000                             ##  the number of observations in each group
    b0=2                                       ##  intercept
    b1=0.2                                    ##  slope of X
    u = rnorm(grp,0,0.05)
    e = rnorm(nindiv*grp,1)
    M.x=rnorm(grp,3,0.05)
    sigma<-diag(grp)                   ## identity matrix
    require(MASS)
    x=mvrnorm(nindiv,M.x, sigma)
    X=as.vector(x)
    data.set = data.frame(group.id = rep(1:grp, each=nindiv),
                  u = rep(u, each=nindiv),
                  M.x = rep(M.x, each=nindiv),
                  indiv.id = 1:(nindiv*grp),
                  X = X,
                  e = e)
    data.set$Y = b0 + b1*data.set$X + b2*data.set$mean.X+ data.set$u + data.set$e

    indiv<- data.frame(data.set$group.id, data.set$X,data.set$Y)          names(indiv)<-c("group.id","X" ,"Y")
my problem is i do not know how can i write the command of the independence between X and the error terms of the two levels of the model (e and u ) in my program code.
and i appreciate your help
thanks mervat  

On Wednesday, July 12, 2017, 9:46:10 AM GMT+2, Thierry Onkelinx <thierry.onkelinx at inbo.be> wrote:

Dear Mervat,
This the r-sig-mixed-models mailing list. Not the write-code-for-me mailing list. Please show us what you have tried. This can easily be solved using expand.grid() and rnorm().
Best regards,
ir. Thierry Onkelinx
Instituut voor natuur- en bosonderzoek / Research Institute for Nature and Forest 
team Biometrie & Kwaliteitszorg / team Biometrics & Quality Assurance 
Kliniekstraat 25
1070 Anderlecht
Belgium

To call in the statistician after the experiment is done may be no more than asking him to perform a post-mortem examination: he may be able to say what the experiment died of. ~ Sir Ronald Aylmer Fisher
The plural of anecdote is not data. ~ Roger Brinner 
The combination of some data and an aching desire for an answer does not ensure that a reasonable answer can be extracted from a given body of data. ~ John Tukey
2017-07-11 13:28 GMT+02:00 mervat mohamed via R-sig-mixed-models <r-sig-mixed-models at r-project. org>:




----- Forwarded Message -----From: mervat mohamed <mervat_moh2006 at yahoo.com>To: mailman-owner at r-project.org <mailman-owner at r-project.org> Sent: Sunday, July 2, 2017, 5:39:42 AM GMT+2Subject: Fw: [R-sig-ME] Fw: simulating a linear random intercept model with exogeneity assumption


     On Wednesday, June 21, 2017 8:27 PM, mervat mohamed via R-sig-mixed-models <r-sig-mixed-models at r-project. org> wrote:




  Show original message    On Monday, June 19, 2017 1:36 PM, mervat mohamed via R-sig-mixed-models <r-sig-mixed-models at r-project. org> wrote:


 Hi r-sig group
I want to know how to write the code of simulating a linear random intercept model with following assumptions using R programing:
the model:
yij=b.+b1 xij + u.j + eij,  where j refer to the group number and i refer to the observation number in the j group
model assumptions:
  1- xij ~ N (3,1.5)  2- u.j ~ N (0,1)  3 - eij ~ N (0,1)  4 - cov (xij , u.j)=0  5 - cov (xij , eij)=0
  6- cov (u.j , eij)= 0
thnx for help
mervat
    [[alternative HTML version deleted]]

______________________________ _________________
R-sig-mixed-models at r-project. org mailing list
https://stat.ethz.ch/mailman/ listinfo/r-sig-mixed-models

  ______________________________ _________________
R-sig-mixed-models at r-project. org mailing list
https://stat.ethz.ch/mailman/ listinfo/r-sig-mixed-models

   
______________________________ _________________
R-sig-mixed-models at r-project. org mailing list
https://stat.ethz.ch/mailman/ listinfo/r-sig-mixed-models





	[[alternative HTML version deleted]]



More information about the R-sig-mixed-models mailing list