Content Description

Characterizing patterns of integration and modularity are hugely important in morphometrics. At their core, these methods seek to summarize patterns of trait covariation, as found in the trait covariance matrix: \(\hat{\mathbf\Sigma} = \mathbf{Z^T}\mathbf{Z}/ (n-1)\). There are three general classes of methods:

`Geomorph implements approaches addressing all three of these (See lecture for detailed description of methods).

library(geomorph)
## Loading required package: RRPP
## Loading required package: rgl
## Loading required package: Matrix

1: Overall Integration

Here we use an effect size based on \(V_{rel}\) to quantify and compare overall levels of integration in a set of traits.

data("plethodon")
Y.gpa <- gpagen(plethodon$land, print.progress = FALSE)
#Separate data by species
coords.gp <- coords.subset(Y.gpa$coords, plethodon$species)

#Z_Vrel by species
Vrel.gp <- Map(function(x) integration.Vrel(x), coords.gp) 
compare.ZVrel(Vrel.gp$Jord, Vrel.gp$Teyah)
## 
## Effect sizes
## 
##  Vrel.gp$Jord Vrel.gp$Teyah 
##    -0.2978931    -0.2642648 
## 
## Effect sizes for pairwise differences in rel.eig effect size
## 
##               Vrel.gp$Jord Vrel.gp$Teyah
## Vrel.gp$Jord    0.00000000    0.09804249
## Vrel.gp$Teyah   0.09804249    0.00000000
## 
## P-values
## 
##               Vrel.gp$Jord Vrel.gp$Teyah
## Vrel.gp$Jord     1.0000000     0.9218986
## Vrel.gp$Teyah    0.9218986     1.0000000

1B: Overall Integration Across Spatial Scales

A method for evaluating integration across spatial scales is also implemented

globalIntegration(Y.gpa$coords) #data are not spatially integrated

##      BEval 
## -0.6369198

2: Integration between Subsets of Traits

Here, one may use PLS to evaluate the covariation between subsets of traits:

data(pupfish) # GPA previously performed
group <- factor(paste(pupfish$Pop, pupfish$Sex, sep = "."))

# Subset 3D array by group, returning a list of 3D arrays
tail.LM <- c(1:3, 5:9, 18:38)
head.LM <- (1:56)[-tail.LM]
tail.coords <- pupfish$coords[tail.LM,,]
head.coords <- pupfish$coords[head.LM,,]

IT <- integration.test(tail.coords, head.coords, print.progress = F)
## 
## No names for A2.  Data are assumed to be ordered as in A.
summary(IT)
## 
## Call:
## integration.test(A = tail.coords, A2 = head.coords, print.progress = F) 
## 
## 
## 
## r-PLS: 0.917
## 
## Effect Size (Z): 5.4983
## 
## P-value: 0.001
## 
## Based on 1000 random permutations
plot(IT)

One could also use a partition of landmarks (NOTE identical \(r_{PLS}\):

land.gp <- rep(1,56); land.gp[tail.LM] <- 2
integration.test(pupfish$coords, partition.gp=land.gp, print.progress = FALSE)
## 
## Call:
## integration.test(A = pupfish$coords, partition.gp = land.gp,  
##     print.progress = FALSE) 
## 
## 
## 
## r-PLS: 0.917
## 
## Effect Size (Z): 5.4019
## 
## P-value: 0.001
## 
## Based on 1000 random permutations

And the statistical results are the same as when using two.b.pls

two.b.pls(tail.coords, head.coords, print.progress = FALSE)
## Data in either A1 or A2 do not have names.  It is assumed data in both A1 and A2 are ordered the same.
## 
## Call:
## two.b.pls(A1 = tail.coords, A2 = head.coords, print.progress = FALSE) 
## 
## 
## 
## r-PLS: 0.917
## 
## Effect Size (Z): 5.4983
## 
## P-value: 0.001
## 
## Based on 1000 random permutations

picknplot.shape example is provided in lab script.

2B: Comparing the Strength of Integration

One may also evaluate whether the strength of integration similar in two or more groups. As above this is accomplished using effect sizes, this time based on \(r_{PLS}\).

tail.coords.gp <- coords.subset(tail.coords, group)
head.coords.gp <- coords.subset(head.coords, group)

# Obtain Integration for groups
integ.tests <- Map(function(x,y) integration.test(x, y, iter=499, 
                print.progress = FALSE), head.coords.gp, tail.coords.gp)
compare.pls(integ.tests)
## 
## Effect sizes
## 
##    Marsh.F    Marsh.M Sinkhole.F Sinkhole.M 
##  3.1410345  1.5130737  0.7429882  2.7156556 
## 
## Effect sizes for pairwise differences in PLS effect size
## 
##              Marsh.F   Marsh.M Sinkhole.F Sinkhole.M
## Marsh.F    0.0000000 0.6586628  1.9027247  0.7518767
## Marsh.M    0.6586628 0.0000000  0.8778808  1.1927859
## Sinkhole.F 1.9027247 0.8778808  0.0000000  2.0927712
## Sinkhole.M 0.7518767 1.1927859  2.0927712  0.0000000
## 
## P-values
## 
##               Marsh.F   Marsh.M Sinkhole.F Sinkhole.M
## Marsh.F    1.00000000 0.5101123 0.05707647 0.45212519
## Marsh.M    0.51011230 1.0000000 0.38000838 0.23295324
## Sinkhole.F 0.05707647 0.3800084 1.00000000 0.03636958
## Sinkhole.M 0.45212519 0.2329532 0.03636958 1.00000000

3: Modularity among Subsets of Traits

Tests of modularity can also be performed in geomorph:

MT <- modularity.test(pupfish$coords,land.gp,CI=FALSE,print.progress = FALSE)
summary(MT)
## 
## Call:
## modularity.test(A = pupfish$coords, partition.gp = land.gp, CI = FALSE,  
##     print.progress = FALSE) 
## 
## 
## 
## CR: 0.8771
## 
## P-value: 0.003
## 
## Effect Size: -4.5427
## 
## Based on 1000 random permutations
plot(MT)

3B: Comparing the Strength of Modularity

One may also evaluate whether the strength of modularity similar in two or more groups. As above this is accomplished using effect sizes based on \(CR\).

coords.gp <- coords.subset(pupfish$coords, group)
modul.tests <- Map(function(x) modularity.test(x, land.gp,print.progress = FALSE), coords.gp) 
compare.CR(modul.tests, CR.null = FALSE)
## 
##  NOTE: more negative effects represent stronger modular signal! 
## 
## 
## Effect sizes
## 
##    Marsh.F    Marsh.M Sinkhole.F Sinkhole.M 
##  -2.816188  -5.778451  -3.363544  -4.534488 
## 
## Effect sizes for pairwise differences in CR effect size
## 
##             Marsh.F   Marsh.M Sinkhole.F Sinkhole.M
## Marsh.F    0.000000 2.5466647  1.8303913   1.017783
## Marsh.M    2.546665 0.0000000  0.1496676   1.693290
## Sinkhole.F 1.830391 0.1496676  0.0000000   1.251317
## Sinkhole.M 1.017783 1.6932905  1.2513171   0.000000
## 
## P-values
## 
##               Marsh.F    Marsh.M Sinkhole.F Sinkhole.M
## Marsh.F    1.00000000 0.01087579 0.06719144 0.30878104
## Marsh.M    0.01087579 1.00000000 0.88102690 0.09040019
## Sinkhole.F 0.06719144 0.88102690 1.00000000 0.21081881
## Sinkhole.M 0.30878104 0.09040019 0.21081881 1.00000000

3B-2: Comparing Alternative Modular Hypotheses

One may also wish to determine whether there is support in a given dataset for a 2-module, 3-module, 4-module, or other hypothesis. This is also accomplished via effect sizes of \(CR\).

land.gps3 <- rep('a',56); land.gps3[39:48]<-'b'; land.gps3[c(6:9,28:38)] <- 'c' 
#3 module hypothesis (tail now a module)
land.gps4 <- rep('a',56); land.gps4[39:48]<-'b'; land.gps4[c(6:9,28:38)] <- 'c'; 
land.gps4[c(10,49:56)] <- 'd'  #4 module hypothesis (eye now a module)

m3.test <- modularity.test(coords.gp$Marsh.F,land.gps3, iter = 499, print.progress = FALSE)
m4.test <- modularity.test(coords.gp$Marsh.F,land.gps4, iter = 499, print.progress = FALSE)

model.Z <- compare.CR(m3.test,m4.test, CR.null = TRUE)
model.Z 
## 
##  NOTE: more negative effects represent stronger modular signal! 
## 
## 
## Effect sizes
## 
## No_Modules    m3.test    m4.test 
##   0.000000  -3.055068  -3.873177 
## 
## Effect sizes for pairwise differences in CR effect size
## 
##            No_Modules  m3.test  m4.test
## No_Modules   0.000000 3.055068 3.873177
## m3.test      3.055068 0.000000 1.030587
## m4.test      3.873177 1.030587 0.000000
## 
## P-values
## 
##              No_Modules     m3.test      m4.test
## No_Modules 1.0000000000 0.002250096 0.0001074256
## m3.test    0.0022500958 1.000000000 0.3027346572
## m4.test    0.0001074256 0.302734657 1.0000000000