[R-sig-ME] within subjects, 3 conditions, 3 lines

[Stanislav Aggerwal]

Thierry, thanks for your help, but your model is completely different from
what I want, so far as I can tell.

The model I want to fit is this. Each subject produces 3 lines, one for
each condition
y = b1 + b2*x + error  #condition a
y = b3 + b4*x + error # condition b
y = b5 + b6*x + error #condition c
Within each subject, these three lines have differing slopes and
intercepts. For each subject, each parameter (b1 to pick an example) is not
identical to that of other subjects but is similar to that from other
subjects.

The intercepts, b1, b3, and b5 vary between subjects. Each subject's b1,
b3, and b5 is randomly sampled from N(p_i,sd_i) where p_i is the population
value for each of the 3 parameters

The intercepts, b2, b4, and b6 vary between subjects. Each subject's b2,
b4, and b6 is randomly sampled from N(p_i,sd_i) where p_i is the population
value for each of the 3 parameters.

In my simulations I used params a1-a6 which are the R parametrisation used
in fitting linear models, instead of b1-b6. E.g. a2 = b1+b2

So far as I can tell (and I looked hard), my code is producing the right
sort of simulated data. And these simulated data look similar to my real
data.

If we pretend that subjects is a fixed effect, we have
lm(y ~ x*cond*subj)
A different slope and intercept for each condition and each subject.

Thanks again,
Bill



On Sat, Jan 24, 2015 at 10:35 AM, Stanislav Aggerwal <
stan.aggerwal at gmail.com> wrote:

> Thanks Thierry!
>
> Stan
>
> On Thu, Jan 22, 2015 at 4:16 PM, ONKELINX, Thierry <
> Thierry.ONKELINX at inbo.be> wrote:
>
>> Dear Stan,
>>
>> Please keep the mailing in cc.  You use rnorm(nsubj,mean=0,sd=.1) as the
>> random intercept. However rnorm() will be evaluate again in a1, a2, ..., a6
>> and thus yielding different values of the random intercept in a1, a2, ...,
>> a6. You want
>>
>> Rf <- rnorm(nsubj,mean=0,sd=.1)
>> a1<- rep(1+ rf,each=5)  #inta=1
>> a2<- rep(1 + rf,each=5)  #intb-inta=2-1=1
>>
>> Note that is more clear to define the fixed, random effect and noise
>> separately and add those components in a separate step.
>>
>> Best regards,
>>
>> Thierry
>>
>> 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
>> + 32 2 525 02 51
>> + 32 54 43 61 85
>> Thierry.Onkelinx at inbo.be
>> www.inbo.be
>>
>> 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
>>
>> Van: Stanislav Aggerwal [mailto:stan.aggerwal at gmail.com]
>> Verzonden: donderdag 22 januari 2015 12:25
>> Aan: ONKELINX, Thierry
>> Onderwerp: Re: [R-sig-ME] within subjects, 3 conditions, 3 lines
>>
>> Thanks Thierry for your help. Could I get some clarification please?
>>
>>
>> You get the right answer for the wrong question. Because you have a bug
>> in the code that simulates the data. You added a second noise term instead
>> of a random intercept. Hint: set all parameters not related to the random
>> intercept to 0 in a1 to a6. Then a1 to a6 should be equal.
>>
>> There are 30*3 (subjects * conditions) random intercepts. Also 30 random
>> slopes. There are 6 params because there are 3 conditions: (slope and
>> intercept)*3
>>
>> a1-a6 have the population values: 1, 1, 2, 3, -1, -2
>> Then each subject has a1-a6 values that randomly vary about these
>> population values.
>> It is the model I meant, and it describes my data.
>> Are you saying I implemented it wrong? I don't understand how.
>>
>> a1<- rep(1+rnorm(nsubj,mean=0,sd=.1),each=5)  #inta=1
>> a2<- rep(1+rnorm(nsubj,mean=0,sd=.1),each=5)  #intb-inta=2-1=1
>> a3<- rep(2+rnorm(nsubj,mean=0,sd=.1),each=5)  #intc-inta=3-1=2
>> a4<- rep(3+rnorm(nsubj,mean=0,sd=.1),each=5)  #slopea=3
>> a5<- rep(-1+rnorm(nsubj,mean=0,sd=.1),each=5) #slopeb-slopea
>> a6<- rep(-2+rnorm(nsubj,mean=0,sd=.1),each=5) #slopec-slopea
>>
>> y[cond=="a"]<-a1    + a4*x[cond=="a"]
>> y[cond=="b"]<-a1+a2 + (a4+a5)*x[cond=="b"]
>> y[cond=="c"]<-a1+a3 + (a4+a6)*x[cond=="c"]
>>
>> Then I add errors. For each subject the data are not perfectly described
>> by a line, due to measurement error. If I did not add this, the only source
>> of noise in the model is due to the variation in slopes and intercepts
>> between subjects.
>> Are you saying I should eliminate these measurement errors?
>> This is the model I intend:
>> For each subject in each condition:
>> y = int + slope*x + error
>> each subject has 3 lines (for a, b, c) with params a1-a6
>>
>> #autocorrelated errors
>> for(i in 1:nsubj)
>>   {
>>   y[subj==i]<-y[subj==i] +
>> as.numeric(filter(rnorm(5*3,mean=0,sd=.5),filter=0.5,method="recursive"))
>>   }
>>
>> If you look at the plot I think it looks correct.
>> Thanks very much!
>> Cheers
>> Stan
>>
>> Disclaimer<https://www.inbo.be/nl/disclaimer-mailberichten-van-het-inbo>
>>
>
>

	[[alternative HTML version deleted]]