Dear Thierry,
Thank you very much for the demonstration. Three questions:
(1) Are two different sets of boxplots needed?
(2) From these boxplots, how can we exactly determine the the failure of "
m2" in capturing the switch of a school from control to treatment (what
are the criteria)?
(3) Suppose "m2" failed, what is generally an alternative to account for
the switch of a school from control to treatment?
Thank you, Simon
On Tue, May 19, 2020 at 2:27 AM Thierry Onkelinx <thierry.onkelinx at inbo.be>
wrote:
Dear Simon,
Something like this
library(lme4)
dat <- read.csv('https://raw.githubusercontent.com/hkil/m/master/z.csv')
dat$year <- factor(dat$year)
m2 <- lmer(y~ year*group + (1|scid/stid), data =dat)
dat$resid <- resid(m2)
library(ggplot2)
ggplot(dat, aes(x = year, y = resid, colour = group)) +
geom_boxplot() +
facet_wrap(~scid)
ggplot(dat, aes(x = scid, y = resid, colour = year, fill = group)) +
geom_boxplot()
Best regards,
ir. Thierry Onkelinx
Statisticus / Statistician
Vlaamse Overheid / Government of Flanders
INSTITUUT VOOR NATUUR- EN BOSONDERZOEK / RESEARCH INSTITUTE FOR NATURE
AND FOREST
Team Biometrie & Kwaliteitszorg / Team Biometrics & Quality Assurance
thierry.onkelinx at inbo.be
Havenlaan 88 bus 73, 1000 Brussel
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
///////////////////////////////////////////////////////////////////////////////////////////
<https://www.inbo.be>
Op ma 18 mei 2020 om 19:35 schreef Simon Harmel <sim.harmel at gmail.com>:
Dear Thierry,
Would you possibly demonstrate what you exactly mean perhaps in R?
Thank you, Simon
On Mon, May 18, 2020, 12:13 PM Thierry Onkelinx <
thierry.onkelinx at inbo.be> wrote:
Dear Simon,
You are wondering if the model captures the change in treatment for a
school. Therefore you need to plot the residuals vs every combination of
school and year. A boxplot for every combination would be useful. If the
change in treatment triggers a shift in residuals, then the current model
fails.
Best regards,
ir. Thierry Onkelinx
Statisticus / Statistician
Vlaamse Overheid / Government of Flanders
INSTITUUT VOOR NATUUR- EN BOSONDERZOEK / RESEARCH INSTITUTE FOR NATURE
AND FOREST
Team Biometrie & Kwaliteitszorg / Team Biometrics & Quality Assurance
thierry.onkelinx at inbo.be
Havenlaan 88 bus 73, 1000 Brussel
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
///////////////////////////////////////////////////////////////////////////////////////////
<https://www.inbo.be>
Op ma 18 mei 2020 om 17:45 schreef Simon Harmel <sim.harmel at gmail.com>:
Dear Thierry,
By "Have a look at the residuals" you mean something like the
following (below)? So no other adjustment is required for the switching
that occurred?
plot(m1, type = c("p","smooth"), col.line = 2)
plot(m1, sqrt(abs(resid(.)))~fitted(.), type = c("p","smooth"),
col.line = 2)
On Mon, May 18, 2020 at 2:01 AM Thierry Onkelinx <
thierry.onkelinx at inbo.be> wrote:
Dear Simon,
The question is rather if the model is able to capture this change.
Have a look at the residuals. If they look OK, then the model handles the
change in treatment.
Best regards,
ir. Thierry Onkelinx
Statisticus / Statistician
Vlaamse Overheid / Government of Flanders
INSTITUUT VOOR NATUUR- EN BOSONDERZOEK / RESEARCH INSTITUTE FOR
NATURE AND FOREST
Team Biometrie & Kwaliteitszorg / Team Biometrics & Quality Assurance
thierry.onkelinx at inbo.be
Havenlaan 88 bus 73, 1000 Brussel
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
///////////////////////////////////////////////////////////////////////////////////////////
<https://www.inbo.be>
Op zo 17 mei 2020 om 01:09 schreef Simon Harmel <sim.harmel at gmail.com
Hello All,
I have a 3-year longitudinal dataset (*see link below the table*).
Up to
year 2 (coded "1"), 8 schools (4 in Treatment, 4 in Control)
cooperated
with the study. But in year 3 (coded "2"), one of the Treatment
schools
(named "good") dropped out.
Also in year 3 (coded "2"), we were made to move one of the *Control
*schools
(named "*orange*") to the *Treatment *group. The full design of the
study
is shown in the Table below.
I want to regress "year" and "group" on "y" (a continuous response)
in lme4
package in R. But is there a way to capture the switch of one of the
control schools to the treatment group?
Thank you very much, Simon
? *Switched from control to treatment*
? *Out as of year coded 2*
*SCHOOL NAMES*
*Year*
*Codes*
*Control*
*Treatment*
0
har
john
bus
orange
caro
good
bla
carm
1
har
john
bus
*orange*
caro
good
bla
carm
2
har
john
bus
X
caro
*orange*
bla
carm
*library(lme4)*
*dat <- read.csv('
https://raw.githubusercontent.com/hkil/m/master/z.csv
<https://raw.githubusercontent.com/hkil/m/master/z.csv>')*
*m1 <- lmer(y~ year*group + (1|stid), data = dat) #### 'stid' =
student id m2 <- lmer(y~ year*group + (1|scid/stid),
data =
dat) #### 'scid' = school id*
[[alternative HTML version deleted]]