From 91449c343271cb2370edb19639cf108ede9ba501 Mon Sep 17 00:00:00 2001 From: "Ruby E. Stephens" Date: Thu, 18 Apr 2024 14:20:22 +1000 Subject: [PATCH] just latitude analysis of interest --- scripts/analysis/PGLS_with_latitude.R | 180 ++++---------------------- 1 file changed, 23 insertions(+), 157 deletions(-) diff --git a/scripts/analysis/PGLS_with_latitude.R b/scripts/analysis/PGLS_with_latitude.R index 98c7892..51fb573 100644 --- a/scripts/analysis/PGLS_with_latitude.R +++ b/scripts/analysis/PGLS_with_latitude.R @@ -124,156 +124,24 @@ pgls_allotb %>% ylab("Floral longevity (log # days)") ggsave("figures/symmetry_loglongevity_abs_latitude_allotb.pdf", width = 9, height = 5) -#** ttest & PGLS ---- -# t-test of longevity by symmetry -ttests <- list() -ttests$allotblong <- t.test(pgls_allotb$spmean_long_days[pgls_allotb$sym_species == "1"], - pgls_allotb$spmean_long_days[pgls_allotb$sym_species == "0"]) -ttests$allotblong -# Welch Two Sample t-test -# -# data: pgls_allotb$spmean_long_days[pgls_allotb$sym_species == "1"] and pgls_allotb$spmean_long_days[pgls_allotb$sym_species == "0"] -# t = 4.3647, df = 710.24, p-value = 1.462e-05 -# alternative hypothesis: true difference in means is not equal to 0 -# 95 percent confidence interval: -# 0.6167554 1.6252529 -# sample estimates: -# mean of x mean of y -# 4.682592 3.561587 - -# t-test of latitude by symmetry -ttests$allotblat <- t.test(pgls_allotb$spmeanlat[pgls_allotb$sym_species == "1"], - pgls_allotb$spmeanlat[pgls_allotb$sym_species == "0"]) -ttests$allotblat -# Welch Two Sample t-test -# -# data: pgls_allotb$spmeanlat[pgls_allotb$sym_species == "1"] and pgls_allotb$spmeanlat[pgls_allotb$sym_species == "0"] -# t = 1.1786, df = 921.11, p-value = 0.2389 -# alternative hypothesis: true difference in means is not equal to 0 -# 95 percent confidence interval: -# -0.5881533 2.3567204 -# sample estimates: -# mean of x mean of y -# 26.37643 25.49215 +#** PGLS ---- -# run PGLS without latitude first, to check symmetry by itself +# list to store model output pgls_models <- list() -pgls_models$allotbsymlong <- nlme::gls(log(spmean_long_days) ~ - scale(as.numeric(sym_species)), - correlation = ape::corBrownian(phy = allotb, - form = ~spp), - data = pgls_allotb, method = "ML") -summary(pgls_models$allotbsymlong) -# Generalized least squares fit by maximum likelihood -# Model: log(spmean_long_days) ~ scale(as.numeric(sym_species)) -# Data: pgls_allotb -# AIC BIC logLik -# 4524.62 4540.402 -2259.31 -# -# Correlation Structure: corBrownian -# Formula: ~spp -# Parameter estimate(s): -# numeric(0) -# -# Coefficients: -# Value Std.Error t-value p-value -# (Intercept) 0.9563604 0.7457652 1.282388 0.1999 -# scale(as.numeric(sym_species)) 0.1426238 0.0413410 3.449937 0.0006 -# -# Correlation: -# (Intr) -# scale(as.numeric(sym_species)) 0.017 -# -# Standardized residuals: -# Min Q1 Med Q3 Max -# -1.10113651 -0.28198308 -0.01443645 0.23146808 0.81545011 -# -# Residual standard error: 3.047146 -# Degrees of freedom: 1423 total; 1421 residual - -# effect of (absolute) latitude on longevity? -pgls_models$allotblonglat <- nlme::gls(log(spmean_long_days) ~ scale(spmeanlatabs), - correlation = ape::corBrownian(phy = allotb, - form = ~spp), - data = pgls_allotb, method = "ML") -summary(pgls_models$allotblonglat) -# Generalized least squares fit by maximum likelihood -# Model: log(spmean_long_days) ~ scale(spmeanlatabs) -# Data: pgls_allotb -# AIC BIC logLik -# 4496.767 4512.548 -2245.383 -# -# Correlation Structure: corBrownian -# Formula: ~spp -# Parameter estimate(s): -# numeric(0) -# -# Coefficients: -# Value Std.Error t-value p-value -# (Intercept) 0.8665400 0.7384333 1.173484 0.2408 -# scale(spmeanlatabs) -0.2110542 0.0332768 -6.342376 0.0000 -# -# Correlation: -# (Intr) -# scale(spmeanlatabs) 0.01 -# -# Standardized residuals: -# Min Q1 Med Q3 Max -# -1.0204262 -0.3280829 0.0318071 0.3029884 0.8935845 -# -# Residual standard error: 3.017469 -# Degrees of freedom: 1423 total; 1421 residual - - -# yes looks like latitude has an effect on longevity -# can't do effect of latitude on symmetry bc binary variable (would need phylo logistic regression) -# but can dummy check the reverse?? -pgls_models$allotbsymlat <- nlme::gls(spmeanlatabs ~ sym_species, - correlation = ape::corBrownian(phy = allotb, - form = ~spp), - data = pgls_allotb, method = "ML") -summary(pgls_models$allotbsymlat) -# Generalized least squares fit by maximum likelihood -# Model: spmeanlatabs ~ sym_species -# Data: pgls_allotb -# AIC BIC logLik -# 11242 11257.78 -5617.998 -# -# Correlation Structure: corBrownian -# Formula: ~spp -# Parameter estimate(s): -# numeric(0) -# -# Coefficients: -# Value Std.Error t-value p-value -# (Intercept) 23.076481 7.901225 2.920620 0.0035 -# sym_species1 -1.478082 0.940960 -1.570824 0.1164 -# -# Correlation: -# (Intr) -# sym_species1 -0.021 -# -# Standardized residuals: -# Min Q1 Med Q3 Max -# -0.7060763 -0.2233425 0.1016666 0.4622944 1.4485092 -# -# Residual standard error: 32.28142 -# Degrees of freedom: 1423 total; 1421 residual -# latitude does not vary by symmetry, which agrees with the scatter plot - -# how do symmetry and latitude interact as fixed effects on longevity +# what is the effect of symmetry and latitude on longevity ? pgls_models$allotbsymlatlong <- nlme::gls(log(spmean_long_days) ~ - scale(as.numeric(sym_species))*scale(spmeanlatabs), + scale(as.numeric(sym_species)) + + scale(spmeanlatabs), correlation = ape::corBrownian(phy = allotb, form = ~spp), data = pgls_allotb, method = "ML") summary(pgls_models$allotbsymlatlong) # Generalized least squares fit by maximum likelihood -# Model: log(spmean_long_days) ~ scale(as.numeric(sym_species)) * scale(spmeanlatabs) +# Model: log(spmean_long_days) ~ scale(as.numeric(sym_species)) + scale(spmeanlatabs) # Data: pgls_allotb # AIC BIC logLik -# 4458.418 4484.721 -2224.209 +# 4491.906 4512.948 -2241.953 # # Correlation Structure: corBrownian # Formula: ~spp @@ -282,29 +150,27 @@ summary(pgls_models$allotbsymlatlong) # # Coefficients: # Value Std.Error t-value p-value -# (Intercept) 0.9536605 0.7281932 1.309626 0.1905 -# scale(as.numeric(sym_species)) 0.1691634 0.0409237 4.133630 0.0000 -# scale(spmeanlatabs) -0.1239595 0.0359223 -3.450767 0.0006 -# scale(as.numeric(sym_species)):scale(spmeanlatabs) -0.1628957 0.0287153 -5.672776 0.0000 +# (Intercept) 0.9090730 0.7370253 1.233435 0.2176 +# scale(as.numeric(sym_species)) 0.1367322 0.0408635 3.346076 0.0008 +# scale(spmeanlatabs) -0.2087726 0.0332832 -6.272614 0.0000 # # Correlation: -# (Intr) sc(.(_)) scl(s) -# scale(as.numeric(sym_species)) 0.019 -# scale(spmeanlatabs) 0.014 0.102 -# scale(as.numeric(sym_species)):scale(spmeanlatabs) -0.011 -0.160 -0.405 -# -# Standardized residuals: -# Min Q1 Med Q3 Max -# -1.115757251 -0.285272896 0.003389668 0.264009912 0.857956007 -# -# Residual standard error: 2.972902 -# Degrees of freedom: 1423 total; 1419 residual +# (Intr) s(.(_) +# scale(as.numeric(sym_species)) 0.017 +# scale(spmeanlatabs) 0.011 0.044 +# +# Standardized residuals: +# Min Q1 Med Q3 Max +# -1.007123509 -0.319472645 0.009280592 0.279001215 0.838333798 +# +# Residual standard error: 3.010204 +# Degrees of freedom: 1423 total; 1420 residual # with all predictors scaled (and therefore numeric) all come out -# as significant, comparing effect sizes symmetry has strongest (+ve) effect, -# closely followed by interaction (-ve) and then latitude (-ve) +# as significant, comparing effect sizes latitude has strongest (-ve) effect, +# followed by symmetry (+ve) -rm(spp, pgls_allotb, allotb, sym_long_lat, ttests, pgls_models) +rm(spp, pgls_allotb, allotb, sym_long_lat, pgls_models) # tried Bayesian models as per Song et al. (2022) paper, using brms package # and running full interaction model symmetry still has significant +ve effect,