The phytools package contains a small function, multi.mantel, for multiple matrix
regression. The idea of this function is merely to fit a model in which distance or correlation matrix
Y ~ X1 + X2 and so on, but then obtain our p-values for both
the full model & the model coefficients via (Mantel) permutations of the rows & columns together in
Y.
A colleague today contacted me with some questions about the function & so I decided it needed a small
re-boot. None of the internal operation of the function has been changed; however, I have now
updated
it with a new object class, as well as print, residuals, and fitted
S3 methods.
For fun, this is what the print method looks like:
print.multi.mantel<-function(x,...){
if(hasArg(digits)) digits<-list(...)$digits
else digits<-6
star<-function(p){
obj<-if(p>0.1) "" else if(p<=0.1&&p>0.05) "." else
if(p<=0.05&&p>0.01) "*" else if(p<=0.01&&p>0.001) "**" else
if(p<=0.001) "***"
obj
}
cat("\nResults from a (multiple) Mantel regression using \"multi.mantel\":\n\n")
cat("Coefficients:\n")
object<-data.frame(x$coefficients,
x$tstatistic,x$probt,
sapply(x$probt,star))
rownames(object)<-names(x$coefficients)
colnames(object)<-c("Estimate","t value","Pr(>|t|)","")
print(object)
cat("---\n")
cat(paste("Signif. codes: 0 \u2018***\u2019 0.001 \u2018**\u2019 0.01",
"\u2018*\u2019 0.05 \u2018.\u2019 0.1 \u2018 \u2019 1\n"))
cat(paste("Pr(>|t|) based on",x$nperm,
"(Mantel) permutations of rows & columns together in Y.\n\n"))
cat(paste("Multiple R-squared:",round(x$r.squared,digits),"\n"))
cat(paste("F-statistic: ",round(x$fstatistic,digits),
", p-value (based on ",x$nperm," permutations): ",
x$probF,"\n\n",sep=""))
}
Here's a small example. It takes as arguments a single dependent matrix, Y, and either a
single independent matrix or a list of such matrices, X. Either matrices or class
"dist" objects are acceptable. Here I use the latter.
library(phytools)
packageVersion("phytools")
## [1] '0.6.48'
dY
## 1 2 3 4 5 6
## 2 0.54603980
## 3 1.69212862 2.23816843
## 4 0.16256209 0.70860190 1.52956653
## 5 0.58192166 1.12796146 1.11020697 0.41935957
## 6 0.83933577 0.29329596 2.53146439 1.00189786 1.42125743
## 7 0.58679271 0.04075291 2.27892134 0.74935480 1.16871437 0.25254306
## 8 1.75571440 2.30175420 0.06358577 1.59315231 1.17379274 2.59505017
## 9 0.93129804 1.47733784 0.76083058 0.76873595 0.34937638 1.77063381
## 10 0.97782473 1.52386454 0.71430389 0.81526264 0.39590307 1.81716050
## 7 8 9
## 2
## 3
## 4
## 5
## 6
## 7
## 8 2.34250711
## 9 1.51809075 0.82441636
## 10 1.56461744 0.77788967 0.04652669
dX
## [[1]]
## 1 2 3 4 5 6
## 2 1.75416715
## 3 0.63607244 1.11809471
## 4 3.05970185 1.30553469 2.42362940
## 5 0.77190528 2.52607244 1.40797773 3.83160713
## 6 0.58102311 1.17314404 0.05504933 2.47867874 1.35292839
## 7 2.09863486 0.34446771 1.46256242 0.96106698 2.87054015 1.51761176
## 8 1.45771945 0.29644770 0.82164701 1.60198239 2.22962474 0.87669634
## 9 0.46813108 1.28603607 0.16794136 2.59157077 1.24003636 0.11289203
## 10 1.23194515 0.52222200 0.59587271 1.82775670 2.00385043 0.65092204
## 7 8 9
## 2
## 3
## 4
## 5
## 6
## 7
## 8 0.64091541
## 9 1.63050378 0.98958837
## 10 0.86668972 0.22577430 0.76381407
##
## [[2]]
## 1 2 3 4 5 6
## 2 0.54790843
## 3 1.16512317 1.71303160
## 4 0.26507653 0.28283189 1.43019970
## 5 0.95864760 1.50655603 0.20647557 1.22372414
## 6 0.36213942 0.18576901 1.52726259 0.09706289 1.32078703
## 7 0.03566574 0.51224269 1.20078891 0.22941080 0.99431334 0.32647369
## 8 1.29414610 1.84205453 0.12902294 1.55922264 0.33549850 1.65628553
## 9 1.73415785 2.28206628 0.56903468 1.99923439 0.77551025 2.09629728
## 10 0.30273318 0.24517525 1.46785635 0.03765665 1.26138078 0.05940624
## 7 8 9
## 2
## 3
## 4
## 5
## 6
## 7
## 8 1.32981184
## 9 1.76982359 0.44001175
## 10 0.26706744 1.59687928 2.03689103
fit<-multi.mantel(dY,dX)
fit
##
## Results from a (multiple) Mantel regression using "multi.mantel":
##
## Coefficients:
## Estimate t value Pr(>|t|)
## (intercept) 0.9667831 4.196181 0.618
## X1 -0.1771598 -1.579220 0.103
## X2 0.3950495 2.802433 0.023 *
## ---
## Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
## Pr(>|t|) based on 1000 (Mantel) permutations of rows & columns together in Y.
##
## Multiple R-squared: 0.213525
## F-statistic: 5.701406, p-value (based on 1000 permutations): 0.028
residuals(fit)
## 1 2 3 4 5
## 2 -0.326426449
## 3 0.377750562 0.792734343
## 4 -0.366883315 -0.138625616 0.427153292
## 5 -0.626824208 0.013532473 0.311292771 -0.352048579
## 6 -0.167576456 -0.539041195 0.971089391 0.435892244 0.172382466
## 7 -0.022286474 -1.067365644 1.096874319 -0.137794548 0.317672448
## 8 0.535928683 0.659786867 -0.808605010 0.294205815 0.469470892
## 9 -0.637629371 -0.163140682 -0.400996988 -0.528721733 -0.704087161
## 10 0.109698119 0.552741765 -0.726790552 0.157408232 -0.714186278
## 6 7 8 9
## 2
## 3
## 4
## 5
## 6
## 7 -0.574353601
## 8 1.129267512 0.963926842
## 9 -0.004290689 0.140999293 -0.140877962
## 10 0.942226161 0.645871995 -0.779741779 -1.589612195
fitted(fit)
## 1 2 3 4 5 6 7
## 2 0.8724663
## 3 1.3143781 1.4454341
## 4 0.5294454 0.8472275 1.1024132
## 5 1.2087459 1.1144290 0.7989142 0.7714081
## 6 1.0069122 0.8323372 1.5603750 0.5660056 1.2488750
## 7 0.6090792 1.1081186 1.1820470 0.8871494 0.8510419 0.8268967
## 8 1.2197857 1.6419673 0.8721908 1.2989465 0.7043218 1.4657827 1.3785803
## 9 1.5689274 1.6404785 1.1618276 1.2974577 1.0534635 1.7749245 1.3770915
## 10 0.8681266 0.9711228 1.4410944 0.6578544 1.1100894 0.8749343 0.9187454
## 8 9
## 2
## 3
## 4
## 5
## 6
## 7
## 8
## 9 0.9652943
## 10 1.5576314 1.6361389
You get the idea.
I simulated the distance matrices for the preceding as follows:
dY<-dist(y<-rnorm(10))
dX<-list(dist(rnorm(10)),dist(sqrt(0.5)*y+rnorm(10,sd=0.5)))
That's it.
Liam, this is a really useful function. Thanks! Do you have any recommendations for how to use it for variation partitioning if you provide it with a list of independent X matrices?
ReplyDelete