Background

Response Variable = dewlappermin
Predictor Variables = dewsize and visibility
Covariate = svl
Other variables = pushbobpermin and movepermin

Question Do dewlap display behaviors of male anoles align with principles of signal efficacy?

Data Analysis

Load Data
Also create separate dataframes for each species to make analysis later a little easier

load(here::here('data/jm_anoles.Rda'))

g <- jmanoles[which(jmanoles$species == "grahami"),]
l <- jmanoles[which(jmanoles$species == "lineatopus"),]
v <- jmanoles[which(jmanoles$species == "valencienni"),]


knitr::kable(head(jmanoles))
species ID dewsize1 dewsize2 dewsize3 dewsizemax sex svl height diameter inclination visibilitygroup visibility minobserved dewlappermin pushbobpermin movepermin X95.MCP
grahami G123 0.092 0.098 NA 0.098 m 43 1.010000 1.500000 0.00000 bottom 0.6122927 0 NA NA NA NA
grahami G113 0.081 0.102 NA 0.102 m 45 0.175000 6.750000 59.33333 bottom 0.6200000 40 0.075 0.000 0.450000 0.9668063
grahami G337 0.093 0.097 NA 0.097 m 46 1.046667 8.750000 25.00000 bottom 0.6697648 12 0.000 0.000 1.166667 NA
grahami G300 0.204 0.259 NA 0.259 m 56 1.201429 2.500000 71.25000 bottom 0.6970845 52 0.000 0.000 1.346154 13.7153600
grahami G140 0.331 0.337 NA 0.337 m 62 2.207400 2.666667 45.00000 bottom 0.7289377 40 1.250 0.825 0.850000 4.7038580
grahami G128 0.213 0.193 NA 0.213 m 56 0.620000 1.000000 20.00000 bottom 0.7297416 0 NA NA NA NA

Correlation Matrix Creating a correlation matrix that compares the correlation among svl, dewsizemax, visibility, dewlappermin, pushbobpermin, and movepermin for all the lizard species together and each species separately.

Nubmers in the matrix represent the correlation coefficient. The color and size of the circle represnet the magnitude and direciton of the correlation. And an x over the correlation coefficient signifies that it is not a statistically signficant correlation.

# correlation matrix function 

crmtrx <- function(df, speciesname){
  
  # subsetting the continuous variables for the matrix
  cordata <- df[,c('svl', 'dewsizemax', 'visibility', 'dewlappermin', 'pushbobpermin', 'movepermin')]
  
  M <- cor(cordata, method = "pearson", use = "complete.obs")

  testRes <- corrplot::cor.mtest(cordata, conf.level = 0.95)

  corrplot::corrplot(M, p.mat = testRes$p, method = 'circle', type = 'lower',
         addCoef.col ='black', number.cex = 0.8, order = 'AOE', diag=FALSE, title = speciesname)
}


crmtrx(jmanoles, 'All')

crmtrx(g, 'grahami')

crmtrx(l, 'lineatopus')

crmtrx(v, 'valencienni')

Regression Models

Before continuing, for valencienni, we have both males and females. All the regression analyses hae males and females clumped together, but we probably can’t do that because when running t tests there are differences in dewlap behavior, dewlap size, and svl between the sexes (however, there is not a difference in the visibility of the habitats they were found displaying at).

t.test(data = v, dewlappermin ~ sex)
## 
##  Welch Two Sample t-test
## 
## data:  dewlappermin by sex
## t = -3.0385, df = 6.0217, p-value = 0.02274
## alternative hypothesis: true difference in means between group f and group m is not equal to 0
## 95 percent confidence interval:
##  -1.1191889 -0.1211846
## sample estimates:
## mean in group f mean in group m 
##      0.02519452      0.64538127
t.test(data = v, dewsizemax ~ sex)
## 
##  Welch Two Sample t-test
## 
## data:  dewsizemax by sex
## t = -5.0268, df = 15.942, p-value = 0.0001253
## alternative hypothesis: true difference in means between group f and group m is not equal to 0
## 95 percent confidence interval:
##  -0.3621594 -0.1472628
## sample estimates:
## mean in group f mean in group m 
##       0.2878889       0.5426000
t.test(data = v, svl ~ sex)
## 
##  Welch Two Sample t-test
## 
## data:  svl by sex
## t = -2.5302, df = 18.385, p-value = 0.02073
## alternative hypothesis: true difference in means between group f and group m is not equal to 0
## 95 percent confidence interval:
##  -8.6171199 -0.8051024
## sample estimates:
## mean in group f mean in group m 
##        58.55556        63.26667
t.test(data = v, visibility ~ sex)
## 
##  Welch Two Sample t-test
## 
## data:  visibility by sex
## t = -1.6153, df = 19.723, p-value = 0.1221
## alternative hypothesis: true difference in means between group f and group m is not equal to 0
## 95 percent confidence interval:
##  -0.1366439  0.0174375
## sample estimates:
## mean in group f mean in group m 
##       0.8264410       0.8860442

Function to run regression models

library(tidyverse)
## -- Attaching packages --------------------------------------- tidyverse 1.3.1 --
## v ggplot2 3.3.5     v purrr   0.3.4
## v tibble  3.1.3     v dplyr   1.0.7
## v tidyr   1.1.3     v stringr 1.4.0
## v readr   2.0.1     v forcats 0.5.1
## -- Conflicts ------------------------------------------ tidyverse_conflicts() --
## x dplyr::filter() masks stats::filter()
## x dplyr::lag()    masks stats::lag()
library(knitr)
library(broom)

# creating a funciton to run the regression analyses 
lmfxn <- function(response, predictors, speciesname){
  
   fixedformula <- paste(response, "~", predictors)
   
   lm(as.formula(fixedformula), data = subset(jmanoles, species == speciesname)) %>%
     tidy() %>%
     kable()

}

Full Model (absolute dewlap size)
dewsizemax + visibility

lmfxn('dewlappermin', 'dewsizemax + visibility', 'grahami')
term estimate std.error statistic p.value
(Intercept) -3.278308 1.206925 -2.716249 0.0103052
dewsizemax 1.367933 1.006390 1.359248 0.1830160
visibility 4.258064 1.448834 2.938960 0.0058808 
lmfxn('dewlappermin', 'dewsizemax + visibility', 'lineatopus')
term estimate std.error statistic p.value
(Intercept) -0.2791601 0.3780039 -0.7385110 0.4645153
dewsizemax 0.3508507 0.1881182 1.8650547 0.0695216
visibility 0.3886007 0.4657376 0.8343769 0.4090241 
lmfxn('dewlappermin', 'dewsizemax + visibility', 'valencienni')
term estimate std.error statistic p.value
(Intercept) 0.1521760 1.0744628 0.1416299 0.8893902
dewsizemax 1.6633053 0.5050969 3.2930421 0.0053356
visibility -0.6502239 1.3826710 -0.4702665 0.6454067

Full Model with Intereaction

lmfxn('dewlappermin', 'dewsizemax * visibility', 'grahami')
term estimate std.error statistic p.value
(Intercept) -2.8712681 2.560654 -1.1213027 0.2702526
dewsizemax -0.4573428 10.136302 -0.0451193 0.9642842
visibility 3.7758139 3.043014 1.2408140 0.2234235
dewsizemax:visibility 2.1431313 11.840918 0.1809937 0.8574796 
lmfxn('dewlappermin', 'dewsizemax * visibility', 'lineatopus')
term estimate std.error statistic p.value
(Intercept) -0.3945753 0.9115879 -0.4328440 0.6675125
dewsizemax 0.6645241 2.2566199 0.2944776 0.7699547
visibility 0.5222358 1.0677339 0.4891067 0.6275057
dewsizemax:visibility -0.3578142 2.5649908 -0.1394992 0.8897736

Full Model (relative dewlap size)
svl + dewsizemax + visibility

lmfxn('dewlappermin', 'svl + dewsizemax + visibility', 'grahami')
term estimate std.error statistic p.value
(Intercept) -3.2834290 1.7674411 -1.8577303 0.0721493
svl 0.0001527 0.0379887 0.0040197 0.9968169
dewsizemax 1.3593061 2.3768273 0.5718994 0.5712656
visibility 4.2564468 1.5246810 2.7916967 0.0086523 
lmfxn('dewlappermin', 'svl + dewsizemax + visibility', 'lineatopus')
term estimate std.error statistic p.value
(Intercept) -0.2777665 0.5502190 -0.5048290 0.6165198
svl -0.0000390 0.0110683 -0.0035262 0.9972045
dewsizemax 0.3516676 0.2999381 1.1724673 0.2481217
visibility 0.3889894 0.4843790 0.8030682 0.4268026 
lmfxn('dewlappermin', 'svl + dewsizemax + visibility', 'valencienni')
term estimate std.error statistic p.value
(Intercept) 0.3208546 2.0319314 0.1579062 0.8769571
svl -0.0034788 0.0350391 -0.0992843 0.9224272
dewsizemax 1.7330502 0.8763628 1.9775488 0.0695777
visibility -0.6323120 1.4456234 -0.4373975 0.6690016

Signal Size Model
Do males with larger dewlaps (absolute size) display their dewlaps more than lizards with smaller dewlaps?
dewsizemax

Do males with relatively larger dewlaps display their dewlaps more than lizards with smaller dewlaps? dewsizemax + svl

lmfxn('dewlappermin', 'dewsizemax', 'grahami')
term estimate std.error statistic p.value
(Intercept) 0.1606349 0.3264514 0.4920637 0.6257466
dewsizemax 1.9761692 1.0870390 1.8179376 0.0776416 
lmfxn('dewlappermin', 'dewsizemax + svl', 'grahami')
term estimate std.error statistic p.value
(Intercept) -1.0239071 1.7210609 -0.5949279 0.5558307
dewsizemax 0.3437192 2.5728032 0.1335972 0.8945088
svl 0.0281425 0.0401359 0.7011813 0.4879615 
lmfxn('dewlappermin', 'svl', 'grahami')
term estimate std.error statistic p.value
(Intercept) -1.2109157 0.9871632 -1.226662 0.2281379
svl 0.0329946 0.0168410 1.959187 0.0580947 
lmfxn('dewlappermin', 'dewsizemax', 'lineatopus')
term estimate std.error statistic p.value
(Intercept) 0.0275355 0.0878549 0.3134199 0.7555507
dewsizemax 0.4144739 0.1713327 2.4191166 0.0200810 
lmfxn('dewlappermin', 'dewsizemax + svl', 'lineatopus')
term estimate std.error statistic p.value
(Intercept) -0.0588854 0.4758426 -0.1237499 0.9021331
dewsizemax 0.3697209 0.2977640 1.2416577 0.2215946
svl 0.0019837 0.0107299 0.1848723 0.8542634 
lmfxn('dewlappermin', 'svl', 'lineatopus')
term estimate std.error statistic p.value
(Intercept) -0.4709745 0.3432405 -1.372141 0.1774806
svl 0.0128149 0.0062892 2.037615 0.0480766 
lmfxn('dewlappermin', 'dewsizemax', 'valencienni')
term estimate std.error statistic p.value
(Intercept) -0.345359 0.1825188 -1.892183 0.0779307
dewsizemax 1.523585 0.3977257 3.830743 0.0016371 
lmfxn('dewlappermin', 'dewsizemax + svl', 'valencienni')
term estimate std.error statistic p.value
(Intercept) -0.0627007 1.7792454 -0.0352400 0.9723859
dewsizemax 1.6376401 0.8239031 1.9876612 0.0667717
svl -0.0053915 0.0337462 -0.1597657 0.8753478 
lmfxn('dewlappermin', 'svl', 'valencienni')
term estimate std.error statistic p.value
(Intercept) -2.9420546 1.1301226 -2.603306 0.0199672
svl 0.0527281 0.0184295 2.861066 0.0118971 
# library(ppcor)
# 
# cor <- data.frame(g$dewlappermin, g$dewsizemax, g$svl)
# spcor.test(g$dewlappermin, g$dewsizemax, g$svl, method = c('pearson'))
library(ggplot2)

ggplt <- ggplot(jmanoles,aes(x = svl, y = dewsizemax, shape = sex, color = species))+
         geom_point()+
         theme_classic() 

ggplt+geom_smooth(method=lm,se=FALSE,fullrange=TRUE,
                  aes(color=species))
## `geom_smooth()` using formula 'y ~ x'

library(ggplot2)

ggplt <- ggplot(jmanoles,aes(x = svl, y = dewlappermin, shape = species, color = species))+
         geom_point()+
         theme_classic() + ylim(0,3.2)

ggplt+geom_smooth(method=lm,se=FALSE,fullrange=TRUE,
                  aes(color=species))
## `geom_smooth()` using formula 'y ~ x'
## Warning: Removed 44 rows containing non-finite values (stat_smooth).
## Warning: Removed 44 rows containing missing values (geom_point).
## Warning: Removed 38 rows containing missing values (geom_smooth).

library(ggplot2)

ggplt <- ggplot(jmanoles,aes(x = dewsizemax, y = dewlappermin, shape = species, color = species))+
         geom_point()+
         theme_classic() + ylim(0,3.2)

ggplt+geom_smooth(method=lm,se=FALSE,fullrange=TRUE,
                  aes(color=species))
## `geom_smooth()` using formula 'y ~ x'
## Warning: Removed 44 rows containing non-finite values (stat_smooth).
## Warning: Removed 44 rows containing missing values (geom_point).
## Warning: Removed 16 rows containing missing values (geom_smooth).

library(ggplot2)

ggplt <- ggplot(jmanoles,aes(x = dewsizemax/svl, y = dewlappermin, shape = species, color = species))+
         geom_point()+
         theme_classic() + ylim(0,3.2)

ggplt+geom_smooth(method=lm,se=FALSE,fullrange=TRUE,
                  aes(color=species))
## `geom_smooth()` using formula 'y ~ x'
## Warning: Removed 44 rows containing non-finite values (stat_smooth).
## Warning: Removed 44 rows containing missing values (geom_point).
## Warning: Removed 16 rows containing missing values (geom_smooth).

Visibility Model
Do males with in more visible habitats display their dewlap more than lizards in more cluttered habitats?
visibility

lmfxn('dewlappermin', 'visibility', 'grahami')
term estimate std.error statistic p.value
(Intercept) -3.243017 1.221168 -2.655668 0.0118357
visibility 4.663041 1.434933 3.249658 0.0025546 
lmfxn('dewlappermin', 'visibility', 'lineatopus')
term estimate std.error statistic p.value
(Intercept) -0.4172093 0.3817256 -1.092956 0.2807935
visibility 0.7406915 0.4384407 1.689377 0.0987391 
lmfxn('dewlappermin', 'visibility', 'valencienni')
term estimate std.error statistic p.value
(Intercept) -1.450292 1.232847 -1.176376 0.2577767
visibility 2.028065 1.439038 1.409320 0.1791349 
library(ggplot2)

#ggplot(data = jmanoles, aes(x = visibility, y = dewlappermin, color = species)) + geom_point()

ggplt <- ggplot(jmanoles,aes(x = visibility, y = dewlappermin, shape = species, color = species))+
         geom_point()+
         theme_classic() + ylim(0,3.2) + xlim(0.59,1)

ggplt+geom_smooth(method=lm,se=FALSE,fullrange=TRUE,
                  aes(color=species))
## `geom_smooth()` using formula 'y ~ x'
## Warning: Removed 44 rows containing non-finite values (stat_smooth).
## Warning: Removed 44 rows containing missing values (geom_point).
## Warning: Removed 46 rows containing missing values (geom_smooth).

Note: this relationship, especially for grahami seems more curvilinear than linear