Winter 2020 Flight Trials: PCA

Flight Trials Winter 2020 Dataset was conducted from 2/17/2020 - 3/10/2020. Soapberry bugs were flight tested twice for multiple hours in the flight mill and observed from 8 AM to (5-8 PM) each day.

Size-adjusted CPCA (notes)

“Morphological relationships change with overall body size and body size often varies among populations. Therefore, quantitative analyses of indi- vidual traits from organisms in different populations or environments (e.g., in studies of phenotypic plasticity) often adjust for differences in body size to isolate changes in allometry.” (McCoy et al, 2006)

“For most studies of morphological plasticity and size correction, it is only the first CPC (CPC1) that is of interest because it describes the general scaling of traits with body size.” (McCoy et al, 2006)

“If the first principal components are not shared (Fig. 1c), then the patterns of morphological variation are so fundamentally different that there can be no globally applied size correction because ‘size’ does not have a common meaning across groups. On the other hand, if all principal components are common to all groups (Fig. 1a), then the groups show identical patterns of within-group covariation, indicating the same all- ometries of traits with size in all groups.” (McCoy et al, 2006)

Looking at the figures….I would say differences in size but not shape.

“As long as the groups to be compared share CPC1 (Fig. 1a, b) and all the traits have strong loadings on CPC1, then it can be interpreted as a common body size dimension:i.e., a ‘size-axis’ that can be used to evaluate between-group differences in morphology that go beyond simple (isometric) changes in size” (McCoy et al, 2006)

“…where the main goal is to correct for the effects of a growth axis that explains a large fraction of the total variation (e.g., >90%)” (McCoy et al, 2006)

No Mass

source("src/clean_flight_data.R") # Script that loads and cleans up the data
source("src/center_flight_data.R")

data <- read_flight_data("data/all_flight_data-Winter2020.csv")
data_all <- data[[1]]
data_tested <- data[[2]]
#Standardize variables
data_tested$thorax_s <- (data_tested$thorax-mean(data_tested$thorax, na.rm=TRUE)) / sd(data_tested$thorax, na.rm=TRUE)
data_tested$body_s <- (data_tested$body-mean(data_tested$body, na.rm=TRUE)) / sd(data_tested$body, na.rm=TRUE)
data_tested$wing_s <- (data_tested$wing-mean(data_tested$wing, na.rm=TRUE)) / sd(data_tested$wing, na.rm=TRUE)
data_tested$beak_s <- (data_tested$beak-mean(data_tested$beak, na.rm=TRUE)) / sd(data_tested$beak, na.rm=TRUE)
data_tested$mass_s <- (data_tested$mass-mean(data_tested$mass, na.rm=TRUE)) / sd(data_tested$mass, na.rm=TRUE)

PCA applies to data tables where rows are considered as individuals and columns as quantitative variables. Let’s filter out the dupliates:

data_tested <- data_tested[data_tested$trial_type == "T1",]
#The output is a PCA of soapberry bug morphology

d <- data_tested %>%
  select(thorax_s, body_s, wing_s, beak_s) %>%
  filter(!is.na(body_s))
colnames(d) <- c("thorax", "body", "wing", "beak")
MorphPCA <- PCA(d, scale.unit=TRUE, graph=TRUE)

#Investiating variables
eig.val <- get_eigenvalue(MorphPCA)
eig.val
##       eigenvalue variance.percent cumulative.variance.percent
## Dim.1 3.50143036        87.535759                    87.53576
## Dim.2 0.37328713         9.332178                    96.86794
## Dim.3 0.10977100         2.744275                    99.61221
## Dim.4 0.01551152         0.387788                   100.00000
fviz_eig(MorphPCA, addlabels = TRUE, ylim = c(0, 50))

#An eigenvalue > 1 indicates that PCs account for more variance than accounted by one of the original variables in standardized data. This is commonly used as a cutoff point for which PCs are retained. This holds true only when the data are standardized, which we have done.
#You can also limit the number of component to that number that accounts for a certain fraction of the total variance. (No formal cut-off)

variables <- get_pca_var(MorphPCA)

corrplot(variables$cos2, is.corr=FALSE)

corrplot(variables$contrib, is.corr=FALSE)

GFdim <- dimdesc(MorphPCA, axes = c(1,2), proba = 0.05)
GFdim$Dim.1
## $quanti
##        correlation       p.value
## body     0.9826965 9.758854e-244
## thorax   0.9593090 2.657695e-183
## wing     0.9451532 2.013968e-162
## beak     0.8497939  8.736348e-94
## 
## attr(,"class")
## [1] "condes" "list "
GFdim$Dim.2
## $quanti
##      correlation      p.value
## beak   0.5212173 1.608037e-24
## body  -0.1428160 9.166353e-03
## wing  -0.2825846 1.635350e-07
## 
## attr(,"class")
## [1] "condes" "list "
vars <- MorphPCA$ind

par(mfrow=c(2,3))
#Correlations

m <- lm(d$thorax ~ vars$coord[,1])
#pval <- paste("p = ", summary(m)$coefficients[8])
vars <- MorphPCA$ind
plot(vars$coord[,1], d$thorax, xlab="PC1", ylab="thorax (mm)")
cor_val <- paste("R = ", round(cor(vars$coord[,1], d$thorax),3))
text(-2,2, cor_val)
#text(-2,1.5, pval)
abline(m, lty=2)

m <- lm(d$body ~ vars$coord[,1])
plot(vars$coord[,1], d$body, xlab="PC1", ylab="body (mm)")
cor_val <- paste("R = ", round(cor(vars$coord[,1], d$body),3))
text(-2,2, cor_val)
abline(m, lty=2)

m <- lm(d$wing ~ vars$coord[,1])
plot(vars$coord[,1], d$wing, xlab="PC1", ylab="wing (mm)")
cor_val <- paste("R = ", round(cor(vars$coord[,1], d$wing),3))
text(-2,2, cor_val)
abline(m, lty=2)

m <- lm(d$beak ~ vars$coord[,1])
plot(vars$coord[,1], d$beak, xlab="PC1", ylab="beak (mm)")
cor_val <- paste("R = ", round(cor(vars$coord[,1], d$beak),3))
text(-2,2, cor_val)
abline(m, lty=2)

m <- lm(d$beak ~ vars$coord[,2])
plot(vars$coord[,2], d$beak, xlab="PC2", ylab="beak (mm)")
cor_val <- paste("R = ", round(cor(vars$coord[,2], d$beak), 3))
text(-0.5,3, cor_val)
abline(m, lty=2)

#Plot by host
host <- data_tested %>%
  select(host_plant)

d$host_plant <- host$host_plant
fviz_pca_ind(MorphPCA,
             geom.ind = "point",
             col.ind = d$host_plant,
             legend.title = "Host Plant",
             addEllipses=TRUE)

# Plot by sex
sex <- data_tested %>%
  select(sex) 

d$sex <- sex$sex
fviz_pca_ind(MorphPCA,
             geom.ind = "point",
             col.ind = d$sex,
             legend.title = "Sex",
             addEllipses=TRUE)

# Plot by distance from the sympatric zone
sym_dist <- data_tested %>%
  select(sym_dist)

d$sym_dist <- sym_dist$sym_dist
fviz_pca_ind(MorphPCA,
             geom.ind = "point",
             col.ind = d$sym_dist,
             legend.title = "Sympatric Zone Distance")

# Plot by population
pop <- data_tested %>%
  select(population)

d$population <- pop$pop
fviz_pca_ind(MorphPCA,
             geom.ind = "point",
             col.ind = d$population,
             legend.title = "Population",
             addEllipses=TRUE)

Mass

d <- data_tested %>%
  select(thorax_s, body_s, wing_s, beak_s, mass_s) %>%
  filter(!is.na(mass_s))
colnames(d) <- c("thorax", "body", "wing", "beak", "mass")
MorphMPCA <- PCA(d, scale.unit=TRUE, graph=TRUE)

Wing2Body

# No mass
d <- data_tested %>%
  select(thorax_s, wing2body_s, beak_s) %>%
  filter(!is.na(beak_s))
colnames(d) <- c("thorax", "wing2body", "beak")
MorphMPCA <- PCA(d, scale.unit=TRUE, graph=TRUE)

# Mass
d <- data_tested %>%
  select(thorax_s, wing2body_s, beak_s, mass_s) %>%
  filter(!is.na(mass_s))
colnames(d) <- c("thorax", "wing2body", "beak", "mass")
MorphMPCA <- PCA(d, scale.unit=TRUE, graph=TRUE)