19
$\begingroup$

The following question is one of those holy grails for me for some time now, I hope someone might be able to offer a good advice.

I wish to perform a non-parametric repeated measures multiway anova using R.

I have been doing some online searching and reading for some time, and so far was able to find solutions for only some of the cases: friedman test for one way nonparametric repeated measures anova, ordinal regression with {car} Anova function for multi way nonparametric anova, and so on. The partial solutions is NOT what I am looking for in this question thread. I have summarized my findings so far in a post I published some time ago (titled: Repeated measures ANOVA with R (functions and tutorials), in case it would help anyone)

If what I read online is true, this task might be achieved using a mixed Ordinal Regression model (a.k.a: Proportional Odds Model).

I found two packages that seems relevant, but couldn't find any vignette on the subject:

So being new to the subject matter, I was hoping for some directions from people here.

Are there any tutorials/suggested-reading on the subject? Even better, can someone suggest a simple example code for how to run and analyse this in R (e.g: "non-parametric repeated measures multiway anova") ?

Improve this question
$\endgroup$
3
  • $\begingroup$ Tal, may I ask whether you found a solution for this? I am having the same problem, and the replies below might be helpful in assisting to find an answer, but don't really provide a definitive answer. I have 9 ordinal DVs and 2 time points, and am looking for the same test you wanted to perform. $\endgroup$
    – Torvon
    Commented Jul 5, 2013 at 23:36
  • 1
    $\begingroup$ Hi Torvon. I have never came by a solution. I guess permutation tests will be the safest route, but I never got the time to sit and make it work. If you do - please come again to post your answer. Cheers, T $\endgroup$
    – Tal Galili
    Commented Jul 6, 2013 at 18:53
  • $\begingroup$ Thank you for your quick response. I'll have to work this out and will let you know. $\endgroup$
    – Torvon
    Commented Jul 6, 2013 at 23:46

3 Answers 3

Reset to default
8
$\begingroup$

The ez package, of which I am the author, has a function called ezPerm() which computes a permutation test, but probably doesn't do interactions properly (the documentation admits as much). The latest version has a function called ezBoot(), which lets you do bootstrap resampling that takes into account repeated measures (by resampling subjects, then within subjects), either using traditional cell means as the prediction statistic or using mixed effects modelling to make predictions for each cell in the design. I'm still not sure how "non-parametric" the bootstrap CIs from mixed effects model predictions are; my intuition is that they might reasonably be considered non-parametric, but my confidence in this area is low given that I'm still learning about mixed effects models.

Improve this answer
$\endgroup$
3
  • $\begingroup$ Hello Mike. Thank you for the answer, and for your package - it is really great! $\endgroup$
    – Tal Galili
    Commented Sep 3, 2010 at 11:05
  • $\begingroup$ @Mike, Your package seems to be the only working for mixed multiple factor designs. The aovp alternative - from the orphaned lmperm package - produces huge variations for the p-values, see this. I have a few questions: Where can I found a bibliographical references for the implementation of ezPerm? How can I interpret that the function probably doesn't do interactions properly? What could be a post-hoc test in this case? Thanks! $\endgroup$
    – toto_tico
    Commented Nov 11, 2015 at 1:10
  • $\begingroup$ @Mike, would something like ezPerm( data = DATA, dv = DV, wid = WID, within = interaction(A,B), perms = 1e3) make sense to double check if the interaction is significative? $\endgroup$
    – toto_tico
    Commented Nov 11, 2015 at 1:14
6
$\begingroup$

When in doubt, bootstrap! Really, I don't know of a canned procedure to handle such a scenario.

Bootstrapping is a generally applicable way of generating some error parameters from the data at hand. Rather than relying on the typical parametric assumptions, bootstrap procedures capitalize on the characteristics of the sample to generate an empirical distribution against which your sample estimates can be compared.

Google scholar is gold...it's been done before...at least once.

Lunneborg, Clifford E.; Tousignant, James P.; 1985 "Efron's Bootstrap with Application to the Repeated Measures Design." Multivariate Behavioral Research; Apr85, Vol. 20 Issue 2, p161, 18p

Improve this answer
$\endgroup$
3
  • 1
    $\begingroup$ Thank you for the lead Brett! I wonder if someone got to implement it by now in R (I'd guess not). $\endgroup$
    – Tal Galili
    Commented Aug 5, 2010 at 1:43
  • 1
    $\begingroup$ Right. R has lots of routines for assisting with bootstrap and other randomization methods, but I don't know that you'll find anything specific to this problem. $\endgroup$
    – Brett
    Commented Aug 5, 2010 at 3:15
  • $\begingroup$ Very nice first sentence. I hope it is not copyrighted, because I plan to use it :D $\endgroup$
    – gui11aume
    Commented Jun 28, 2012 at 15:55
1
$\begingroup$

There is a "trick" mentioned in some forums and mailing lists - I also found it mentioned in Joop Hox' book "Multilevel Analysis" (second edition, 2010), pp. 189.

The idea is: you reformat your long data into a long long dataset in which you create a new DV that includes all your DV responses, and use an index variable that holds information about the nature of the DVs to predict this outcome.

Let's assume you have 9 depression symptoms (ordinal), 2 measurement points, and 300 subjects. So while you have 300 rows in your normal dataset, and in your long dataset you'll have 600 rows, this new dataset will have 9 (symptoms) x 2 (time) x 300 (subjects) rows.

The new DV variable "symptoms" now contains the symptom severity of participants on 9 symptoms, the variables "index" contains the information about the nature of the symptom (1 to 9), and then there are the two variables "time" and "UserID".

You can now use the ordinal package to run this. 

data<-read.csv("data_long_long.csv", head=T)

data$symptoms <- factor(data$symptoms)
data$time <- factor(data$time)
data$index <-factor(data$index)

m1<-clmm2(symptoms ~ index+time, random=UserID, data = data, Hess=TRUE, nAGQ=10)

In my specific case, I was interested whether there was a significant interaction between index and time, so I ran one additional model and compared them:

m2<-clmm2(symptoms ~ index+time, random=UserID, data = data, Hess=TRUE, nAGQ=10)
anova(m1,m2)

CLMM2 uses a random intercept model (to the best of my knowledge, the package ordinal does not do random slopes), if you do not with for a random intercept model you can run the models instead using CLM, e.g.: 

m3<-clm(symptoms ~ index+time, data = data)
Improve this answer
$\endgroup$
1
  • 1
    $\begingroup$ aren't m1 and m2 exactly the same, maybe you meant clmm2(symptoms ~ index*time, random=UserID, data = data, Hess=TRUE, nAGQ=10) anova(m1,m2). Do you know by any chance how to express this with the new clmm function? I am not sure about the notation. $\endgroup$
    – toto_tico
    Commented Nov 15, 2015 at 20:41

Not the answer you're looking for? Browse other questions tagged or ask your own question.