You have 38 women who have done 17 profiles each. Each profile can be of three different types (a, b, or c), randomly assigned to each woman. The women can either be ovulating or not (binary variable) and X is the mean of all 17 profile ratings for each subject.
You now want to assess the impact of ovulation on score for profiles a and b. You can do this with a mixed model approach and it requires that you restructure your data so that each profile rating is one row. If there are 17*3 profiles it means that each woman will have 51 rows, but 34 of them will have missing data for the score variable because each woman only rated 17 profiles. The other variables, in this case only ovulation, will be the same for all of the 51 rows. This is called converting data from wide form to long form, and you can do it in this way:
library(reshape)
DF.long <- melt(DF, id.vars=c("subject", "ovul", value.name="rating", variable.name="profile")
# id.vars are the variables that you want to keep per subject, the others will be converted.
# rating is the rating score for each trial
# profile is the name of the profile for each trial
Now, you will need to convert all names in the rating
variable, "profile_1a", "profile_2a" etc to "a", "profile_1b", "profile_2b" to b and so on. When you're done, you're ready to analyze. The first method is to analyze a and b separately. You will need to use subject as a random effect to take the repeated measurements into account:
mod_a <- lmer(rating ~ ovul + (1|subject), data=subset(DF_long, profile=="a")
mod_b <- lmer(rating ~ ovul + (1|subject), data=subset(DF_long, profile=="b")
You now have the effect of ovul on a and b for the two models. If you want to see if ovul has a different effect on a than on b, you can specify an interaction model:
mod_a <- lmer(rating ~ profile * ovul + (1|subject), data=subset(DF_long, profile=="a" | profile=="b")
You will now get one value for profile, that tells you the difference between the profiles, and another value for ovul, that tells you the average effect of ovul. The interaction term profile:ovul will tell you the difference between the profles in the effect of ovul on score.
Best Answer
@Stefan is probably right. If you have a single measurement per tick then you should leave the explicit
mouse.id:tick.id
grouping variable out of the model, i.e. uselog10.load ~ treat + (1|mouse.id)
. Alternatively, if this is a nested design (each mouse gets one treatment, each tick is fed on a single mouse, multiple ticks per mice) you could also follow Murtaugh 2007 "Simplicity and complexity in ecological data analysis" Ecology 88 and simply compute the average load per mouse, then use a simplelm(log10.load~treat,data=aggregated_data)
. (If you insist on fitting the model the way it is you can usecontrol=lmerControl(check.nobs.vs.nlev="ignore")
, but I don't advise it.)